Dynamics of dissipative topological defects in coupled phase oscillators

TL;DR

This paper investigates how dissipative topological defects evolve in coupled phase oscillators, revealing that defect dynamics depend on coupling rates and influence system coherence, with implications for understanding phase transitions.

Contribution

It introduces a numerical study of defect dynamics in coupled oscillators, highlighting the impact of coupling rate on defect number and system order, connecting to the Kibble-Zurek mechanism.

Findings

Faster coupling rates lead to fewer topological defects.

Defect number scales with the coupling rate, similar to KZM.

Reducing defects enhances long-term coherence and system order.

Abstract

The dynamics of dissipative topological defects in a system of coupled phase oscillators, arranged in one and two-dimensional arrays, is numerically investigated using the Kuramoto model. After an initial rapid decay of the number of topological defects, due to vortex-anti-vortex annihilation, we identify a long-time (quasi) steady state where the number of defects is nearly constant. We find that the number of topological defects at long times is significantly smaller when the coupling between the oscillators is increased at a finite rate rather than suddenly turned on. Moreover, the number of topological defects scales with the coupling rate, analogous to the cooling rate in KibbleZurek mechanism (KZM). Similar to the KZM, the dynamics of topological defects is governed by two competing time scales: the dissipation rate and the coupling rate. Reducing the number of topological defects improves the long time coherence and order parameter of the system and enhances its probability to reach a global minimal loss state that can be mapped to the ground state of a classical XY spin Hamiltonian.

Equations

$$\frac{d{\varphi }_m}{d(t/{\tau }_p)}={\Omega }_m{\tau }_p+\sum _{n=(m{)}_{NN}}{K}_{mn}\sin ({\varphi }_n-{\varphi }_m),$$

$$r=\frac{1}{M}\left|\sum _{m=1}^M{e}^{i{\varphi }_m}\right|,$$

$$K(T)=\{\begin{array}{ll}CT,& \text{if}CT\le K_{max}\\ K_{max},& \text{otherwise}\end{array}$$

$${\Omega }_{max}(T)=\{\begin{array}{ll}\left[1-CT\right]{\tau }_p{\Omega }_{max},& \text{if}\left[1-CT\right]\ge 0\\ 0,& \text{otherwise.}\end{array}$$

$$D=\frac{1}{2\pi }\sum _{m=1}^M\left\{{\varphi }_{m+1}-{\varphi }_m\right\},$$

$$D_{{}_{ij}}=$$

Full text

I Introduction

Topological defects occur in many fields including atomic and molecular physics, condensed matter, cosmology, cold atoms, fluids mechanics, spin systems and optics [1, 2, 3, 4, 5, 6]. They appear across a phase transition and limit the order of systems. Their origin and behavior were first characterized by the Kibble-Zurek mechanism (KZM) as a system quenched across a phase transition into an ordered state via competing time scales [7, 8, 9]. In the KZM, the number of topological defects obeys a power law as a function of the quench rate. This power law was related to the critical exponent associated to the correlation length and the relaxation time of the system in a universal manner [8, 10].

The KZM was investigated in diverse environments: atomic gases [11, 12, 13], condensed matter systems [9, 14] and nonlinear optics [15]. Observing the KZM experimentally is subject to strong limitations, the cooling rate must be slow enough to keep the system homogeneous but not too slow to prevent finite size effect [11]. In particular for systems cooled from the outside this implies limitations into the interesting regime of fast cooling, where the system fails to maintain cooling uniformity [12, 14, 16, 17, 11].

Recently, dissipative topological defects were numerically and experimentally investigated using a one-dimensional ring array of phased-locked lasers, where their formation was shown to be related to the KZM in a universal manner by two competing time scales [18]. The ratio between these two time scales depends on the system parameters, with which it is possible to enable the system to dissipate to a fully ordered, defect-free state that can be exploited for solving computational problems in various fields [19]. The dissipative topological defects were investigated with a fixed coupling between the lasers but intensity fluctuations acted as an effective heat bath with controllable cooling rate. The lasers platform is essentially a system of coupled oscillators where the formation of topological defects occurs without external cooling. They reach a stable phase-locked state due to a dissipative coupling that minimizes the loss of the system, and can be mapped to the ground state of a classical XY spin Hamiltonian [19, 20, 21]. With a complex landscape of many local minima, the system does not reach the ground state solution and remains stuck in local minima. For the ring geometry, these local minima have a nonzero topological charge (winding number) and were defined as dissipative topological defects [18, 22, 23]. Under the assumption of constant field amplitudes, such coupled lasers are well approximated as Kuramoto phase oscillators [24, 21].

In this paper, we investigated the dynamics of topological defects in one and two-dimensional arrays of coupled oscillators, where we varied the coupling strenght in time to correspond to temperature variation in the KZM. In our investigation, we used the Kuramoto model [25], where we assumed that all oscillators have the same amplitude and interact with nearest neighbors only. We swept the coupling in time at different rates to observe the dynamics of dissipative topological defects, and then established their connections with KZM where uncoupled oscillators are equivalent to the infinite temperature state and strongly coupled phase-locked oscillators are equivalent to a zero temperature state. In addition to the variation of the coupling, we showed that introducing disorder at the initial state leads to a more ordered stable state.

II Array of coupled oscillators and the Kuramoto model

The Kuramoto model describes the temporal evolution of the phase of coupled oscillators [25, 26]. It was successfully applied in many areas such as neural networks, complex systems, chemical and biological oscillators [25, 24] and recently in coupled lasers for observing topological defects [18, 27]. The model is valid when all oscillators are nearly identical and the coupling between them is weak.

Our oscillators are coupled to their nearest neighbors ($NN$) only [18]. Their phases evolve in time as a function of the detuning and the coupling strength between the oscillators:

[TABLE]

where $\phi_{m}$ is the phase of the $m^{th}$ oscillator, $\Omega_{m}$ is the frequency detuning from a common frequency $\Omega_{0}=0$, $\tau_{p}$ is the lifetime between two iterations, $K_{mn}$ is the coupling between oscillators $m$ and $n$ and the sum occurs only on the nearest neighbors of the $m^{th}$ oscillator.

Equation (1) is solved using the fourth-order Runge-Kutta method. The time is expressed in unit of iterations $T=t/\tau_{p}$ and the initial phases of the oscillators were chosen randomly between $\left[-\pi~{}\mathrm{to}~{}\pi\right]$. Using Eq. (1), we determined the dynamics of dissipative topological defects in a square and ring array of coupled phase oscillators, we show the relation between time dependent coupling strength and KZM, and we introduce disorder into the system. The results are presented in section III to VI.

III Dissipative Topological Defects in a 2-D square array

We investigated the dynamics of dissipative topological defects in a square array of $M=900$ oscillators, using Eq. (1) to calculate the phases evolution. A topological defect is defined as either a vortex or an anti-vortex, identified by using a defect number $D$, where $D=1$ corresponds to a vortex and $D=-1$ to an anti-vortex (see Appendix A). We assumed that there were no detuning between the oscillators ($\Omega_{m}=0$), and that the coupling between the oscillators was kept constant $K_{mn}=0.1$. The results are presented in Fig. 1.

Figure 1(a) shows, for periodic boundary conditions, the simulated time evolution (iteration unit $T$) of the phases of the oscillators (second column) and defects (third column). As evident, the phases become smoother with time due to the dissipative dynamics and the number of phase singularities which connect the smooth regions decrease. The dissipative dynamics involve vortex$-$vortex, anti-vortex$-$anti-vortex, and vortex$-$anti-vortex interactions, where same charges repel each other, and opposite charges attract each other. Accordingly, vortex$-$anti-vortex attraction forms coupled pairs (surrounded by the green dashed circles in the third column of Fig. 1 that eventually annihilate. This annihilation process follows a well-defined potential, with logarithmic scaling [28]. After a long time ($T>95000$), few topological defects remain and the large distance between them critically slows down the anihilation process, leading to a quasi steady-state (stable state). The average distance between the topological defects can be then related to a correlation length thereby indicating the coherence of the system [11]. As the number of topological defects decrease, the distance between a pair and the correlation length of the system increase. Figure 1(b) shows the simulated results for open boundary conditions. The same dissipative dynamics as for the periodic boundary conditions is observed. However, some topological defects dissipate into the boundaries. These are shown by the yellow star circles in the last column of Fig. 1. As the system size is larger, the difference between periodic and open boundary conditions becomes smaller. In the thermodynamic limit (and also for a square array of $40000$ oscillators that we study below), there is no difference between the periodic and open boundary conditions.

Figure 2 shows a quantitative time evolution of the number of topological defects for both periodic and open boundary conditions. As evident, the number of topological defects is large in the initial state ($T=1$), and decreases in time due to the dissipative dynamics, until only few topological defects remain in the stable state ($T=95000$). For the periodic boundary conditions, Fig. 2(a), both the number of vortices and anti-vortices curves overlap in time, indicating that they can only dissipate by annihilating each other. For the open boundary conditions, Fig. 2(b), the number of vortices differs slightly from the number of anti-vortices indicating that topological defects can also dissipate into the boundaries.

To quantify the ordering in the system, we calculated the order parameter [25, 26]:

[TABLE]

where the sum occurs over all $M$ oscillators. Equation (2) represents the analog of thermal averaging for a classical spin-system [21]. The order parameter $r$ is in the range $\left[0~{}\mathrm{to}~{}1\right]$. When $r=1$, all the oscillators have the same phase, i.e. a fully ordered state. Topological defect reduces the order parameter of the system.

Using Eq. (2), we calculated the evolution in time of the order parameter in the square array of coupled oscillators shown in Fig. 1 with periodic and open boundary conditions. The results are presented in Fig. 3. As shown, the order parameter increases from a low value (initial state) to a finite value (stable state). The higher stable state order parameter for the open boundaries ($0.4$ versus $0.15$) is consistent with the smaller number of topological defect in Fig. 2, due to the boundary dissipation.

IV Time-dependent coupling strength

Next, we introduced time varying coupling (coupling rate) as an analogue of cooling rate in KZM [29]. For non-zero detuning, below a certain coupling strength the oscillators remain unsynchronized [30, 31, 32], in a disorder state. For coupling higher than a characteristic critical coupling, the oscillators phase-lock, to form an ordered state. Thus, varying coupling with time across the transition from disordered state to an ordered state, provides a more direct analogy to KZM [9]. Time varying coupling can be given as:

[TABLE]

where $C$ is the quench rate, $K_{max}=2.5$ is the maximal coupling (arbitrarily chosen). An infinite quench rate corresponds to a sudden jump in the coupling strength.

Using Eq.(1), we investigated the dynamics of dissipative topological defects and KZM for a square array of $M=40000$ oscillators with nearest neighbors coupling. The detuning between the oscillators was randomly chosen in the range $[0~{}\mathrm{to}~{}\Omega_{max}\approx\tau_{p}^{-1}\pi/8]$ and the coupling was increased in time according to Eq. (3). For this large system size, periodic and open boundary conditions give the same results.

The results for the periodic boundary conditions are presented in Fig.4. Figure  4(a) shows the number of topological defects as a function of time. At very high quench rate $C=80\times 10^{-4}$, the dynamics of the number of topological defects is similar to those of constant coupling (Fig. 2). Such high quench rate corresponds to an early sudden jump of the coupling (see Fig. 10(a) of Appendix A). As the quench rate is lowered, the number of topological defects starts to decrease later but eventually reaches a lower value at the stable state. This long time crossing is also evident in the corresponding time evolution of order parameter, as shown in Fig. 4(b).

In the stable state, the observation of increased topological defects with the increase in coupling rate is analogous to the KZM, where there are more defects at higher cooling rate [9]. We verified these observations for different system size (number of oscillators), and found a similar behaviour (see Appendix B), where the initial number of topological defects and the number of topological defects at the stable state increase with the system size (see Appendix C).

We also determined the density of defects at the stable state as a function of the quench rate, where the density of defects is defined as the number of topological defects divided by the total number of oscillators. The results are presented in Fig. 5. As evident, the density of defects at the stable state saturates beyond $C\approx 5\times 10^{-4}$. As the quench rate decreases, the density of defects is reduced.

Next, we calculated the density of defects as a function of time for different maximal detuning $\Omega_{max}$. The results are presented in Fig. 6. As evident, when $\Omega_{max}$ is very small, the system converges to the same number of topological defects for all quench rates. As $\Omega_{max}$ increases, there is more inversion of the two competing time scales and less topological defects at long times, indicating an increase of the ordering of the system in stable state. Such effect was demonstrated in a disorder-induced ordering phenomena [33, 34, 35].

V Time-dependent detuning

In this section, we exchange the role of the coupling strength and the detuning between the oscillators in the KZM. Specifically, we use the maximal detuning $\Omega_{max}$ as the control parameter (linearly varied in time) and let the coupling strength between the oscillators be constant in time. Specifically:

[TABLE]

After incorporating Eq. (4) into Eq. (1), we calculated the density of defects as a function of time for different quench rates of $\Omega_{max}$ and (constant) coupling strengths $K$. The results are presented in Fig. 7. As evident, when $K$ is small, the system converges to the same number of topological defects for all quench rates. As $K$ increases, there is more inversion of the two competing time scales and less topological defects at long times, indicating an increase of the ordering of the system in stable state. Comparing the results of Fig. 6 and Fig. 7 confirm that quenching the coupling strength in a fixed spread of detuning is equivalent to quenching the spread of detuning in a fixed coupling strength. This reflects the close relation between maximal detuning and coupling strength in the Kuramoto model.

VI Dissipative topological defects in a one-dimensional ring

In this section, we consider a one-dimensional ring array geometry of $M=50$ coupled oscillators, i.e. a chain with periodic boundary conditions. The detuning between the oscillators is randomly initialized in the range of $[0~{}\mathrm{to}~{}\Omega_{max}\approx\tau_{p}^{-1}\frac{\pi}{34}]$. The coupling strength between the oscillators varies linearly in time, Eq. (3), and the coupling occurs over nearest neighbors oscillators only.

For a ring array geometry, the defect number $D$ is globally defined for the whole array as [18, 22]:

[TABLE]

where $\{\}$ wrap angle in $\left[-\pi~{}\mathrm{to}~{}\pi\right]$. $D>0$ corresponds to a vortex and $D<0$ to an anti-vortex. After incorporating Eq. (3) into Eq. (1), we calculated the phases of the $50$ oscillators at the stable state ($T=150000$) for two different realizations, each with a different initial phase distribution, and both with a quench rate of $C=3\times 10^{-6}$. The results are presented in Fig. 8. Two different stable states are observed: a vortex in Fig. 8(a) with $D=1$ and $r\approx 0$; and an in-phase-locked state in Fig. 8(b) with $D=0$ and $r\approx 1$.

The statistics on the number of topological defects is obtained by averaging over $1500$ realizations. The probability of having a topological defect is the ratio of the number of realizations that lead to a topological defect at the stable state over the total number of realizations. The results are presented in Fig. 9. As evident, the probability of having topological defects at the stable state saturates to about $70\%$ at quench rates $C=4\times 10^{-6}$ and above. For lower quench rates, the probability of having topological defects monotonically decreases, in analogy to KZM in one dimension [18]. However, as the quench rate goes to zero, the probability of having a topological defect approaches a finite non-zero value ($\approx 0.15$ for our system). This value grows up for lower values of $\Omega_{max}$ indicating the need of finite detuning spread to suppress topological defects for slow quench rates.

VII Concluding Remarks

We numerically investigated the dynamics of dissipative topological defects in one and two-dimensional arrays of coupled oscillators, by using the Kuramoto model. An analogy with the Kibble-Zurek mechanism was shown, where the dynamics of the topological defects was governed by two competing time scales and their density scaled with the coupling rate. In the short time scale regime, the number of topological defects rapidly decreased in time with increasing coupling strength. However, in the stable state, the number of topological defects reduced with the coupling rate. This reduction depends both on the detuning between the oscillators and coupling rate. For a certain coupling rate, introducing detuning between the oscillators initially in the system lead to a more ordered stable state with less topological defects (disorder-induced ordering).

In the Kibble-Zurek mechanism, the number of topological defects obeys a power law as a function of the quench rate. This power law was related to the critical exponents associated to the correlation length and the relaxation time of the system. Accordingly, we plan to extend our investigations to determine whether the scaling of our system will exhibit any power-law behavior, and so as directly resemble the Kibble-Zurek mechanism. Then, we could analyze the critical exponents of scaling behavior, and find the universality class of coupled phase oscillators with dissipative topological defects.

Acknowledgements.

We wish to acknowledge the Israel Science Foundation and the Israel-U.S. Binational Science foundation for their support.

Appendix A Defect number D for a square array and time evolution of coupling and detuning.

For a square array of coupled oscillators, the defect number $D_{ij}$ is defined for each oscillator of the array, as:

[TABLE]

where the brackets $\{\}$ wrap phases in the range of $\left[-\pi~{}\mathrm{to}~{}\pi\right]$, $i$ the row number and $j$ the column number ($i$ and $j$ in the range of $\left[1~{}\mathrm{to}~{}M\right]$). For a ring array of coupled oscillators, the defect number $D$ is defined by Eq. (5) for the whole array.

We calculated the coupling strength and detuning as a function in time using Eq. (3) and Eq. (4). The results are presented in Fig. 10. Figure 10(a) shows the variation of coupling in time for the different quench rates that were used to calculate the results of Figs. 4, 5 and  6. Figure  10(b) shows the variation of detuning in time for the quench rate that was used to calculate the results of Fig. 7.

Appendix B Dissipation of topological defects as a function of array size

We calculated the number of topological defects as a function of time for four different quench rates and for different sizes of the system (number of oscillators). Aside from the system size, all other parameters were the same as those used for obtaining the results of Fig. 4. The results are presented in Fig. 11, where the number of oscillators was relatively large and varied from $2000$ to $60000$. As depicted in Fig. 11, the smaller the system size, the smaller the number of topological defects. Yet, the inversion of the two competing time scales occurs at the same time for the four different system sizes. Surprisingly, the system size does not govern the inversion of the two time scales.

Appendix C Distribution of the topological defects at T=1

We determined the probability distribution of the number of number of vortices, anti-vortices and topological defects at the initial state (at $T=1$). Since the definition of an anti-vortex is the opposite of a vortex, the probability distribution of the number of vortices must be identical to that of the number of anti-vortices, and proportional to the distribution of the number of topological defects. For this reason, only the distribution of the number of vortices at the initial state was considered. The probability distribution was calculated using $10000$ realizations, where in each realization the number of vortices was calculated for a square array of $M$ oscillators with random phase initialized in the range $\left[-\pi~{}\mathrm{to}~{}\pi\right]$.

We calculated the probability distribution of the number of vortices at the initial state for a square array of $40000$ oscillators and the mean and variance of the Normal fit as a function of the number of oscillators. The results are presented in Fig. 12. Figure  12(a) shows the probability distribution as a function of the number of vortices at the initial state. As evident by the red curve, the probability distribution fit well the Normal distribution. Using the Normal fit, we calculated the corresponding mean $\mu$ and variance $\sigma^{2}$ as a function of the number of oscillators. The results are shown in Fig. 12(b). The mean and the variance of the Normal fit are a linear function of the number of oscillators $M$ [23].