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