Quantifier.hpp 15 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377
  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 = 1
  45. };
  46. enum outputs {
  47. OUT = 1
  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_port({IN, "in"});
  65. this->output_port({OUT, "out"});
  66. this->observables({{UP, "up"},
  67. {DOWN, "down"},
  68. {VALUE, "value"}});
  69. _adaptive = context.parameters().allow_offsets;
  70. _adaptive_state = _adaptive ? POSSIBLE : IMPOSSIBLE;
  71. _zero_init_offset = context.parameters().zero_init_offset;
  72. _step_size = context.parameters().quantum;
  73. assert(_step_size > 0);
  74. _past_length = context.parameters().archive_length;
  75. assert(_past_length > 2);
  76. }
  77. virtual ~Quantifier() { }
  78. virtual void dconf(typename Time::type t, typename Time::type e,
  79. const common::Bag<Time>& bag)
  80. {
  81. dint(t);
  82. dext(t, e, bag);
  83. }
  84. virtual void dint(typename Time::type /* t */)
  85. {
  86. switch (_state) {
  87. case INIT:
  88. break;
  89. case IDLE:
  90. break;
  91. case RESPONSE:
  92. _state = IDLE;
  93. break;
  94. }
  95. }
  96. virtual void dext(typename Time::type t, typename Time::type e,
  97. const common::Bag<Time>& bag)
  98. {
  99. std::for_each(bag.begin(), bag.end(),
  100. [this, t, e](const common::ExternalEvent<Time>& event) {
  101. IntegratorData data;
  102. double shifting_factor;
  103. double value;
  104. int cnt;
  105. event.data()(data);
  106. value = data.value;
  107. if (_state == INIT) {
  108. init_step_number_and_offset(value);
  109. update_thresholds();
  110. _state = RESPONSE;
  111. } else {
  112. cnt = 0;
  113. while (value >= _up_threshold or value <= _down_threshold) {
  114. cnt++;
  115. if (value >= _up_threshold) {
  116. _step_number++;
  117. } else {
  118. _step_number--;
  119. }
  120. switch (_adaptive_state) {
  121. case IMPOSSIBLE:
  122. update_thresholds();
  123. break;
  124. case POSSIBLE:
  125. if (value >= _up_threshold) {
  126. store_change(_step_size, t);
  127. } else {
  128. store_change(-_step_size, t);
  129. }
  130. shifting_factor = shift_quanta();
  131. assert(shifting_factor >= 0
  132. and shifting_factor <= 1);
  133. if (shifting_factor != 0 and shifting_factor != 1) {
  134. if (value >= _up_threshold) {
  135. update_thresholds(shifting_factor,
  136. DIRECTION_DOWN);
  137. } else {
  138. update_thresholds(shifting_factor,
  139. DIRECTION_UP);
  140. }
  141. _adaptive_state = DONE;
  142. } else {
  143. update_thresholds();
  144. }
  145. break;
  146. case DONE:
  147. init_step_number_and_offset(value);
  148. _adaptive_state = POSSIBLE;
  149. update_thresholds();
  150. break;
  151. }
  152. }
  153. }
  154. });
  155. _state = RESPONSE;
  156. }
  157. virtual void start(typename Time::type /* time */)
  158. {
  159. _offset = 0;
  160. _state = INIT;
  161. }
  162. virtual typename Time::type ta(typename Time::type /* time */)
  163. {
  164. switch (_state) {
  165. case INIT:
  166. case IDLE:
  167. return Time::infinity;
  168. case RESPONSE:
  169. return 0.0;
  170. }
  171. return Time::infinity;
  172. }
  173. virtual common::Bag<Time> lambda(typename Time::type /* time */) const
  174. {
  175. common::Bag<Time> msgs;
  176. const QuantifierData data = {_up_threshold, _down_threshold};
  177. msgs.push_back(common::ExternalEvent<Time>(OUT, data));
  178. return msgs;
  179. }
  180. virtual common::Value observe(const typename Time::type& /* t */,
  181. unsigned int index) const
  182. {
  183. switch (index) {
  184. case UP:
  185. return (double) _up_threshold;
  186. case DOWN:
  187. return (double) _down_threshold;
  188. case VALUE:
  189. return (double) (_up_threshold - _down_threshold);
  190. default:
  191. return common::Value();
  192. }
  193. }
  194. private:
  195. typedef enum {
  196. DIRECTION_UP, DIRECTION_DOWN
  197. } Direction;
  198. void init_step_number_and_offset(double value)
  199. {
  200. _step_number = static_cast<long int>(std::floor(value / _step_size));
  201. if (_zero_init_offset) {
  202. _offset = 0;
  203. } else {
  204. _offset = value - static_cast<double>(_step_number) * _step_size;
  205. }
  206. }
  207. bool monotonous(unsigned int range)
  208. {
  209. if ((range + 1) > _archive.size()) {
  210. return false;
  211. }
  212. for (size_t i = 0; i < range; i++) {
  213. if (_archive[i].value * _archive[i + 1].value < 0) {
  214. return false;
  215. }
  216. }
  217. return true;
  218. }
  219. bool oscillating(unsigned int range)
  220. {
  221. if ((range + 1) > _archive.size()) {
  222. return false;
  223. }
  224. for (size_t i = _archive.size() - range; i < _archive.size() - 1; i++) {
  225. if (_archive[i].value * _archive[i + 1].value > 0) {
  226. return false;
  227. }
  228. }
  229. return true;
  230. }
  231. double shift_quanta()
  232. {
  233. double factor = 0;
  234. if (oscillating(_past_length - 1) and
  235. _archive.back().date - _archive.front().date != 0) {
  236. double acc;
  237. double local_estim;
  238. int cnt;
  239. acc = 0;
  240. cnt = 0;
  241. for (size_t i = 0; i < _archive.size() - 2; ++i) {
  242. if (0 != (_archive[i + 2].date - _archive[i].date)) {
  243. if ((_archive.back().value * _archive[i + 1].value) > 0) {
  244. local_estim =
  245. 1 - (_archive[i + 1].date - _archive[i].date) /
  246. (_archive[i + 2].date - _archive[i].date);
  247. } else {
  248. local_estim = (_archive[i + 1].date - _archive[i].date) /
  249. (_archive[i + 2].date - _archive[i].date);
  250. }
  251. acc += local_estim;
  252. cnt++;
  253. }
  254. }
  255. acc = acc / cnt;
  256. factor = acc;
  257. _archive.resize(0);
  258. }
  259. return factor;
  260. }
  261. void store_change(double val, const typename Time::type& time)
  262. {
  263. record_t record;
  264. record.date = time;
  265. record.value = val;
  266. _archive.push_back(record);
  267. while (_archive.size() > _past_length) {
  268. _archive.pop_front();
  269. }
  270. }
  271. void update_thresholds()
  272. {
  273. auto step_number = static_cast<double>(_step_number);
  274. _up_threshold = _offset + _step_size * (step_number + 1);
  275. _down_threshold = _offset + _step_size * (step_number - 1);
  276. }
  277. void update_thresholds(double factor)
  278. {
  279. auto step_number = static_cast<double>(_step_number);
  280. _up_threshold = _offset + _step_size * (step_number + (1 - factor));
  281. _down_threshold = _offset + _step_size * (step_number - (1 - factor));
  282. }
  283. void update_thresholds(double factor, Direction d)
  284. {
  285. auto step_number = static_cast<double>(_step_number);
  286. if (d == DIRECTION_UP) {
  287. _up_threshold = _offset + _step_size * (step_number + (1 - factor));
  288. _down_threshold = _offset + _step_size * (step_number - 1);
  289. } else {
  290. _up_threshold = _offset + _step_size * (step_number + 1);
  291. _down_threshold = _offset + _step_size * (step_number - (1 - factor));
  292. }
  293. }
  294. typedef enum vars {
  295. UP, DOWN, VALUE
  296. } Observable;
  297. typedef enum {
  298. INIT, IDLE, RESPONSE
  299. } State;
  300. typedef enum {
  301. IMPOSSIBLE, POSSIBLE, DONE
  302. } AdaptiveState;
  303. struct record_t {
  304. double value;
  305. typename Time::type date;
  306. };
  307. // parameters
  308. bool _adaptive;
  309. bool _zero_init_offset;
  310. unsigned int _past_length;
  311. double _step_size;
  312. // state
  313. int _state;
  314. int _adaptive_state;
  315. unsigned int _step_number; // long int
  316. double _offset;
  317. double _up_threshold;
  318. double _down_threshold;
  319. std::deque<record_t> _archive;
  320. };
  321. }
  322. }
  323. }
  324. #endif