pyclustering  0.10.1
pyclustring is a Python, C++ data mining library.
fsync.py
1 """!
2 
3 @brief Oscillatory Neural Network based on Kuramoto model in frequency domain.
4 @details Implementation based on paper @cite book::chemical_oscillatorions_waves.
5 
6 @authors Andrei Novikov (pyclustering@yandex.ru)
7 @date 2014-2020
8 @copyright BSD-3-Clause
9 
10 """
11 
12 
13 import numpy
14 import random
15 import pyclustering.utils
16 
17 from scipy.integrate import odeint
18 
19 from pyclustering.nnet import network, conn_type, conn_represent
20 
21 
23  """!
24  @brief Represents output dynamic of Sync in frequency domain.
25 
26  """
27 
28 
29  def __init__(self, amplitude, time):
30  """!
31  @brief Constructor of Sync dynamic in frequency domain.
32 
33  @param[in] amplitude (list): Dynamic of oscillators on each step of simulation.
34  @param[in] time (list): Simulation time where each time-point corresponds to amplitude-point.
35 
36  """
37 
38  self.__amplitude = amplitude
39  self.__time = time
40 
41 
42  @property
43  def output(self):
44  """!
45  @brief (list) Returns output dynamic of the Sync network (amplitudes of each oscillator in the network) during simulation.
46 
47  """
48 
49  return self.__amplitude
50 
51 
52  @property
53  def time(self):
54  """!
55  @brief (list) Returns time-points corresponds to dynamic-points points.
56 
57  """
58 
59  return self.__time
60 
61 
62  def __len__(self):
63  """!
64  @brief (uint) Returns number of simulation steps that are stored in dynamic.
65 
66  """
67 
68  return len(self.__amplitude)
69 
70 
71  def __getitem__(self, index):
72  """!
73  @brief Indexing of the dynamic.
74 
75  """
76  if index == 0:
77  return self.__time
78 
79  elif index == 1:
80  return self.__amplitude
81 
82  else:
83  raise NameError('Out of range ' + index + ': only indexes 0 and 1 are supported.')
84 
85 
86  def allocate_sync_ensembles(self, tolerance=0.1):
87  """!
88  @brief Allocate clusters in line with ensembles of synchronous oscillators where each synchronous ensemble corresponds to only one cluster.
89 
90  @param[in] tolerance (double): Maximum error for allocation of synchronous ensemble oscillators.
91 
92  @return (list) Grours of indexes of synchronous oscillators, for example, [ [index_osc1, index_osc3], [index_osc2], [index_osc4, index_osc5] ].
93 
94  """
95 
96  return pyclustering.utils.allocate_sync_ensembles(self.__amplitude, tolerance, 0.0)
97 
98 
99  def extract_number_oscillations(self, index, amplitude_threshold):
100  """!
101  @brief Extracts number of oscillations of specified oscillator.
102 
103  @param[in] index (uint): Index of oscillator whose dynamic is considered.
104  @param[in] amplitude_threshold (double): Amplitude threshold when oscillation is taken into account, for example,
105  when oscillator amplitude is greater than threshold then oscillation is incremented.
106 
107  @return (uint) Number of oscillations of specified oscillator.
108 
109  """
110 
111  return pyclustering.utils.extract_number_oscillations(self.__amplitude, index, amplitude_threshold);
112 
113 
114 
116  """!
117  @brief Visualizer of output dynamic of sync network in frequency domain.
118 
119  """
120 
121  @staticmethod
122  def show_output_dynamic(fsync_output_dynamic):
123  """!
124  @brief Shows output dynamic (output of each oscillator) during simulation.
125 
126  @param[in] fsync_output_dynamic (fsync_dynamic): Output dynamic of the fSync network.
127 
128  @see show_output_dynamics
129 
130  """
131 
132  pyclustering.utils.draw_dynamics(fsync_output_dynamic.time, fsync_output_dynamic.output, x_title = "t", y_title = "amplitude");
133 
134 
135  @staticmethod
136  def show_output_dynamics(fsync_output_dynamics):
137  """!
138  @brief Shows several output dynamics (output of each oscillator) during simulation.
139  @details Each dynamic is presented on separate plot.
140 
141  @param[in] fsync_output_dynamics (list): list of output dynamics 'fsync_dynamic' of the fSync network.
142 
143  @see show_output_dynamic
144 
145  """
146 
147  pyclustering.utils.draw_dynamics_set(fsync_output_dynamics, "t", "amplitude", None, None, False, False);
148 
149 
150 
152  """!
153  @brief Model of oscillatory network that uses Landau-Stuart oscillator and Kuramoto model as a synchronization mechanism.
154  @details Dynamic of each oscillator in the network is described by following differential Landau-Stuart equation with feedback:
155 
156  \f[
157  \dot{z}_{i} = (i\omega_{i} + \rho^{2}_{i} - |z_{i}|^{2} )z_{i} + \frac{1}{N}\sum_{j=0}^{N}k_{ij}(z_{j} - z_{i});
158  \f]
159 
160  Where left part of the equation is Landau-Stuart equation and the right is a Kuramoto model for synchronization.
161  For solving this equation Runge-Kutta 4 method is used by default.
162 
163  Example:
164  @code
165  # Prepare oscillatory network parameters.
166  amount_oscillators = 3;
167  frequency = 1.0;
168  radiuses = [1.0, 2.0, 3.0];
169  coupling_strength = 1.0;
170 
171  # Create oscillatory network
172  oscillatory_network = fsync_network(amount_oscillators, frequency, radiuses, coupling_strength);
173 
174  # Simulate network during 200 steps on 10 time-units of time-axis.
175  output_dynamic = oscillatory_network.simulate(200, 10, True); # True is to collect whole output dynamic.
176 
177  # Visualize output result
178  fsync_visualizer.show_output_dynamic(output_dynamic);
179  @endcode
180 
181  Example of output dynamic of the network:
182  @image html fsync_sync_examples.png
183 
184  """
185 
186  __DEFAULT_FREQUENCY_VALUE = 1.0;
187  __DEFAULT_RADIUS_VALUE = 1.0;
188  __DEFAULT_COUPLING_STRENGTH = 1.0;
189 
190 
191  def __init__(self, num_osc, factor_frequency = 1.0, factor_radius = 1.0, factor_coupling = 1.0, type_conn = conn_type.ALL_TO_ALL, representation = conn_represent.MATRIX):
192  """!
193  @brief Constructor of oscillatory network based on synchronization Kuramoto model and Landau-Stuart oscillator.
194 
195  @param[in] num_osc (uint): Amount oscillators in the network.
196  @param[in] factor_frequency (double|list): Frequency of oscillators, it can be specified as common value for all oscillators by
197  single double value and for each separately by list.
198  @param[in] factor_radius (double|list): Radius of oscillators that affects amplitude, it can be specified as common value for all oscillators by
199  single double value and for each separately by list.
200  @param[in] factor_coupling (double): Coupling strength between oscillators.
201  @param[in] type_conn (conn_type): Type of connection between oscillators in the network (all-to-all, grid, bidirectional list, etc.).
202  @param[in] representation (conn_represent): Internal representation of connection in the network: matrix or list.
203 
204  """
205 
206  super().__init__(num_osc, type_conn, representation);
207 
208  self.__frequency = factor_frequency if isinstance(factor_frequency, list) else [ fsync_network.__DEFAULT_FREQUENCY_VALUE * factor_frequency for _ in range(num_osc) ];
209  self.__radius = factor_radius if isinstance(factor_radius, list) else [ fsync_network.__DEFAULT_RADIUS_VALUE * factor_radius for _ in range(num_osc) ];
210  self.__coupling_strength = fsync_network.__DEFAULT_COUPLING_STRENGTH * factor_coupling;
211  self.__properties = [ self.__oscillator_property(index) for index in range(self._num_osc) ];
212 
213  random.seed();
214  self.__amplitude = [ random.random() for _ in range(num_osc) ];
215 
216 
217  def simulate(self, steps, time, collect_dynamic = False):
218  """!
219  @brief Performs static simulation of oscillatory network.
220 
221  @param[in] steps (uint): Number simulation steps.
222  @param[in] time (double): Time of simulation.
223  @param[in] collect_dynamic (bool): If True - returns whole dynamic of oscillatory network, otherwise returns only last values of dynamics.
224 
225  @return (list) Dynamic of oscillatory network. If argument 'collect_dynamic' is True, than return dynamic for the whole simulation time,
226  otherwise returns only last values (last step of simulation) of output dynamic.
227 
228  @see simulate()
229  @see simulate_dynamic()
230 
231  """
232 
233  dynamic_amplitude, dynamic_time = ([], []) if collect_dynamic is False else ([self.__amplitude], [0]);
234 
235  step = time / steps;
236  int_step = step / 10.0;
237 
238  for t in numpy.arange(step, time + step, step):
239  self.__amplitude = self.__calculate(t, step, int_step);
240 
241  if collect_dynamic is True:
242  dynamic_amplitude.append([ numpy.real(amplitude)[0] for amplitude in self.__amplitude ]);
243  dynamic_time.append(t);
244 
245  if collect_dynamic is False:
246  dynamic_amplitude.append([ numpy.real(amplitude)[0] for amplitude in self.__amplitude ]);
247  dynamic_time.append(time);
248 
249  output_sync_dynamic = fsync_dynamic(dynamic_amplitude, dynamic_time);
250  return output_sync_dynamic;
251 
252 
253  def __calculate(self, t, step, int_step):
254  """!
255  @brief Calculates new amplitudes for oscillators in the network in line with current step.
256 
257  @param[in] t (double): Time of simulation.
258  @param[in] step (double): Step of solution at the end of which states of oscillators should be calculated.
259  @param[in] int_step (double): Step differentiation that is used for solving differential equation.
260 
261  @return (list) New states (phases) for oscillators.
262 
263  """
264 
265  next_amplitudes = [0.0] * self._num_osc;
266 
267  for index in range (0, self._num_osc, 1):
268  z = numpy.array(self.__amplitude[index], dtype = numpy.complex128, ndmin = 1);
269  result = odeint(self.__calculate_amplitude, z.view(numpy.float64), numpy.arange(t - step, t, int_step), (index , ));
270  next_amplitudes[index] = (result[len(result) - 1]).view(numpy.complex128);
271 
272  return next_amplitudes;
273 
274 
275  def __oscillator_property(self, index):
276  """!
277  @brief Calculate Landau-Stuart oscillator constant property that is based on frequency and radius.
278 
279  @param[in] index (uint): Oscillator index whose property is calculated.
280 
281  @return (double) Oscillator property.
282 
283  """
284 
285  return numpy.array(1j * self.__frequency[index] + self.__radius[index]**2, dtype = numpy.complex128, ndmin = 1);
286 
287 
288  def __landau_stuart(self, amplitude, index):
289  """!
290  @brief Calculate Landau-Stuart state.
291 
292  @param[in] amplitude (double): Current amplitude of oscillator.
293  @param[in] index (uint): Oscillator index whose state is calculated.
294 
295  @return (double) Landau-Stuart state.
296 
297  """
298 
299  return (self.__properties[index] - numpy.absolute(amplitude) ** 2) * amplitude;
300 
301 
302  def __synchronization_mechanism(self, amplitude, index):
303  """!
304  @brief Calculate synchronization part using Kuramoto synchronization mechanism.
305 
306  @param[in] amplitude (double): Current amplitude of oscillator.
307  @param[in] index (uint): Oscillator index whose synchronization influence is calculated.
308 
309  @return (double) Synchronization influence for the specified oscillator.
310 
311  """
312 
313  sync_influence = 0.0;
314 
315  for k in range(self._num_osc):
316  if self.has_connection(index, k) is True:
317  amplitude_neighbor = numpy.array(self.__amplitude[k], dtype = numpy.complex128, ndmin = 1);
318  sync_influence += amplitude_neighbor - amplitude;
319 
320  return sync_influence * self.__coupling_strength / self._num_osc;
321 
322 
323  def __calculate_amplitude(self, amplitude, t, argv):
324  """!
325  @brief Returns new amplitude value for particular oscillator that is defined by index that is in 'argv' argument.
326  @details The method is used for differential calculation.
327 
328  @param[in] amplitude (double): Current amplitude of oscillator.
329  @param[in] t (double): Current time of simulation.
330  @param[in] argv (uint): Index of the current oscillator.
331 
332  @return (double) New amplitude of the oscillator.
333 
334  """
335 
336  z = amplitude.view(numpy.complex);
337  dzdt = self.__landau_stuart(z, argv) + self.__synchronization_mechanism(z, argv);
338 
339  return dzdt.view(numpy.float64);
pyclustering.nnet.fsync.fsync_network.__calculate_amplitude
def __calculate_amplitude(self, amplitude, t, argv)
Returns new amplitude value for particular oscillator that is defined by index that is in 'argv' argu...
Definition: fsync.py:323
pyclustering.nnet.fsync.fsync_network.__radius
__radius
Definition: fsync.py:209
pyclustering.nnet.fsync.fsync_dynamic.__time
__time
Definition: fsync.py:39
pyclustering.nnet.fsync.fsync_network.__init__
def __init__(self, num_osc, factor_frequency=1.0, factor_radius=1.0, factor_coupling=1.0, type_conn=conn_type.ALL_TO_ALL, representation=conn_represent.MATRIX)
Constructor of oscillatory network based on synchronization Kuramoto model and Landau-Stuart oscillat...
Definition: fsync.py:191
pyclustering.nnet.fsync.fsync_network.__properties
__properties
Definition: fsync.py:211
pyclustering.nnet.fsync.fsync_dynamic.extract_number_oscillations
def extract_number_oscillations(self, index, amplitude_threshold)
Extracts number of oscillations of specified oscillator.
Definition: fsync.py:99
pyclustering.nnet.fsync.fsync_network.__landau_stuart
def __landau_stuart(self, amplitude, index)
Calculate Landau-Stuart state.
Definition: fsync.py:288
pyclustering.utils.draw_dynamics_set
def draw_dynamics_set(dynamics, xtitle=None, ytitle=None, xlim=None, ylim=None, xlabels=False, ylabels=False)
Draw lists of dynamics of neurons (oscillators) in the network.
Definition: __init__.py:957
pyclustering.nnet.fsync.fsync_network.__synchronization_mechanism
def __synchronization_mechanism(self, amplitude, index)
Calculate synchronization part using Kuramoto synchronization mechanism.
Definition: fsync.py:302
pyclustering.nnet.network._num_osc
int _num_osc
Definition: __init__.py:88
pyclustering.nnet.fsync.fsync_dynamic.__amplitude
__amplitude
Definition: fsync.py:38
pyclustering.nnet.fsync.fsync_visualizer
Visualizer of output dynamic of sync network in frequency domain.
Definition: fsync.py:115
pyclustering.utils.extract_number_oscillations
def extract_number_oscillations(osc_dyn, index=0, amplitude_threshold=1.0)
Extracts number of oscillations of specified oscillator.
Definition: __init__.py:592
pyclustering.utils.allocate_sync_ensembles
def allocate_sync_ensembles(dynamic, tolerance=0.1, threshold=1.0, ignore=None)
Allocate clusters in line with ensembles of synchronous oscillators where each synchronous ensemble c...
Definition: __init__.py:631
pyclustering.nnet.fsync.fsync_visualizer.show_output_dynamic
def show_output_dynamic(fsync_output_dynamic)
Shows output dynamic (output of each oscillator) during simulation.
Definition: fsync.py:122
pyclustering.utils.draw_dynamics
def draw_dynamics(t, dyn, x_title=None, y_title=None, x_lim=None, y_lim=None, x_labels=True, y_labels=True, separate=False, axes=None)
Draw dynamics of neurons (oscillators) in the network.
Definition: __init__.py:829
pyclustering.nnet.fsync.fsync_dynamic
Represents output dynamic of Sync in frequency domain.
Definition: fsync.py:22
pyclustering.nnet.fsync.fsync_network.__coupling_strength
__coupling_strength
Definition: fsync.py:210
pyclustering.nnet.fsync.fsync_dynamic.__len__
def __len__(self)
(uint) Returns number of simulation steps that are stored in dynamic.
Definition: fsync.py:62
pyclustering.nnet.fsync.fsync_network.__oscillator_property
def __oscillator_property(self, index)
Calculate Landau-Stuart oscillator constant property that is based on frequency and radius.
Definition: fsync.py:275
pyclustering.nnet.network
Common network description that consists of information about oscillators and connection between them...
Definition: __init__.py:82
pyclustering.nnet.network.has_connection
def has_connection(self, i, j)
Returns True if there is connection between i and j oscillators and False - if connection doesn't exi...
Definition: __init__.py:351
pyclustering.nnet.fsync.fsync_dynamic.time
def time(self)
(list) Returns time-points corresponds to dynamic-points points.
Definition: fsync.py:53
pyclustering.nnet.fsync.fsync_dynamic.__getitem__
def __getitem__(self, index)
Indexing of the dynamic.
Definition: fsync.py:71
pyclustering.nnet.fsync.fsync_network
Model of oscillatory network that uses Landau-Stuart oscillator and Kuramoto model as a synchronizati...
Definition: fsync.py:151
pyclustering.nnet
Neural and oscillatory network module. Consists of models of bio-inspired networks.
Definition: __init__.py:1
pyclustering.nnet.fsync.fsync_dynamic.allocate_sync_ensembles
def allocate_sync_ensembles(self, tolerance=0.1)
Allocate clusters in line with ensembles of synchronous oscillators where each synchronous ensemble c...
Definition: fsync.py:86
pyclustering.nnet.fsync.fsync_network.__calculate
def __calculate(self, t, step, int_step)
Calculates new amplitudes for oscillators in the network in line with current step.
Definition: fsync.py:253
pyclustering.nnet.fsync.fsync_visualizer.show_output_dynamics
def show_output_dynamics(fsync_output_dynamics)
Shows several output dynamics (output of each oscillator) during simulation.
Definition: fsync.py:136
pyclustering.nnet.fsync.fsync_network.__amplitude
__amplitude
Definition: fsync.py:214
pyclustering.nnet.fsync.fsync_network.simulate
def simulate(self, steps, time, collect_dynamic=False)
Performs static simulation of oscillatory network.
Definition: fsync.py:217
pyclustering.nnet.fsync.fsync_dynamic.__init__
def __init__(self, amplitude, time)
Constructor of Sync dynamic in frequency domain.
Definition: fsync.py:29
pyclustering.utils
Utils that are used by modules of pyclustering.
Definition: __init__.py:1
pyclustering.nnet.fsync.fsync_network.__frequency
__frequency
Definition: fsync.py:208
pyclustering.nnet.fsync.fsync_dynamic.output
def output(self)
(list) Returns output dynamic of the Sync network (amplitudes of each oscillator in the network) duri...
Definition: fsync.py:43