Quantifier.hpp 15 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383
  1. /**
  2. * @file kernel/pdevs/qss/Quantifier.hpp
  3. * @author The ARTIS Development Team
  4. * See the AUTHORS or Authors.txt file
  5. */
  6. /*
  7. * ARTIS - the multimodeling and simulation environment
  8. * This file is a part of the ARTIS environment
  9. *
  10. * Copyright (C) 2013-2019 ULCO http://www.univ-littoral.fr
  11. *
  12. * This program is free software: you can redistribute it and/or modify
  13. * it under the terms of the GNU General Public License as published by
  14. * the Free Software Foundation, either version 3 of the License, or
  15. * (at your option) any later version.
  16. *
  17. * This program is distributed in the hope that it will be useful,
  18. * but WITHOUT ANY WARRANTY; without even the implied warranty of
  19. * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  20. * GNU General Public License for more details.
  21. *
  22. * You should have received a copy of the GNU General Public License
  23. * along with this program. If not, see <http://www.gnu.org/licenses/>.
  24. */
  25. #ifndef QSS_QUANTIFIER
  26. #define QSS_QUANTIFIER
  27. #include <artis-star/kernel/pdevs/Dynamics.hpp>
  28. #include <artis-star/kernel/pdevs/qss/Data.hpp>
  29. #include <cmath>
  30. namespace artis {
  31. namespace pdevs {
  32. namespace qss {
  33. struct QuantifierParameters {
  34. bool allow_offsets;
  35. bool zero_init_offset;
  36. double quantum;
  37. unsigned int archive_length;
  38. };
  39. template<class Time>
  40. class Quantifier
  41. : public artis::pdevs::Dynamics<Time, Quantifier<Time>, QuantifierParameters> {
  42. public:
  43. enum inputs {
  44. IN, RESET
  45. };
  46. enum outputs {
  47. OUT
  48. };
  49. enum states {
  50. STATE, ADAPTIVE_STATE, STEP_NUMBER, OFFSET, UP_THRESHOLD, DOWN_THRESHOLD
  51. };
  52. Quantifier(const std::string& name,
  53. const Context<Time, Quantifier<Time>, QuantifierParameters>& context)
  54. :
  55. artis::pdevs::Dynamics<Time, Quantifier<Time>, QuantifierParameters>(name,
  56. context)
  57. {
  58. DECLARE_STATES(int,
  59. ((STATE, &Quantifier<Time>::_state), (ADAPTIVE_STATE, &Quantifier<Time>::_adaptive_state)));
  60. DECLARE_STATES(unsigned int,
  61. ((STEP_NUMBER, &Quantifier<Time>::_step_number)));
  62. DECLARE_STATES(double,
  63. ((OFFSET, &Quantifier<Time>::_offset), (UP_THRESHOLD, &Quantifier<Time>::_up_threshold), (DOWN_THRESHOLD, &Quantifier<Time>::_down_threshold)));
  64. this->input_ports({{IN, "in"},
  65. {RESET, "reset"}});
  66. this->output_port({OUT, "out"});
  67. this->observables({{UP, "up"},
  68. {DOWN, "down"},
  69. {VALUE, "value"}});
  70. _adaptive = context.parameters().allow_offsets;
  71. _adaptive_state = _adaptive ? POSSIBLE : IMPOSSIBLE;
  72. _zero_init_offset = context.parameters().zero_init_offset;
  73. _step_size = context.parameters().quantum;
  74. assert(_step_size > 0);
  75. _past_length = context.parameters().archive_length;
  76. assert(_past_length > 2);
  77. }
  78. virtual ~Quantifier() { }
  79. virtual void dconf(typename Time::type t, typename Time::type e,
  80. const common::Bag<Time>& bag)
  81. {
  82. dint(t);
  83. dext(t, e, bag);
  84. }
  85. virtual void dint(typename Time::type /* t */)
  86. {
  87. switch (_state) {
  88. case INIT:
  89. break;
  90. case IDLE:
  91. break;
  92. case RESPONSE:
  93. _state = IDLE;
  94. break;
  95. }
  96. }
  97. virtual void dext(typename Time::type t, typename Time::type e,
  98. const common::Bag<Time>& bag)
  99. {
  100. std::for_each(bag.begin(), bag.end(),
  101. [this, t, e](const common::ExternalEvent<Time>& event) {
  102. if (event.on_port(IN)) {
  103. IntegratorData data;
  104. double shifting_factor;
  105. double value;
  106. int cnt;
  107. event.data()(data);
  108. value = data.value;
  109. if (_state == INIT) {
  110. init_step_number_and_offset(value);
  111. update_thresholds();
  112. _state = RESPONSE;
  113. } else {
  114. cnt = 0;
  115. while (value >= _up_threshold or value <= _down_threshold) {
  116. cnt++;
  117. if (value >= _up_threshold) {
  118. _step_number++;
  119. } else {
  120. _step_number--;
  121. }
  122. switch (_adaptive_state) {
  123. case IMPOSSIBLE:
  124. update_thresholds();
  125. break;
  126. case POSSIBLE:
  127. if (value >= _up_threshold) {
  128. store_change(_step_size, t);
  129. } else {
  130. store_change(-_step_size, t);
  131. }
  132. shifting_factor = shift_quanta();
  133. assert(shifting_factor >= 0
  134. and shifting_factor <= 1);
  135. if (shifting_factor != 0 and shifting_factor != 1) {
  136. if (value >= _up_threshold) {
  137. update_thresholds(shifting_factor,
  138. DIRECTION_DOWN);
  139. } else {
  140. update_thresholds(shifting_factor,
  141. DIRECTION_UP);
  142. }
  143. _adaptive_state = DONE;
  144. } else {
  145. update_thresholds();
  146. }
  147. break;
  148. case DONE:
  149. init_step_number_and_offset(value);
  150. _adaptive_state = POSSIBLE;
  151. update_thresholds();
  152. break;
  153. }
  154. }
  155. }
  156. } else if (event.on_port(RESET)) {
  157. _offset = 0;
  158. _state = INIT;
  159. }
  160. });
  161. _state = RESPONSE;
  162. }
  163. virtual void start(typename Time::type /* time */)
  164. {
  165. _offset = 0;
  166. _state = INIT;
  167. }
  168. virtual typename Time::type ta(typename Time::type /* time */)
  169. {
  170. switch (_state) {
  171. case INIT:
  172. case IDLE:
  173. return Time::infinity;
  174. case RESPONSE:
  175. return 0.0;
  176. }
  177. return Time::infinity;
  178. }
  179. virtual common::Bag<Time> lambda(typename Time::type /* time */) const
  180. {
  181. common::Bag<Time> msgs;
  182. const QuantifierData data = {_up_threshold, _down_threshold};
  183. msgs.push_back(common::ExternalEvent<Time>(OUT, data));
  184. return msgs;
  185. }
  186. virtual common::Value observe(const typename Time::type& /* t */,
  187. unsigned int index) const
  188. {
  189. switch (index) {
  190. case UP:
  191. return (double) _up_threshold;
  192. case DOWN:
  193. return (double) _down_threshold;
  194. case VALUE:
  195. return (double) (_up_threshold - _down_threshold);
  196. default:
  197. return common::Value();
  198. }
  199. }
  200. private:
  201. typedef enum {
  202. DIRECTION_UP, DIRECTION_DOWN
  203. } Direction;
  204. void init_step_number_and_offset(double value)
  205. {
  206. _step_number = static_cast<long int>(std::floor(value / _step_size));
  207. if (_zero_init_offset) {
  208. _offset = 0;
  209. } else {
  210. _offset = value - static_cast<double>(_step_number) * _step_size;
  211. }
  212. }
  213. bool monotonous(unsigned int range)
  214. {
  215. if ((range + 1) > _archive.size()) {
  216. return false;
  217. }
  218. for (size_t i = 0; i < range; i++) {
  219. if (_archive[i].value * _archive[i + 1].value < 0) {
  220. return false;
  221. }
  222. }
  223. return true;
  224. }
  225. bool oscillating(unsigned int range)
  226. {
  227. if ((range + 1) > _archive.size()) {
  228. return false;
  229. }
  230. for (size_t i = _archive.size() - range; i < _archive.size() - 1; i++) {
  231. if (_archive[i].value * _archive[i + 1].value > 0) {
  232. return false;
  233. }
  234. }
  235. return true;
  236. }
  237. double shift_quanta()
  238. {
  239. double factor = 0;
  240. if (oscillating(_past_length - 1) and
  241. _archive.back().date - _archive.front().date != 0) {
  242. double acc;
  243. double local_estim;
  244. int cnt;
  245. acc = 0;
  246. cnt = 0;
  247. for (size_t i = 0; i < _archive.size() - 2; ++i) {
  248. if (0 != (_archive[i + 2].date - _archive[i].date)) {
  249. if ((_archive.back().value * _archive[i + 1].value) > 0) {
  250. local_estim =
  251. 1 - (_archive[i + 1].date - _archive[i].date) /
  252. (_archive[i + 2].date - _archive[i].date);
  253. } else {
  254. local_estim = (_archive[i + 1].date - _archive[i].date) /
  255. (_archive[i + 2].date - _archive[i].date);
  256. }
  257. acc += local_estim;
  258. cnt++;
  259. }
  260. }
  261. acc = acc / cnt;
  262. factor = acc;
  263. _archive.resize(0);
  264. }
  265. return factor;
  266. }
  267. void store_change(double val, const typename Time::type& time)
  268. {
  269. record_t record;
  270. record.date = time;
  271. record.value = val;
  272. _archive.push_back(record);
  273. while (_archive.size() > _past_length) {
  274. _archive.pop_front();
  275. }
  276. }
  277. void update_thresholds()
  278. {
  279. auto step_number = static_cast<double>(_step_number);
  280. _up_threshold = _offset + _step_size * (step_number + 1);
  281. _down_threshold = _offset + _step_size * (step_number - 1);
  282. }
  283. void update_thresholds(double factor)
  284. {
  285. auto step_number = static_cast<double>(_step_number);
  286. _up_threshold = _offset + _step_size * (step_number + (1 - factor));
  287. _down_threshold = _offset + _step_size * (step_number - (1 - factor));
  288. }
  289. void update_thresholds(double factor, Direction d)
  290. {
  291. auto step_number = static_cast<double>(_step_number);
  292. if (d == DIRECTION_UP) {
  293. _up_threshold = _offset + _step_size * (step_number + (1 - factor));
  294. _down_threshold = _offset + _step_size * (step_number - 1);
  295. } else {
  296. _up_threshold = _offset + _step_size * (step_number + 1);
  297. _down_threshold = _offset + _step_size * (step_number - (1 - factor));
  298. }
  299. }
  300. typedef enum vars {
  301. UP, DOWN, VALUE
  302. } Observable;
  303. typedef enum {
  304. INIT, IDLE, RESPONSE
  305. } State;
  306. typedef enum {
  307. IMPOSSIBLE, POSSIBLE, DONE
  308. } AdaptiveState;
  309. struct record_t {
  310. double value;
  311. typename Time::type date;
  312. };
  313. // parameters
  314. bool _adaptive;
  315. bool _zero_init_offset;
  316. unsigned int _past_length;
  317. double _step_size;
  318. // state
  319. int _state;
  320. int _adaptive_state;
  321. unsigned int _step_number; // long int
  322. double _offset;
  323. double _up_threshold;
  324. double _down_threshold;
  325. std::deque<record_t> _archive;
  326. };
  327. }
  328. }
  329. }
  330. #endif