Frustration Induced Chimeras and Motion in Two Dimensional Swarmalators

R. Senthamizhana, R. Gopala , V. K. Chandrasekara

a
Department of Physics, Centre for Nonlinear Science and Engineering, School of Electrical and Electronics Engineering,SASTRA Deemed
University, Thanjavur, 613 401, Tamil Nadu, India

Abstract
Swarmalators are phase oscillators capable of simultaneous swarming and synchronization, making them potential
candidates for replicating complex dynamical states. In this work, we explore the effects of a frustration parameter
in the phase interaction functions of a two-dimensional swarmalator model inspired by the solvable Sakaguchi-
arXiv:2501.05135v1 [nlin.AO] 9 Jan 2025

swarmalators that move in a one-dimensional ring. The impact of the frustration parameter in these models has
been a topic of great interest. Real-world coupled systems with frustration exhibit remarkable collective dynamical
states, underscoring the relevance of this study. The frustration parameter induces various states exhibiting non-
stationarity, chimeric clustering, and global translational motion, where swarmalators move spontaneously in
two-dimensional space. We investigate the characteristics of these states and their responses to changes in the
frustration parameter. Notably, the emergence of chimeric states suggests the crucial role of non-stationarity in
phase interactions for spontaneous population clustering. Additionally, we examine how phase non-stationarity
influences the spatial positions of swarmalators and provide a classification of these states based on different order
parameters.
Keywords: Swarmalator model, phase-lag, frustration parameter, phase frustration

1. Introduction base for understanding collective motion and behav-

ior in systems of self-propelled particles, serving as a
Synchronization is a captivating phenomenon fundamental framework in studying swarm dynamics,
prevalent in nature and various artificial systems [1]. [12] and active matter physics [13]. Also, this model
It manifests in diverse forms, from the firing of neu- has been instrumental in understanding behaviors in
rons [2] and the flashing of fireflies [3] to the os- natural systems like bird flocking, fish schooling, etc.
cillations of metronomes on a shared platform [4]. [14, 15]. Despite being studied independently, efforts
Arthur Winfree introduced the first significant math- have been made to unify synchronization and swarm-
ematical model to study synchronization in 1967 [5]. ing. For instance, models such as mobile oscillators
Winfree’s model was groundbreaking and paved the and chemotactic oscillators have emerged, combining
way for further research. Building on this, Yoshiki these phenomena [16, 17, 18].
Kuramoto presented a simplified and analytically
tractable model in 1975 [6]. The Kuramoto model A foundational work to combine these phenom-
has since become a cornerstone in understanding syn- ena was conducted by Tanaka and Iwasa, who intro-
chronization, extensively studied in systems such as duced a generalized model of chemotactic oscillators
power systems analysis [7], ecological systems [8], capable of exhibiting diverse dynamical behaviors
and chemical systems [9], providing profound insights [19, 20]. Building on this, a significant advancement
into the underlying mechanisms of collective behav- in the field is the “swarmalator” model introduced by
ior. O’Keeffe, Hong, and Strogatz in 2017 [21], which
Building on the understanding of synchroniza- integrates synchronization and swarming behaviors
tion, another intriguing phenomenon, swarming, in- within a unified framework, allowing the study of sys-
volves the spatial aggregation of multiple agents [10]. tems where individuals exhibit both spatial cluster-
While synchronization demonstrates self-organization ing and phase coherence. This model applies to var-
in time, swarming shows self-organization in space. A ious natural systems, including Japanese tree frogs,
notable model for swarming was introduced by Vic- sperms, and vinegar eels [22, 23, 24]. Since its in-
sek in 1995 [11]. The Vicsek model provides the troduction, the swarmalator model has been explored
in different configurations, considering the effects of
noise, delay, pinning, forcing, higher-order interac-
Email addresses: gopalphysics@gmail.com (R. Gopal), tions, Higher-harmonic interactions, and time-varying
chandru25nld@gmail.com (V. K. Chandrasekar) interactions [25, 26, 27, 28, 29, 30, 31, 32, 33].
Preprint submitted to Chaos, Solitons & Fractals January 10, 2025
 Since the two-dimensional swarmalator model in- coupled phase and spatial dynamics proposed by
troduced by O’Keeffe et al. exhibits a chaotic na- O’Keeffe et. al. [21] to investigate various collective
ture [34], a solvable model of swarmalators, which dynamical states. In this study, we consider a general-
includes a pair of Kuramoto models, has also been ized swarmalator model, incorporating frustration pa-
developed [35] and studied under various settings rameters as given in the Sakaguchi-Kuramoto model
[36, 37]. An exciting modification involves including [38], described by the equations (ẋ, θ̇) as follows,
a frustration parameter or phase lag in the phase inter-
action function. Sakaguchi and Kuramoto extensively  N
X x j − xi
1
studied the impact of this frustration parameter in the ẋi = vi + ǫ1 (A + J cos(θ j − θi + α1 ))

Kuramoto model in 1986 [38]. Recently, the study of N 
j,i x j − xi
frustration parameters has also been incorporated into 
the solvable model of swarmalators by Lizárraga et x j − xi 
−B  ,
ǫ2 
al., [39], identifying several interesting states, such as x j − xi
the active asynchronous and turbulent states. N
K X sin(θ j − θi + α2 )
The effect of the frustration parameter in the Ku- θ̇i = ωi + ǫ3 . (1)
ramoto model has been studied extensively under var- N j,i x j − xi
ious settings [40, 41], revealing rich dynamics in-
cluding the emergence of frequency synchronization Where, xi = (xi , yi ) denotes the spatial coordinate
[42], complex phase locking patterns [43], and the vector of the ith particle, and θi denotes the internal
presence of chimera states [44] where coherent and phase of the ith swarmalator. N is the total number
incoherent oscillations coexist within the same sys- of swarmalators. The parameters vi and ωi , represent-
tem.One notable outcome is the induction of frustra- ing the velocity and natural frequency of the ith swar-
tion, which prevents the oscillators from achieving a malator, respectively, can be set to zero by selecting an
perfect synchronization [45, 46]. Our research ex- appropriate reference frame without any loss of gen-
plores the introduction of phase frustration in the two- erality. The value of A is set to unity to ensure that
dimensional swarmalator model by adding a phase lag the spatial attraction function in ẋi is positive definite,
which act as a frustration parameter in the phase in- while the repulsion strength B is assigned as 1. The
teraction functions. Since the phase is coupled to the interaction strengths K and J quantify the influence
spatial position in our model, this non-stationarity re- of surrounding swarmalators on the phase and spatial
sults in unique states ranging from radially symmet- coordinates of the ith swarmalator, respectively. To
ric chimera to states with global translational motion. achieve long-range attraction and short-range repul-
Although the two-dimensional model of swarmala- sion among the swarmalators, the exponents are set
tors is not exactly solvable, it yields a plethora of as ǫ1 = 1 and ǫ2 = 2. Additionally, ǫ3 is set to 1
states exhibiting interesting dynamics, paving the way to ensure that the influence of spatial distance in the
for a deeper understanding of several naturally occur- phase equation θ̇i remains long-range, consistent with
ring swarming phenomena, such as chimeric flash- the phase-dependent spatial attraction.
ing in Photuris frontalis fireflies where both synchro- In this model, the parameters α2 and α1 (both < π2 )
nized cluster and independent flashers exist in a single represent phase shifts. The interplay between α1 ,
swarm [47] and swarm motility in several species of α2 , and the coupling constants K and J induces non-
bacteria like Bacillus subtilis, in which a multicellular stationarity in the system. Extending frustration pa-
group of bacteria executes a coordinated movement rameters to both phase interaction functions in the
through the synchronization of flagella as a result of xy model reveals multiple non-stationary states. The
hydrodynamic interactions [48, 49]. emerging dynamical states and the order parameters
This paper is organized as follows: the introduction used are explored in the subsequent sections.
is provided in Sec. I, followed by a description of the
model in Sec. II. Sec. III further explores the charac- 3. Results and Discussion
teristics of chimeric, translational, and other existing
states and the effect of changes in the frustration pa- All the simulations were performed using the RK-
rameter on the emerging translational states. Sec. IV 4 numerical algorithm with a population size of N =
discusses the order parameters used and the existence 500 and a step size of dt = 0.01. The swarmalators
of different states with respect to a range of interaction were initialized with positions x and y uniformly dis-
strengths. Finally, Sec. V concludes the paper. tributed in the range (−1, 1). Their phases θ were ran-
domly initialized following a uniform distribution in
2. Model the range [0, 2π). For various pairs of J and K with
α2 = α1 < π2 , the swarmalator system exhibits several
We examine an N-particle model known as swar- unique dynamical states. These include static states,
malators in a two dimensional space (x ∈ R2 ) with such as the static asynchronous state (SAS) and static
2
chimera (SC), as well as active states, including the have also analyzed the frequency of each swarmala-
active phase wave (APW), active synchronized state tor with the help of the frequency histogram, shown
(AS), and active chimera (AC). Additionally, the sys- in Figure 2, and it reveals the existence of distinct
tem demonstrates globally dynamic states like global synchronized and desynchronized clusters within the
translational motion (GTM) and synchronized global population in SC state.
translational motion (SGTM). It is worth noting that
varying α1 and α2 independently near π2 does not in-
troduce new states in the JK parameter space. These
states are discussed in detail in the following subsec-
tions.

3.1. Static states

By introducing the frustration parameters α2 and α1
in Eq.(1), two distinct static states are observed: the
static asynchronous state and the static chimera state.
Static asynchronous states (SAS) arise for small Figure 2: Histogram of phase frequency ν shows the frequency dis-
values of J and K, such as J = 0.05 , K = 0.05, tribution of (a) static asynchronous state (SAS) with near zero fre-
and α1 = α2 = 1.568 as depicted in Figure 1 (a). quencies and (b) static chimera state (SC) with multiple frequen-
In this state, the swarmalators are fixed in position, cies.
forming a circular disk-like pattern. A random dis-
From figure 2 (a); it is evident that the frequency
tribution of phases characterizes static asynchronous
of all the swarmalators is identical in the static asyn-
state. This randomness in phase distribution indicates
chronous state, indicating they are frequency-locked
a lack of global synchronization despite the spatial or-
despite their randomly distributed phases. Even
ganization.
though the SAS shows a near-zero frequency, as
shown in figure 2 (a), it still exhibits non-stationarity,
and the low values of the interaction strengths J and
K explain the near zero frequency. In contrast, the
static chimera state shows that the synchronized clus-
ter is concentrated around the frequency ν = −1 and
a few randomly distributed frequencies as seen in the
figure 2 (b). This confirms that the SC state exhibits
chimeric properties.
In static chimera state, although oscillators are clus-
tered into synchronous and asynchronous groups, the
Figure 1: Static states obtained for different pairs of interaction spatial distribution of these clusters exhibits radial
strengths J and K , with α1 = α2 = 1.568 symmetry. Notably, the synchronous cluster ceases to
: (a) Static asynchronous state (SAS) for J = 0.05, exist beyond a certain radius from the centroid of the
K = 0.05, (b) Static chimera (SC) showing the population. The radius of the synchronized population
radially symmetric clustering behavior but without Rc is calculated from the time-averaged order param-
any within-population motion for J = 0.05, K = 1.0. eter from Equation 2 which calculates the degree of
synchronization of swarmalators inside the specified
In the static chimera (SC) state, the swarmalators radius R,
show the co-existence of synchronized and desyn-
chronized population. Here, the synchronized popu- * +
lation occupies the central position in the radial ar- 1 X
hr(R)it = eiθi (R) (2)
rangement. In contrast, the desynchronized popula- NR i∈{i:ri ≤R} t
tion forms a ring-like structure around the synchro-
nized group, as shown in Figure 1 (b). This formation Where, ri is the distance of the ith swarmalator from
of distinct populations occurs spontaneously, without the spatial center of the population, R is the arbi-
any conditional separation of the swarmalator popu- trary radius from the spatial center of the population,
lations. Similar to the static asynchronous state, the NR is the number of swarmalators within the radius
phases continuously evolve due to the frustration in R. Here, the order parameter value remains near one
the phase interaction functions. One may also note when the radius includes only the synchronized popu-
that the introduction of frustration in phase oscillators lation. However, as the value of that radius increases
results in chimera states with clusters having locked and begins to encompass the desynchronized popula-
and nonlocked frequencies in Ref.[50]. Similarly, we tion, the respective order parameter value decreases,
3
as shown in Figure 3 (a). The radius at which the or-
der parameter starts to fall below 0.98 is marked as the
critical radius Rc for the synchronized population.

Figure 4: Snapshot of the active phase wave (APW) obtained for

J = 1.0 , K = 0.05 and α1 = α2 = 1.568. (b) Spatial angle φi vs. θ
of the APW, depicting the correlation between θ and φ.

Figure 3: (a) Time averaged order parameter vs. arbitrary radius R

shows the order parameter from Eq. 2, calculated for the population
NR as the radius R varies. The critical radius Rc , beyond which
would not be possible in the absence of frustration
the synchronized population ceases to exist, is indicated by the red parameters. Also, the swarmalators execute a radial
dashed line. (b) Frequency vs. radius confirms that the population oscillation, moving from the center of the population
beyond Rc is desynchronized. to the outer edge and then back to the center. This
movement repeats over time. The change in the ra-
Figure 3(b) indicates that the frequency ν of swar- dial distance of a swarmalator from the center of the
malators within the critical radius Rc takes a sin- population is correlated with the rate of change of θ
gle value, whereas the desynchronized population ex- (frequency ν). This correlation is confirmed by cal-
hibits multiple frequencies. This observation confirms culating the power spectral density of the time series
the chimeric characteristics of the SC in the spatial do- of both the radial displacement and the frequency of a
main. random swarmalator, and we can confirm from figure
5 (c) that the spectral density shows the existence of
3.2. Active states
same frequency spectrum for both the time series.
So far, we have discussed scenarios where the swar-
malators exhibit a static spatial property. However,
the swarmalators exhibit active behavior under certain
parametric conditions, where they do not remain sta-
tionary and change positions over time. The resulting
active states are discussed as follows.
Due to the inclusion of the frustration parameters
α1 and α2 , five different active states are observed
for various pairs of J and K. One of these states is
the active phase wave (APW), documented in vari-
ous literature sources [21, 33]. In the APW state, the
swarmalators form an annular structure and execute
spatial counter-rotation within the ring structure, as
shown in Figure 4 (a). As the term ‘phase wave’ sug-
gests, the swarmalators’ phases θi are correlated
 with
yi −cy
their spatial angles φi , defined as φi = tan−1 xi −c x ,
where (c x , cy ) represents the centroid of the popula-
tion, as depicted in Figure 4(b). Unlike the previously
reported APW state, due to the inclusion of frustra- Figure 5: (a) Active synchronized (AS) state obtained for J = 1.0
tion parameters, Figure 4 demonstrates the APW state ,K = −0.25, and (b) active chimera (AC) state obtained for J =
with few number of swarmalators occupies the central 0.55, K = 0.8, (c) The power spectral density obtained from the
time series of the radial displacement of the i = 100th swarmalator
region. and it’s phase frequency shows the correlation between them and
For J = 1.0, K = −0.25, and α1 = α2 = 1.568, the other fixed parameters are α1 = α2 = 1.568
the active synchronized (AS) state is observed. In
the active synchronized state (Figure 5 (a)), the swar- Like the static chimera state, the active chimera de-
malators achieve fully synchronized phases despite picted in Figure 5 (b) also possesses synchronized and
the negative phase coupling K. This phase synchro- desynchronized populations. The main difference be-
nization is particularly interesting because the sys- tween these states is that the swarmalators do not re-
tem’s non-stationarity favors phase synchronization main in either of the population constantly. Due to
mediated by the spatial attraction strength J, which the radial oscillation of the swarmalators in the active
4
chimera, similar to the active synchronized state, the grows and oscillates randomly, confirming desynchro-
swarmalators from the synchronized population in the nization.
central region of the disc move to the outer desynchro-
nized ring and then fall back to the central synchro- 3.3. Active states with translational motion
nized population with the evolution of time. The num-
In this section, we discuss states with translational
ber of swarmalators in the synchronized population is
motion, which is characterized by the uniform collec-
more significant in the active chimera state compared
tive motion of the swarm [54], even though individual
to the static chimera state at any instant of time. This
swarmalators execute motion within the population.
can be verified using the frequency histogram of the
This collective motion in xy space is due to the in-
AC state shown in Figure 6 (a).
clusion of frustration parameters α1 and α2 (see equa-
tion 1), which induce non-stationarity in both phase
and spatial position of each swarmalator. Here, we
observed two states with non-stationary spatial dy-
namics: the global translational motion (GTM) state
(J = 1.0 and K = 0.32) with wide distribution of
phases and the synchronized global translational mo-
tion (SGTM) state (J = −1.0 and K = 1.0), in which
the swarmalators show synchronized phases, which
are shown in Figure 7(a) and 7(b), respectively.

Figure 7: Snapshot of (a) global translational motion (GTM) for

Figure 6: (a)The frequency ν versus swarmalator index i histogram J = 1.0, K = 0.32, (b)synchronized global translational motion
for the active chimera state shows the disproportionate existence of (SGTM) for J = −1.0 and K = 1.0, and α1 = α2 = 1.568 for both
frequencies of synchronized and desynchronized populations, (b) GTM and SGTM
shows the critical radius (Rc ) at which the order parameter < r >t
starts to decrease, (c) radial displacement time series shows the ran- Here, in GTM, a subgroup in the central region of
dom sample swarmalators (i = 220, 470) shows near zero phase dif- the spatial structure which has the highest degree of
ference ∆θi when the swarmalators move inside the critical radius
Rc
directional synchrony aligns towards the direction of
the entire population’s global collective motion and
~ of
is identified by calculating the velocity vector v(t)
The critical radius Rc calculated using the equa-
each swarmalator which is given by
tion 2 is marked at the radius where the order pa-
rameter value falls below the threshold value of 0.98,
as depicted in Figure 6(b). Since there is an active ~ ~
exchange of swarmalators between synchronized and ~ = ∆r = (x(t + ∆t) − x(t))i + (y(t + ∆t) − y(t)) j
v(t)
∆t ∆t
desynchronized populations, this active chimera state (3)
resembles a breathing chimera state observed in phase
oscillators, where the synchronized and desynchro- To analyze the directional properties of these veloc-
nized domains oscillates in time [51, 52, 53] lead- ~ is translated to a common
ity vectors, each vector v(t)
ing to the continuous switching of oscillators between origin as seen in figure 8 (a). Translation operation
both domains. This behavior is confirmed by Fig- ensures that the analysis is centered on the directional
ure 6(c), where the time series of radial displacement characteristics rather than the positions of the vectors.
of two random swarmalators i = 220, 470 shows the After translation, the vectors are normalized, resulting
switching between the synchronized and desynchro- in vectors with their heads placed on a unit circle as
nized regions in space by moving inside and out of given in figure 8 (b). This transformation allows for a
the critical radius Rc . When both swarmalators fall be- precise analysis of the predominant orientation of the
low Rc , they exhibit a near-zero phase difference ∆θi , vectors. The heads of the normalized velocity vec-
emphasizing synchronization and when both oscilla- tors are analogous to the phase oscillators in unit cir-
tors move above the radius Rc , their phase difference cle, therefore we used the Kuramoto order parameter
5
given in the equation 4 to determine at which direc- Here, ∆c x (t), ∆cy (t) is the velocity of centroid in x
tion most of the vectors are aligned, the angle of the and y coordinates respectively at time t. The compar-
order parameter shows the direction of alignment of ison reveals that the directional evolution of the pre-
the predominant velocity vectors and the magnitude, dominant group calculated from the order parameter
represented by the length of the arrow, indicates the aligns with that of the centroid, underscoring the in-
degree of synchrony. fluence of the directionally synchronized swarmala-
tors on the overall motion of the population as shown
in figure 9 (b).
The above study clearly emphasizes that the global
translational motion of the swarmalators depends on
the within-population directional synchrony of the
swarmalator group in the central region of the spa-
tial structure and the population in the central region
is actively replaced by other swarmalators, preserving
the number of directionally synchronized swarmala-
tors. In addition to the above dynamics, unlike all
other states, the spatial distribution of swarmalators
Figure 8: (a) Velocity vectors of swarmalators in GTM translated
to a common origin. (b) Normalized velocity vectors with the Ku- in the circular disc arrangement is not uniform at any
ramoto order parameter (equation 4) depicted as an arrow. instant, as is clearly evident from Figure 7(a), where
swarmalators form two less dense regions adjacent to
the central area.
In the case of the SGTM state, the swarmalator
exhibits high degree of phase synchronization when
compared to GTM. A small population at the cen-
ter of the radial arrangement exhibits desynchronized
phases. SGTM doesn’t show within population direc-
tional synchrony, because of it’s radial oscillatory na-
ture similar to the active chimera state and it is shown
in Figure 10.

Figure 9: (a) Time series of the magnitude of the Kuramoto order SGTM
π
parameter obtained from the velocity vectors. (b) Time evolution −3.50
of the angle Θ1 of the Kuramoto order parameter. Θ2 is the time
evolution of the velocity vector obtained from the centroid of the
entire population of swarmalators. −4.75
vy

θ
0

The order parameter for these velocity vectors are

calculated over time and the time series of the mag- −6.00
nitude r of the Kuramoto order parameter in figure 9 −2.5 −1.0 0.5
−π

(a) demonstrates that r ≈ 0.5 over time, indicating vx

consistent directional synchrony within the population
regardless of individual movements of the swarmala- Figure 10: The velocity vectors of swarmalators in SGTM state
tors. Additionally, the angular evolution of the Ku- shows the radial alignment with vectors pointing towards and away
ramoto order parameter, Θ1 (t), given by from the center. The vector’s bases are drawn from the position of
it’s corresponding swarmalator and scaled by the factor of 0.01 for
better visualization.
 N

 1 X 
Θ1 (t) = arg  eiψ j (t)  , (4)
N j=1 We can observe that the properties of GTM and
SGTM exhibit distinct characteristics and possess dif-
v̂y (t)i
! fering magnitudes of translational velocity. However,
−1
ψi (t) = tan (5) to compare the centroid velocity of GTM and SGTM,
v̂ x (t)i we plotted the displacement of the centroid in xy space
Where v̂ x (t)i and v̂y (t)i are the velocity vectors corre- (see figure 11), in which the arrows represents the di-
sponding to ith swarmalator at time t. Futher, Θ1 (t) rection of movement of the centroid, and it can be ob-
compared with the angular evolution of the centroid’s served from the figure 11, that the translational mo-
velocity vector, Θ2 (t) from equation 6. tion of SGTM has a relatively low centroid velocity,
! results in the small displacement (See Fig. 11(b)),
∆cy (t)
Θ2 (t) = tan−1 (6) compared with the GTM state, where we observe a
∆c x (t) relatively large displacement (See Fig. 11(a)). This
6
emphasizes the crucial role of directional synchrony malators exhibit GTM only when α2 ≥ π4 , regardless
in the global spatial displacement of GTM states. of the value of α1 , highlighting the importance of α2
for the initiation of GTM. Furthermore, as α1 and α2
approaches π2 in both the cases, the range of GTM de-
creases. Therefore, one can confirm that the frustra-
tion parameters are crucial and the interplay between
interaction strengths J, K and parameters α1 , α2 gives
rise to the states with translational motion. It is worth
noting that selecting values of α1 and α2 near π2 en-
hances the effect of non-stationarity in the system.

4. Order Parameter-Based State Classification

Figure 11: Evolution of the centroid over time for (a) GTM and (b)
SGTM. The blue square marker indicates the initial position, and In our study, distinct order parameters are exam-
the red arrows show the direction of the centroid’s evolution. ined to investigate the existence of various dynam-
ical states. To differentiate states based on the de-
gree of synchronization, the Kuramoto order param-
eter is commonly utilized [55]. The time-averaged
Kuramoto order parameter magnitude hrit is defined
as:

Z T N
1 1 X iθ j (t)
hrit = e dt (7)
T 0 N j=1

Here, θ j is the phase of the jth swarmalator. States

such as SAS, GTM, and APW, which lack global
phase synchronization, yield values of hrit ≈ 0. In
contrast, chimera states like AC and SC exhibit values
in the range 0 < hrit < 1, while AS and SGTM states,
characterized by global synchronization, have values
of hrit ≈ 1.
Furthermore, states involving translational motion,
such as GTM and SGTM, can be characterized by cal-
culating the maximum displacement of the centroid
(c x(t) , cy(t) ) of the swarmalators at time t, expressed as
Figure 12: The change in total centroid displacement δ, calculated follows:
from Eq. 9, with respect to the frustration parameter α1 (with α2
fixed) in (a) and with varying α2 (with α1 fixed) in (b) demonstrates  
the transition of swarmalators from a non-GTM to a GTM state.  1 XN N
1 X 
(c x(t) , cy(t) ) =  x j (t), y j (t) , (8)
N j=1 N j=1
Both GTM and SGTM states emerge due to the in-
clusion of frustration parameters α1 and α2 . Vary-
ing the frustration parameters affects the displacement q
range and ability of these states to execute global δ= (c x(0) − c x(T ) )2 + (cy(0) − cy(T ) )2 (9)
translational motion. Since GTM shows enhanced
spatial displacement than SGTM, we can see the where δ represents the resultant magnitude of the
impact of frustration parameter in GTM state more centroid displacement between the initial time and the
clearly. To show the effect of change in frustration pa- final time T .
rameters in GTM state, we have illustrated the change For GTM and SGTM, δ > 0, while for all other
in total centroid displacement δ (see equation 9) of states, δ = 0.
GTM with respect to the various values of the frustra- Finally, to distinguish between active and static
tion parameters α1 and α2 in Figure 12. Figure 12(a) states, the displacement δ is calculated between two
shows that GTM occurs only after α1 ≈ 3π 8 , with a de- arbitrary times t1 and t2 such that t2 > t1 . How-
layed onset of GTM as α2 increases and minimal vari- ever, static states exhibit small displacement value δ
ation in GTM range across different α2 values. Simi- that fluctuates due to the non-stationary behavior in
larly, the results from Figure 12 (b) indicate that swar- phase and spatial-phase coupling, and active states
7
Figure 13: Variation in order parameter values with respect to change in interaction strengths (a) For J = 0.8, varying K reveals the existence
of SGTM and AS through changes in δ. (b) For J = 0.8, varying K indicates the presence of APW, GTM, and AC via < r >t and δ. (c) For
J = 0.1, varying K demonstrates the existence of SAS and AC through H, δ and < r >t . (d) varying J for K = 0.5, H reveals the presence of
SC among other states.

placement calculated from Eq. 9, and δm is the thresh-

old for the order parameter H(δ), set at 0.07 + 0.02.
The upper limit of error range 0.02 arises from slight
variations in observed values due to the system’s
chaotic nature [34]. Here the threshold value δm is
determined by averaging out the maximum spatial
displacement of each swarmalator from the popula-
tion’s centroid in the static chimera state for differ-
ent neighboring parameter values, where each swar-
malator exhibits small vibrational motion due to phase
non-stationarity. Value of δ greater than δm gives
H = 1 representing active state, while H = 0 which
Figure 14: The J vs K plot demonstrates the existence of various corresponds to static states.
dynamical states across different pairs of interaction strengths with
α2 = α1 = 1.568.
The defined order parameters successfully delineate
all the identified states, and the corresponding one-
parameter plots for different interaction strengths are
show large variations in the displacement. To address shown in Figure 13.
these variations in the displacement values in the both
In Figure 13(a), for J = 0.8 and −1 < K < 0, the
static and active states, a Heaviside step function is
translational motion of the SGTM leads to larger dis-
applied, as described by the following expression
placement values, resulting in higher δ values. Con-
 versely, for AS, δ fluctuates around zero due to its
1

 if δ > δm non-translational nature. In Figure 13(b), for J = 0.8
H(δ) =  (10) and 0 < K < 1, the AC state exhibits greater phase
0
 if δ ≤ δm
synchronization, as evidenced by higher values of the
where H(δ) is the transformed value, δ is the dis- Kuramoto order parameter in the range 0 < hrit < 1.0.
8
 States hrit δ H(δ)
Static Asynchronous State (SAS) ∼0 ∼0 0
Static Chimera (SC) 0 < hrit < 1 ∼0 0
Active Chimera (AC) 0 < hrit < 1 ∼0 1
Active Sync (AS) ∼1 ∼0 1
Global Translational Motion (GTM) ∼0 >0 1
Synchronized Global Translational Motion (SGTM) ∼1 >0 1
Active Phase Wave (APW) ∼0 ∼0 1

Table 1: Summary of the order parameters hrit , δ, and H(δ) for various states.

However, for GTM and APW, δ serves as the distin- symmetry is preserved, resulting in the observed sym-
guishing factor due to the translational nature of GTM metry in the J − K parameter space. This also explains
and the non-translational nature of APW. Since all the emergence of phase synchronization behaviors in
these states are active, H shows unitary value. the K < 0 regime.
In Figure 13(c), for J = 0.1 and 0 < K < 1, both
H and hrit yield zero value for SAS, while for AC, 5. conclusion
these parameters are non-zero, with hrit ≈ 0.5, re-
flecting its chimeric nature. δ fluctuates around zero We investigated the effect of the frustration parame-
as both states are non-translational. In Figure 13(d), ters in a two-dimensional swarmalator model. Includ-
for a smoothly varying spatial interaction strength J ing the frustration parameter induces non-stationarity
from −1 to 1 with fixed K = 0.5, H = 0 indicates the in the phases, leading to several intriguing dynami-
presence of SC around J = 0, whereas H = 1 appears cal states, such as chimera and global translational
elsewhere. The higher δ values for J < 0 suggest the motion (GTM) states. In chimeric states, we iden-
existence of SGTM, with hrit ≈ 1 for AS in J < 0, tified both static and active chimeras. In the static
and 0 < hrit < 1 for AC in J > 0. The Ranges of chimera, the swarmalators form clusters in phases
order parameter values that correspond to each state with synchronized and desynchronized properties,
is illustrated in Table 1. where the synchronized population occupies the cen-
Figure 14 presents a two-parameter plot that tral region, and the desynchronized population forms
sketches the existence of all the collective dynamical a ring around the central region. The active chimera
states mentioned above for different pairs of J and K. exhibits a breathing-type property, with the swarmala-
These regions are obtained from the order parameters tors switching back and forth between synchronized
defined above, utilizing the characteristic variations and desynchronized states while oscillating radially in
depicted in Figure 13. States such as global trans- spatial coordinates.
lational motion, active phase wave, and static asyn- Based on the coupling between phase and spatial
chronous state appear in the (+J, +K) and (−J, −K) movement in swarmalators, non-stationarity in phase
quadrants. In contrast, synchronized global transla- dynamics gives rise to complex states such as global
tional motion and active synchronized states are found translational motion (GTM) and synchronized global
in the (+J, −K) and (−J, +K) quadrants. The active translational motion (SGTM). In the GTM state, the
chimera state is observed in all quadrants. From Fig- entire population of swarmalators moves cohesively
ure 14, it is evident that introducing frustration pa- in the xy space while simultaneously exhibiting active
rameters α1 and α2 alters the landscape of swarmala- within-population movements. Compared to SGTM,
tor dynamics to a greater extent. Notably, we found GTM demonstrates more pronounced motion with
that the state distribution exhibits inverted symmetry higher translational velocity. Additionally, the system
about J = 0, which is a direct consequence of the exhibits states such as active phase wave (APW) and
system’s inherent symmetry under the transformation static asynchronous state (SAS), highlighting the rich
(θi , K, J) → (−θi , −K, −J). For α1 = α2 = π/2, this dynamical behaviors of swarmalators under the influ-
symmetry ensures that any dynamical state observed ence of frustration parameters. Our findings will be
at a point (J, K) in the parameter space will have a cor- helpful for significant implications for understanding
responding state at (−J, −K). In our study, with α1 and decentralized control mechanisms in artificial swarm
α2 taken near π/2 (specifically, α1 , α2 = 1.568), this systems such as robotic swarms. By manipulating the
9
frustration parameters, it might be possible to achieve [26] N. Blum et al., Physical Review. E 109, 014205, (2024).
desired collective behaviors in swarms of autonomous [27] G. K. Sar, D. Ghosh, and K. O’Keeffe, Physical Review. E
107, 024215, (2023).
robots, which can be applied to various fields such [28] G. K. Sar, K. O’Keeffe, and D. Ghosh, Chaos an Interdisci-
as search and rescue operations [56], environmen- plinary Journal of Nonlinear Science 33, 111103, (2023).
tal monitoring [57], and agricultural automation [58]. [29] J. U. F. Lizarraga and M. A. M. De Aguiar, Chaos an Interdis-
The future scope of this work lies in the ability to con- ciplinary Journal of Nonlinear Science 30, 053112, (2020).
[30] M. S. Anwar et al., Communications Physics 7, 59, (2024).
trol specific dynamical states like GTM’s directional [31] Smith, L.D., SIAM Journal on Applied Dynamical Systems,
evolution and spatial distribution of chimera states, 23(2), pp. 1133–1158, (2024).
which might enhance the efficiency and adaptability [32] Senthamizhan, R., Gopal, R. and Chandrasekar, V.K., Physi-
of swarm systems in complex and dynamic environ- cal Review E, 109(6), 064303, (2024).
[33] Sar, G.K. et al., New Journal of Physics, 24(4), p. 043004,
ments. (2022).
[34] S. Ansarinasab, F. Nazarimehr, F. Ghassemi, D. Ghosh, and
S. Jafari, Appl. Math. Comput. 468, 128508 (2024).
Acknowledgements [35] S. Yoon et al., Physical Review Letters 129, 208002, (2022).
[36] K. O’Keeffe, S. Ceron, and K. Petersen, Physical Review. E
The work of R.G. and V.K.C. forms part of 105, 014211, (2022).
a research project sponsored by ANRF-DST-CRG [37] G. K. Sar, K. O’Keeffe, and D. Ghosh, Chaos an Interdisci-
plinary Journal of Nonlinear Science 33, 111103, (2023).
Project Grant No. C.R.G./2023/003505. R.G. and [38] Sakaguchi, H. and Kuramoto, Y., Progress of Theoretical
V.K.C. thanks DST, New Delhi, for computational Physics, 76(3), pp. 576–581, (1986).
facilities under the DST-FIST programme (Grant [39] J. U. F. Lizárraga and M. a. M. De Aguiar, Physical Review.
No. SR/FST/PS-1/2020/135) to the Department of E 108, 024212, (2023).
[40] M. Manoranjani, S. Gupta, and V. K. Chandrasekar, Chaos
Physics. an Interdisciplinary Journal of Nonlinear Science 31, 083130,
(2021).
[41] M. Manoranjani et al., Physical Review. E 105, 034307,
References (2022).
[42] C.-H. Hsia et al., Journal of Differential Equations 268, 7897
[1] Bayani, A., Nazarimehr, F., Jafari, S. et al., Nat Commun 15, (2020).
4955 (2024). [43] B. Moyal et al., Physical Review. E 109, 034211, (2024).
[2] Penn, Y., Segal, M. and Moses, E., Proceedings of the Na- [44] E. A. Martens, C. Bick, and M. J. Panaggio, Chaos an Inter-
tional Academy of Sciences, 113(12), pp. 3341–3346, (2016). disciplinary Journal of Nonlinear Science 26, 094819, (2016).
[3] M. McCrea, B. Ermentrout, and J. E. Rubin, Journal of the [45] S.Y. Ha, Y. Kim, and Z. Li, SIAM J. Appl. Dyn. Syst. 13, 466
Royal Society Interface 19, (2022). (2014).
[4] J. Pantaleone, American Journal of Physics 70, 992 (2002). [46] P. S. Skardal and A. Arenas, Sci. Adv. 1, e1500339 (2015).
[5] A. T. Winfree, Journal of Theoretical Biology 16, 15 (1967). [47] R. Sarfati and O. Peleg, Science Advances 8, 46, (2022).
[6] Y. Kuramoto, in Springer eBooks, pp. 420–422, (2005). [48] D. B. Kearns, Nature Reviews Microbiology 8, 634 (2010).
[7] Y. Guo et al., International Journal of Electrical Power & En- [49] D. B. Kearns and R. Losick, Molecular Microbiology 49, 581
ergy Systems 129, 106804 (2021). (2003).
[8] J. Vandermeer et al., Royal Society Open Science 8(3), [50] J. Rode et al., Frontiers in Applied Mathematics and Statistics
(2021). 5, (2019).
[9] Kuramoto, Y., Berlin: Springer, (1984). [51] M. I. Bolotov et al., JETP Letters 106, 393 (2017).
[10] M. Schranz et al., Frontiers in Robotics and AI 7, (2020). [52] O. E. Omel’chenko, Journal of Nonlinear Science 32m, 22,
[11] T. Vicsek et al., Physical Review Letters 75, 1226 (1995). (2022).
[12] X. Wang, H. Zhao, and L. Li, Applied Sciences 13, 11513 [53] J. H. Sheeba, V. K. Chandrasekar, and M. Lakshmanan, Phys-
(2023). ical Review E 81, 046203, (2010).
[13] M. L. R. Puzzo et al., Journal of Physics Condensed Matter [54] A. Cavagna et al., Nature Physics 13, 914 (2017).
34, 314001 (2022). [55] J. A. Acebrón et al., Reviews of Modern Physics 77, 137
[14] H. Christodoulidi et al., in WORLD SCIENTIFIC eBooks, (2005).
pp. 383–398, (2014). [56] R. D. Arnold, H. Yamaguchi, and T. Tanaka, Journal of Inter-
[15] N. Phuoc et al., Communications in Physics 23, 121 (2013). national Humanitarian Action 3, 18, (2018).
[16] Meli, V.N. et al., Chaos, Solitons &; Fractals, 177, p. 114278, [57] M. Duarte et al., OCEANS 2016 - Shanghai (2016).
(2023). [58] D. Albiero et al., Computers and Electronics in Agriculture
[17] Fujiwara, N., Kurths, J. and Dı́az-Guilera, A., Physical Re- 193, 106608 (2022).
view E, 83(2), 025101, (2011).
[18] Sawai, S. and Aizawa, Y., Journal of the Physical Society of
Japan, 67(8), pp. 2557–2560, (1998).
[19] Tanaka, D., Physical Review Letters, 99(13), 134103, (2007).
[20] Iwasa, M. and Tanaka, D., Physics Letters A, 381(36), pp.
3054–3061, (2017).
[21] K. P. O’Keeffe, H. Hong, and S. H. Strogatz, Nature Commu-
nications 8, 1504, (2017).
[22] I. Aihara et al., Scientific Reports 4, 3891, (2014).
[23] I. H. Riedel, K. Kruse, and J. Howard, Science 309, 300
(2005).
[24] A. Peshkov, S. McGaffigan, and A. C. Quillen, Soft Matter
18, 1174 (2022).
[25] Hong, H. et al., Physical Review Research, 5(2), 023105,
(2023).

10