"""Kuramoto oscillator network with empirical coupling strengths."""
from typing import List, Optional, Tuple
import numpy as np
from mercurial.params.empirical import KuramotoParams
[docs]
class KuramotoOscillator:
def __init__(self, natural_freq: float, initial_phase: float = 0.0):
self.omega = natural_freq
self.phi = initial_phase % (2 * np.pi)
[docs]
def update(
self, dt: float, coupling: float, mean_field: complex, noise_amplitude: float = 0.05
) -> None:
coupling_term = coupling * np.imag(mean_field * np.exp(-1j * self.phi))
noise = noise_amplitude * np.sqrt(dt) * np.random.normal()
self.phi += (self.omega + coupling_term) * dt + noise
self.phi %= 2 * np.pi
[docs]
class KuramotoNetwork:
def __init__(
self,
n_oscillators: int,
natural_freqs: Optional[List[float]] = None,
params: Optional[KuramotoParams] = None,
coupling_strength: Optional[float] = None,
noise_amplitude: Optional[float] = None,
):
if params is None:
params = KuramotoParams()
self.K = coupling_strength if coupling_strength is not None else params.coupling_moderate
self.sigma = noise_amplitude if noise_amplitude is not None else params.phase_noise
self.N = n_oscillators
if natural_freqs is None:
mean_omega = 2 * np.pi * 40.0
std_omega = 2 * np.pi * 5.0
self.omegas = np.random.normal(mean_omega, std_omega, n_oscillators)
else:
self.omegas = np.array(natural_freqs)
self.oscillators = [KuramotoOscillator(omega) for omega in self.omegas]
[docs]
def mean_field(self) -> complex:
phases = np.array([osc.phi for osc in self.oscillators])
return np.mean(np.exp(1j * phases))
[docs]
def update(self, dt: float) -> None:
mf = self.mean_field()
for osc in self.oscillators:
osc.update(dt, self.K, mf, self.sigma)
[docs]
def order_parameter(self) -> Tuple[float, float]:
phases = np.array([osc.phi for osc in self.oscillators])
z = np.mean(np.exp(1j * phases))
return np.abs(z), np.angle(z)