Motion of spiral waves in the Complex Ginzburg-Landau equation M. Aguareles S.J. Chapman January 14, 2009 9 Abstract 0 0 Solutions of the general cubic complex Ginzburg-Landau equation comprising 2 multiple spiral waves are considered. For parameters close to the vortex limit, and n a forasystemofspiralwaveswithwell-separatedcentres,lawsofmotionofthecentres J are found which vary depending on the order of magnitude of the separation of the 4 centres. Inparticular, thedirection oftheinteraction changes fromalongthelineof 1 centres to perpendicular to the line of centres as the separation increases, with the ] strengthoftheinteraction algebraicatsmallseparationsandexponentially smallat S large separations. The corresponding asymptotic wavenumber and frequency are P . determined. These depend on the positions of the centres of the spirals, and so n i evolve slowly as the spirals move. l n [ 1 Introduction 1 v 0 The complex Ginzburg-Landau equation is one of the most-studied nonlinear models in 7 9 physics. It describes on a qualitative level, and in many important cases on a quanti- 1 tative level, a great number of phenomena, from nonlinear waves to second-order phase . 1 transitions, including superconductivity, superfluidity, Bose-Einstein condensation, liquid 0 9 crystals, and string theory [3]. 0 The equation arises as the amplitude equation in the vicinity of a Hopf bifurcation in : v spatially extended systems, and is therefore generic for active media displaying wave pat- i X terns. The simplest examples of such media are chemical oscillations such as the famous r Belousov-Zhabotinsky reaction. More complex examples include thermal convection of a binary fluids [14] and transverse patterns of high intensity light [8]. The general cubic complex Ginzburg-Landau equation is given by ∂Ψ = Ψ−(1+ia) |Ψ|2Ψ+(1+ib)∇2Ψ, (1) ∂t where a and b are real parameters and the complex field Ψ represents the amplitude and phase of the modulations of the oscillatory pattern. Of particular interest are “defect” solutions of (1). These are topologically stable solutions in which Ψ has a single zero, around which the phase of Ψ varies by a non-zero integer multiple of 2π. When a = b these solutions are known as “vortices”, and the 1 constant phase lines are rays emanating from the zero. When a 6= b the defect solutions are known as “spirals”, with the constant phase lines behaving as rotating Archimedean spirals except in the immediate vicinity of the core. It is often convenient to factor out the rotation of the spiral, by writing 1+ωb t′ 1+b2 Ψ = e−iωt ψ, t = , x = x′. r1+ab 1+ωb r1+bω This gives, on dropping the primes, ∂ψ (1−ib) = ∇2ψ +(1−|ψ|2)ψ +iqψ(1−k2 −|ψ|2), (2) ∂t where a−b ω −b q = , q(1−k2) = ; 1+ba 1+bω rotatingsinglespiralwavesarenowstationarysolutionsof (2). Theconstantk isknownas the asymptotic wavenumber, since it is easily shown that at infinity arg(ψ) ∼ nφ±kr. An important property of spiral wave solutions is that the asymptotic wavenumber k is not a freeparameter, butisuniquely determinedbyq [7]. Physical systems corresponding to(2) generally containnot onebut many defects. Such complex patterns may beunderstoodin termsofthepositionofthesedefects. Thusifthemotionofthedefectscanbedetermined, much of the dynamics of the solution can be understood. Defect solutions of (2) behave very differently depending on whether q = 0 (corre- sponding to a = b) or q 6= 0 (corresponding to a 6= b). When q = 0 the wavenumber k = 0, and a great amount is known about the solutions to (2). In particular, in a seminal work, Neu [9] analysed a system of n vortices asymptotically in the limit in which the separation of vortices is much greater than the core radius, using the theory of matched asymptotic expansions. By approximating the solution using near-field or “inner” ex- pansions in the vicinity of each vortex core and matching these to a far-field or “outer” expansion away from vortex cores, Neu derived a law of motion for each vortex in terms of the positions of the others, thus reducing (2) to the solution of 2n ordinary differential equations (for the x- and y-coordinates of each vortex). The interaction between defects in this case is long-range, essentially decaying like r−1 for large r. Neu’s analysis has become the template for the analysis of the motion of a system of singularities in many equations, including more detailed models of superconductivity [10], [4], [5]. As we shall show, the key property of (2) that facilitated Neu’s analysis is the linearity of the far- field equations. Thus the contribution from many vortices in the far field can be obtained by a simple linear superposition, and the motion of vortices is determined through the interaction of this far-field with the individual core solutions. When q > 0 the wavenumber k > 0, and the situation is much more complicated, even for a single defect. Hagan [7] studied single spiral wave solutions of (2) and demonstrated that the asymptotic wavenumber k is uniquely determined by q. Using perturbation techniques he found asymptotic expressions for ψ and the asymptotic wavenumber k as a function of q and the winding number of the spiral n. For small values of q a single defect has a multi-layer structure, with the solution comprising inner, outer and far-field regions in which different approximations hold. The transition from the outer region to the far-field occurs exponentially far (in q) from the centre (at what we shall term the 2 outer core radius), and the asymptotic wavenumber k is correspondingly exponentially small in q. This outer core radius is the radius at which the level phase lines switch from being essentially radial to essentially azimuthal. Thus for non-zero q there is a new lengthscale in the problem, with each spiral core having two lengthscales. In studying the motion of spirals it is no longer enough to say that they are well-separated by comparison to the core radius; now it must be determined whether the separation is large compared to the inner core, but small compared to the outer, or whether the separation is large compared to the outer core radius so that the interaction is truly far-field. When the separation lies between the inner and outer core radii the interaction is algebraic, but when the separation is large compared to the outer core radius the interaction of defects decays exponentially. The fact that the outer equation for the phase of ψ is nonlinear, so that the con- tributions from multiple defects may not simply be added, along with the exponential scaling of the outer variable, explains the difficulty in applying Neu’s techniques to the general case of non-zero q. Thus, despite much work and some partial results [2, 12], the interaction of defects in the case of non-zero q was not completely understood. However, recently in [1] a set of laws of motion for N spirals with unit winding number was derived systematically in the limit 0 ≤ q ≪ 1. The aim of the present work is to give the details of that calculation. We start in Section 2 with reviewing the general asymptotic scheme that determines the spirals’ mobility when q = 0. In Section 3 we inspect the equilibrium solutions to (2) for q > 0 formed by a single spiral, highlighting the existence of two distinguished outer regions where either the azimuthal or the radial components of the phase function dominate. These two sections serve as a template for the analysis of multiple-spirals patterns, which is performed in Sections 4 and 5. In Section 4 we derive a law of motion forspiralswhich areseparatedby distances comparabletotheouter coreradius; since this is a distinguished limit we refer to it as the canonical separation. The interaction at the canonical separation is found to be exponentially small and it takes place in the direction perpendicular to the line of centres of the spirals. In Section 5 spirals are assumed to be separated by distances lying between the inner and outer core radii; we denote this as the near-field separation. In this case the interaction becomes algebraic with a component along the line of the centres. Changingtheseparationofthespiralswhilekeepingtheparameterq fixedisequivalent to varying q at fixed separation. We will see that the near-field separation is needed to interpolate between the canonical separation and the case q = 0. 2 Interaction of vortices in the Ginzburg-Landau equa- tion with real coefficients Without lost of generality and to simplify the calculations we consider equation (2) with b = 0. We return briefly to the case of general b in the appendix. For q = 0 this reads Ψ = Ψ(1−|Ψ|2)+∇2Ψ, (3) t 3 which by writing Ψ = feiχ with f and χ real and separating real and imaginary parts we may write as f = ∇2f −f|∇χ|2 +(1−f2)f, (4) t fχ = 2∇χ∇f +f∇2χ. (5) t We wish to determine the law of motion for well separated vortices following [9]. We assume that the separation is O(ǫ−1), with ǫ ≪ 1. This leads to an “outer region”, in which x is scaled with ǫ−1, and an “inner region” in the vicinity of each vortex. Matching the asymptotic expansions of the solutions in each of these regions leads to the law of motion for each vortex centre. 2.1 Outer region In the outer region we rescale x and t by setting X = ǫx and T = ǫ2µt; here µ is a small parameter (the timescale for vortex motion) which will be determined later. We will find that µ is logarithmic in ǫ. With this rescaling (4), (5) read ǫ2µf = ǫ2(∇2f −f|∇χ|2)+(1−f2)f, (6) T µfχ = 2∇χ·∇f +f∇2χ. (7) T Expanding in powers of ǫ as f ∼ f +ǫ2f +··· , 0 1 χ ∼ χ +ǫ2χ +··· 0 1 we find f = 1, (8) 0 ∂χ µ 0 = ∇2χ (9) 0 ∂T Now expanding χ for small µ as 0 χ ∼ χ +µχ +··· 0 00 01 and substituting into (9) we obtain at leading order ∇2χ = 0, 00 with solution N χ = n φ , 00 j j Xj=1 where φ is the polar angle measured from the centre of the jth vortex X , and n is the j j j winding number (or degree) of the jth vortex. At the next order in µ we find ∂χ ∇2χ = 00, 01 ∂T 4 with solution N 1 dX j χ = − n R logR e · , 01 j j j φj 2 dT Xj=1 where R is the distance from the jth vortex. Continuing to O(µ2) we find j N 1 dX dX χ = n R2logR e · j e · j . 02 8 j j j (cid:18) φj dT (cid:19)(cid:18) rj dT (cid:19) Xj=1 In general we find that χ is O(RmlogR) as R → 0. 0m 2.2 Inner Region We rescale near the ℓth vortex by setting X = X +ǫx to give ℓ dX ǫµ ǫf − ℓ ·∇f = ∇2f −f|∇χ|2 +(1−f2)f, (10) T (cid:18) dT (cid:19) dX ǫµ ǫfχ −f ℓ ·∇χ = 2∇χ∇f +f∇2χ, (11) T (cid:18) dT (cid:19) or equivalently dX ǫµ ǫΨ − ℓ ·∇Ψ = Ψ(1−|Ψ|2)+∇2Ψ. (12) T (cid:18) dT (cid:19) Expanding in powers of ǫ as Ψ ∼ Ψ +ǫΨ +... we find at leading order 0 1 ∇2Ψ +Ψ (1−|Ψ |2) = 0. (13) 0 0 0 This is just the equation for a single static vortex, with solution Ψ = f eiχ0 = f (r)ei(nℓφ+C(T)), (14) 0 0 0 where d2f 1df n2 0 + 0 −f ℓ +f −f3 = 0, dr2 r dr 0r2 0 0 f(0) = 0, f → 1 as r → ∞, and C(T) is determined by matching. It is well known that this equation has a unique increasing monotone solution [6]. Continuing with the expansion we find at first order in ǫ that dX −µ∇Ψ · ℓ = ∇2Ψ +Ψ (1−|Ψ |2)−Ψ (Ψ Ψ∗ +Ψ∗Ψ ). (15) 0 dT 1 1 0 0 0 1 0 1 5 2.2.1 Solvability Condition Since the linear operator on the right-hand side of equation (15) is self-adjoint, and the homogeneous version of (15) is satisfied by the partial derivatives of Ψ (as can be seen 0 by differentiating (13)), by the Fredholm Alternative there is a solvability condition on (15) which we can write as dX ∂Ψ∗ ∂(∇Ψ ·d) − ℜ µ ∇Ψ · ℓ (∇Ψ ·d)∗ dD = ℜ (∇Ψ ·d) 1 −Ψ∗ 0 dl Z (cid:26) (cid:18) 0 dT (cid:19) 0 (cid:27) Z (cid:26) 0 ∂r 1 ∂r (cid:27) D ∂D where d is an arbitrary constant vector and D is an arbitrary region in the plane. Taking D to be a ball of radius r we find, after some calculations, dX r f2 2π ∂χ χ −µπ ℓ ·d s(f′)2 +n2 0 ds = n 1 + 1 e ·ddφ (16) (cid:18) dT (cid:19)Z (cid:18) 0 k s (cid:19) kZ (cid:18) ∂r r (cid:19) φ 0 0 where χ ∼ χ +ǫχ +···. After matching with the outer region to determine χ , equation 0 1 1 (16) will determine vortex velocity dX /dT. ℓ 2.3 Asymptotic matching 2.3.1 Inner limit of the outer We express the leading-order (in ǫ) outer χ in terms of the inner variable x given by 0 X = X +ǫx so that R = ǫr, φ = φ. Then, Taylor-expanding χ , ℓ ℓ ℓ 0 χ ∼ χ +µχ +... 0 00 01 ǫr dX ∼ n φ+G(X )−n logǫre · ℓµ+ǫ∇G(X )·x+O(ǫ2), (17) ℓ ℓ ℓ φ ℓ 2 dT where n dX j j G(X) = n φ −µ |X−X |log|X−X | ·e j j j j φj 2 dT X X j6=ℓ j6=k n dX dX +µ2 j|X−X |2log|X−X | e · j e · j +O(µ3). j j φj rj 8 (cid:18) dT (cid:19)(cid:18) dT (cid:19) X j6=ℓ 2.3.2 Outer limit of the inner The leading-order phase in the inner region is χ = n φ+C(T). 0 ℓ This matches with the first term in (17) providing C(T) = G(X ). ℓ The two-term inner expansion for the phase is χ +ǫχ . This should match with the 0 1 two-term inner expansion of the outer (17). Since µ = O(1/log1/ǫ) and all logarithmic terms need to be matched at the same time, to perform this matching we need to take the full µ expansion of both terms. Fortunately only the expansion of G(X ) and ∇G(X ) ℓ ℓ involve infinitely many terms in µ, and these are evaluated at X and are therefore ℓ independent of x. We see that the expansions match if r dX ℓ χ ∼ −n logǫre · µ+∇G(X )·x (18) 1 ℓ φ ℓ 2 dT as r → ∞, which is good to all orders in µ. 6 2.4 Law of motion We can now use the matching condition (18) in the solvability condition (16) to find a law of motion for each vortex. As r → ∞ the left-hand side of (16) is dX −µπ(n2logr +a)d· ℓ (19) ℓ dT where r f2 a = lim s(f′)2 +n2 0 ds−logr , r→∞(cid:20)Z (cid:18) 0 ℓ s (cid:19) (cid:21) 0 is a constant independent of µ and ǫ. Using (18) in the right-hand side of (16) gives dX dX −µπ(1/2+logr)n2d· ℓ −µn2logǫπd· ℓ +2πd·∇G⊥(X ), (20) ℓ dT ℓ dT ℓ where ⊥ represents rotation of a vector by π/2. Since µ is a small constant, we find that the only way to make the left and right hand sides balance is by taking µ = 1/log(1/ǫ), giving dX dX n2 dX −aµd· ℓ = n2d· ℓ −µ ℓ d· ℓ +2n d·∇G⊥(X ). (21) dT ℓ dT 2 dT ℓ ℓ Since d is arbitrary this can be rearranged to give the law of motion dX 2n ∇G⊥(X ) ℓ ℓ ℓ = − +O(ǫ) (22) dT n2 +µ(a−n2/2) ℓ ℓ where n e µ dX ∇G(X ) = j φj +n log|X −X | j ·e e ℓ |X −X | j2 ℓ j (cid:18) dT rj(cid:19) φj X ℓ j j6=ℓ µ dX −n (1+log|X −X |) j ·e e +O(µ2). j2 ℓ j (cid:18) dT φj(cid:19) rj Notethat the expression (22) is accurate to allorders inµ; theexpansion in µis necessary only to evaluate ∇G. To leading order in µ the law of motion reads dX 2 n e ℓ ∼ j rj . (23) dT n |X −X | ℓ X ℓ j j6=ℓ In particular, for two vortices at positions (X ,0) and (X ,0) with X < X the laws 1 2 1 2 of motion are given by dX n 2 1 2 = ,0 +O(µ), (24) dT n (cid:18)X −X (cid:19) 1 1 2 dX n 2 2 1 = ,0 +O(µ). (25) dT n (cid:18)X −X (cid:19) 2 2 1 The direction of motion is always along the line of centres, with like vortices repelling and opposites attracting. 7 2.5 An alternative matching procedure The analysis above is based follows the method that was used in [9]. However, when q > 0 there is the added complication of three small parameters, µ, ǫ and q, rather than just two, µ and ǫ. This has several implications on the way we compute the asymptotic expansions, and affects the way the matching procedure must be carried out. In the analysis above the key step which determines the law of motion is matching the first-order inner solution χ with the correction to the outer phase due to the other 1 vortices. When q = 0 we are fortunate that we do not need to expand the inner χ in 1 powers of µ, so that we in effect retain all orders of µ while matching in ǫ. For q > 0 we are not so lucky, and we need to expand both the inner and outer solutions in powers of µ to make any progress. Since we then no longer have the full µ expansion of either region, we cannot match them together. We can get around this deficiency by writing down and solving equations for the outer limit of the inner solution. By solving these equations we are able to resum the infinite series in µ that is present in the inner region, and this is exactly what we need when matching. To illustrate this new procedure we apply it here to the q = 0 case, where we check the results against the known solution above. 2.5.1 Outer limit of the leading-order inner Rather than solving the inner equations and writing the solution in terms of the outer variable and expanding, we rewriting the leading-order inner equations in terms of the outer variable R = ǫr to obtain 0 = ǫ2(∇2f −f |∇χ |2)+(1−f2)f , (26) 0 0 0 0 0 0 = ǫ2∇·(f2∇χ ). (27) 0 0 We now expand in powers of ǫ to give the outer limit of the leading-order inner solution as χ ∼ χ +ǫ2χ +··· 0 00 01 f ∼ f +ǫ2f +··· . 0 00 01 b b b b Substituting these expansions into (26), (27) gives 1 f = 1, f = − |∇χ |2, (28) 00 01 00 2 b0 = ∇2χ00,b b (29) with solution χ = n φ. b 00 k b 2.5.2 Outer limit of the first-order inner We write the first-order inner equation in terms of the outer variable to give dX −ǫµ k ·∇f = ǫ2∇2f −ǫ2f |∇χ |2 −2ǫ2f ∇χ ·∇χ +f −3f2f , (30) dT 0 1 1 0 0 0 1 1 0 1 dX −µǫf2 k ·∇χ = ǫ2∇·(f2∇χ )+ǫ2∇·(2f f ∇χ ). (31) 0 dT 0 0 1 0 1 0 8 We now expand in powers of ǫ as χ 10 χ ∼ +χ +··· , 1 11 ǫ b f ∼ f +ǫfb +··· , 1 10 11 to give b b f = 0, 10 f = −∇χ ·∇χ , b11 00 10 dX dX 1 µ ∇2χb10 = −µ b k ·∇bχ00 = −nkµ k · eφ = nk (V1sinφ−V2cosφ), dT dT R R where we hbave written b dX ℓ = (V ,V ). 1 2 dT Thus µRlogR µRlogRdX k χ = n (V sinφ−V cosφ) = −n ·e 10 k 1 2 k φ 2 2 dT plus a homogbeneous solution which comes from matching with the outer. We see that this homogeneous solution is X·∇G(X ), k giving µRlogRdX k χ = −n ·e +X·∇G(X ). (32) 10 k φ k 2 dT This expression is the obuter limit of the full µ expansion of the first-order inner solution. Rewriting in terms of the inner variable gives µrlogǫdX k χ ∼ −n ·e +x·∇G(X )+··· , 1 k φ k 2 dT which is (18) as expected. Note that first term in (32), which is the particular integral, was obtained previously from the O(µ) terms in the outer solution; this time we have used only the leading-order outer solution. Thus this method allows us to match with the inner solution when we only know a few terms in the µ-expansion of the outer solution. 3 Equilibrium spiral wave solutions We now consider equation (2) and analyse equilibrium spiral wave solutions in the limit where the parameter q is small. These equilibria, which correspond to single spirals with arbitrary winding numbers, were studied by Hagan [7] who showed that the asymptotic wavenumber k, and thus the frequency of the corresponding periodic solution to (1), is uniquely determined by the parameter q. We shall recast Hagan’s results in a more systematic asymptotic framework in a way which generalises to the many-spiral solutions considered in 4 and 5. § § With ψ = feiχ and f and χ real, we seek solutions of the form f = f(r) and χ = nφ+ϕ(r) with f(r) and ϕ (r) bounded as r → ∞ and r ϕ (0) = 0, f(r) ∼ Crn as r → 0, r 9 where C is some positive constant, and a subscript denotes partial differentiation. In fact it can be shown [7] that bounded solutions of (2) satisfy f(r) → (1−k2)1/2 and ϕ → −k r as r → ∞. We start by introducing an auxiliary parameter ǫ. For the multi-spiral case ǫ will represent the inverse of the spiral separation; here it represents the inverse of the outer core radius, which will be defined shortly. Rescaling (2) onto this new (outer) lengthscale by setting X = ǫx gives iǫ2α2 0 = ǫ2∇2ψ +(1+iq)(1−|ψ|2)ψ − ψ, (33) q where we have introduced the new parameter α = qk/ǫ. This can now be seen as an eigenvalue problem for α(q) which provides the relationship between k and q. The outer core radius of the spiral is the value of ǫ which makes α of order one as q → 0. 3.1 Outer Region With f = f(r) and χ = nφ+ϕ(r) equation (33) becomes f′ n2 0 = ǫ2 f′′ + −ǫ2f +(ϕ′)2 +(1−f2)f, (34) (cid:18) R(cid:19) (cid:18)R2 (cid:19) ǫ2 ǫ2α2 0 = (Rf2ϕ′)′ +qf2(1−f2)− f2. (35) R q Expanding in powers of ǫ as f ∼ f (X;q)+ǫ2f (X;q)+··· , 0 1 ϕ ∼ ϕ (X;q)+ǫ2ϕ (X;q)+··· , 0 1 we find 1 n2 f = 1, f = − +(ϕ′)2 , 0 1 2 (cid:18)R2 (cid:19) ϕ′ n2 α2 ϕ′′ + 0 +q +(ϕ′)2 − = 0. (36) 0 R (cid:18)R2 0 (cid:19) q Equation (36) is a Riccati equation and can be linearised through the transformation ϕ = (1/q)logH to give 0 0 H′ q2n2 H′′ + 0 +H −α2 = 0, (37) 0 R 0(cid:18) R2 (cid:19) with the general solution H (R) = K (αR) + λI (αR) where λ is an arbitrary real 0 inq inq number, and K and I are the modified Bessel functions of the first and second kind. inq inq Only when λ = 0 is the function ϕ monotone [7], so that the iso-phase contours are spirals. Then 1 χ ∼ nφ+ log(K (αR)). (38) 0 inq q 10