Skip to content

Synchronization

Synchronization: coupled phase oscillators, pulse-coupled oscillators, swarmalators, chaos, and Kuramoto order.

AggregateAveragingInteraction(label, *, source, target, r=4.0, kappa=0.4)

Bases: MapInteraction

Logistic map of an aggregate, x_t(t + 1) = r a (1 - a) with a = (1 - kappa) x_t + kappa x_s.

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
source Annotated[str, RoleSpec()]

The path of the unit whose value the aggregate reads.

required
target Annotated[str, RoleSpec()]

The path of the unit whose value the map sets.

required
r Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The growth rate of the logistic map.

4.0
kappa Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The weight of the source in the aggregate.

0.4

Examples:

>>> from types import SimpleNamespace
>>> from sage.all__sagemath_symbolics import var
>>> from skcomplex.models.physics.synchronization import AggregateAveragingInteraction
>>> averaging = AggregateAveragingInteraction('averaging', source='a', target='b').check().prepare()
>>> source, target = SimpleNamespace(x=var('u')), SimpleNamespace(x=var('v'))
>>> namespace = SimpleNamespace(source=source, target=target, r=4, kappa=var('kappa'))
>>> [law.value.factor() for law in averaging.run(namespace).laws_]
[-4*(kappa*u - kappa*v + v)*(kappa*u - kappa*v + v - 1)]
Source code in src/skcomplex/models/physics/synchronization/_aggregate_averaging.py
46
47
48
49
50
51
52
53
54
55
56
57
58
def __init__(
    self,
    label: str,
    *,
    source: Annotated[str, RoleSpec()],
    target: Annotated[str, RoleSpec()],
    r: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 4.0,
    kappa: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.4,
) -> None:
    super().__init__(label, target=target, state='x', function=self._average)
    self.source = source
    self.r = r
    self.kappa = kappa

DelayedCouplingInteraction(label, *, source, target, K=1.0, tau=1.0)

Bases: RateInteraction

Delayed sine coupling of one phase oscillator to another, theta_t' += K sin(theta_s(t - tau) - theta_t).

The target reads the phase of the source as it was a delay tau earlier.

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
source Annotated[str, RoleSpec()]

The path of the oscillator whose delayed phase the coupling reads.

required
target Annotated[str, RoleSpec()]

The path of the oscillator whose phase the coupling turns.

required
K Annotated[float, ParameterSpec()]

The coupling strength.

1.0
tau Annotated[float, ParameterSpec(domain=open(0, oo))]

The delay.

1.0

Examples:

>>> from types import SimpleNamespace
>>> from sage.all__sagemath_symbolics import var
>>> from skcomplex.models.physics.synchronization import DelayedCouplingInteraction
>>> coupling = DelayedCouplingInteraction('coupling', source='a', target='b').check().prepare()
>>> namespace = SimpleNamespace(source=SimpleNamespace(theta=var('a')), target=SimpleNamespace(theta=var('b')))
>>> [law.value for law in coupling.run(SimpleNamespace(**vars(namespace), K=var('K'), tau=1)).laws_]
[K*sin(-b + delayed(a, 1))]
Source code in src/skcomplex/models/physics/synchronization/_kuramoto.py
61
62
63
64
65
66
67
68
69
70
71
72
73
def __init__(
    self,
    label: str,
    *,
    source: Annotated[str, RoleSpec()],
    target: Annotated[str, RoleSpec()],
    K: Annotated[float, ParameterSpec()] = 1.0,
    tau: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, target=target, state='theta')
    self.source = source
    self.K = K
    self.tau = tau

DiffusiveCouplingInteraction(label, *, source, target, state='x2', k=0.1)

Bases: RateInteraction

Diffusive coupling of a state of one element to the same state of another, x_t' += k (x_s - x_t).

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
source Annotated[str, RoleSpec()]

The path of the element whose state the coupling reads.

required
target Annotated[str, RoleSpec()]

The path of the element whose state the coupling changes.

required
state str

The name of the coupled state.

'x2'
k Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The coupling strength.

0.1

Examples:

>>> from types import SimpleNamespace
>>> from sage.all__sagemath_symbolics import var
>>> from skcomplex.models.physics.synchronization import DiffusiveCouplingInteraction
>>> coupling = DiffusiveCouplingInteraction('coupling', source='a', target='b').check().prepare()
>>> source, target = SimpleNamespace(x2=var('u')), SimpleNamespace(x2=var('v'))
>>> [law.value for law in coupling.run(SimpleNamespace(source=source, target=target, k=var('k'))).laws_]
[k*(u - v)]
Source code in src/skcomplex/models/physics/synchronization/_diffusive_coupling.py
45
46
47
48
49
50
51
52
53
54
55
56
def __init__(
    self,
    label: str,
    *,
    source: Annotated[str, RoleSpec()],
    target: Annotated[str, RoleSpec()],
    state: str = 'x2',
    k: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
) -> None:
    super().__init__(label, target=target, state=state)
    self.source = source
    self.k = k

EdgeCouplingInteraction(label, *, target, links, K=1.0, alpha=0.0)

Bases: KuramotoCouplingInteraction

Kuramoto coupling along the weighted edges of a network, theta_i' += K w_ij sin(theta_j - theta_i).

The coupling runs once per edge of the network of its oscillators and reads the weight w of the edge from the group on the edges.

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
target Annotated[str, RoleSpec()]

The path of the group of oscillators on the nodes.

required
links Annotated[str, RoleSpec()]

The path of the group on the edges, whose parameter w weighs every edge.

required
K float

The coupling strength.

1.0
alpha float

The phase lag.

0.0

Examples:

>>> from skcomplex.models.physics.synchronization import EdgeCouplingInteraction
>>> EdgeCouplingInteraction('coupling', target='nodes', links='edges', K=0.5).weights
'w'
Source code in src/skcomplex/models/physics/synchronization/_edge_coupling.py
40
41
42
43
44
45
46
47
48
49
50
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    links: Annotated[str, RoleSpec()],
    K: float = 1.0,
    alpha: float = 0.0,
) -> None:
    super().__init__(label, target=target, K=K, alpha=alpha, weights='w')
    self.links = links

GlobalCouplingInteraction(label, *, target, state='x', K=1.0)

Bases: CouplingInteraction

Diffusive coupling of one state of every member to the mean of the group, x_i' += K (mean_j x_j - x_i).

The mean runs over every member, the member itself included.

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
target Annotated[str, RoleSpec()]

The path of the group.

required
state str

The name of the coupled state.

'x'
K Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The coupling strength.

1.0

Examples:

>>> from types import SimpleNamespace
>>> from sage.all__sagemath_symbolics import var
>>> from skcomplex.models.physics.synchronization import GlobalCouplingInteraction
>>> coupling = GlobalCouplingInteraction('coupling', target='oscillators', state='y').check().prepare()
>>> member, other = SimpleNamespace(y=var('y')), SimpleNamespace(y=var('z'))
>>> coupling.function(SimpleNamespace(K=var('K'), member=member, other=other))
-K*(y - z)
Source code in src/skcomplex/models/physics/synchronization/_amplitude_death.py
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    state: str = 'x',
    K: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(
        label,
        target=target,
        read=state,
        write=state,
        function=self._pull,
        reduce='mean',
        mode='group',
        include_self=True,
    )
    self.state = state
    self.K = K

IntegrateAndFireOscillatorElement(label, visual=None, *, x=0.0, S=2.0, gamma=1.0, epsilon=0.1)

Bases: Element

Integrate-and-fire oscillator that charges as x' = S - gamma x, and sends pulses of strength epsilon.

The oscillator charges towards S / gamma. Its firing at 1 is a crossing of the group it belongs to, which resets the most charged oscillators to 0 and lifts every other one by epsilon.

Parameters:

Name Type Description Default
label str

The label of the oscillator or the group.

required
visual VisualSpec | None

How the oscillator is drawn, or None for default_visual_.

None
x Annotated[float, StateSpec()]

The initial charge.

0.0
S Annotated[float, ParameterSpec(domain=open(0, oo))]

The charging input.

2.0
gamma Annotated[float, ParameterSpec(domain=open(0, oo))]

The leak rate.

1.0
epsilon Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The strength of a pulse.

0.1

Attributes:

Name Type Description
default_visual_

The visual of the oscillator when it sets none.

Examples:

>>> from skcomplex.models.physics.synchronization import IntegrateAndFireOscillatorElement
>>> IntegrateAndFireOscillatorElement('oscillator').epsilon
0.1
Source code in src/skcomplex/models/physics/synchronization/_integrate_and_fire_oscillator.py
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x: Annotated[float, StateSpec()] = 0.0,
    S: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 2.0,
    gamma: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    epsilon: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x
    self.S = S
    self.gamma = gamma
    self.epsilon = epsilon

LinkElement(label, visual=None, *, w=1.0)

Bases: Element

Edge of a network with a weight w, which couplings along the edge read.

Parameters:

Name Type Description Default
label str

The label of the edge or the group on the edges.

required
visual VisualSpec | None

How the edge is drawn, or None for default_visual_.

None
w Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The weight of the edge.

1.0

Attributes:

Name Type Description
default_visual_

The visual of the edge when it sets none.

Examples:

>>> from skcomplex.models.physics.synchronization import LinkElement
>>> LinkElement('link', w=0.5).w
0.5
Source code in src/skcomplex/models/physics/synchronization/_link.py
38
39
40
41
42
43
44
45
46
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    w: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.w = w

LogisticCouplingInteraction(label, *, target, state='x')

Bases: CouplingInteraction

Coupled logistic maps, x_i(t + 1) = (1 - epsilon) f(x_i) + epsilon mean_j f(x_j) with f(x) = r x (1 - x).

The map of every node mixes with the mean of the maps of its neighbours j on its lattice or network, with the parameters r and epsilon of the node. A node with no neighbours keeps (1 - epsilon) f(x_i). It is the map form of coupling, with the mean of the maps of the neighbours combined with the map of the node.

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
target Annotated[str, RoleSpec()]

The path of the group of maps, on the nodes of a lattice or network.

required
state str

The name of the state that the maps iterate.

'x'

Examples:

>>> import networkx as nx
>>> from skcomplex.core.elements import LatticeElement
>>> from skcomplex.core.spaces import DiscreteTime, NetworkSpace
>>> from skcomplex.core.system import System
>>> from skcomplex.engine.runtime import Simulation
>>> from skcomplex.models.physics.synchronization import LogisticCouplingInteraction, LogisticNodeElement
>>> maps = LatticeElement(
...     'maps',
...     element=LogisticNodeElement,
...     space=NetworkSpace(nx.path_graph(2)),
...     values={'x': [0.25, 0.5], 'epsilon': 0.5},
... )
>>> coupling = LogisticCouplingInteraction('coupling', target='maps')
>>> system = System(elements=[maps], interactions=[coupling], time=DiscreteTime())
>>> Simulation(system).check().prepare().run(until=1).trajectory_.values('maps.x')[-1].tolist()
[0.875, 0.875]
Source code in src/skcomplex/models/physics/synchronization/_logistic_coupling.py
52
53
54
55
56
57
58
59
60
61
62
63
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()], state: str = 'x') -> None:
    super().__init__(
        label,
        target=target,
        read=state,
        write=state,
        function=self._get_map,
        reduce='mean',
        form='map',
        combine=self._get_mix,
    )
    self.state = state

LogisticNodeElement(label, visual=None, *, x=0.5, r=4.0, epsilon=0.5)

Bases: Element

Logistic map on a node of a network, with its value, growth rate, and weight of the neighbours.

LogisticCouplingInteraction iterates x_i(t + 1) = (1 - epsilon) f(x_i) + epsilon mean_j f(x_j), where the map f(x) = r x (1 - x) of the node mixes with the mean of the maps of its neighbours j on the network.

Parameters:

Name Type Description Default
label str

The label of the node or the group on the nodes.

required
visual VisualSpec | None

How the node is drawn, or None for default_visual_.

None
x Annotated[float, StateSpec()]

The initial value.

0.5
r Annotated[float, ParameterSpec()]

The growth rate of the logistic map.

4.0
epsilon Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The weight of the neighbours.

0.5

Attributes:

Name Type Description
default_visual_

The visual of the node when it sets none.

Examples:

>>> from skcomplex.models.physics.synchronization import LogisticNodeElement
>>> LogisticNodeElement('node').epsilon
0.5
Source code in src/skcomplex/models/physics/synchronization/_logistic_node.py
47
48
49
50
51
52
53
54
55
56
57
58
59
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x: Annotated[float, StateSpec()] = 0.5,
    r: Annotated[float, ParameterSpec()] = 4.0,
    epsilon: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.5,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x
    self.r = r
    self.epsilon = epsilon

LogisticUnitElement(label, visual=None, *, x=0.3)

Bases: Element

Value x in [0, 1] of a logistic map, which maps of the coupling iterate.

Parameters:

Name Type Description Default
label str

The label of the unit.

required
visual VisualSpec | None

How the unit is drawn, or None for default_visual_.

None
x Annotated[float, StateSpec()]

The initial value.

0.3

Attributes:

Name Type Description
default_visual_

The visual of the unit when it sets none.

Examples:

>>> from skcomplex.models.physics.synchronization import LogisticUnitElement
>>> LogisticUnitElement('unit').x
0.3
Source code in src/skcomplex/models/physics/synchronization/_logistic_unit.py
35
36
37
38
39
40
41
42
43
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x: Annotated[float, StateSpec()] = 0.3,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x

RosslerOscillatorElement(label, visual=None, *, x1=1.0, x2=0.0, x3=0.0, omega=1.0, a=0.2, b=0.2, c=5.7)

Bases: Element

Rössler oscillator x1' = -omega x2 - x3, x2' = omega x1 + a x2, x3' = b + x3 (x1 - c).

It is chaotic for a = b = 0.2 and c = 5.7.

Parameters:

Name Type Description Default
label str

The label of the oscillator.

required
visual VisualSpec | None

How the oscillator is drawn, or None for default_visual_.

None
x1 Annotated[float, StateSpec()]

The initial first variable.

1.0
x2 Annotated[float, StateSpec()]

The initial second variable.

0.0
x3 Annotated[float, StateSpec()]

The initial third variable.

0.0
omega Annotated[float, ParameterSpec()]

The natural frequency.

1.0
a Annotated[float, ParameterSpec()]

The first Rössler parameter.

0.2
b Annotated[float, ParameterSpec()]

The second Rössler parameter.

0.2
c Annotated[float, ParameterSpec()]

The third Rössler parameter.

5.7

Attributes:

Name Type Description
default_visual_

The visual of the oscillator when it sets none.

Examples:

>>> from skcomplex.models.physics.synchronization import RosslerOscillatorElement
>>> RosslerOscillatorElement('oscillator').c
5.7
Source code in src/skcomplex/models/physics/synchronization/_rossler_oscillator.py
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x1: Annotated[float, StateSpec()] = 1.0,
    x2: Annotated[float, StateSpec()] = 0.0,
    x3: Annotated[float, StateSpec()] = 0.0,
    omega: Annotated[float, ParameterSpec()] = 1.0,
    a: Annotated[float, ParameterSpec()] = 0.2,
    b: Annotated[float, ParameterSpec()] = 0.2,
    c: Annotated[float, ParameterSpec()] = 5.7,
) -> None:
    super().__init__(label, visual=visual)
    self.x1 = x1
    self.x2 = x2
    self.x3 = x3
    self.omega = omega
    self.a = a
    self.b = b
    self.c = c

SimplexCouplingInteraction(label, *, target, K1=0.0, K2=1.0)

Bases: HyperedgeInteraction

Kuramoto coupling through the simplices of a hypergraph, pairs and triangles.

On a link, the focal oscillator i gets K1 sin(theta_j - theta_i), and on a triangle K2 sin(theta_j + theta_k - 2 theta_i) from its two other members j and k.

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
target Annotated[str, RoleSpec()]

The path of the group of oscillators on the nodes.

required
K1 Annotated[float, ParameterSpec()]

The coupling of every link.

0.0
K2 Annotated[float, ParameterSpec()]

The coupling of every triangle.

1.0

Examples:

>>> from skcomplex.models.physics.synchronization import SimplexCouplingInteraction
>>> SimplexCouplingInteraction('coupling', target='oscillators', K2=0.5).K2
0.5
Source code in src/skcomplex/models/physics/synchronization/_hypergraphs.py
55
56
57
58
59
60
61
62
63
64
65
66
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    K1: Annotated[float, ParameterSpec()] = 0.0,
    K2: Annotated[float, ParameterSpec()] = 1.0,
) -> None:
    super().__init__(label)
    self.target = target
    self.K1 = K1
    self.K2 = K2

SwarmalatorCouplingInteraction(label, first, second, cutoff=1000000.0, position='x', *, A=1.0, B=1.0, J=1.0, K=0.0, period=None, n_units=1)

Bases: PairCouplingInteraction

Coupling of two swarmalators through their separation r and their phase difference dtheta.

Each swarmalator moves toward the other at the speed (A + J cos(dtheta)) / N, away from it at B / (N |r|), and turns its phase toward the other's at the rate K sin(dtheta) / (N |r|), as in O'Keeffe, Ceron, and Petersen. It is a pair coupling of the positions and the phases jointly.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
first Annotated[str, RoleSpec()]

The path of the group of the first members.

required
second Annotated[str, RoleSpec()]

The path of the group of the second members.

required
cutoff Annotated[float, ParameterSpec(domain=open(0, oo))]

The distance beyond which pairs do not interact, large for the all-to-all model.

1000000.0
position str

The name of the vector state that holds the positions of the members.

'x'
A Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The strength of the attraction.

1.0
B Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The strength of the repulsion.

1.0
J Annotated[float, ParameterSpec()]

The influence of the phase similarity on the attraction.

1.0
K Annotated[float, ParameterSpec()]

The phase coupling, synchronizing when positive.

0.0
period float | None

The side of a periodic box, for minimum image distances, or None for an open space.

None
n_units int

The number of swarmalators N that normalizes the coupling.

1

Examples:

>>> from skcomplex.models.physics.synchronization import SwarmalatorCouplingInteraction
>>> SwarmalatorCouplingInteraction('coupling', first='units', second='units', J=0.5).J
0.5
Source code in src/skcomplex/models/physics/synchronization/_swarmalators.py
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
def __init__(
    self,
    label: str,
    first: Annotated[str, RoleSpec()],
    second: Annotated[str, RoleSpec()],
    cutoff: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1000000.0,
    position: str = 'x',
    *,
    A: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    B: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    J: Annotated[float, ParameterSpec()] = 1.0,
    K: Annotated[float, ParameterSpec()] = 0.0,
    period: float | None = None,
    n_units: int = 1,
) -> None:
    super().__init__(
        label,
        first,
        second,
        cutoff=cutoff,
        position=position,
        write=(position, 'theta'),
        function=self._get_terms,
        period=period,
    )
    self.A = A
    self.B = B
    self.J = J
    self.K = K
    self.n_units = n_units

SwarmalatorElement(label, visual=None, *, theta=0.0, omega=1.0, x=0.0)

Bases: SitedOscillatorElement

Phase oscillator whose place x is a position that moves by the rates of interactions.

Parameters:

Name Type Description Default
label str

The label of the swarmalator or the group.

required
visual VisualSpec | None

How the swarmalator is drawn, or None for default_visual_.

None
theta float

The initial phase.

0.0
omega float

The natural frequency.

1.0
x Annotated[ArrayLike, StateSpec(vector=True, role=position)]

The initial position. A number applies to every component.

0.0

Examples:

>>> from skcomplex.models.physics.synchronization import SwarmalatorElement
>>> SwarmalatorElement('unit', x=[0.0, 1.0]).x
[0.0, 1.0]
Source code in src/skcomplex/models/physics/synchronization/_swarmalator.py
39
40
41
42
43
44
45
46
47
48
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    theta: float = 0.0,
    omega: float = 1.0,
    x: Annotated[npt.ArrayLike, StateSpec(vector=True, role='position')] = 0.0,
) -> None:
    super().__init__(label, visual=visual, theta=theta, omega=omega, x=x)

WinfreeMeanFieldInteraction(label, *, target, K=1.0)

Bases: CouplingInteraction

Winfree coupling theta_i' += K mean_j P(theta_j) R(theta_i) over every member j of the group.

The influence P(theta) = 1 + cos(theta) that an oscillator sends and the sensitivity R(theta) = -sin(theta) with which it responds make every oscillator respond to the mean influence of the group.

Parameters:

Name Type Description Default
label str

The label of the interaction, unique among the parts of its system.

required
target Annotated[str, RoleSpec()]

The path of the group of oscillators.

required
K Annotated[float, ParameterSpec()]

The coupling strength.

1.0

Examples:

>>> from types import SimpleNamespace
>>> from sage.all__sagemath_symbolics import var
>>> from skcomplex.models.physics.synchronization import WinfreeMeanFieldInteraction
>>> coupling = WinfreeMeanFieldInteraction('coupling', target='oscillators')
>>> member, other = SimpleNamespace(theta=var('a')), SimpleNamespace(theta=var('b'))
>>> coupling.function(SimpleNamespace(K=var('K'), member=member, other=other))
-K*(cos(b) + 1)*sin(a)
Source code in src/skcomplex/models/physics/synchronization/_winfree_mean_field.py
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    K: Annotated[float, ParameterSpec()] = 1.0,
) -> None:
    super().__init__(
        label,
        target=target,
        read='theta',
        write='theta',
        function=self._response,
        reduce='mean',
        mode='group',
        include_self=True,
    )
    self.K = K

aggregate_averaging_synchronization(x1=0.3, x2=0.6, *, r=4.0, kappa=0.4)

Return two logistic maps coupled by aggregate averaging, which synchronize their chaos for strong coupling.

Each map iterates on the aggregate of its own value and the other's, x_i(t + 1) = r a_i (1 - a_i) with a_1 = (1 - kappa) x_1 + kappa x_2 and a_2 = (1 - kappa) x_2 + kappa x_1. The synchronized orbit is stable for kappa > (r - 1) / (2 r), which is 3 / 8 for r = 4, and the maps agree in one step for kappa = 1 / 2.

Parameters:

Name Type Description Default
x1 float

The initial value of the first map.

0.3
x2 float

The initial value of the second map.

0.6
r float

The growth rate of the logistic map.

4.0
kappa float

The coupling strength.

0.4

Returns:

Name Type Description
system System

The model, with the values first.x and second.x.

Examples:

>>> from skcomplex.models.physics.synchronization import aggregate_averaging_synchronization
>>> aggregate_averaging_synchronization(kappa=0.1).interactions[0].kappa
0.1
Source code in src/skcomplex/models/physics/synchronization/_chaos.py
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'chaos', 'map', 'logistic map'),
    references=('Gros 2015',),
    defaults={'until': 50, 'dt': 1},
)
def aggregate_averaging_synchronization(
    x1: float = 0.3,
    x2: float = 0.6,
    *,
    r: float = 4.0,
    kappa: float = 0.4,
) -> System:
    """Return two logistic maps coupled by aggregate averaging, which synchronize their chaos for strong coupling.

    Each map iterates on the aggregate of its own value and the other's, `x_i(t + 1) = r a_i (1 - a_i)` with
    `a_1 = (1 - kappa) x_1 + kappa x_2` and `a_2 = (1 - kappa) x_2 + kappa x_1`. The synchronized orbit is stable for
    `kappa > (r - 1) / (2 r)`, which is `3 / 8` for `r = 4`, and the maps agree in one step for `kappa = 1 / 2`.

    Args:
        x1:
            The initial value of the first map.

        x2:
            The initial value of the second map.

        r:
            The growth rate of the logistic map.

        kappa:
            The coupling strength.

    Returns:
        system:
            The model, with the values `first.x` and `second.x`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import aggregate_averaging_synchronization
        >>> aggregate_averaging_synchronization(kappa=0.1).interactions[0].kappa
        0.1
    """
    return System(
        elements=[LogisticUnitElement('first', x=x1), LogisticUnitElement('second', x=x2)],
        interactions=[
            AggregateAveragingInteraction(f'averaging_{target}', source=source, target=target, r=r, kappa=kappa)
            for source, target in (('second', 'first'), ('first', 'second'))
        ],
    )

coupled_map_lattices_on_networks(n_nodes=100, *, network='erdos_renyi', degree=20, rewiring=0.1, r=4.0, epsilon=0.9, seed=0)

Return chaotic logistic maps on the nodes of a network, each pulled towards the mean of its neighbours.

Node i follows x_i(t + 1) = (1 - epsilon) f(x_i) + (epsilon / k_i) sum_j A_ij f(x_j) with the logistic map f(x) = r x (1 - x), a mean over the neighbours of the network. For a coupling epsilon large enough against the Lyapunov exponent of the map and the spectral gap of the network, the maps synchronize on one chaotic orbit, and the variance of the states over the nodes falls to zero. Below it they stay incoherent.

Parameters:

Name Type Description Default
n_nodes int

The number of maps, the nodes of the network.

100
network str

The generator, 'watts_strogatz', 'erdos_renyi', or 'barabasi_albert'.

'erdos_renyi'
degree int

The mean degree of the network.

20
rewiring float

The probability of rewiring every edge of the Watts-Strogatz ring.

0.1
r float

The parameter of the logistic map, chaotic at 4.

4.0
epsilon float

The coupling between 0 and 1.

0.9
seed int

The seed of the network and the initial states.

0

Returns:

Name Type Description
system System

The model, with the states maps.x and their variance over the nodes, spread.

Raises:

Type Description
ValueError

The network is unknown.

Examples:

>>> from skcomplex.models.physics.synchronization import coupled_map_lattices_on_networks
>>> coupled_map_lattices_on_networks(20).elements[0].check().prepare().size_
20
Source code in src/skcomplex/models/physics/synchronization/_networks.py
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'map', 'chaos', 'network'),
    references=('Arenas et al. 2008', 'Kaneko 1989'),
    defaults={'until': 100.0, 'dt': 1.0},
    observables=lambda system: [ObservableSpec('spread', group='maps', variable='x', reduce='var')],
)
def coupled_map_lattices_on_networks(
    n_nodes: int = 100,
    *,
    network: str = 'erdos_renyi',
    degree: int = 20,
    rewiring: float = 0.1,
    r: float = 4.0,
    epsilon: float = 0.9,
    seed: int = 0,
) -> System:
    """Return chaotic logistic maps on the nodes of a network, each pulled towards the mean of its neighbours.

    Node `i` follows `x_i(t + 1) = (1 - epsilon) f(x_i) + (epsilon / k_i) sum_j A_ij f(x_j)` with the logistic map
    `f(x) = r x (1 - x)`, a mean over the neighbours of the network. For a coupling `epsilon` large enough against the
    Lyapunov exponent of the map and the spectral gap of the network, the maps synchronize on one chaotic orbit, and
    the variance of the states over the nodes falls to zero. Below it they stay incoherent.

    Args:
        n_nodes:
            The number of maps, the nodes of the network.

        network:
            The generator, `'watts_strogatz'`, `'erdos_renyi'`, or `'barabasi_albert'`.

        degree:
            The mean degree of the network.

        rewiring:
            The probability of rewiring every edge of the Watts-Strogatz ring.

        r:
            The parameter of the logistic map, chaotic at 4.

        epsilon:
            The coupling between 0 and 1.

        seed:
            The seed of the network and the initial states.

    Returns:
        system:
            The model, with the states `maps.x` and their variance over the nodes, `spread`.

    Raises:
        ValueError:
            The network is unknown.

    Examples:
        >>> from skcomplex.models.physics.synchronization import coupled_map_lattices_on_networks
        >>> coupled_map_lattices_on_networks(20).elements[0].check().prepare().size_
        20
    """
    space = _network(network, n_nodes, degree, rewiring, seed)
    rng = np.random.default_rng(seed + 1)
    maps = LatticeElement(
        'maps',
        element=LogisticNodeElement,
        space=space,
        values={'x': rng.uniform(0.0, 1.0, n_nodes), 'r': r, 'epsilon': epsilon},
        visual=VisualSpec('sites', {'states': ['x'], 'low': 0.0, 'high': 1.0}),
    )
    coupling = LogisticCouplingInteraction('coupling', target='maps')
    return System(elements=[maps], interactions=[coupling], time=DiscreteTime())

coupled_rossler(n_oscillators=2, *, k=0.04, mismatch=0.02, a=0.2, b=0.2, c=5.7, seed=0)

Return chaotic Rössler oscillators on a ring, diffusively coupled through their second variable.

Oscillator i follows x1' = -ω_i x2 - x3, x2' = ω_i x1 + a x2 + k sum_j (x2_j - x2_i), and x3' = b + x3 (x1 - c), summed over its neighbours on the ring. Two oscillators with frequencies 1 ± mismatch pass from no synchronization to chaotic phase synchronization as k grows. Identical oscillators, with no mismatch, synchronize completely when every nonzero eigenvalue of the ring's Laplacian times k falls where the master stability function is negative.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators on the ring.

2
k float

The coupling strength.

0.04
mismatch float

The half width of the interval of natural frequencies around 1.

0.02
a float

The first Rössler parameter.

0.2
b float

The second Rössler parameter.

0.2
c float

The third Rössler parameter.

5.7
seed int

The seed of the initial states.

0

Returns:

Name Type Description
system System

The model, with the states oscillator_<i>.x1, oscillator_<i>.x2, and oscillator_<i>.x3.

Examples:

>>> from skcomplex.models.physics.synchronization import coupled_rossler
>>> [element.omega for element in coupled_rossler().elements]
[0.98, 1.02]
Source code in src/skcomplex/models/physics/synchronization/_chaos.py
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'chaos', 'ode', 'master stability function', 'networks'),
    references=('Rössler 1976', 'Pecora and Carroll 1998', 'Datseris and Parlitz 2022', 'Arenas et al. 2008'),
    defaults={'until': 100.0, 'dt': 0.1},
)
def coupled_rossler(
    n_oscillators: int = 2,
    *,
    k: float = 0.04,
    mismatch: float = 0.02,
    a: float = 0.2,
    b: float = 0.2,
    c: float = 5.7,
    seed: int = 0,
) -> System:
    """Return chaotic Rössler oscillators on a ring, diffusively coupled through their second variable.

    Oscillator `i` follows `x1' = -ω_i x2 - x3`, `x2' = ω_i x1 + a x2 + k sum_j (x2_j - x2_i)`, and
    `x3' = b + x3 (x1 - c)`, summed over its neighbours on the ring. Two oscillators with frequencies `1 ± mismatch`
    pass from no synchronization to chaotic phase synchronization as `k` grows. Identical oscillators, with no mismatch,
    synchronize completely when every nonzero eigenvalue of the ring's Laplacian times `k` falls where the master
    stability function is negative.

    Args:
        n_oscillators:
            The number of oscillators on the ring.

        k:
            The coupling strength.

        mismatch:
            The half width of the interval of natural frequencies around 1.

        a:
            The first Rössler parameter.

        b:
            The second Rössler parameter.

        c:
            The third Rössler parameter.

        seed:
            The seed of the initial states.

    Returns:
        system:
            The model, with the states `oscillator_<i>.x1`, `oscillator_<i>.x2`, and `oscillator_<i>.x3`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import coupled_rossler
        >>> [element.omega for element in coupled_rossler().elements]
        [0.98, 1.02]
    """
    omegas = np.linspace(1 - mismatch, 1 + mismatch, n_oscillators)
    starts = np.random.default_rng(seed).uniform(-1.0, 1.0, (n_oscillators, 2))
    oscillators = [
        RosslerOscillatorElement(
            f'oscillator_{index}',
            x1=float(starts[index, 0]),
            x2=float(starts[index, 1]),
            x3=0.0,
            omega=round(float(omegas[index]), 12),
            a=a,
            b=b,
            c=c,
        )
        for index in range(n_oscillators)
    ]
    links = [(index, index + 1) for index in range(n_oscillators - 1)]
    if n_oscillators > PAIR:
        links.append((n_oscillators - 1, 0))
    couplings = [
        DiffusiveCouplingInteraction(
            f'coupling_{source}_{target}',
            source=f'oscillator_{source}',
            target=f'oscillator_{target}',
            k=k,
        )
        for first, second in links
        for source, target in ((first, second), (second, first))
    ]
    return System(elements=oscillators, interactions=couplings)

delayed_kuramoto(n_oscillators=5, *, K=2.0, tau=1.0, spread=0.2, seed=0)

Return Kuramoto phase oscillators coupled all to all through signals that arrive after a delay.

Oscillator i follows theta_i' = omega_i + (K / N) sum_j sin(theta_j(t - tau) - theta_i), with natural frequencies spread evenly around 1 and initial phases spread around the circle.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators.

5
K float

The coupling strength.

2.0
tau float

The delay of the coupling.

1.0
spread float

The width of the interval of natural frequencies around 1.

0.2
seed int

The seed of the initial phases.

0

Returns:

Name Type Description
system System

The model, with the phases oscillator_<i>.theta.

Examples:

>>> from skcomplex.models.physics.synchronization import delayed_kuramoto
>>> len(delayed_kuramoto(n_oscillators=3).interactions)
6
Source code in src/skcomplex/models/physics/synchronization/_kuramoto.py
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'delay', 'oscillators'),
    references=('Yeung and Strogatz 1999', 'Kuramoto 1984'),
    defaults={'until': 50.0, 'dt': 0.1},
)
def delayed_kuramoto(
    n_oscillators: int = 5,
    *,
    K: float = 2.0,
    tau: float = 1.0,
    spread: float = 0.2,
    seed: int = 0,
) -> System:
    """Return Kuramoto phase oscillators coupled all to all through signals that arrive after a delay.

    Oscillator `i` follows `theta_i' = omega_i + (K / N) sum_j sin(theta_j(t - tau) - theta_i)`, with natural
    frequencies spread evenly around 1 and initial phases spread around the circle.

    Args:
        n_oscillators:
            The number of oscillators.

        K:
            The coupling strength.

        tau:
            The delay of the coupling.

        spread:
            The width of the interval of natural frequencies around 1.

        seed:
            The seed of the initial phases.

    Returns:
        system:
            The model, with the phases `oscillator_<i>.theta`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import delayed_kuramoto
        >>> len(delayed_kuramoto(n_oscillators=3).interactions)
        6
    """
    phases = [2 * math.pi * ((seed * 0.618 + index / n_oscillators) % 1.0) for index in range(n_oscillators)]
    omegas = [1.0 + spread * (index / max(n_oscillators - 1, 1) - 0.5) for index in range(n_oscillators)]
    elements = [
        PhaseOscillatorElement(f'oscillator_{index}', theta=phases[index], omega=omegas[index])
        for index in range(n_oscillators)
    ]
    interactions = [
        DelayedCouplingInteraction(
            f'coupling_{source}_{target}',
            source=f'oscillator_{source}',
            target=f'oscillator_{target}',
            K=K / n_oscillators,
            tau=tau,
        )
        for source in range(n_oscillators)
        for target in range(n_oscillators)
        if source != target
    ]
    return System(elements=elements, interactions=interactions)

heterogeneous_global_coupling(n_systems=100, *, k=2.0, beta=0.1, gamma=2.0, heterogeneity=0.8, random_state=0)

Return globally coupled oscillators with heterogeneous time scales, whose fixed point the coupling restabilizes.

System j follows x_j' = τ_j g(x_j, y_j) + k (X - x_j) and y_j' = τ_j h(x_j, y_j) + k (Y - y_j), where X and Y are the population means and g + ih = (beta + i gamma) z - |z|² z is the Hopf normal form of z = x + iy. Its fixed point 0 is an unstable focus with eigenvalues beta ± i gamma. The time scales τ_j are spread evenly over [1 - heterogeneity, 1 + heterogeneity], with mean 1 and variance sigma². The order parameters X, Y, V = <δτ δx>, and W = <δτ δy> close on four equations. The global fixed point is stable when their eigenvalues nu - k / 2 ± √(k² / 4 + sigma² nu²), for nu = beta ± i gamma, have negative real parts, which holds roughly for 2 beta (1 + sigma) < k < (gamma² - beta²) sigma² / beta.

Parameters:

Name Type Description Default
n_systems int

The number of coupled systems.

100
k float

The coupling constant.

2.0
beta float

The growth rate β of the unstable focus.

0.1
gamma float

The rotation rate gamma of the unstable focus.

2.0
heterogeneity float

The half width of the interval of time scales, below 1.

0.8
random_state int | Generator | None

The seed or generator of the initial phases.

0

Returns:

Name Type Description
system System

The model, with the states systems.x and systems.y and the observables amplitude, the mean amplitude, and the order parameters X, Y, V, and W.

Examples:

>>> from skcomplex.models.physics.synchronization import heterogeneous_global_coupling
>>> heterogeneous_global_coupling(4).elements[0].values['tau']
array([0.4, 0.8, 1.2, 1.6])
Source code in src/skcomplex/models/physics/synchronization/_amplitude_death.py
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'amplitude death', 'mean field', 'heterogeneity', 'ode'),
    references=("De Monte, d'Ovidio, and Mosekilde 2003", 'Bertin 2021'),
    defaults={'until': 40.0, 'dt': 0.2},
    observables=lambda system: [
        ObservableSpec('amplitude', _amplitude, group='systems'),
        ObservableSpec('X', group='systems', variable='x'),
        ObservableSpec('Y', group='systems', variable='y'),
        ObservableSpec('V', _time_scale_deviation_x, group='systems'),
        ObservableSpec('W', _time_scale_deviation_y, group='systems'),
    ],
)
def heterogeneous_global_coupling(
    n_systems: int = 100,
    *,
    k: float = 2.0,
    beta: float = 0.1,
    gamma: float = 2.0,
    heterogeneity: float = 0.8,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return globally coupled oscillators with heterogeneous time scales, whose fixed point the coupling restabilizes.

    System `j` follows `x_j' = τ_j g(x_j, y_j) + k (X - x_j)` and `y_j' = τ_j h(x_j, y_j) + k (Y - y_j)`, where `X`
    and `Y` are the population means and `g + ih = (beta + i gamma) z - |z|² z` is the Hopf normal form of
    `z = x + iy`. Its fixed point 0 is an unstable focus with eigenvalues `beta ± i gamma`. The time scales `τ_j` are
    spread evenly over `[1 - heterogeneity, 1 + heterogeneity]`, with mean 1 and variance `sigma²`. The order
    parameters `X`, `Y`, `V = <δτ δx>`, and `W = <δτ δy>` close on four equations. The global fixed point is stable
    when their eigenvalues `nu - k / 2 ± √(k² / 4 + sigma² nu²)`, for `nu = beta ± i gamma`, have negative real parts,
    which holds roughly for `2 beta (1 + sigma) < k < (gamma² - beta²) sigma² / beta`.

    Args:
        n_systems:
            The number of coupled systems.

        k:
            The coupling constant.

        beta:
            The growth rate `β` of the unstable focus.

        gamma:
            The rotation rate `gamma` of the unstable focus.

        heterogeneity:
            The half width of the interval of time scales, below 1.

        random_state:
            The seed or generator of the initial phases.

    Returns:
        system:
            The model, with the states `systems.x` and `systems.y` and the observables `amplitude`, the mean amplitude,
            and the order parameters `X`, `Y`, `V`, and `W`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import heterogeneous_global_coupling
        >>> heterogeneous_global_coupling(4).elements[0].values['tau']
        array([0.4, 0.8, 1.2, 1.6])
    """
    taus = 1 + heterogeneity * (2 * (np.arange(n_systems) + 1 / 2) / n_systems - 1)
    x, y = _spread_on_circle(n_systems, random_state)
    radius = math.sqrt(abs(beta))
    systems = GroupElement(
        'systems',
        element=StuartLandauElement,
        size=n_systems,
        values={'x': radius * x, 'y': radius * y, 'mu': beta, 'omega': gamma, 'tau': taus},
    )
    return System(
        elements=[systems],
        interactions=[
            GlobalCouplingInteraction(f'coupling_{state}', target='systems', state=state, K=k) for state in 'xy'
        ],
    )

hypergraph_kuramoto(n_oscillators=60, *, K1=0.0, K2=8.0, width=0.5, links=1.0, triangles=1.0, phases=None, seed=0)

Return Kuramoto oscillators coupled in pairs and by three-body interactions on the triangles of a hypergraph.

Oscillator i follows theta_i' = omega_i + (K1 / <k1>) sum_j sin(theta_j - theta_i) + (K2 / <k2>) sum_jk sin(theta_j + theta_k - 2 theta_i) over the links ij and the triangles ijk of a random hypergraph, where <k1> and <k2> are the mean numbers of links and of triangles of a node. All-to-all, the triangle term is the three-body coupling of Skardal and Arenas, and the Ott-Antonsen reduction of Lorentzian frequencies of half width width gives r' = -width r + K1 (r - r^3) / 2 + K2 (r^3 - r^5) / 2. Without links, the incoherent state stays stable for every K2, and for K2 > 8 width a synchronized state of r^2 = (1 + sqrt(1 - 8 width / K2)) / 2 coexists with it, so the order parameter depends on the initial phases and sweeps of the coupling show hysteresis.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators, the nodes of the hypergraph.

60
K1 float

The pairwise coupling, divided by the mean number of links of a node.

0.0
K2 float

The three-body coupling, divided by the mean number of triangles of a node.

8.0
width float

The half width of the Lorentzian distribution of the natural frequencies.

0.5
links float

The probability that a pair of oscillators is a link.

1.0
triangles float

The probability that a triple of oscillators is a triangle.

1.0
phases Sequence[float] | None

The initial phases, or None for phases drawn uniformly from the circle.

None
seed int

The seed of the hypergraph, the frequencies, and the phases.

0

Returns:

Name Type Description
system System

The model, with the phases oscillators.theta and the order parameter order.

Examples:

>>> from skcomplex.models.physics.synchronization import hypergraph_kuramoto
>>> system = hypergraph_kuramoto(6)
>>> len(system.elements[0].space.check().prepare().hyperedges_)
35
Source code in src/skcomplex/models/physics/synchronization/_hypergraphs.py
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'phase', 'hypergraph'),
    references=('Skardal and Arenas 2019', 'Majhi et al. 2022'),
    defaults={'until': 40.0, 'dt': 0.1},
    observables=lambda system: _order_parameter('order', 'oscillators'),
)
def hypergraph_kuramoto(
    n_oscillators: int = 60,
    *,
    K1: float = 0.0,
    K2: float = 8.0,
    width: float = 0.5,
    links: float = 1.0,
    triangles: float = 1.0,
    phases: Sequence[float] | None = None,
    seed: int = 0,
) -> System:
    """Return Kuramoto oscillators coupled in pairs and by three-body interactions on the triangles of a hypergraph.

    Oscillator `i` follows `theta_i' = omega_i + (K1 / <k1>) sum_j sin(theta_j - theta_i) + (K2 / <k2>) sum_jk
    sin(theta_j + theta_k - 2 theta_i)` over the links `ij` and the triangles `ijk` of a random hypergraph, where `<k1>`
    and `<k2>` are the mean numbers of links and of triangles of a node. All-to-all, the triangle term is the
    three-body coupling of Skardal and Arenas, and the Ott-Antonsen reduction of Lorentzian frequencies of half width
    `width` gives `r' = -width r + K1 (r - r^3) / 2 + K2 (r^3 - r^5) / 2`. Without links, the incoherent state stays
    stable for every `K2`, and for `K2 > 8 width` a synchronized state of `r^2 = (1 + sqrt(1 - 8 width / K2)) / 2`
    coexists with it, so the order parameter depends on the initial phases and sweeps of the coupling show hysteresis.

    Args:
        n_oscillators:
            The number of oscillators, the nodes of the hypergraph.

        K1:
            The pairwise coupling, divided by the mean number of links of a node.

        K2:
            The three-body coupling, divided by the mean number of triangles of a node.

        width:
            The half width of the Lorentzian distribution of the natural frequencies.

        links:
            The probability that a pair of oscillators is a link.

        triangles:
            The probability that a triple of oscillators is a triangle.

        phases:
            The initial phases, or `None` for phases drawn uniformly from the circle.

        seed:
            The seed of the hypergraph, the frequencies, and the phases.

    Returns:
        system:
            The model, with the phases `oscillators.theta` and the order parameter `order`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import hypergraph_kuramoto
        >>> system = hypergraph_kuramoto(6)
        >>> len(system.elements[0].space.check().prepare().hyperedges_)
        35
    """
    hypergraph = _simplices(n_oscillators, links, triangles, seed)
    space = NetworkSpace(nx.complete_graph(n_oscillators), hyperedges=hypergraph)
    sizes = [len(members) for members in hypergraph.edges.members(dtype=list)]
    mean_links = max(PAIR * sizes.count(PAIR) / n_oscillators, 1.0)
    mean_triangles = max(3 * sizes.count(3) / n_oscillators, 1.0)
    rng = np.random.default_rng(seed + 1)
    omega = width * stats.cauchy.ppf(rng.uniform(size=n_oscillators))
    theta = rng.uniform(0, 2 * math.pi, n_oscillators) if phases is None else np.asarray(phases, dtype=float)
    oscillators = LatticeElement(
        'oscillators',
        element=PhaseOscillatorElement,
        space=space,
        values={'theta': theta, 'omega': omega},
        visual=VisualSpec('sites', {'states': ['theta'], 'low': 0.0, 'high': 2 * math.pi}),
    )
    coupling = SimplexCouplingInteraction(
        'coupling',
        target='oscillators',
        K1=K1 / mean_links,
        K2=K2 / mean_triangles,
    )
    return System(elements=[oscillators], interactions=[coupling])

kuramoto(n_oscillators=50, *, K=4.0, sigma=1.0, seed=0)

Return Kuramoto oscillators coupled all to all, theta_i' = omega_i + (K / N) sum_j sin(theta_j - theta_i).

Every oscillator reads the group means of cos theta and sin theta, so the coupling costs time linear in the number of oscillators. The natural frequencies are drawn from a normal distribution of width sigma, and the phases uniformly from the circle. Above the critical coupling K_c = 2 sigma sqrt(2 / pi) a synchronized cluster forms, and the order parameter r, the modulus of the mean of exp(i theta), grows from about 1 / sqrt(N) towards 1.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators.

50
K float

The coupling strength.

4.0
sigma float

The standard deviation of the natural frequencies.

1.0
seed int

The seed of the natural frequencies and the initial phases.

0

Returns:

Name Type Description
system System

The model, with the phases oscillators.theta and the order parameter order.

Examples:

>>> from skcomplex.models.physics.synchronization import kuramoto
>>> system = kuramoto(10)
>>> system.elements[0].size, system.space
(10, None)
Source code in src/skcomplex/models/physics/synchronization/_kuramoto.py
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'phase', 'mean field'),
    references=('Kuramoto 1984', 'Datseris and Parlitz 2022', 'Gros 2015', 'Jensen 2023'),
    defaults={'until': 50.0, 'dt': 0.5},
    observables=lambda system: _order_parameter('order', 'oscillators'),
)
def kuramoto(n_oscillators: int = 50, *, K: float = 4.0, sigma: float = 1.0, seed: int = 0) -> System:
    """Return Kuramoto oscillators coupled all to all, `theta_i' = omega_i + (K / N) sum_j sin(theta_j - theta_i)`.

    Every oscillator reads the group means of `cos theta` and `sin theta`, so the coupling costs time linear in the
    number of oscillators. The natural frequencies are drawn from a normal distribution of width `sigma`, and the
    phases uniformly from the circle. Above the critical coupling `K_c = 2 sigma sqrt(2 / pi)` a synchronized cluster
    forms, and the order parameter `r`, the modulus of the mean of `exp(i theta)`, grows from about `1 / sqrt(N)`
    towards 1.

    Args:
        n_oscillators:
            The number of oscillators.

        K:
            The coupling strength.

        sigma:
            The standard deviation of the natural frequencies.

        seed:
            The seed of the natural frequencies and the initial phases.

    Returns:
        system:
            The model, with the phases `oscillators.theta` and the order parameter `order`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import kuramoto
        >>> system = kuramoto(10)
        >>> system.elements[0].size, system.space
        (10, None)
    """
    rng = np.random.default_rng(seed)
    oscillators = GroupElement(
        'oscillators',
        element=PhaseOscillatorElement,
        size=n_oscillators,
        values={'theta': _spread_phases(n_oscillators, seed), 'omega': rng.normal(0.0, sigma, n_oscillators)},
    )
    coupling = KuramotoCouplingInteraction(
        'coupling',
        target='oscillators',
        K=K,
        mode='group',
        reduce='mean',
        include_self=True,
    )
    return System(
        elements=[oscillators],
        interactions=[coupling],
    )

kuramoto_on_networks(n_oscillators=200, *, network='watts_strogatz', degree=6, rewiring=0.1, K=1.0, sigma=0.5, seed=0)

Return Kuramoto oscillators coupled along the edges of a network, theta_i' = omega_i + K sum_j A_ij sin(...).

Oscillator i follows theta_i' = omega_i + (K / <k>) sum_j w_ij sin(theta_j - theta_i) over its neighbours j, where the weights w_ij = 1 live on a group on the edges and <k> is the mean degree. The coupling runs once per edge. The natural frequencies are drawn from a normal distribution of width sigma, and the phases uniformly from the circle. Above a critical coupling of the order of the all-to-all value 2 sigma sqrt(2 / pi), the order parameter r = |<exp(i theta)>| grows from about 1 / sqrt(N) towards 1, and small-world shortcuts lower it.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators, the nodes of the network.

200
network str

The generator, 'watts_strogatz', 'erdos_renyi', or 'barabasi_albert'.

'watts_strogatz'
degree int

The mean degree of the network, even for the Watts-Strogatz ring.

6
rewiring float

The probability of rewiring every edge of the Watts-Strogatz ring.

0.1
K float

The coupling strength, divided by the mean degree.

1.0
sigma float

The standard deviation of the natural frequencies.

0.5
seed int

The seed of the network, the frequencies, and the phases.

0

Returns:

Name Type Description
system System

The model, with the phases oscillators.theta, the weights links.w, and the order parameter order.

Raises:

Type Description
ValueError

The network is unknown.

Examples:

>>> from skcomplex.models.physics.synchronization import kuramoto_on_networks
>>> system = kuramoto_on_networks(20, degree=4)
>>> system.elements[1].check().prepare().size_
40
Source code in src/skcomplex/models/physics/synchronization/_networks.py
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'phase', 'network'),
    references=('Arenas et al. 2008', 'Jensen 2023'),
    defaults={'until': 40.0, 'dt': 0.2},
    observables=lambda system: _order_parameter('order', 'oscillators'),
)
def kuramoto_on_networks(
    n_oscillators: int = 200,
    *,
    network: str = 'watts_strogatz',
    degree: int = 6,
    rewiring: float = 0.1,
    K: float = 1.0,
    sigma: float = 0.5,
    seed: int = 0,
) -> System:
    """Return Kuramoto oscillators coupled along the edges of a network, `theta_i' = omega_i + K sum_j A_ij sin(...)`.

    Oscillator `i` follows `theta_i' = omega_i + (K / <k>) sum_j w_ij sin(theta_j - theta_i)` over its neighbours `j`,
    where the weights `w_ij = 1` live on a group on the edges and `<k>` is the mean degree. The coupling runs once per
    edge. The natural frequencies are drawn from a normal distribution of width `sigma`, and the phases uniformly from
    the circle. Above a critical coupling of the order of the all-to-all value `2 sigma sqrt(2 / pi)`, the order
    parameter `r = |<exp(i theta)>|` grows from about `1 / sqrt(N)` towards 1, and small-world shortcuts lower it.

    Args:
        n_oscillators:
            The number of oscillators, the nodes of the network.

        network:
            The generator, `'watts_strogatz'`, `'erdos_renyi'`, or `'barabasi_albert'`.

        degree:
            The mean degree of the network, even for the Watts-Strogatz ring.

        rewiring:
            The probability of rewiring every edge of the Watts-Strogatz ring.

        K:
            The coupling strength, divided by the mean degree.

        sigma:
            The standard deviation of the natural frequencies.

        seed:
            The seed of the network, the frequencies, and the phases.

    Returns:
        system:
            The model, with the phases `oscillators.theta`, the weights `links.w`, and the order parameter `order`.

    Raises:
        ValueError:
            The network is unknown.

    Examples:
        >>> from skcomplex.models.physics.synchronization import kuramoto_on_networks
        >>> system = kuramoto_on_networks(20, degree=4)
        >>> system.elements[1].check().prepare().size_
        40
    """
    space = _network(network, n_oscillators, degree, rewiring, seed)
    rng = np.random.default_rng(seed + 1)
    mean_degree = 2 * space.graph.number_of_edges() / n_oscillators
    oscillators = LatticeElement(
        'oscillators',
        element=PhaseOscillatorElement,
        space=space,
        values={'theta': rng.uniform(0, 2 * math.pi, n_oscillators), 'omega': rng.normal(0.0, sigma, n_oscillators)},
        visual=VisualSpec('sites', {'states': ['theta'], 'low': 0.0, 'high': 2 * math.pi}),
    )
    links = GroupElement('links', element=LinkElement, space=space, on='edges')
    coupling = EdgeCouplingInteraction('coupling', target='oscillators', links='links', K=K / mean_degree)
    return System(elements=[oscillators, links], interactions=[coupling])

kuramoto_sakaguchi_chimera(n_oscillators=32, *, A=0.2, beta=0.1, spread=0.1, seed=0)

Return two populations of identical phase oscillators with phase lag, where one synchronizes and one drifts.

Oscillator i of population p follows theta_i' = sum_q (K_pq / N) sum_j sin(theta_j - theta_i - alpha), with the stronger coupling mu = (1 + A) / 2 inside a population, the weaker nu = (1 - A) / 2 between them, and the lag alpha = pi / 2 - beta. Every population reads the means of cos theta and sin theta of both populations. The first population starts nearly synchronized and the second spread over the circle. The populations settle in a chimera, where one is fully synchronized and the other keeps drifting with a partial order parameter, which breathes for A = 0.2. Which population ends up coherent depends on the initial phases.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators of every population.

32
A float

The disparity of the couplings inside and between the populations.

0.2
beta float

The distance of the phase lag from pi / 2.

0.1
spread float

The width of the initial phases of the first population.

0.1
seed int

The seed of the initial phases.

0

Returns:

Name Type Description
system System

The model, with the phases first.theta and second.theta and the order parameters first_order and second_order.

Examples:

>>> from skcomplex.models.physics.synchronization import kuramoto_sakaguchi_chimera
>>> len(kuramoto_sakaguchi_chimera(4).interactions)
4
Source code in src/skcomplex/models/physics/synchronization/_kuramoto.py
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'phase', 'chimera'),
    references=('Abrams et al. 2008', 'Jensen 2023'),
    defaults={'until': 100.0, 'dt': 1.0},
    observables=lambda system: [*_order_parameter('first_order', 'first'), *_order_parameter('second_order', 'second')],
)
def kuramoto_sakaguchi_chimera(
    n_oscillators: int = 32,
    *,
    A: float = 0.2,
    beta: float = 0.1,
    spread: float = 0.1,
    seed: int = 0,
) -> System:
    """Return two populations of identical phase oscillators with phase lag, where one synchronizes and one drifts.

    Oscillator `i` of population `p` follows `theta_i' = sum_q (K_pq / N) sum_j sin(theta_j - theta_i - alpha)`, with
    the stronger coupling `mu = (1 + A) / 2` inside a population, the weaker `nu = (1 - A) / 2` between them, and the
    lag `alpha = pi / 2 - beta`. Every population reads the means of `cos theta` and `sin theta` of both populations.
    The first population starts nearly synchronized and the second spread over the circle. The populations settle in a
    chimera, where one is fully synchronized and the other keeps drifting with a partial order parameter, which
    breathes for `A = 0.2`. Which population ends up coherent depends on the initial phases.

    Args:
        n_oscillators:
            The number of oscillators of every population.

        A:
            The disparity of the couplings inside and between the populations.

        beta:
            The distance of the phase lag from `pi / 2`.

        spread:
            The width of the initial phases of the first population.

        seed:
            The seed of the initial phases.

    Returns:
        system:
            The model, with the phases `first.theta` and `second.theta` and the order parameters `first_order` and
            `second_order`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import kuramoto_sakaguchi_chimera
        >>> len(kuramoto_sakaguchi_chimera(4).interactions)
        4
    """
    rng = np.random.default_rng(seed)
    alpha, mu, nu = math.pi / 2 - beta, (1 + A) / 2, (1 - A) / 2
    populations = [
        GroupElement(
            'first',
            element=PhaseOscillatorElement,
            size=n_oscillators,
            values={'theta': rng.normal(0.0, spread, n_oscillators), 'omega': 0.0},
            visual=VisualSpec('phase_ring', {'phase': 'theta'}),
        ),
        GroupElement(
            'second',
            element=PhaseOscillatorElement,
            size=n_oscillators,
            values={'theta': _spread_phases(n_oscillators, seed + 1), 'omega': 0.0},
            visual=VisualSpec('phase_ring', {'phase': 'theta'}),
        ),
    ]
    interactions: list[BaseInteraction] = [
        KuramotoCouplingInteraction(
            f'{target}_from_{source}',
            target=target,
            source=source,
            K=mu if target == source else nu,
            alpha=alpha,
            reduce='mean',
            mode='group',
            include_self=True,
        )
        for target in ('first', 'second')
        for source in ('first', 'second')
    ]
    return System(
        elements=populations,
        interactions=interactions,
    )

mirollo_strogatz(n_oscillators=5, *, S=2.0, gamma=1.0, epsilon=0.1, seed=0)

Return integrate-and-fire oscillators coupled by pulses, which end up firing in unison.

Every oscillator of the group charges as x' = S - gamma x, a concave rise, and fires when it reaches 1. Firing is an event whose guard is the group maximum of the charges, and it resets the most charged oscillators to 0 and lifts every other one by epsilon. An oscillator lifted to 1 fires at the same instant and is absorbed into the group that fired, so the groups merge until all oscillators fire together.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators.

5
S float

The charging input.

2.0
gamma float

The leak rate, below S so that the oscillators fire.

1.0
epsilon float

The strength of a pulse.

0.1
seed int

The seed of the initial charges.

0

Returns:

Name Type Description
system System

The model, with the charges oscillators.x of the group members and the firings as fire events.

Examples:

>>> from skcomplex.models.physics.synchronization import mirollo_strogatz
>>> mirollo_strogatz(3).elements[0].check().prepare().size_
3
Source code in src/skcomplex/models/physics/synchronization/_pulses.py
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'events', 'pulse coupling', 'integrate and fire'),
    references=('Mirollo and Strogatz 1990', 'Arenas et al. 2008'),
    defaults={'until': 20.0, 'dt': 0.1},
)
def mirollo_strogatz(
    n_oscillators: int = 5,
    *,
    S: float = 2.0,
    gamma: float = 1.0,
    epsilon: float = 0.1,
    seed: int = 0,
) -> System:
    """Return integrate-and-fire oscillators coupled by pulses, which end up firing in unison.

    Every oscillator of the group charges as `x' = S - gamma x`, a concave rise, and fires when it reaches 1. Firing is
    an event whose guard is the group maximum of the charges, and it resets the most charged oscillators to 0 and lifts
    every other one by `epsilon`. An oscillator lifted to 1 fires at the same instant and is absorbed into the group
    that fired, so the groups merge until all oscillators fire together.

    Args:
        n_oscillators:
            The number of oscillators.

        S:
            The charging input.

        gamma:
            The leak rate, below `S` so that the oscillators fire.

        epsilon:
            The strength of a pulse.

        seed:
            The seed of the initial charges.

    Returns:
        system:
            The model, with the charges `oscillators.x` of the group members and the firings as `fire` events.

    Examples:
        >>> from skcomplex.models.physics.synchronization import mirollo_strogatz
        >>> mirollo_strogatz(3).elements[0].check().prepare().size_
        3
    """
    charges = np.random.default_rng(seed).uniform(0.0, 1.0, n_oscillators)
    oscillators = GroupElement(
        'oscillators',
        element=IntegrateAndFireOscillatorElement,
        size=n_oscillators,
        values={'x': charges, 'S': S, 'gamma': gamma, 'epsilon': epsilon},
    )
    fire = CrossingInteraction(
        'fire',
        target='oscillators',
        guard=lambda s: maximum(s.target.x) - 1,
        reset={
            'x': lambda s: where(
                s.target.x >= maximum(s.target.x),
                0,
                where(s.target.x + s.target.epsilon >= 1, 0, s.target.x + s.target.epsilon),
            ),
        },
        direction=1,
    )
    return System(elements=[oscillators], interactions=[fire])

nonlocal_ring_chimera(n_oscillators=100, *, R=14, alpha=1.46, omega=0.0, seed=0)

Return identical phase oscillators on a ring, each coupled to its R nearest neighbours on both sides.

Oscillator i follows phi_i' = omega - (1 / 2R) sum_{j = i - R}^{i + R} sin(phi_i - phi_j + alpha), where the sum includes the oscillator itself. Near alpha = pi / 2, coherent and incoherent domains coexist on the ring as a chimera state, which is a long transient on small rings.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators on the ring.

100
R int

The range of the coupling.

14
alpha float

The phase lag.

1.46
omega float

The natural frequency.

0.0
seed int

The seed of the initial phases.

0

Returns:

Name Type Description
system System

The model, with the phases ring.theta and the places of the oscillators ring.x.

Examples:

>>> from skcomplex.models.physics.synchronization import nonlocal_ring_chimera
>>> nonlocal_ring_chimera(20, R=3).interactions[0].neighborhood.radius
3.5
Source code in src/skcomplex/models/physics/synchronization/_kuramoto.py
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'phase', 'chimera', 'nonlocal coupling'),
    references=('Kuramoto and Battogtokh 2002', 'Datseris and Parlitz 2022'),
    defaults={'until': 100.0, 'dt': 1.0},
    observables=lambda system: _order_parameter('order', 'ring'),
)
def nonlocal_ring_chimera(
    n_oscillators: int = 100,
    *,
    R: int = 14,
    alpha: float = 1.46,
    omega: float = 0.0,
    seed: int = 0,
) -> System:
    """Return identical phase oscillators on a ring, each coupled to its `R` nearest neighbours on both sides.

    Oscillator `i` follows `phi_i' = omega - (1 / 2R) sum_{j = i - R}^{i + R} sin(phi_i - phi_j + alpha)`, where the sum
    includes the oscillator itself. Near `alpha = pi / 2`, coherent and incoherent domains coexist on the ring as a
    chimera state, which is a long transient on small rings.

    Args:
        n_oscillators:
            The number of oscillators on the ring.

        R:
            The range of the coupling.

        alpha:
            The phase lag.

        omega:
            The natural frequency.

        seed:
            The seed of the initial phases.

    Returns:
        system:
            The model, with the phases `ring.theta` and the places of the oscillators `ring.x`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import nonlocal_ring_chimera
        >>> nonlocal_ring_chimera(20, R=3).interactions[0].neighborhood.radius
        3.5
    """
    weight = 1 / (2 * R)
    ring = GroupElement(
        'ring',
        element=SitedOscillatorElement,
        size=n_oscillators,
        values={
            'theta': _spread_phases(n_oscillators, seed),
            'omega': omega,
            'x': np.arange(n_oscillators, dtype=float).reshape(-1, 1),
        },
        visual=VisualSpec('phase_ring', {'phase': 'theta', 'position': 'x'}),
    )
    return System(
        elements=[ring],
        interactions=[
            KuramotoCouplingInteraction(
                'coupling',
                target='ring',
                K=weight,
                alpha=alpha,
                include_self=True,
                neighborhood=RadiusNeighborhood(R + 0.5),
            ),
        ],
        space=ContinuousSpace((float(n_oscillators),), boundary='periodic'),
    )

oscillator_death(n_oscillators=100, *, K=3.0, width=1.5, distribution='lorentzian', random_state=0)

Return globally coupled Stuart-Landau oscillators, whose oscillations all die for strong coupling and spread.

Oscillator j follows z_j' = (1 + iω_j) z_j - |z_j|² z_j + (K / N) Σ_k (z_k - z_j) for z_j = x_j + i y_j. Alone it settles on the unit circle. The frequencies are the midpoint quantiles of a density of width ω0, the Lorentzian ω0 / (π (ω² + ω0²)) or the quartic g0 / (ω⁴ + ω0⁴). For many oscillators, the common fixed point z = 0 is stable for K > 1 and ω0 above a critical width, 1 for the Lorentzian and (1 + √(2K - 1)) / √2 for the quartic density, and every amplitude then decays to zero. For the Lorentzian, the slowest perturbation decays at the rate ω0 - 1.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators.

100
K float

The coupling strength.

3.0
width float

The width ω0 of the density of frequencies.

1.5
distribution str

The density of frequencies, 'lorentzian' or 'quartic'.

'lorentzian'
random_state int | Generator | None

The seed or generator of the initial phases.

0

Returns:

Name Type Description
system System

The model, with the states oscillators.x and oscillators.y and the observables amplitude, the mean amplitude, X and Y, the mean field, and order, its modulus.

Raises:

Type Description
ValidationError

The distribution is neither 'lorentzian' nor 'quartic'.

Examples:

>>> from skcomplex.models.physics.synchronization import oscillator_death
>>> oscillator_death(4, distribution='quartic').elements[0].check().prepare().size_
4
Source code in src/skcomplex/models/physics/synchronization/_amplitude_death.py
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'amplitude death', 'mean field', 'ode'),
    references=('Ermentrout 1990', 'Mirollo and Strogatz 1990', 'Bertin 2021'),
    defaults={'until': 20.0, 'dt': 0.2},
    observables=lambda system: [
        ObservableSpec('amplitude', _amplitude, group='oscillators'),
        ObservableSpec('X', group='oscillators', variable='x'),
        ObservableSpec('Y', group='oscillators', variable='y'),
        ObservableSpec('order', _order),
    ],
)
def oscillator_death(
    n_oscillators: int = 100,
    *,
    K: float = 3.0,
    width: float = 1.5,
    distribution: str = 'lorentzian',
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return globally coupled Stuart-Landau oscillators, whose oscillations all die for strong coupling and spread.

    Oscillator `j` follows `z_j' = (1 + iω_j) z_j - |z_j|² z_j + (K / N) Σ_k (z_k - z_j)` for `z_j = x_j + i y_j`.
    Alone it settles on the unit circle. The frequencies are the midpoint quantiles of a density of width `ω0`, the
    Lorentzian `ω0 / (π (ω² + ω0²))` or the quartic `g0 / (ω⁴ + ω0⁴)`. For many oscillators, the common fixed point
    `z = 0` is stable for `K > 1` and `ω0` above a critical width, 1 for the Lorentzian and
    `(1 + √(2K - 1)) / √2` for the quartic density, and every amplitude then decays to zero. For the Lorentzian, the
    slowest perturbation decays at the rate `ω0 - 1`.

    Args:
        n_oscillators:
            The number of oscillators.

        K:
            The coupling strength.

        width:
            The width `ω0` of the density of frequencies.

        distribution:
            The density of frequencies, `'lorentzian'` or `'quartic'`.

        random_state:
            The seed or generator of the initial phases.

    Returns:
        system:
            The model, with the states `oscillators.x` and `oscillators.y` and the observables `amplitude`, the mean
            amplitude, `X` and `Y`, the mean field, and `order`, its modulus.

    Raises:
        ValidationError:
            The distribution is neither `'lorentzian'` nor `'quartic'`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import oscillator_death
        >>> oscillator_death(4, distribution='quartic').elements[0].check().prepare().size_
        4
    """
    omegas = _frequency_quantiles(n_oscillators, width, distribution)
    x, y = _spread_on_circle(n_oscillators, random_state)
    oscillators = GroupElement(
        'oscillators',
        element=StuartLandauElement,
        size=n_oscillators,
        values={'x': x, 'y': y, 'mu': 1.0, 'omega': omegas, 'tau': 1.0},
    )
    return System(
        elements=[oscillators],
        interactions=[
            GlobalCouplingInteraction(f'coupling_{state}', target='oscillators', state=state, K=K) for state in 'xy'
        ],
    )

ott_antonsen(r=0.1, *, K=2.0, Delta=0.5, omega=0.0)

Return the Ott-Antonsen equation for the order parameter of infinitely many Kuramoto oscillators.

For natural frequencies with a Lorentzian distribution of center omega and half width Delta, the complex order parameter z = x + i y follows exactly z' = (i omega - Delta) z + (K / 2)(1 - |z|²) z, a Stuart-Landau oscillator with the growth rate K / 2 - Delta and the saturation K / 2. Its modulus settles at sqrt(1 - 2 Delta / K) above the critical coupling K_c = 2 Delta, and at 0 below it.

Parameters:

Name Type Description Default
r float

The initial modulus of the order parameter.

0.1
K float

The coupling strength.

2.0
Delta float

The half width of the distribution of natural frequencies.

0.5
omega float

The center of the distribution of natural frequencies.

0.0

Returns:

Name Type Description
system System

The model, with the real part field.x and the imaginary part field.y of the order parameter.

Examples:

>>> from skcomplex.models.physics.synchronization import ott_antonsen
>>> ott_antonsen(K=3.0).elements[0].kappa
1.5
Source code in src/skcomplex/models/physics/synchronization/_kuramoto.py
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'ode', 'mean field', 'dimension reduction'),
    references=('Ott and Antonsen 2008', 'Bick et al. 2020'),
    defaults={'until': 30.0, 'dt': 0.1},
)
def ott_antonsen(
    r: float = 0.1,
    *,
    K: float = 2.0,
    Delta: float = 0.5,
    omega: float = 0.0,
) -> System:
    """Return the Ott-Antonsen equation for the order parameter of infinitely many Kuramoto oscillators.

    For natural frequencies with a Lorentzian distribution of center `omega` and half width `Delta`, the complex order
    parameter `z = x + i y` follows exactly `z' = (i omega - Delta) z + (K / 2)(1 - |z|²) z`, a Stuart-Landau
    oscillator with the growth rate `K / 2 - Delta` and the saturation `K / 2`. Its modulus settles at
    `sqrt(1 - 2 Delta / K)` above the critical coupling `K_c = 2 Delta`, and at 0 below it.

    Args:
        r:
            The initial modulus of the order parameter.

        K:
            The coupling strength.

        Delta:
            The half width of the distribution of natural frequencies.

        omega:
            The center of the distribution of natural frequencies.

    Returns:
        system:
            The model, with the real part `field.x` and the imaginary part `field.y` of the order parameter.

    Examples:
        >>> from skcomplex.models.physics.synchronization import ott_antonsen
        >>> ott_antonsen(K=3.0).elements[0].kappa
        1.5
    """
    return System(elements=[StuartLandauElement('field', x=r, y=0.0, mu=K / 2 - Delta, omega=omega, kappa=K / 2)])

swarmalators(n_units=100, *, J=1.0, K=0.0, A=1.0, B=1.0, random_state=0)

Return swarmalators, identical oscillators whose phases set their attraction and whose positions their coupling.

Every swarmalator i moves by (1 / N) sum_j [(x_j - x_i) / |x_j - x_i| (A + J cos(theta_j - theta_i)) - B (x_j - x_i) / |x_j - x_i|^2] and turns by (K / N) sum_j sin(theta_j - theta_i) / |x_j - x_i|. With J > 0 like phases attract. With K = 0 and J = 1 the swarm settles into a static phase wave, a disk whose phase winds once around its center, so one of the order parameters S± = |mean(exp(i (phi ± theta)))| of the spatial angle phi and the phase approaches 1. With K > 0 the phases synchronize into a uniform disk, and with K < 0 and small J they stay asynchronous.

Parameters:

Name Type Description Default
n_units int

The number of swarmalators N.

100
J float

The influence of the phase similarity on the attraction.

1.0
K float

The phase coupling.

0.0
A float

The strength of the attraction.

1.0
B float

The strength of the repulsion.

1.0
random_state int | Generator | None

The seed or generator of the initial positions in the square of side 2 and of the initial phases.

0

Returns:

Name Type Description
system System

The model, with the positions units.x and the phases units.theta.

Examples:

>>> from skcomplex.models.physics.synchronization import swarmalators
>>> swarmalators(10, K=-0.1).interactions[0].K
-0.1
Source code in src/skcomplex/models/physics/synchronization/_swarmalators.py
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
@register_model(
    field='physics/synchronization',
    tags=('swarmalators', 'synchronization', 'swarming', 'pair interaction', 'ode'),
    references=("O'Keeffe, Ceron, and Petersen 2017", 'te Vrugt and Wittkowski 2024'),
    defaults={'until': 100.0, 'dt': 1.0, 'solver': RK4Solver(dt=0.1)},
)
def swarmalators(
    n_units: int = 100,
    *,
    J: float = 1.0,
    K: float = 0.0,
    A: float = 1.0,
    B: float = 1.0,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return swarmalators, identical oscillators whose phases set their attraction and whose positions their coupling.

    Every swarmalator `i` moves by `(1 / N) sum_j [(x_j - x_i) / |x_j - x_i| (A + J cos(theta_j - theta_i)) -
    B (x_j - x_i) / |x_j - x_i|^2]` and turns by `(K / N) sum_j sin(theta_j - theta_i) / |x_j - x_i|`. With `J > 0`
    like phases attract. With `K = 0` and `J = 1` the swarm settles into a static phase wave, a disk whose phase winds
    once around its center, so one of the order parameters `S± = |mean(exp(i (phi ± theta)))|` of the spatial angle
    `phi` and the phase approaches 1. With `K > 0` the phases synchronize into a uniform disk, and with `K < 0` and
    small `J` they stay asynchronous.

    Args:
        n_units:
            The number of swarmalators `N`.

        J:
            The influence of the phase similarity on the attraction.

        K:
            The phase coupling.

        A:
            The strength of the attraction.

        B:
            The strength of the repulsion.

        random_state:
            The seed or generator of the initial positions in the square of side 2 and of the initial phases.

    Returns:
        system:
            The model, with the positions `units.x` and the phases `units.theta`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import swarmalators
        >>> swarmalators(10, K=-0.1).interactions[0].K
        -0.1
    """
    rng = np.random.default_rng(random_state)
    units = GroupElement(
        'units',
        element=SwarmalatorElement,
        size=n_units,
        values={'x': rng.uniform(-1.0, 1.0, (n_units, 2)), 'theta': rng.uniform(-np.pi, np.pi, n_units), 'omega': 0.0},
        visual=VisualSpec('points', {'position': 'x', 'color_by': 'theta', 'radius': 0.03}),
    )
    coupling = SwarmalatorCouplingInteraction(
        'coupling',
        first='units',
        second='units',
        A=A,
        B=B,
        J=J,
        K=K,
        n_units=n_units,
    )
    return System(elements=[units], interactions=[coupling], space=ContinuousSpace(n_dim=2))

winfree(n_oscillators=50, *, K=0.5, spread=0.1, seed=0)

Return the Winfree model, theta_i' = omega_i + (K / N) sum_j P(theta_j) R(theta_i), with the sum over all j.

The influence P(theta) = 1 + cos(theta) and the sensitivity R(theta) = -sin(theta) make every oscillator respond to the mean influence of the population, which every oscillator reads as a group mean. The natural frequencies are spread evenly over [1 - spread, 1 + spread]. Weak coupling leaves the oscillators incoherent, moderate coupling locks them to one frequency, and strong coupling stops them, which is oscillator death.

Parameters:

Name Type Description Default
n_oscillators int

The number of oscillators.

50
K float

The coupling strength.

0.5
spread float

The half width of the interval of natural frequencies around 1.

0.1
seed int

The seed of the initial phases.

0

Returns:

Name Type Description
system System

The model, with the phases oscillators.theta and the order parameter order.

Examples:

>>> from skcomplex.models.physics.synchronization import winfree
>>> len(winfree(5).interactions)
1
Source code in src/skcomplex/models/physics/synchronization/_pulses.py
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
@register_model(
    field='physics/synchronization',
    tags=('synchronization', 'oscillators', 'phase', 'oscillator death'),
    references=('Winfree 1967', 'Ariaratnam and Strogatz 2001', 'Bick et al. 2020'),
    defaults={'until': 50.0, 'dt': 0.5},
    observables=lambda system: _order_parameter('order', 'oscillators'),
)
def winfree(n_oscillators: int = 50, *, K: float = 0.5, spread: float = 0.1, seed: int = 0) -> System:
    """Return the Winfree model, `theta_i' = omega_i + (K / N) sum_j P(theta_j) R(theta_i)`, with the sum over all `j`.

    The influence `P(theta) = 1 + cos(theta)` and the sensitivity `R(theta) = -sin(theta)` make every oscillator respond
    to the mean influence of the population, which every oscillator reads as a group mean. The natural frequencies are
    spread evenly over `[1 - spread, 1 + spread]`. Weak coupling leaves the oscillators incoherent, moderate coupling
    locks them to one frequency, and strong coupling stops them, which is oscillator death.

    Args:
        n_oscillators:
            The number of oscillators.

        K:
            The coupling strength.

        spread:
            The half width of the interval of natural frequencies around 1.

        seed:
            The seed of the initial phases.

    Returns:
        system:
            The model, with the phases `oscillators.theta` and the order parameter `order`.

    Examples:
        >>> from skcomplex.models.physics.synchronization import winfree
        >>> len(winfree(5).interactions)
        1
    """
    oscillators = GroupElement(
        'oscillators',
        element=PhaseOscillatorElement,
        size=n_oscillators,
        values={
            'theta': np.random.default_rng(seed).uniform(0, 2 * np.pi, n_oscillators),
            'omega': np.linspace(1 - spread, 1 + spread, n_oscillators),
        },
    )
    return System(
        elements=[oscillators],
        interactions=[WinfreeMeanFieldInteraction('coupling', target='oscillators', K=K)],
    )