sync.py
1 """!
2 
3 @brief Neural Network: Oscillatory Neural Network based on Kuramoto model
4 @details Implementation based on paper @cite article::syncnet::1, @cite article::nnet::sync::1, @cite inproceedings::net::sync::1.
5 
6 @authors Andrei Novikov (pyclustering@yandex.ru)
7 @date 2014-2020
8 @copyright GNU Public License
9 
10 @cond GNU_PUBLIC_LICENSE
11  PyClustering is free software: you can redistribute it and/or modify
12  it under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  PyClustering is distributed in the hope that it will be useful,
17  but WITHOUT ANY WARRANTY; without even the implied warranty of
18  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19  GNU General Public License for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with this program. If not, see <http://www.gnu.org/licenses/>.
23 @endcond
24 
25 """
26 
27 import math
28 import numpy
29 import random
30 import warnings
31 
32 try:
33  import matplotlib.pyplot as plt
34  import matplotlib.animation as animation
35 except Exception as error_instance:
36  warnings.warn("Impossible to import matplotlib (please, install 'matplotlib'), pyclustering's visualization "
37  "functionality is not available (details: '%s')." % str(error_instance))
38 
39 import pyclustering.core.sync_wrapper as wrapper
40 
41 from pyclustering.core.wrapper import ccore_library
42 
43 from scipy.integrate import odeint
44 
45 from pyclustering.nnet import network, conn_represent, conn_type, initial_type, solve_type
46 from pyclustering.utils import pi, draw_dynamics, draw_dynamics_set, set_ax_param
47 
48 
50  """!
51  @brief Provides services to calculate order parameter and local order parameter that are used
52  for synchronization level estimation.
53 
54  """
55  @staticmethod
56  def calculate_sync_order(oscillator_phases):
57  """!
58  @brief Calculates level of global synchronization (order parameter) for input phases.
59  @details This parameter is tend 1.0 when the oscillatory network close to global synchronization and it tend to 0.0 when
60  desynchronization is observed in the network.
61 
62  @param[in] oscillator_phases (list): List of oscillator phases that are used for level of global synchronization.
63 
64  @return (double) Level of global synchronization (order parameter).
65 
66  @see calculate_order_parameter()
67 
68  """
69 
70  exp_amount = 0.0
71  average_phase = 0.0
72 
73  for phase in oscillator_phases:
74  exp_amount += math.expm1(abs(1j * phase))
75  average_phase += phase
76 
77  exp_amount /= len(oscillator_phases)
78  average_phase = math.expm1(abs(1j * (average_phase / len(oscillator_phases))))
79 
80  return abs(average_phase) / abs(exp_amount)
81 
82 
83  @staticmethod
84  def calculate_local_sync_order(oscillator_phases, oscillatory_network):
85  """!
86  @brief Calculates level of local synchorization (local order parameter) for input phases for the specified network.
87  @details This parameter is tend 1.0 when the oscillatory network close to local synchronization and it tend to 0.0 when
88  desynchronization is observed in the network.
89 
90  @param[in] oscillator_phases (list): List of oscillator phases that are used for level of local (partial) synchronization.
91  @param[in] oscillatory_network (sync): Instance of oscillatory network whose connections are required for calculation.
92 
93  @return (double) Level of local synchronization (local order parameter).
94 
95  """
96 
97  exp_amount = 0.0
98  num_neigh = 0.0
99 
100  for i in range(0, len(oscillatory_network), 1):
101  for j in range(0, len(oscillatory_network), 1):
102  if oscillatory_network.has_connection(i, j) is True:
103  exp_amount += math.exp(-abs(oscillator_phases[j] - oscillator_phases[i]))
104  num_neigh += 1.0
105 
106  if num_neigh == 0:
107  num_neigh = 1.0
108 
109  return exp_amount / num_neigh
110 
111 
112 
114  """!
115  @brief Represents output dynamic of Sync.
116 
117  """
118 
119  @property
120  def output(self):
121  """!
122  @brief (list) Returns output dynamic of the Sync network (phase coordinates of each oscillator in the network) during simulation.
123 
124  """
125  if (self._ccore_sync_dynamic_pointer is not None) and ((self._dynamic is None) or (len(self._dynamic) == 0)):
126  self._dynamic = wrapper.sync_dynamic_get_output(self._ccore_sync_dynamic_pointer)
127 
128  return self._dynamic
129 
130 
131  @property
132  def time(self):
133  """!
134  @brief (list) Returns sampling times when dynamic is measured during simulation.
135 
136  """
137  if (self._ccore_sync_dynamic_pointer is not None) and ((self._time is None) or (len(self._time) == 0)):
138  self._time = wrapper.sync_dynamic_get_time(self._ccore_sync_dynamic_pointer)
139 
140  return self._time
141 
142 
143  def __init__(self, phase, time, ccore=None):
144  """!
145  @brief Constructor of Sync dynamic.
146 
147  @param[in] phase (list): Dynamic of oscillators on each step of simulation. If ccore pointer is specified than it can be ignored.
148  @param[in] time (list): Simulation time.
149  @param[in] ccore (ctypes.pointer): Pointer to CCORE sync_dynamic instance in memory.
150 
151  """
152 
153  self._dynamic = phase
154  self._time = time
155  self._ccore_sync_dynamic_pointer = ccore
156 
157 
158  def __del__(self):
159  """!
160  @brief Default destructor of Sync dynamic.
161 
162  """
163  if self._ccore_sync_dynamic_pointer is not None:
164  wrapper.sync_dynamic_destroy(self._ccore_sync_dynamic_pointer)
165 
166 
167  def __len__(self):
168  """!
169  @brief Returns number of simulation steps that are stored in dynamic.
170  @return (uint) Number of simulation steps that are stored in dynamic.
171 
172  """
173  if (self._ccore_sync_dynamic_pointer is not None):
174  return wrapper.sync_dynamic_get_size(self._ccore_sync_dynamic_pointer)
175 
176  return len(self._dynamic)
177 
178 
179  def __getitem__(self, index):
180  """!
181  @brief Indexing of the dynamic.
182 
183  """
184  if index == 0:
185  return self.time
186 
187  elif index == 1:
188  return self.output
189 
190  else:
191  raise NameError('Out of range ' + index + ': only indexes 0 and 1 are supported.')
192 
193 
194  def allocate_sync_ensembles(self, tolerance = 0.01, indexes = None, iteration = None):
195  """!
196  @brief Allocate clusters in line with ensembles of synchronous oscillators where each synchronous ensemble corresponds to only one cluster.
197 
198  @param[in] tolerance (double): Maximum error for allocation of synchronous ensemble oscillators.
199  @param[in] indexes (list): List of real object indexes and it should be equal to amount of oscillators (in case of 'None' - indexes are in range [0; amount_oscillators]).
200  @param[in] iteration (uint): Iteration of simulation that should be used for allocation.
201 
202  @return (list) Groups (lists) of indexes of synchronous oscillators.
203  For example [ [index_osc1, index_osc3], [index_osc2], [index_osc4, index_osc5] ].
204 
205  """
206 
207  if self._ccore_sync_dynamic_pointer is not None:
208  ensembles = wrapper.sync_dynamic_allocate_sync_ensembles(self._ccore_sync_dynamic_pointer, tolerance, iteration)
209 
210  if indexes is not None:
211  for ensemble in ensembles:
212  for index in range(len(ensemble)):
213  ensemble[index] = indexes[ensemble[index]]
214 
215  return ensembles
216 
217  if (self._dynamic is None) or (len(self._dynamic) == 0):
218  return []
219 
220  number_oscillators = len(self._dynamic[0])
221  last_state = None
222 
223  if iteration is None:
224  last_state = self._dynamic[len(self._dynamic) - 1]
225  else:
226  last_state = self._dynamic[iteration]
227 
228  clusters = []
229  if number_oscillators > 0:
230  clusters.append([0])
231 
232  for i in range(1, number_oscillators, 1):
233  cluster_allocated = False
234  for cluster in clusters:
235  for neuron_index in cluster:
236  last_state_shifted = abs(last_state[i] - 2 * pi)
237 
238  if ( ( (last_state[i] < (last_state[neuron_index] + tolerance)) and (last_state[i] > (last_state[neuron_index] - tolerance)) ) or
239  ( (last_state_shifted < (last_state[neuron_index] + tolerance)) and (last_state_shifted > (last_state[neuron_index] - tolerance)) ) ):
240  cluster_allocated = True
241 
242  real_index = i
243  if indexes is not None:
244  real_index = indexes[i]
245 
246  cluster.append(real_index)
247  break
248 
249  if cluster_allocated is True:
250  break
251 
252  if cluster_allocated is False:
253  clusters.append([i])
254 
255  return clusters
256 
257 
258  def allocate_phase_matrix(self, grid_width = None, grid_height = None, iteration = None):
259  """!
260  @brief Returns 2D matrix of phase values of oscillators at the specified iteration of simulation.
261  @details User should ensure correct matrix sizes in line with following expression grid_width x grid_height that should be equal to
262  amount of oscillators otherwise exception is thrown. If grid_width or grid_height are not specified than phase matrix size
263  will by calculated automatically by square root.
264 
265  @param[in] grid_width (uint): Width of the allocated matrix.
266  @param[in] grid_height (uint): Height of the allocated matrix.
267  @param[in] iteration (uint): Number of iteration of simulation for which correlation matrix should be allocated.
268  If iternation number is not specified, the last step of simulation is used for the matrix allocation.
269 
270  @return (list) Phase value matrix of oscillators with size [number_oscillators x number_oscillators].
271 
272  """
273 
274  output_dynamic = self.output
275 
276  if (output_dynamic is None) or (len(output_dynamic) == 0):
277  return []
278 
279  current_dynamic = output_dynamic[len(output_dynamic) - 1]
280  if iteration is not None:
281  current_dynamic = output_dynamic[iteration]
282 
283  width_matrix = grid_width
284  height_matrix = grid_height
285  number_oscillators = len(current_dynamic)
286  if (width_matrix is None) or (height_matrix is None):
287  width_matrix = int(math.ceil(math.sqrt(number_oscillators)))
288  height_matrix = width_matrix
289 
290  if (number_oscillators != width_matrix * height_matrix):
291  raise NameError("Impossible to allocate phase matrix with specified sizes, amout of neurons should be equal to grid_width * grid_height.");
292 
293  phase_matrix = [[0.0 for _ in range(width_matrix)] for _ in range(height_matrix)]
294  for i in range(height_matrix):
295  for j in range(width_matrix):
296  phase_matrix[i][j] = current_dynamic[j + i * width_matrix]
297 
298  return phase_matrix
299 
300 
301  def allocate_correlation_matrix(self, iteration = None):
302  """!
303  @brief Allocate correlation matrix between oscillators at the specified step of simulation.
304 
305  @param[in] iteration (uint): Number of iteration of simulation for which correlation matrix should be allocated.
306  If iternation number is not specified, the last step of simulation is used for the matrix allocation.
307 
308  @return (list) Correlation matrix between oscillators with size [number_oscillators x number_oscillators].
309 
310  """
311 
312  if self._ccore_sync_dynamic_pointer is not None:
313  return wrapper.sync_dynamic_allocate_correlation_matrix(self._ccore_sync_dynamic_pointer, iteration)
314 
315  if (self._dynamic is None) or (len(self._dynamic) == 0):
316  return []
317 
318  dynamic = self._dynamic
319  current_dynamic = dynamic[len(dynamic) - 1]
320 
321  if iteration is not None:
322  current_dynamic = dynamic[iteration]
323 
324  number_oscillators = len(dynamic[0])
325  affinity_matrix = [[0.0 for i in range(number_oscillators)] for j in range(number_oscillators)]
326 
327  for i in range(number_oscillators):
328  for j in range(number_oscillators):
329  phase1 = current_dynamic[i]
330  phase2 = current_dynamic[j]
331 
332  affinity_matrix[i][j] = abs(math.sin(phase1 - phase2))
333 
334  return affinity_matrix
335 
336 
337  def calculate_order_parameter(self, start_iteration=None, stop_iteration=None):
338  """!
339  @brief Calculates level of global synchorization (order parameter).
340  @details This parameter is tend 1.0 when the oscillatory network close to global synchronization and it tend to 0.0 when
341  desynchronization is observed in the network. Order parameter is calculated using following equation:
342 
343  \f[
344  r_{c}=\frac{1}{Ne^{i\varphi }}\sum_{j=0}^{N}e^{i\theta_{j}};
345  \f]
346 
347  where \f$\varphi\f$ is a average phase coordinate in the network, \f$N\f$ is an amount of oscillators in the network.
348 
349  @param[in] start_iteration (uint): The first iteration that is used for calculation, if 'None' then the last iteration is used.
350  @param[in] stop_iteration (uint): The last iteration that is used for calculation, if 'None' then 'start_iteration' + 1 is used.
351 
352  Example:
353  @code
354  oscillatory_network = sync(16, type_conn = conn_type.ALL_TO_ALL);
355  output_dynamic = oscillatory_network.simulate_static(100, 10);
356 
357  print("Order parameter at the last step: ", output_dynamic.calculate_order_parameter());
358  print("Order parameter at the first step:", output_dynamic.calculate_order_parameter(0));
359  print("Order parameter evolution between 40 and 50 steps:", output_dynamic.calculate_order_parameter(40, 50));
360  @endcode
361 
362  @return (list) List of levels of global synchronization (order parameter evolution).
363 
364  @see order_estimator
365 
366  """
367 
368  (start_iteration, stop_iteration) = self.__get_start_stop_iterations(start_iteration, stop_iteration)
369 
370  if self._ccore_sync_dynamic_pointer is not None:
371  return wrapper.sync_dynamic_calculate_order(self._ccore_sync_dynamic_pointer, start_iteration, stop_iteration)
372 
373  sequence_order = []
374  for index in range(start_iteration, stop_iteration):
375  sequence_order.append(order_estimator.calculate_sync_order(self.output[index]))
376 
377  return sequence_order
378 
379 
380  def calculate_local_order_parameter(self, oscillatory_network, start_iteration = None, stop_iteration = None):
381  """!
382  @brief Calculates local order parameter.
383  @details Local order parameter or so-called level of local or partial synchronization is calculated by following expression:
384 
385  \f[
386  r_{c}=\left | \sum_{i=0}^{N} \frac{1}{N_{i}} \sum_{j=0}e^{ \theta_{j} - \theta_{i} } \right |;
387  \f]
388 
389  where N - total amount of oscillators in the network and \f$N_{i}\f$ - amount of neighbors of oscillator with index \f$i\f$.
390 
391  @param[in] oscillatory_network (sync): Sync oscillatory network whose structure of connections is required for calculation.
392  @param[in] start_iteration (uint): The first iteration that is used for calculation, if 'None' then the last iteration is used.
393  @param[in] stop_iteration (uint): The last iteration that is used for calculation, if 'None' then 'start_iteration' + 1 is used.
394 
395  @return (list) List of levels of local (partial) synchronization (local order parameter evolution).
396 
397  """
398 
399  (start_iteration, stop_iteration) = self.__get_start_stop_iterations(start_iteration, stop_iteration)
400 
401  if self._ccore_sync_dynamic_pointer is not None:
402  network_pointer = oscillatory_network._ccore_network_pointer
403  return wrapper.sync_dynamic_calculate_local_order(self._ccore_sync_dynamic_pointer, network_pointer, start_iteration, stop_iteration)
404 
405  sequence_local_order = []
406  for index in range(start_iteration, stop_iteration):
407  sequence_local_order.append(order_estimator.calculate_local_sync_order(self.output[index], oscillatory_network))
408 
409  return sequence_local_order
410 
411 
412  def __get_start_stop_iterations(self, start_iteration, stop_iteration):
413  """!
414  @brief Aplly rules for start_iteration and stop_iteration parameters.
415 
416  @param[in] start_iteration (uint): The first iteration that is used for calculation.
417  @param[in] stop_iteration (uint): The last iteration that is used for calculation.
418 
419  @return (tuple) New the first iteration and the last.
420 
421  """
422  if start_iteration is None:
423  start_iteration = len(self) - 1
424 
425  if stop_iteration is None:
426  stop_iteration = start_iteration + 1
427 
428  return start_iteration, stop_iteration
429 
430 
431 
433  """!
434  @brief Visualizer of output dynamic of sync network (Sync).
435 
436  """
437 
438  @staticmethod
439  def show_output_dynamic(sync_output_dynamic):
440  """!
441  @brief Shows output dynamic (output of each oscillator) during simulation.
442 
443  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network.
444 
445  @see show_output_dynamics
446 
447  """
448 
449  draw_dynamics(sync_output_dynamic.time, sync_output_dynamic.output, x_title="t", y_title="phase", y_lim=[0, 2 * 3.14])
450 
451 
452  @staticmethod
453  def show_output_dynamics(sync_output_dynamics):
454  """!
455  @brief Shows several output dynamics (output of each oscillator) during simulation.
456  @details Each dynamic is presented on separate plot.
457 
458  @param[in] sync_output_dynamics (list): list of output dynamics 'sync_dynamic' of the Sync network.
459 
460  @see show_output_dynamic
461 
462  """
463 
464  draw_dynamics_set(sync_output_dynamics, "t", "phase", None, [0, 2 * 3.14], False, False)
465 
466 
467  @staticmethod
468  def show_correlation_matrix(sync_output_dynamic, iteration=None):
469  """!
470  @brief Shows correlation matrix between oscillators at the specified iteration.
471 
472  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network.
473  @param[in] iteration (uint): Number of iteration of simulation for which correlation matrix should be
474  allocated. If iteration number is not specified, the last step of simulation is used for the matrix
475  allocation.
476 
477  """
478 
479  _ = plt.figure()
480  correlation_matrix = sync_output_dynamic.allocate_correlation_matrix(iteration)
481 
482  plt.imshow(correlation_matrix, cmap = plt.get_cmap('cool'), interpolation='kaiser', vmin=0.0, vmax=1.0)
483  plt.show()
484 
485 
486  @staticmethod
487  def show_phase_matrix(sync_output_dynamic, grid_width=None, grid_height=None, iteration=None):
488  """!
489  @brief Shows 2D matrix of phase values of oscillators at the specified iteration.
490  @details User should ensure correct matrix sizes in line with following expression grid_width x grid_height that should be equal to
491  amount of oscillators otherwise exception is thrown. If grid_width or grid_height are not specified than phase matrix size
492  will by calculated automatically by square root.
493 
494  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network whose phase matrix should be shown.
495  @param[in] grid_width (uint): Width of the phase matrix.
496  @param[in] grid_height (uint): Height of the phase matrix.
497  @param[in] iteration (uint): Number of iteration of simulation for which correlation matrix should be allocated.
498  If iternation number is not specified, the last step of simulation is used for the matrix allocation.
499 
500  """
501 
502  _ = plt.figure()
503  phase_matrix = sync_output_dynamic.allocate_phase_matrix(grid_width, grid_height, iteration)
504 
505  plt.imshow(phase_matrix, cmap = plt.get_cmap('jet'), interpolation='kaiser', vmin=0.0, vmax=2.0 * math.pi)
506  plt.show()
507 
508 
509  @staticmethod
510  def show_order_parameter(sync_output_dynamic, start_iteration=None, stop_iteration=None):
511  """!
512  @brief Shows evolution of order parameter (level of global synchronization in the network).
513 
514  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network whose evolution of global synchronization should be visualized.
515  @param[in] start_iteration (uint): The first iteration that is used for calculation, if 'None' then the first is used
516  @param[in] stop_iteration (uint): The last iteration that is used for calculation, if 'None' then the last is used.
517 
518  """
519 
520  (start_iteration, stop_iteration) = sync_visualizer.__get_start_stop_iterations(sync_output_dynamic, start_iteration, stop_iteration)
521 
522  order_parameter = sync_output_dynamic.calculate_order_parameter(start_iteration, stop_iteration)
523  axis = plt.subplot(111)
524  plt.plot(sync_output_dynamic.time[start_iteration:stop_iteration], order_parameter, 'b-', linewidth=2.0)
525  set_ax_param(axis, "t", "R (order parameter)", None, [0.0, 1.05])
526 
527  plt.show()
528 
529 
530  @staticmethod
531  def show_local_order_parameter(sync_output_dynamic, oscillatory_network, start_iteration=None, stop_iteration=None):
532  """!
533  @brief Shows evolution of local order parameter (level of local synchronization in the network).
534 
535  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network whose evolution of global synchronization should be visualized.
536  @param[in] oscillatory_network (sync): Sync oscillatory network whose structure of connections is required for calculation.
537  @param[in] start_iteration (uint): The first iteration that is used for calculation, if 'None' then the first is used
538  @param[in] stop_iteration (uint): The last iteration that is used for calculation, if 'None' then the last is used.
539 
540  """
541  (start_iteration, stop_iteration) = sync_visualizer.__get_start_stop_iterations(sync_output_dynamic, start_iteration, stop_iteration)
542 
543  order_parameter = sync_output_dynamic.calculate_local_order_parameter(oscillatory_network, start_iteration, stop_iteration)
544  axis = plt.subplot(111)
545  plt.plot(sync_output_dynamic.time[start_iteration:stop_iteration], order_parameter, 'b-', linewidth=2.0)
546  set_ax_param(axis, "t", "R (local order parameter)", None, [0.0, 1.05])
547 
548  plt.show()
549 
550 
551  @staticmethod
552  def animate_output_dynamic(sync_output_dynamic, animation_velocity = 75, save_movie = None):
553  """!
554  @brief Shows animation of output dynamic (output of each oscillator) during simulation on a circle from [0; 2pi].
555 
556  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network.
557  @param[in] animation_velocity (uint): Interval between frames in milliseconds.
558  @param[in] save_movie (string): If it is specified then animation will be stored to file that is specified in this parameter.
559 
560  """
561 
562  figure = plt.figure()
563 
564  dynamic = sync_output_dynamic.output[0]
565  artist, = plt.polar(dynamic, [1.0] * len(dynamic), 'o', color='blue')
566 
567  def init_frame():
568  return [artist]
569 
570  def frame_generation(index_dynamic):
571  dynamic = sync_output_dynamic.output[index_dynamic]
572  artist.set_data(dynamic, [1.0] * len(dynamic))
573 
574  return [artist]
575 
576  phase_animation = animation.FuncAnimation(figure, frame_generation, len(sync_output_dynamic), interval = animation_velocity, init_func = init_frame, repeat_delay = 5000);
577 
578  if save_movie is not None:
579  phase_animation.save(save_movie, writer='ffmpeg', fps=15, bitrate=1500)
580  else:
581  plt.show()
582 
583 
584  @staticmethod
585  def animate_correlation_matrix(sync_output_dynamic, animation_velocity = 75, colormap = 'cool', save_movie = None):
586  """!
587  @brief Shows animation of correlation matrix between oscillators during simulation.
588 
589  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network.
590  @param[in] animation_velocity (uint): Interval between frames in milliseconds.
591  @param[in] colormap (string): Name of colormap that is used by matplotlib ('gray', 'pink', 'cool', spring', etc.).
592  @param[in] save_movie (string): If it is specified then animation will be stored to file that is specified in this parameter.
593 
594  """
595 
596  figure = plt.figure()
597 
598  correlation_matrix = sync_output_dynamic.allocate_correlation_matrix(0)
599  artist = plt.imshow(correlation_matrix, cmap = plt.get_cmap(colormap), interpolation='kaiser', vmin = 0.0, vmax = 1.0)
600 
601  def init_frame():
602  return [ artist ]
603 
604  def frame_generation(index_dynamic):
605  correlation_matrix = sync_output_dynamic.allocate_correlation_matrix(index_dynamic)
606  artist.set_data(correlation_matrix)
607 
608  return [artist]
609 
610  correlation_animation = animation.FuncAnimation(figure, frame_generation, len(sync_output_dynamic), init_func = init_frame, interval = animation_velocity , repeat_delay = 1000, blit = True)
611 
612  if save_movie is not None:
613  correlation_animation.save(save_movie, writer='ffmpeg', fps=15, bitrate=1500)
614  else:
615  plt.show()
616 
617 
618  @staticmethod
619  def animate_phase_matrix(sync_output_dynamic, grid_width=None, grid_height=None, animation_velocity=75, colormap='jet', save_movie=None):
620  """!
621  @brief Shows animation of phase matrix between oscillators during simulation on 2D stage.
622  @details If grid_width or grid_height are not specified than phase matrix size will by calculated automatically by square root.
623 
624  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network.
625  @param[in] grid_width (uint): Width of the phase matrix.
626  @param[in] grid_height (uint): Height of the phase matrix.
627  @param[in] animation_velocity (uint): Interval between frames in milliseconds.
628  @param[in] colormap (string): Name of colormap that is used by matplotlib ('gray', 'pink', 'cool', spring', etc.).
629  @param[in] save_movie (string): If it is specified then animation will be stored to file that is specified in this parameter.
630 
631  """
632 
633  figure = plt.figure()
634 
635  def init_frame():
636  return frame_generation(0)
637 
638  def frame_generation(index_dynamic):
639  figure.clf()
640  axis = figure.add_subplot(111)
641 
642  phase_matrix = sync_output_dynamic.allocate_phase_matrix(grid_width, grid_height, index_dynamic)
643  axis.imshow(phase_matrix, cmap=plt.get_cmap(colormap), interpolation='kaiser', vmin=0.0, vmax=2.0 * math.pi)
644  artist = figure.gca()
645 
646  return [artist]
647 
648  phase_animation = animation.FuncAnimation(figure, frame_generation, len(sync_output_dynamic), init_func = init_frame, interval = animation_velocity , repeat_delay = 1000);
649 
650  if save_movie is not None:
651  phase_animation.save(save_movie, writer='ffmpeg', fps=15, bitrate=1500)
652  else:
653  plt.show()
654 
655 
656  @staticmethod
657  def __get_start_stop_iterations(sync_output_dynamic, start_iteration, stop_iteration):
658  """!
659  @brief Apply rule of preparation for start iteration and stop iteration values.
660 
661  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network.
662  @param[in] start_iteration (uint): The first iteration that is used for calculation.
663  @param[in] stop_iteration (uint): The last iteration that is used for calculation.
664 
665  @return (tuple) New values of start and stop iterations.
666 
667  """
668  if start_iteration is None:
669  start_iteration = 0
670 
671  if stop_iteration is None:
672  stop_iteration = len(sync_output_dynamic)
673 
674  return start_iteration, stop_iteration
675 
676 
677  @staticmethod
678  def animate(sync_output_dynamic, title=None, save_movie=None):
679  """!
680  @brief Shows animation of phase coordinates and animation of correlation matrix together for the Sync dynamic output on the same figure.
681 
682  @param[in] sync_output_dynamic (sync_dynamic): Output dynamic of the Sync network.
683  @param[in] title (string): Title of the animation that is displayed on a figure if it is specified.
684  @param[in] save_movie (string): If it is specified then animation will be stored to file that is specified in this parameter.
685 
686  """
687 
688  dynamic = sync_output_dynamic.output[0]
689  correlation_matrix = sync_output_dynamic.allocate_correlation_matrix(0)
690 
691  figure = plt.figure(1)
692  if title is not None:
693  figure.suptitle(title, fontsize = 26, fontweight = 'bold')
694 
695  ax1 = figure.add_subplot(121, projection='polar')
696  ax2 = figure.add_subplot(122)
697 
698  artist1, = ax1.plot(dynamic, [1.0] * len(dynamic), marker='o', color='blue', ls='')
699  artist2 = ax2.imshow(correlation_matrix, cmap = plt.get_cmap('Accent'), interpolation='kaiser')
700 
701  def init_frame():
702  return [artist1, artist2]
703 
704  def frame_generation(index_dynamic):
705  dynamic = sync_output_dynamic.output[index_dynamic]
706  artist1.set_data(dynamic, [1.0] * len(dynamic))
707 
708  correlation_matrix = sync_output_dynamic.allocate_correlation_matrix(index_dynamic)
709  artist2.set_data(correlation_matrix)
710 
711  return [artist1, artist2]
712 
713  dynamic_animation = animation.FuncAnimation(figure, frame_generation, len(sync_output_dynamic), interval=75, init_func=init_frame, repeat_delay=5000)
714 
715  if save_movie is not None:
716  dynamic_animation.save(save_movie, writer='ffmpeg', fps=15, bitrate=1500)
717  else:
718  plt.show()
719 
720 
721 
723  """!
724  @brief Model of oscillatory network that is based on the Kuramoto model of synchronization.
725 
726  @details CCORE option can be used to use the pyclustering core - C/C++ shared library for processing that significantly increases performance.
727 
728  """
729 
730  def __init__(self, num_osc, weight = 1, frequency = 0, type_conn = conn_type.ALL_TO_ALL, representation = conn_represent.MATRIX, initial_phases = initial_type.RANDOM_GAUSSIAN, ccore = True):
731  """!
732  @brief Constructor of oscillatory network is based on Kuramoto model.
733 
734  @param[in] num_osc (uint): Number of oscillators in the network.
735  @param[in] weight (double): Coupling strength of the links between oscillators.
736  @param[in] frequency (double): Multiplier of internal frequency of the oscillators.
737  @param[in] type_conn (conn_type): Type of connection between oscillators in the network (all-to-all, grid, bidirectional list, etc.).
738  @param[in] representation (conn_represent): Internal representation of connection in the network: matrix or list.
739  @param[in] initial_phases (initial_type): Type of initialization of initial phases of oscillators (random, uniformly distributed, etc.).
740  @param[in] ccore (bool): If True simulation is performed by CCORE library (C++ implementation of pyclustering).
741 
742  """
743 
744  self._ccore_network_pointer = None; # Pointer to CCORE Sync implementation of the network.
745 
746  if ( (ccore is True) and ccore_library.workable() ):
747  self._ccore_network_pointer = wrapper.sync_create_network(num_osc, weight, frequency, type_conn, initial_phases);
748  self._num_osc = num_osc;
749  self._conn_represent = conn_represent.MATRIX;
750 
751  else:
752  super().__init__(num_osc, type_conn, representation);
753 
754  self._weight = weight;
755 
756  self._phases = list();
757  self._freq = list();
758 
759  random.seed();
760  for index in range(0, num_osc, 1):
761  if (initial_phases == initial_type.RANDOM_GAUSSIAN):
762  self._phases.append(random.random() * 2 * pi);
763 
764  elif (initial_phases == initial_type.EQUIPARTITION):
765  self._phases.append( pi / num_osc * index);
766 
767  self._freq.append(random.random() * frequency);
768 
769 
770  def __del__(self):
771  """!
772  @brief Destructor of oscillatory network is based on Kuramoto model.
773 
774  """
775 
776  if (self._ccore_network_pointer is not None):
777  wrapper.sync_destroy_network(self._ccore_network_pointer);
778  self._ccore_network_pointer = None;
779 
780 
781  def sync_order(self):
782  """!
783  @brief Calculates current level of global synchorization (order parameter) in the network.
784  @details This parameter is tend 1.0 when the oscillatory network close to global synchronization and it tend to 0.0 when
785  desynchronization is observed in the network. Order parameter is calculated using following equation:
786 
787  \f[
788  r_{c}=\frac{1}{Ne^{i\varphi }}\sum_{j=0}^{N}e^{i\theta_{j}};
789  \f]
790 
791  where \f$\varphi\f$ is a average phase coordinate in the network, \f$N\f$ is an amount of oscillators in the network.
792 
793  Example:
794  @code
795  oscillatory_network = sync(16, type_conn = conn_type.ALL_TO_ALL);
796  output_dynamic = oscillatory_network.simulate_static(100, 10);
797 
798  if (oscillatory_network.sync_order() < 0.9): print("Global synchronization is not reached yet.");
799  else: print("Global synchronization is reached.");
800  @endcode
801 
802  @return (double) Level of global synchronization (order parameter).
803 
804  @see sync_local_order()
805 
806  """
807 
808  if (self._ccore_network_pointer is not None):
809  return wrapper.sync_order(self._ccore_network_pointer);
810 
811  return order_estimator.calculate_sync_order(self._phases);
812 
813 
814  def sync_local_order(self):
815  """!
816  @brief Calculates current level of local (partial) synchronization in the network.
817 
818  @return (double) Level of local (partial) synchronization.
819 
820  @see sync_order()
821 
822  """
823 
824  if (self._ccore_network_pointer is not None):
825  return wrapper.sync_local_order(self._ccore_network_pointer);
826 
827  return order_estimator.calculate_local_sync_order(self._phases, self);
828 
829 
830  def _phase_kuramoto(self, teta, t, argv):
831  """!
832  @brief Returns result of phase calculation for specified oscillator in the network.
833 
834  @param[in] teta (double): Phase of the oscillator that is differentiated.
835  @param[in] t (double): Current time of simulation.
836  @param[in] argv (tuple): Index of the oscillator in the list.
837 
838  @return (double) New phase for specified oscillator (don't assign here).
839 
840  """
841 
842  index = argv;
843  phase = 0;
844  for k in range(0, self._num_osc):
845  if (self.has_connection(index, k) == True):
846  phase += math.sin(self._phases[k] - teta);
847 
848  return ( self._freq[index] + (phase * self._weight / self._num_osc) );
849 
850 
851  def simulate(self, steps, time, solution = solve_type.FAST, collect_dynamic = True):
852  """!
853  @brief Performs static simulation of Sync oscillatory network.
854 
855  @param[in] steps (uint): Number steps of simulations during simulation.
856  @param[in] time (double): Time of simulation.
857  @param[in] solution (solve_type): Type of solution (solving).
858  @param[in] collect_dynamic (bool): If True - returns whole dynamic of oscillatory network, otherwise returns only last values of dynamics.
859 
860  @return (list) Dynamic of oscillatory network. If argument 'collect_dynamic' = True, than return dynamic for the whole simulation time,
861  otherwise returns only last values (last step of simulation) of dynamic.
862 
863  @see simulate_dynamic()
864  @see simulate_static()
865 
866  """
867 
868  return self.simulate_static(steps, time, solution, collect_dynamic);
869 
870 
871  def simulate_dynamic(self, order = 0.998, solution = solve_type.FAST, collect_dynamic = False, step = 0.1, int_step = 0.01, threshold_changes = 0.0000001):
872  """!
873  @brief Performs dynamic simulation of the network until stop condition is not reached. Stop condition is defined by input argument 'order'.
874 
875  @param[in] order (double): Order of process synchronization, distributed 0..1.
876  @param[in] solution (solve_type): Type of solution.
877  @param[in] collect_dynamic (bool): If True - returns whole dynamic of oscillatory network, otherwise returns only last values of dynamics.
878  @param[in] step (double): Time step of one iteration of simulation.
879  @param[in] int_step (double): Integration step, should be less than step.
880  @param[in] threshold_changes (double): Additional stop condition that helps prevent infinite simulation, defines limit of changes of oscillators between current and previous steps.
881 
882  @return (list) Dynamic of oscillatory network. If argument 'collect_dynamic' = True, than return dynamic for the whole simulation time,
883  otherwise returns only last values (last step of simulation) of dynamic.
884 
885  @see simulate()
886  @see simulate_static()
887 
888  """
889 
890  if (self._ccore_network_pointer is not None):
891  ccore_instance_dynamic = wrapper.sync_simulate_dynamic(self._ccore_network_pointer, order, solution, collect_dynamic, step, int_step, threshold_changes);
892  return sync_dynamic(None, None, ccore_instance_dynamic);
893 
894  # For statistics and integration
895  time_counter = 0;
896 
897  # Prevent infinite loop. It's possible when required state cannot be reached.
898  previous_order = 0;
899  current_order = self.sync_local_order();
900 
901  # If requested input dynamics
902  dyn_phase = [];
903  dyn_time = [];
904  if (collect_dynamic == True):
905  dyn_phase.append(self._phases);
906  dyn_time.append(0);
907 
908  # Execute until sync state will be reached
909  while (current_order < order):
910  # update states of oscillators
911  self._phases = self._calculate_phases(solution, time_counter, step, int_step);
912 
913  # update time
914  time_counter += step;
915 
916  # if requested input dynamic
917  if (collect_dynamic == True):
918  dyn_phase.append(self._phases);
919  dyn_time.append(time_counter);
920 
921  # update orders
922  previous_order = current_order;
923  current_order = self.sync_local_order();
924 
925  # hang prevention
926  if (abs(current_order - previous_order) < threshold_changes):
927  # print("Warning: sync_network::simulate_dynamic - simulation is aborted due to low level of convergence rate (order = " + str(current_order) + ").");
928  break;
929 
930  if (collect_dynamic != True):
931  dyn_phase.append(self._phases);
932  dyn_time.append(time_counter);
933 
934  output_sync_dynamic = sync_dynamic(dyn_phase, dyn_time, None);
935  return output_sync_dynamic;
936 
937 
938  def simulate_static(self, steps, time, solution = solve_type.FAST, collect_dynamic = False):
939  """!
940  @brief Performs static simulation of oscillatory network.
941 
942  @param[in] steps (uint): Number steps of simulations during simulation.
943  @param[in] time (double): Time of simulation.
944  @param[in] solution (solve_type): Type of solution.
945  @param[in] collect_dynamic (bool): If True - returns whole dynamic of oscillatory network, otherwise returns only last values of dynamics.
946 
947  @return (list) Dynamic of oscillatory network. If argument 'collect_dynamic' = True, than return dynamic for the whole simulation time,
948  otherwise returns only last values (last step of simulation) of dynamic.
949 
950  @see simulate()
951  @see simulate_dynamic()
952 
953  """
954 
955  if (self._ccore_network_pointer is not None):
956  ccore_instance_dynamic = wrapper.sync_simulate_static(self._ccore_network_pointer, steps, time, solution, collect_dynamic);
957  return sync_dynamic(None, None, ccore_instance_dynamic);
958 
959  dyn_phase = [];
960  dyn_time = [];
961 
962  if (collect_dynamic == True):
963  dyn_phase.append(self._phases);
964  dyn_time.append(0);
965 
966  step = time / steps;
967  int_step = step / 10.0;
968 
969  for t in numpy.arange(step, time + step, step):
970  # update states of oscillators
971  self._phases = self._calculate_phases(solution, t, step, int_step);
972 
973  # update states of oscillators
974  if (collect_dynamic == True):
975  dyn_phase.append(self._phases);
976  dyn_time.append(t);
977 
978  if (collect_dynamic != True):
979  dyn_phase.append(self._phases);
980  dyn_time.append(time);
981 
982  output_sync_dynamic = sync_dynamic(dyn_phase, dyn_time);
983  return output_sync_dynamic;
984 
985 
986  def _calculate_phases(self, solution, t, step, int_step):
987  """!
988  @brief Calculates new phases for oscillators in the network in line with current step.
989 
990  @param[in] solution (solve_type): Type solver of the differential equation.
991  @param[in] t (double): Time of simulation.
992  @param[in] step (double): Step of solution at the end of which states of oscillators should be calculated.
993  @param[in] int_step (double): Step differentiation that is used for solving differential equation.
994 
995  @return (list) New states (phases) for oscillators.
996 
997  """
998 
999  next_phases = [0.0] * self._num_osc; # new oscillator _phases
1000 
1001  for index in range (0, self._num_osc, 1):
1002  if (solution == solve_type.FAST):
1003  result = self._phases[index] + self._phase_kuramoto(self._phases[index], 0, index);
1004  next_phases[index] = self._phase_normalization(result);
1005 
1006  elif ( (solution == solve_type.RK4) or (solution == solve_type.RKF45) ):
1007  result = odeint(self._phase_kuramoto, self._phases[index], numpy.arange(t - step, t, int_step), (index , ));
1008  next_phases[index] = self._phase_normalization(result[len(result) - 1][0]);
1009 
1010  else:
1011  raise NameError("Solver '" + str(solution) + "' is not supported");
1012 
1013  return next_phases;
1014 
1015 
1016  def _phase_normalization(self, teta):
1017  """!
1018  @brief Normalization of phase of oscillator that should be placed between [0; 2 * pi].
1019 
1020  @param[in] teta (double): phase of oscillator.
1021 
1022  @return (double) Normalized phase.
1023 
1024  """
1025 
1026  norm_teta = teta;
1027  while (norm_teta > (2.0 * pi)) or (norm_teta < 0):
1028  if (norm_teta > (2.0 * pi)):
1029  norm_teta -= 2.0 * pi;
1030  else:
1031  norm_teta += 2.0 * pi;
1032 
1033  return norm_teta;
1034 
1035 
1036  def get_neighbors(self, index):
1037  """!
1038  @brief Finds neighbors of the oscillator with specified index.
1039 
1040  @param[in] index (uint): index of oscillator for which neighbors should be found in the network.
1041 
1042  @return (list) Indexes of neighbors of the specified oscillator.
1043 
1044  """
1045 
1046  if ( (self._ccore_network_pointer is not None) and (self._osc_conn is None) ):
1047  self._osc_conn = wrapper.sync_connectivity_matrix(self._ccore_network_pointer);
1048 
1049  return super().get_neighbors(index);
1050 
1051 
1052  def has_connection(self, i, j):
1053  """!
1054  @brief Returns True if there is connection between i and j oscillators and False - if connection doesn't exist.
1055 
1056  @param[in] i (uint): index of an oscillator in the network.
1057  @param[in] j (uint): index of an oscillator in the network.
1058 
1059  """
1060 
1061  if ( (self._ccore_network_pointer is not None) and (self._osc_conn is None) ):
1062  self._osc_conn = wrapper.sync_connectivity_matrix(self._ccore_network_pointer);
1063 
1064  return super().has_connection(i, j);
def simulate_dynamic(self, order=0.998, solution=solve_type.FAST, collect_dynamic=False, step=0.1, int_step=0.01, threshold_changes=0.0000001)
Performs dynamic simulation of the network until stop condition is not reached.
Definition: sync.py:871
def animate_correlation_matrix(sync_output_dynamic, animation_velocity=75, colormap='cool', save_movie=None)
Shows animation of correlation matrix between oscillators during simulation.
Definition: sync.py:585
def __del__(self)
Default destructor of Sync dynamic.
Definition: sync.py:158
def __get_start_stop_iterations(self, start_iteration, stop_iteration)
Aplly rules for start_iteration and stop_iteration parameters.
Definition: sync.py:412
def _calculate_phases(self, solution, t, step, int_step)
Calculates new phases for oscillators in the network in line with current step.
Definition: sync.py:986
def calculate_local_sync_order(oscillator_phases, oscillatory_network)
Calculates level of local synchorization (local order parameter) for input phases for the specified n...
Definition: sync.py:84
def time(self)
(list) Returns sampling times when dynamic is measured during simulation.
Definition: sync.py:132
def animate_output_dynamic(sync_output_dynamic, animation_velocity=75, save_movie=None)
Shows animation of output dynamic (output of each oscillator) during simulation on a circle from [0; ...
Definition: sync.py:552
Represents output dynamic of Sync.
Definition: sync.py:113
def show_phase_matrix(sync_output_dynamic, grid_width=None, grid_height=None, iteration=None)
Shows 2D matrix of phase values of oscillators at the specified iteration.
Definition: sync.py:487
def allocate_correlation_matrix(self, iteration=None)
Allocate correlation matrix between oscillators at the specified step of simulation.
Definition: sync.py:301
def simulate_static(self, steps, time, solution=solve_type.FAST, collect_dynamic=False)
Performs static simulation of oscillatory network.
Definition: sync.py:938
Utils that are used by modules of pyclustering.
Definition: __init__.py:1
def sync_local_order(self)
Calculates current level of local (partial) synchronization in the network.
Definition: sync.py:814
def animate_phase_matrix(sync_output_dynamic, grid_width=None, grid_height=None, animation_velocity=75, colormap='jet', save_movie=None)
Shows animation of phase matrix between oscillators during simulation on 2D stage.
Definition: sync.py:619
def allocate_sync_ensembles(self, tolerance=0.01, indexes=None, iteration=None)
Allocate clusters in line with ensembles of synchronous oscillators where each synchronous ensemble c...
Definition: sync.py:194
def show_output_dynamic(sync_output_dynamic)
Shows output dynamic (output of each oscillator) during simulation.
Definition: sync.py:439
def __len__(self)
Returns number of simulation steps that are stored in dynamic.
Definition: sync.py:167
def has_connection(self, i, j)
Returns True if there is connection between i and j oscillators and False - if connection doesn&#39;t exi...
Definition: __init__.py:366
def show_output_dynamics(sync_output_dynamics)
Shows several output dynamics (output of each oscillator) during simulation.
Definition: sync.py:453
def _phase_kuramoto(self, teta, t, argv)
Returns result of phase calculation for specified oscillator in the network.
Definition: sync.py:830
def __getitem__(self, index)
Indexing of the dynamic.
Definition: sync.py:179
def calculate_local_order_parameter(self, oscillatory_network, start_iteration=None, stop_iteration=None)
Calculates local order parameter.
Definition: sync.py:380
def show_correlation_matrix(sync_output_dynamic, iteration=None)
Shows correlation matrix between oscillators at the specified iteration.
Definition: sync.py:468
def calculate_sync_order(oscillator_phases)
Calculates level of global synchronization (order parameter) for input phases.
Definition: sync.py:56
Common network description that consists of information about oscillators and connection between them...
Definition: __init__.py:97
Provides services to calculate order parameter and local order parameter that are used for synchroniz...
Definition: sync.py:49
def __init__(self, phase, time, ccore=None)
Constructor of Sync dynamic.
Definition: sync.py:143
def show_order_parameter(sync_output_dynamic, start_iteration=None, stop_iteration=None)
Shows evolution of order parameter (level of global synchronization in the network).
Definition: sync.py:510
def __del__(self)
Destructor of oscillatory network is based on Kuramoto model.
Definition: sync.py:770
def get_neighbors(self, index)
Finds neighbors of the oscillator with specified index.
Definition: sync.py:1036
def has_connection(self, i, j)
Returns True if there is connection between i and j oscillators and False - if connection doesn&#39;t exi...
Definition: sync.py:1052
def __init__(self, num_osc, weight=1, frequency=0, type_conn=conn_type.ALL_TO_ALL, representation=conn_represent.MATRIX, initial_phases=initial_type.RANDOM_GAUSSIAN, ccore=True)
Constructor of oscillatory network is based on Kuramoto model.
Definition: sync.py:730
def animate(sync_output_dynamic, title=None, save_movie=None)
Shows animation of phase coordinates and animation of correlation matrix together for the Sync dynami...
Definition: sync.py:678
def _phase_normalization(self, teta)
Normalization of phase of oscillator that should be placed between [0; 2 * pi].
Definition: sync.py:1016
def sync_order(self)
Calculates current level of global synchorization (order parameter) in the network.
Definition: sync.py:781
def allocate_phase_matrix(self, grid_width=None, grid_height=None, iteration=None)
Returns 2D matrix of phase values of oscillators at the specified iteration of simulation.
Definition: sync.py:258
def calculate_order_parameter(self, start_iteration=None, stop_iteration=None)
Calculates level of global synchorization (order parameter).
Definition: sync.py:337
Visualizer of output dynamic of sync network (Sync).
Definition: sync.py:432
def show_local_order_parameter(sync_output_dynamic, oscillatory_network, start_iteration=None, stop_iteration=None)
Shows evolution of local order parameter (level of local synchronization in the network).
Definition: sync.py:531
def simulate(self, steps, time, solution=solve_type.FAST, collect_dynamic=True)
Performs static simulation of Sync oscillatory network.
Definition: sync.py:851
def output(self)
(list) Returns output dynamic of the Sync network (phase coordinates of each oscillator in the networ...
Definition: sync.py:120
Neural and oscillatory network module.
Definition: __init__.py:1
Model of oscillatory network that is based on the Kuramoto model of synchronization.
Definition: sync.py:722