Chimera states (Abrams and Strogatz, 2004; Sethia et al., 2008; Laing, 2009; Hagerstrom et al., 2012; Martens et al., 2013; Omelchenko et al., 2013; Gopal et al., 2014; Schöll, 2016; Hart et al., 2019; Majhi et al., 2019; Parastesh et al., 2020; Wang and Liu, 2020) became a paradigm of collective phenomena in dynamical systems that started with the first report by Kuramoto and Battogtokh (2002) of two coexisting synchronous and asynchronous groups of phase oscillators arranged in a ring and coupled in a non-local fashion. The main question was how the symmetry in a completely synchronous ensemble of identical oscillators breaks into two clusters, one synchronous group and another asynchronous group. How stable is this chimera state? In the beginning, it was apprehended that the chimera state is a transient behavior; the transient time increases with the size of a network (Wolfrum and Omel’chenko, 2011; Rosin et al., 2014). Later, it has been established that chimera states are possible stable states (Pecora et al., 2014; Omel’chenko, 2018; Laing, 2019) in an ensemble of identical oscillators and with symmetry in the connectivity matrix or the topology of a network. By this time, this phenomenon has been widely explored in single-layer networks (Abrams and Strogatz, 2004; Sethia et al., 2008; Laing, 2009; Hagerstrom et al., 2012; Martens et al., 2013; Omelchenko et al., 2013; Gopal et al., 2014; Hart et al., 2019; Majhi et al., 2019; Parastesh et al., 2020; Wang and Liu, 2020), multilayer networks (Ghosh and Jalan, 2016; Maksimenko et al., 2016; Sawicki et al., 2018; Ruzzene et al., 2020), and 3D networks (Maistrenko et al., 2015; Kasimatis et al., 2018; Kundu et al., 2019) with different forms of chimeras such as traveling chimera (Bera et al., 2016a; Omel’chenko, 2019; Dudkowski et al., 2019; Alvarez-Socorro et al., 2021) and spiral chimera (Martens et al., 2010; Gu et al., 2013). Various dynamical models (Hizanidis et al., 2015; Banerjee et al., 2016; Bera et al., 2016b; Saha et al., 2019) with different coupling schemes (Meena et al., 2016; Bera et al., 2017) and global coupling (Sethia and Sen, 2014; Yeldesbay et al., 2014; Hens et al., 2015a; Mishra et al., 2015) have been used for observing chimera patterns. The concepts of amplitude chimera (Hens et al., 2015a; Anna et al., 2016; Banerjee et al., 2018) and amplitude death chimera (Zakharova et al., 2014; Banerjee, 2015; Dutta and Banerjee, 2015) have been introduced. All the examples of chimera states revolve around the symmetry breaking of a complete coherent state into coherent and incoherent groups. This perception has been extended further to the observation of two coexisting subgroups of in-phase and antiphase oscillators (Maistrenko et al., 2017), which was also referred to as a chimera state. The necessary requirement of a large set of oscillators for chimera states to observe has also been relaxed with a smaller size of the network: Chimera patterns emerge in a set of four oscillators (Hart et al., 2016; Meena et al., 2016; Senthilkumar and Chandrasekar, 2019) and even three oscillators (Wojewoda et al., 2016; Maistrenko et al., 2017). In a laser system, delay in coupling has been used (Hart et al., 2016) for the observation of chimera states in four oscillators. A modified Kuramoto phase oscillators with inertia (Maistrenko et al., 2017) were used to demonstrate both in-phase and antiphase chimeras in three oscillators and later confirmed in experiments with three attractively coupled pendula (Wojewoda et al., 2016). In almost all the reported studies, attractive coupling has been used for the observation of chimeras, while a few examples are found to use a combination of both attractive and repulsive coupling (Hens et al., 2015a; Mishra et al., 2015) to originate chimera states in globally coupled oscillators.

We focus here on the role of repulsive interactions in the origin of chimera patterns in a small ensemble of three identical dynamical units. As reported earlier (Mishra et al., 2017; Ray et al., 2020), chimera states may emerge in a large ensemble of globally coupled Josephson junctions under repulsive interactions. Indeed, the repulsive coupling can induce chimera states in a smallest set of three Josephson junctions arranged in a ring. Interestingly, a single Josephson junction shows typical neuron-like spiking and bursting behaviors (Hens et al., 2015b; Hongray et al., 2015; Dana et al., 2006; Mishra et al., 2021), which is one of the reasons that encourage us to investigate the role of repulsive interactions (inhibitory in the sense of neuronal interaction) in a ring of three junctions. A single Josephson junction is represented by a resistance–capacitance shunted junction (RCSJ) circuit (Dana et al., 2001; Dana et al., 2006; Mishra et al., 2021), which remains in an excitable state until an external bias current is applied across the junction. The junction shows spiking limit cycle oscillation for bias current above a critical value. The dynamics of the junction also depends upon a damping-like parameter that is related to the shunted resistance and capacitance of the junction. It shows a bistable region, in the parameter space, where the limit cycle coexists with a steady state for a range of low damping. The RCSJ model can be represented by a second-order phase dynamics of the junction that governs the voltage drop across it. The RCSJ model also represents the dynamics of the simple pendulum model and discrete sine-Gordon equation. The dynamics of the junction is rotational when the trajectory of the junction makes a complete rotation around a cylindrical surface like an inverted pendulum, thereby originating a large-amplitude spiking oscillation. On the other hand, the junction may show, for an appropriate choice of parameters, librational motion when the trajectory of the junction dynamics is restricted to a small region of the cylindrical phase space. This appears like the small-amplitude oscillation of a pendulum.

Three Josephson junctions with identical parameters and by our choice of parameters are kept in rotational motion in an uncoupled state but kept away from the bistable region to avoid further complexity. The ring of three identical oscillators represents an all-to-all (or globally) coupled network and hence perfect symmetry. The junctions in the ring show 2π/3 out-of-phase motion in a cyclic order in time but emerge into a state of complete synchrony (CS) for repulsive interactions above a critical value. When the repulsive interactions in the ring are increased further, the symmetry is broken with the emergence of a sequence of desynchronous states, followed by chimera patterns: 1) two junctions are in complete synchrony (CS) in a state of rotational motion and one junction in incoherent rotational motion when they are all in a chaotic state, 2) two junctions in CS in rotational motion and one junction in incoherent libration when the dynamics of the junctions is periodic; this chimera pattern flips to two coherent junctions in libration and one in incoherent rotation for a change in initial conditions. We numerically delineate the different collective states in a two-parameter plane of the junctions using coherence measures and a libration index and use the master stability function (MSF) to demarcate the stable region of CS in the parameter plane.

A schematic diagram of the ring of three junctions is shown in Figure 1A. The RCSJ circuit is shown at right that represents the dynamics of each node. The dynamics of the ring of junctions is represented by the phase dynamics θ _i and voltage v _i (i = 1, 2, 3) across the ith junction given as follows (Josephson, 1962; Mishra et al., 2017): θ i ̇ = v i , (1) v i ̇ = I − sin θ i − α v i + ε ( v i − 1 − 2 v i + v i + 1 ) , where v _i+1 = v ₁, if (i + 1) > n and v _i−1 = v ₃, if (i − 1) < 1. The initial conditions are randomly chosen from the range θ _i (0) ∈ [ − π, π] and v _i (0) ∈ [0.1, 1]. The damping parameter α = [ h / 2 π e I 2 R C ] 1 / 2 , where h is Planck’s constant and e is the electronic charge, and R and C are the resistance and capacitance of a junction, respectively. I is the external bias current normalized by the critical current I _c = 1.0 at each junction. The interaction between the nodes is established through junction voltages with a coupling strength ɛ, which is considered repulsive (negative) to observe our targeted chimera patterns. Figure 1B shows a phase diagram in an α − I plane for a single junction. The single junction remains excitable (FP) for bias current I < I _c and transits to limit cycle oscillation via saddle-node-invariant circle (SNIC) bifurcation for I > I _c and α > 1.19. In a lower range of α < 1.19 and I < I _c , the junction shows bistability (BS) when the limit cycle coexists with a steady state. The bistable state transits to limit cycle oscillation via fold bifurcation when I > 1.0. The bistable region is bounded by fold bifurcation (dashed horizontal line) and homoclinic (HC) bifurcation (blue curve) in the α − I parameter plane. For our observation of chimeras, we make a choice of parameters I > I _c and α > 1.19 so as to obtain limit cycle oscillation (rotational motion) in all three junctions in uncoupled state and keep them away from the bistable region to avoid further complexity in dynamics.

To check the CS state, we first calculate the error functions between all the pairs of nodes, separately, as follows: e 12 = ( v 1 − v 2 ) 2 ; e 13 = ( v 1 − v 3 ) 3 ; e 23 = ( v 2 − v 3 ) 2 , (2) where . indicates time average. We check all the three error functions simultaneously and detect a chimera state when two of the partial error functions are zero (say, e ₁₂ = 0 and e ₁₃ = 0).

The global error variable is then calculated as follows: e r r = ( e 12 + e 13 + e 23 ) / 3 . (3)

For CS in the network, err = 0 and err ≠ 0 indicate partial synchrony or incoherence. We use a libration index (Mishra et al., 2017) to distinguish librational and rotational motion of the junctions in different collective states: L I = 1 3 ∑ j = 1 3 Θ j , (4) where Θ _j = Θ(δ − m _j ). Θ(.) is the Heaviside step function, where δ is an arbitrarily chosen small threshold and m j = 1 − 0.5 * max ( cos [ θ j ( t ) ] ) − min ( cos [ θ j ( t ) ] ) . The libration index becomes LI = 0 for oscillators in libration and LI = 1 when they are in rotational motion. A value of 0 < LI < 1 indicates the coexistence of librational and rotational motions in the ring of junctions.

For demarcation of the stable CS state in the parameter plane of the junctions in Figure 2, we numerically estimate the MSF (Pecora and Carroll, 1998) of identical junctions by using the variational equation of the ring of junctions. The dynamics of the ith node is x ̇ i = F ( x i ) + ε ∑ j G i j H ( x j ) . Now, let x = (x ₁, x ₂, …, x _n ), F(X) = [F(x ₁ ), F(x ₂ ), …, F(x _n )], H(X) = [H(x ₁ ), H(x ₂ ), …, H(x _n )] is the coupling matrix and G be the adjacency matrix {G _ij }, then x ̇ = F ( x ) + ε G ⊗ H ( x ) , (5) where ⊗ denotes the direct product. The variational equation of Eq. 5 by assuming ξ _i as a small perturbation in the ith node when ξ = (ξ ₁, ξ ₂, …, ξ _n ), (i = 1, 2, ..n) is given as follows: ξ ̇ = [ 1 n ⊗ D F + ε G ⊗ D H ] ξ , (6) when H is the coupling matrix E = (0 0; 1 0), implying coupling in second variable and D H = E. By a block diagonalization of the variational equation, each block will be of the form as follows: ξ k ̇ = [ D F + ε γ k D H ] ξ k , (7) where γ _k is the eigenvalue of G, k = 0, 1, 2, …, n − 1. For k = 0, we have the variational equation for the synchronization manifold γ ₀ = 0. The Jacobian matrices D F and D H are the same for each block, since they are evaluated on the synchronized state. Thus, for each k, the form of each block [Eq. 7] is same with only the scalar multiplier ɛγ _k differing. This leads us to the following formulation of the master stability equation and the associated MSF: We calculate the maximum Lyapunov exponents λ _max for the generic variational equation as follows: ζ ̇ = [ D F + ( α + i β ) D H ] ζ , (8) as a function of α and β. This yields the λ _max as a point on the real axis since we check the real part of the exponent, and G has only the real eigenvalues in our case. The sign of λ _max will reveal the stability of the eigenmodes, and hence, we have the MSF, which is used for delineating the boundary (blue dashed lines) of the CS region (A, yellow) as shown in Figure 2, where λ _max becomes negative.

We draw a phase diagram to delineate the regions of a variety of collective states including the chimera patterns for the three junctions in the ring in an α − ɛ parameter plane in Figure 2. All the junctions are in out-of-phase for attractive coupling (ɛ > 0) when the junction voltages evolve in 2π/3 cyclic motion however become completely synchronized (CS) in the periodic state for a range of repulsive coupling (ɛ < 0). The boundaries of ɛ for CS are demarcated at two ends (by two blue dashed lines), which are indicated by a transition of λ _max to a negative value as estimated by using the MSF and perfectly matching with the numerically simulated boundaries of CS (region A, yellow) in the phase diagram. CS breaks down for attractive coupling at one end (red region) and larger repulsive coupling as well at the other end (regions B and G). CS transits to asynchronous periodic motion (B, red) of the junctions and a mixed state (G, green) for larger repulsive coupling. The chimera pattern (partial synchrony) coexists with complete asynchronous oscillation in region G (chaotic or higher periodic rotational motion). The ring of junctions becomes chaotic in region C (orange) and remains completely asynchronous. They form a chaotic chimera pattern with two coherent and one incoherent junctions in region D when all the junctions maintain rotational motion. Finally, a periodic state re-emerges in region E (cyan) with a chimera pattern. Two junctions become coherent in rotational periodic motion that coexists with a single junction in incoherent periodic libration. These chimera patterns flip with a new pattern when two junctions become coherent in libration and one in incoherent rotation; this flipping occurs due to sensitivity in initial conditions. Region F represents mixed states of coexisting chimera patterns and asynchronous states.

Three identical junctions in the ring evolve out of phase in time in a cyclic order of 2π/3 phase lag for attractive coupling (ϵ > 0), and they evolve into a state of complete identical phase and amplitude for a repulsive coupling (ϵ < 0) above a threshold as decided by the MSF. The time evolution of the three junctions is plotted in Figures 3A,B when the junctions emerge in an out-of-phase state and CS state, for ɛ = 0.02 and ɛ = − 0.01, respectively. Both the out-of-phase and CS states are true for N > 3 number of junctions in the ring (results not shown here but checked for N = 4, 7, 8). For further illustration, the dynamics of three junction voltages v are demonstrated in Figure 4. We have already recognized two reasonably broad parameter regions of chimera patterns in the parameter planes (D and E). Now, we demonstrate their exemplary dynamics using the trajectories of the three junctions in a θ − v cylindrical plane and their corresponding spatiotemporal dynamics and how they evolve against varying repulsive coupling ɛ at a fixed α = 1.45. Figures 4A,B represent the CS state for repulsive interactions as confirmed by a plot of the trajectories (red) of three junctions in a θ − v cylindrical plane (left panel) and its corresponding spatiotemporal plot at the right panel. For stronger repulsive interaction, Figures 4C,D show that CS breaks down, although the junctions maintain rotational periodic motion; however, their trajectories are distinctly different (blue, red, and green) which is confirmed by their incoherent motion in the spatiotemporal plot (right panel). The junctions become chaotic for larger repulsive current yet continue with rotational incoherent motion as depicted by their trajectories and spatiotemporal evolution in Figures 4E,F. The first sign of the chimera pattern appears for further increase in repulsive coupling when all the junctions are in chaotic rotational motion, as shown in Figure 4G. However, two of the junctions (indices 1 and 2) are evolving in the coherence of amplitude and phase, while the third junction (index 3) is incoherent in rotational motion as depicted in the spatiotemporal plot in Figure 4H. The collective dynamics of three junctions then switches to periodic motion for a further increase in repulsive coupling. Interestingly, the junctions split into two groups when a pair of junctions remains in coherent rotation (blue), as shown in Figure 4I, and one junction 3) switches to librational motion (red) and incoherent to the other two junctions (1 and 2). This state of coexisting coherence and incoherence in two groups as typically defined as a chimera pattern is confirmed by their spatiotemporal evolution in Figure 4J.

The role of repulsive interactions in the origin of chimera patterns in a small set of three dynamical units in a ring has been explored here. A ring of three nodes with all-to-all coupling is quite often considered as a network motif (Milo et al., 2002). We used the exemplary model of the superconducting Josephson junction as dynamical units that basically represent a second-order phase dynamics under damping and external bias. Noteworthy that this model also represents the dynamics of a simple pendulum and the discrete sine-Gordon equation. For attractive coupling, three junctions maintain an out-of-phase in 2π/3 alternate rotational motion in a cyclic order, which is counterintuitive to the common notion of synchrony under attractive coupling. Interestingly, three junctions emerge into CS for repulsive coupling only. We are still working on this strange collective behavior but yet to find a complete explanation, although we make some inroads into why and when such a behavior may emerge that is to be reported in the future. However, we have checked this new paradigm of CS using the MSF that establishes its stability in a range of repulsive coupling independent of the number of junctions in the ring. Chimera patterns emerge in three junctions beyond the CS state when it breaks down for larger repulsive coupling. We define chimera states here as two coherent junctions and one junction incoherent to the other two that follow the standard definition of chimera states in the literature. We identified two different patterns of chimera states for repulsive interactions. Three junctions may remain in chaotic rotational motion, yet two of them evolve in coherence, while the third junction also evolves in rotational motion but incoherent to the other two junctions. In another kind of chimera pattern, for larger repulsive coupling, two junctions may evolve in coherent rotational motion, while the third junction switches to small-amplitude librational motion and becomes incoherent to the other two junctions. This chimera pattern is sensitive to initial conditions when the chimera pattern may flip to coherent oscillators in coherence, while the third is in incoherent rotational motion. We cannot designate the chimera patterns as strong chimera (Zhang and Motter, 2021). However, some recent reports consider this latter chimera pattern as a solitary state (Jaros et al., 2018; Semenova et al., 2018; Rybalova et al., 2021) since a single oscillator behaves differently, in the dynamical sense, from the other two coherent oscillators.