text:
# Stationary modulated-amplitude waves in the 1-𝐷 complex Ginzburg-Landau equation

## Introduction

The cubic complex Ginzburg-Landau equation (CGLe) is a generic amplitude equation describing Hopf bifurcation in spatially extended systems, i.e., $`I_o`$ systems , with reflection symmetry . It is of great interest due to its genericity and applications to onset of wave pattern-forming instabilities  in various physical systems such as fluid dynamics, optics, chemistry and biology. It exhibits rich dynamics and has become a paradigm for the transition to spatio-temporal chaos.

We consider the one-dimensional CGLe for the complex amplitude field $`A(x,t)`$:

$`A_t=\mu A+(1+i\alpha )A_{xx}(1+i\beta )|A|^2A`$ (1)
where $`A(x,t):^2`$, and $`\mu ,\alpha ,\beta `$, $`x𝒟`$. $`𝒟`$ is the spatial domain on which the equation is defined. Interesting domains for us are either the whole real axis or a finite box of length $`L`$ with periodic boundary conditions. $`\mu `$ is the control parameter. Only $`\mu >0`$ is considered because we study the supercritical Ginzburg-Landau equation; one could set $`\mu =1`$ by appropriate rescaling of the time, space and amplitude, but we keep it as a parameter for closer connection with experimental results and previous literature. Coefficients $`\alpha `$ and $`\beta `$ parametrize the linear and nonlinear dispersion.

If both $`\alpha `$ and $`\beta `$ are set to $`0`$, we recover the real Ginzburg-Landau equation (RGLe) in which only the diffusion term and the stabilizing cubic term compete with each other and the linear term. A Lyapunov functional exists in that case  and the RGLe behaves like a gradient system. Another limit — the nonlinear Schrödinger equation — results from setting $`\alpha ,\beta \mathrm{}`$; we then have an integrable nonlinear PDE. For parameter values in the intermediate range, long-time behavior of the CGLe can vary from stationary to periodic and to spatiotemporal chaos . In this paper, we concentrate on the stationary solutions of the CGLe in a finite box of length $`L`$ with periodic boundary conditions, and the case $`\alpha \beta `$. Stationary solutions are the simplest non-trivial solutions, related to propagating solutions by an appropriate change of frame of reference $`(x,t)(xvt,t)`$ with fixed $`v`$.

Searching for coherent structures allows one to reduce a partial differential equation into an ordinary one, and such solutions of the CGLe are believed to be extremely important in many regimes, including the spatiotemporal chaos . Recently, numerical integrations of the CGLe have focused on a class of solutions called modulated-amplitude waves (MAWs) and their role in the nonlinear evolution of the Eckhaus instability of initially homogeneous plane waves .

MAWs can bifurcate from the trivial solution $`A=0`$ (case I) or plane wave solutions of zero wavenumber (case II). Analytical aspects of modulated solutions of the CGLe have been addressed by Newton and Sirovich who have applied a perturbation analysis to study the bifurcation in case II , and discussed the secondary bifurcation of those MAWs . Takáč  proved the existence of MAW solutions using a standard bifurcation analysis in the infinite-dimensional phase space of the CGLe, in both cases I and II, together with a stability analysis in case I by means of the center manifold theorem.

In this article we reformulate the CGLe equation assuming a coherent structure form for the solutions, and obtain a fourth-order ordinary differential equation (ODE) with a consistency condition. This form is algebraically convenient, because the deduced system of four first-order ODEs contains only quadratic non-linearity. In the Benjamin-Feir-Newell regime, where plane waves solutions are always unstable, we give a proof of existence of MAWs in both case I and II using our ODE. For weak perturbations in case I or II, we write approximate analytic solutions in the ODE phase space. Coming back to the full CGLe, we then prove the stability of those MAWs in a finite box in case II, and prove that the bifurcation is supercritical, as suggested by recent numerical work .

In the next section, we discuss symmetries and solutions of the CGLe. In section 2 we transform the steady CGLe for MAWs into an equivalent ODE, and give the sufficient condition to identify the solutions of these two equations. In section 3 this ODE is used to construct a 4-$`D`$ dynamical system and prove the existence of symmetric stationary solutions of the CGLe in the two cases I and II. In section 4 the approximate analytic form of the solutions is given and compared to numerical calculations, and the stability of MAWs in case II is proved. Several theorems needed in the proofs are reproduced in appendix B.

## 1 Basic properties of the CGLe

### 1.1 Symmetries

The equation (1) is invariant under temporal and spatial translations. Moreover, it is invariant under a global gauge transformation $`AA\mathrm{exp}(i\varphi )`$, where $`\varphi `$, and under $`xx`$ reflection. As a consequence, it preserves parity of A, i.e., if $`A(x,0)=\pm A(x,0)`$, then $`A(x,t)=\pm A(x,t)`$ for any later time $`t>0`$. If $`A(x,t)`$ has no parity, then $`A(x,t)`$ gives another solution.

### 1.2 Stokes solutions and their stability

The global phase invariance implies that the CGLe has nonlinear plane wave solutions of form

$$A(x,t)=R_0\mathrm{exp}(i(qx\omega t)),$$
(2)
where $`R_0^2=\mu q^2`$ is the amplitude squared, $`\omega =\mu \beta +(\alpha \beta )q^2`$ is the frequency, and $`q,q^2\mu `$ is the wavenumber. They are called Stokes solutions  and are parametrized by the wavenumber $`q`$. The two limit cases of interest to us are highlighted on figure 1: a plane wave of wavenumber $`\mu ^{1/2}`$ and of vanishing amplitude (case I), and the wave with zero wavenumber and maximum amplitude $`\mu `$ (case II). In case II, the solution oscillates uniformly in time; we call it the homogeneously oscillating state (HOS).

For the infinite system, the Benjamin-Feir-Newell  criterion states that all plane wave solutions are unstable with respect to long wavelength perturbations (i.e., of wavenumber $`k0`$) if $`1+\alpha \beta <0`$. If $`1+\alpha \beta >0`$, we have to consider the Eckhaus instability criterion; only a band of wavenumbers are stable against long wavelength perturbations (figure 1):

$$q^2<q_\mathrm{E}^2\frac{(1+\alpha \beta )\mu }{3+\alpha \beta +2\beta ^2}.$$
(3)
For a finite periodic system the wavenumbers for both the original states and the perturbations are quantized. These criteria have been reexamined by Matkowsky and Volpert using linear stability analysis .

### 1.3 Coherent structures and MAWs

Coherent structures play a very important role in the study of pattern formation and dynamical properties of the CGLe . They are uniformly propagating structures of the form

$$A(x,t)=R(xvt)\mathrm{exp}(i\varphi (xvt)e^{i\omega t})$$

which can be expressed as solutions of a 3-$`D`$ nonlinear dynamical system obtained by substituting the above ansatz into the CGLe. There are two free parameters: the frequency $`\omega `$ and the group velocity $`v`$.

The fixed points of the 3-$`D`$ system are the plane waves described in the previous section. The homoclinic  and heteroclinic  connections between the fixed points correspond to localized coherent structures. The Nozaki-Bekki solutions  belong to this category; they connect asymptotic plane waves with different wavenumbers. In numerical simulations in large domains, nearly coherent structures are frequently observed in chaotic regimes, thus suggesting those objects are also relevant to spatiotemporally chaotic dynamics.

Recent numerical studies reveal another kind of coherent structure: modulated amplitude waves (MAWs) for the CGLe . They correspond to limit cycles of the 3-$`D`$ nonlinear system. When $`v=0`$, MAWs are stationary. The formation of MAWs is the first instability encountered when a plane wave state crosses the Eckhaus or Benjamin-Feir stability line. The MAW structure is frequently observed in experiments  and considered as a key to interpretation of patterns and bifurcations exhibited during the system’s transition to spatio-temporal chaos . Traveling MAWs have been observed in numerical simulations of the CGLe in periodic boxes, with parameter $`q`$ between 0 and $`\mu ^{1/2}`$, i.e., in between cases I and II; we are interested here only in stationary MAWs that appear either in case I or case II.

In this paper, we propose a new real-valued ODE to describe steady solutions of the CGLe. A 4-$`D`$ dynamical system derived from this ODE enables us to apply the successive approximation method , to prove the existence of stationary MAWs and to give the analytical form of the approximate solutions in both case I and case II. Numerical integrations of the exact CGLe are then compared to the approximate analytic result. Furthermore, we show non-analyticity at discrete points of solutions in case I, and prove the stability of the MAWs in case II. Some theorems needed in our proof are reproduced in the appendix B. In what follows, $`\mathrm{diag}(\mathrm{})`$ denotes a (block) diagonal matrix and $`\mathrm{col}(\mathrm{})`$ a column vector.

## 2 Stationary case

Since we are only interested in the steady solutions of the CGLe, we substitute the ansatz

$$A(x,t)=R(x)\mathrm{exp}(i\varphi (x)i\omega t),(R,\varphi )^2$$
(4)
into (1). We then have

$`(1+\alpha ^2)G_x`$ $`=K`$ $`(\beta \alpha )R^4(\omega \mu \alpha )R^2`$ (5)
$`(1+\alpha ^2)G^2`$ $`=M`$ $`(1+\alpha ^2)R^3R_{xx}+(\alpha \omega +\mu )R^4(1+\alpha \beta )R^6.`$ (6)
where $`G\varphi _xR^2`$ is reminiscent of “angular momentum”. Note that if $`\alpha =\beta `$, this “angular momentum” is conserved — it is constant in space — provided that $`\omega =\mu \alpha `$. In that case, (6) can be solved in terms of elliptic functions . We will only consider the case $`\alpha \beta `$ in the following. Equations (5) and (6) are invariant under $`(G,x)(G,x)`$. Note that for plane waves, $`K=0`$ and $`G`$ is a constant. If $`K`$ is not always zero, differentiating (6) and dividing the result by (5) gives

$$2G=M_x/K,$$
(7)
and by (6)

$`M={\displaystyle \frac{1+\alpha ^2}{4}}{\displaystyle \frac{M_x^2}{K^2}}.`$ (8)
Furthermore, we can factorize $`R^2`$ from $`M_x`$ and $`K`$ and write $`M_x=R^2N`$ and $`K=R^2P`$, where

$`N`$ $``$ $`(1+\alpha ^2){\displaystyle \frac{1}{2}}(R^2)_{xxx}+(\alpha \omega +\mu )2(R^2)_x(1+\alpha \beta )3R^2(R^2)_x`$
$`P`$ $``$ $`(\beta \alpha )R^2(\omega \mu \alpha ).`$ (9)
The last relation can be used to express $`R^2`$ in terms of $`P`$:

$`R^2={\displaystyle \frac{(\omega \mu \alpha )+P}{\beta \alpha }}=e+dP=R_0^2+{\displaystyle \frac{P}{\beta \alpha }},`$ (10)
where $`d1/(\beta \alpha )`$ and $`e(\omega \mu \alpha )/(\beta \alpha )`$.

Note that $`e=R_0^2`$ is the square of the homogeneous amplitude $`R_0(q,\omega )`$ of the Stokes plane wave solution (2) of frequency $`\omega `$ and wavevector $`q(\omega )`$. $`P`$ then appears as the modulation of the amplitude squared with respect to the Stokes solution, and so it is an appropriate variable to describe a MAW.

Substituting $`K`$ and $`M_x`$ into (8), we have

$$\frac{1+\alpha ^2}{4}\frac{N^2}{P^2}=M.$$
(11)
If $`P0`$ (11) is equivalent to (5) and (6). It is easy to check that if we regard (7) as a definition of $`G`$, and use $`K,M,N,P`$ expressed in terms of $`R`$, equation (5) and (6) will be recovered as a result of (8) and (11). Differentiating both sides of (11) results in

$$\frac{1+\alpha ^2}{2}(PN_xNP_x)=R^2P^3.$$
(12)
In this step we have extended the solution set of (11), because as we integrate (12) back, we get

$$\frac{1+\alpha ^2}{4}\frac{N^2}{P^2}=M+𝒞,$$
(13)
where $`𝒞`$ is an integration constant. Only when $`𝒞=0`$, a solution of (12) is a solution of (11). For this reason, when obtaining solutions of (12), we have to check the consistency condition

$$\frac{1+\alpha ^2}{4}\frac{N^2}{P^2}M=0$$
(14)
to make sure that we have a solution of (11), thus a solution of (5) and (6). Note that if $`K`$ vanishes we have to go back to (5) and (6), since in that case (11) is not well defined. Let us rewrite $`N`$ in terms of $`P`$:

$`N={\displaystyle \frac{2}{1+\alpha ^2}}\left(aP_{xxx}+bP_x+cPP_x\right),`$ (15)
where $`a,b,c`$ are constants

$`a`$ $``$ $`{\displaystyle \frac{(1+\alpha ^2)^2}{4(\beta \alpha )}}`$
$`b`$ $``$ $`{\displaystyle \frac{1+\alpha ^2}{2}}\left({\displaystyle \frac{2(\alpha \omega +\mu )}{\beta \alpha }}{\displaystyle \frac{3(1+\alpha \beta )(\omega \mu \alpha )}{(\beta \alpha )^2}}\right)`$ (16)
$`c`$ $``$ $`{\displaystyle \frac{3(1+\alpha \beta )(1+\alpha ^2)}{2(\beta \alpha )^2}}.`$
After some algebra (here relegated to appendix A), we get an equation for $`P`$ only:

$$\left(\frac{\stackrel{~}{M}_x}{P}\right)_x=\frac{\lambda }{a}\stackrel{~}{M}+kP,\stackrel{~}{M}\lambda P_{xx}+dP^2+\stackrel{~}{e}P.$$
(17)
$`\lambda `$ is a fixed real constant that depends on $`\alpha `$ and $`\beta `$ only, and that takes two different values given in appendix A. $`\lambda `$ is a transient variable used in the proof and derivation but our solutions to the CGLe do not depend on $`\lambda `$ and do not distinguish the two values of $`\lambda `$ (see section 4). $`\stackrel{~}{e}`$ and $`k`$ are real parameters introduced as $`\stackrel{~}{e}+\frac{a}{\lambda }k=e`$. So (17) has two free parameters: $`\omega `$, introduced by the ansatz (4) as the carrier frequency of the solution, and $`k`$. the consistency condition (14) fixes one parameter.

## 3 4-$`D`$ dynamical system and the existence of periodic solutions

Let us take $`\tau `$ as the spatial variable, $`P=P(\tau )`$ in (9), and rewrite (17) as a system of first order equations in $`\tau `$. With $`\stackrel{~}{N}=\stackrel{~}{M}_\tau /P`$ and $`Q=P_\tau `$, from (17) we have

$$\{\begin{array}{ccc}\hfill \dot{\stackrel{~}{M}}& =& \stackrel{~}{N}P\hfill \\ \hfill \dot{\stackrel{~}{N}}& =& \frac{\lambda }{a}\stackrel{~}{M}+kP\hfill \\ \hfill \dot{P}& =& Q\hfill \\ \hfill \dot{Q}& =& \frac{1}{\lambda }(\stackrel{~}{M}dP^2\stackrel{~}{e}P)\hfill \end{array},$$
(18)
where the dot represents the derivation with respect to the spatial variable $`\tau `$.

It is easy to check that $`P=0`$ is a solution of the original equations (5) and (6), corresponding to the plane wave solution of the CGLe with frequency $`\omega `$. We will study the behavior near $`P=0`$ and prove the existence of periodic solutions for small $`P`$. In the CGLe, this corresponds to a weakly modulated amplitude wave which bifurcates from a plane wave solution. If $`Pϵ`$, where $`ϵ`$ is a small parameter, so are $`\stackrel{~}{M},\stackrel{~}{N},Q`$ by their definitions. Write

$$(\stackrel{~}{M},\stackrel{~}{N},P,Q)=(ϵx,ϵy,ϵz,ϵw)$$

and set $`k=k_1+ϵk_2`$. Substituting these into the 4-$`D`$ system, we have

$$\frac{d}{d\tau }\left(\begin{array}{c}x\\ y\\ z\\ w\end{array}\right)=A\left(\begin{array}{c}x\\ y\\ z\\ w\end{array}\right)+ϵ\left(\begin{array}{c}yz\\ k_2z\\ 0\\ \frac{d}{\lambda }z^2\end{array}\right),\text{ where }A=\left(\begin{array}{cccc}0& 0& 0& 0\\ \frac{\lambda }{a}& 0& k_1& 0\\ 0& 0& 0& 1\\ \frac{1}{\lambda }& 0& \frac{\stackrel{~}{e}}{\lambda }& 0\end{array}\right).$$

The linear part $`A`$ describes the behavior of the system in the neighborhood of the trivial fixed point $`(0,0,0,0)`$. Note that the system is invariant under $`(t,y,w)(t,y,w)`$. We use this property to simplify our analysis. Moreover, this system defines an incompressible flow since $`\stackrel{}{X}=0`$, where $`\stackrel{}{X}=(x,y,z,w)`$. It follows from (13) that the system has one integration constant $`𝒞`$. This constant induces a foliation of the phase space into three-dimensional manifolds. Physical solutions, i.e., the solutions of the original CGLe, are restricted to $`𝒞=0`$, the manifold that satisfies the consistency condition (14).

These properties strongly restrict the possible distribution of eigenvalues of $`A`$. We restrict our analysis to the case $`\stackrel{~}{e}/\lambda >0`$, then $`A`$ has eigenvalues $`\{0,0,i\omega _1,i\omega _1\}`$ with $`\omega _1=\sqrt{\stackrel{~}{e}/\lambda }`$. In that case, periodic solutions or MAWs can exist as we will prove in the following. The evolution of the system along either of the two degenerate eigenvalue 0 directions respects the constant $`𝒞`$ foliation: if the solution is on a constant $`𝒞`$ manifold at initial time, it remains there for any later time.

We now discuss the condition $`\stackrel{~}{e}/\lambda >0`$ in terms of an instability of the underlying plane wave. We can rewrite $`\stackrel{~}{e}/\lambda `$ using (16) and (52). Assuming that the solution we are searching for is close to a plane wave, we can use the wavenumber $`q`$ instead of the frequency $`\omega `$, using the dispersion relation (2) for plane waves:

$`{\displaystyle \frac{\stackrel{~}{e}}{\lambda }}={\displaystyle \frac{b}{a}}`$ $`=`$ $`{\displaystyle \frac{2}{1+\alpha ^2}}\left[2(\alpha \omega +\mu ){\displaystyle \frac{3(1+\alpha \beta )(\omega \mu \alpha )}{\beta \alpha }}\right]`$
$`=`$ $`{\displaystyle \frac{2}{1+\alpha ^2}}\left[(3+\alpha \beta +2\alpha ^2)q^2(1+\alpha \beta )\mu \right].`$
If we write

$$q_\mathrm{M}^2\frac{(1+\alpha \beta )\mu }{3+\alpha \beta +2\alpha ^2},$$
(20)
we have

$`{\displaystyle \frac{\stackrel{~}{e}}{\lambda }}>0`$ $``$ $`|\begin{array}{cc}\hfill q^2>q_\mathrm{M}^2& \mathrm{if}(1+\alpha \beta )>0\hfill \\ \hfill q^2<q_\mathrm{M}^2& \mathrm{if}(1+\alpha \beta )<2(1+\alpha ^2)\hfill \\ \hfill q[\sqrt{\mu },\sqrt{\mu }]& \mathrm{if}2(1+\alpha ^2)<1+\alpha \beta <0\hfill \end{array}.`$ (24)
The corresponding regions are illustrated on Fig. 2. Note that $`q_\mathrm{M}(\alpha ,\beta ,\mu )=q_\mathrm{E}(\beta ,\alpha ,\mu )`$. If $`|\alpha |=|\beta |`$, the positivity of $`\stackrel{~}{e}/\lambda `$ is assured when the corresponding plane wave is Eckhaus unstable. If $`|\alpha ||\beta |`$, the positivity does not coincide anymore with the Eckhaus criterion; this is not surprising considering that we do not restrict our analysis to long wavelength perturbations of plane waves, but that the solutions we are seeking may have any wavenumber.

In the following we distinguish two cases. In the first case eigenvalue $`0`$ has a simple elementary divisor, i.e., has two distinct eigenvectors citejhale. This coincides with case I: the MAW solution bifurcates from the $`A=0`$ state, with $`\omega \mu \alpha `$ and hence $`\stackrel{~}{e}/\lambda 4\mu >0`$, for $`\mu >0`$. In the second case, eigenvalue $`0`$ has only one eigenvector. This coincides with case II: the MAW is superimposed over a plane wave with $`\omega \mu \beta `$, so $`q0`$, and

$$\frac{\stackrel{~}{e}}{\lambda }\frac{2\mu (1+\alpha \beta )}{1+\alpha ^2}>\mathrm{\hspace{0.17em}0},$$

The positivity is insured if the system is Benjamin-Feir-Newell unstable, $`(1+\alpha \beta )<0`$.

In terms of $`\stackrel{~}{M},\stackrel{~}{N},P,Q`$, the consistency condition (14) can be written as

$$(1+\alpha ^2)M=\left(\frac{a}{\lambda }\stackrel{~}{N}\lambda Q\right)^2$$
(25)
where in new variables

$`M`$ $`=`$ $`{\displaystyle \frac{d(1+\alpha ^2)}{2\lambda }}(dP+e)(\stackrel{~}{M}dP^2\stackrel{~}{e}P){\displaystyle \frac{d^2(1+\alpha ^2)}{4}}Q^2`$
$`(\alpha \omega +\mu )(dP+e)^2(1+\alpha \beta )(dP+e)^3.`$
Recalling (6), we may express $`G`$ by

$$G=\frac{\frac{a}{\lambda }\stackrel{~}{N}\lambda Q}{1+\alpha ^2}.$$
(26)
Here we are allowed to fix the sign of the right hand side expression because of the $`(G,x)(G,x)`$ reflection symmetry of (5) and (6).

### 3.1 Case I

We want the eigenvalue $`0`$ to have non-degenerate eigenvectors, for this, we set

$$\frac{\frac{\lambda }{a}}{\frac{1}{\lambda }}=\frac{k_1}{\frac{\stackrel{~}{e}}{\lambda }},\text{ i.e., }k_1=\frac{\lambda \stackrel{~}{e}}{a}$$

Consequently, we have

$$e=\stackrel{~}{e}+\frac{a}{\lambda }k=\frac{ϵa}{\lambda }k_2.$$
(27)
Notice that $`e0`$ to the zeroth order, so $`R_00`$ and $`\omega \mu \alpha `$, which means that the solution to be considered bifurcates from the zero solution $`A=0`$, corresponding to a plane wave around the marginal stability curve, with wavenumber $`q\pm \mu ^{1/2}`$. This solution is therefore outside the Eckhaus stability region when $`1+\alpha \beta >0`$.

The four eigenvectors of $`A`$ are:

$$\left(\begin{array}{c}0\\ 1\\ 0\\ 0\end{array}\right)\left(\begin{array}{c}\stackrel{~}{e}\\ 0\\ 1\\ 0\end{array}\right)\left(\begin{array}{c}0\\ ik_1\omega _1^1\\ 1\\ i\omega _1\end{array}\right)\left(\begin{array}{c}0\\ ik_1\omega _1^1\\ 1\\ i\omega _1\end{array}\right).$$

Let

$$D=\left(\begin{array}{cccc}\hfill 0& \stackrel{~}{e}& 0& 0\hfill \\ \hfill 1& 0& 0& \frac{\lambda ^2}{a}\hfill \\ \hfill 0& 1& 1& 0\hfill \\ \hfill 0& 0& 0& 1\hfill \end{array}\right),D^1=\left(\begin{array}{cccc}\hfill 0& 1& 0& a^1\lambda ^2\hfill \\ \hfill \stackrel{~}{e}^1& 0& 0& 0\hfill \\ \hfill \stackrel{~}{e}^1& 0& 1& 0\hfill \\ \hfill 0& 0& 0& 1\hfill \end{array}\right).$$

and $`\stackrel{}{\stackrel{~}{X}}(\stackrel{~}{x},\stackrel{~}{y},\stackrel{~}{z},\stackrel{~}{w})=D^1\stackrel{}{X}`$. The dynamical equations for the new variables become

$$\frac{d}{d\tau }\stackrel{}{\stackrel{~}{X}}=M(\omega _1)\stackrel{}{\stackrel{~}{X}}+ϵ\left(\begin{array}{c}k_2(\stackrel{~}{y}+\stackrel{~}{z})+\frac{\lambda d}{a}(\stackrel{~}{y}+\stackrel{~}{z})^2\\ \frac{1}{\stackrel{~}{e}}(\stackrel{~}{x}+\frac{\lambda ^2}{a}\stackrel{~}{w})(\stackrel{~}{y}+\stackrel{~}{z})\\ \frac{1}{\stackrel{~}{e}}(\stackrel{~}{x}+\frac{\lambda ^2}{a}\stackrel{~}{w})(\stackrel{~}{y}+\stackrel{~}{z})\\ \frac{d}{\lambda }(\stackrel{~}{y}+\stackrel{~}{z})^2\end{array}\right),$$

where

$$M(\omega _1)=D^1AD=\left(\begin{array}{cccc}\hfill 0& 0& 0& 0\hfill \\ \hfill 0& 0& 0& 0\hfill \\ \hfill 0& 0& 0& 1\hfill \\ \hfill 0& 0& \omega _1^2& 0\hfill \end{array}\right).$$

The angular frequency of the solution $`\mathrm{\Omega }`$ should be close to $`\omega _1`$, $`\mathrm{\Omega }^2=\omega _1^2+ϵ\gamma `$, with the shift $`\gamma `$ to be determined later. Next, we change variables to:

$$\{\begin{array}{ccc}\hfill \stackrel{~}{x}& =& x_1\hfill \\ \hfill \stackrel{~}{y}& =& x_2\hfill \\ \hfill \stackrel{~}{z}& =& z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \hfill \\ \hfill \stackrel{~}{w}& =& \mathrm{\Omega }z_1\mathrm{cos}\mathrm{\Omega }\tau \mathrm{\Omega }z_2\mathrm{sin}\mathrm{\Omega }\tau \hfill \end{array}$$
(28)
The 4-$`D`$ system of equations then takes form:

$`\dot{x_1}`$ $`=`$ $`ϵ\left[k_2\left(x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)+{\displaystyle \frac{\lambda d}{a}}\left(x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)^2\right]`$
$`\dot{x_2}`$ $`=`$ $`{\displaystyle \frac{ϵ}{\stackrel{~}{e}}}\left[x_1+{\displaystyle \frac{\mathrm{\Omega }\lambda ^2}{a}}\left(z_1\mathrm{cos}\mathrm{\Omega }\tau z_2\mathrm{sin}\mathrm{\Omega }\tau \right)\right]\left(x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)`$
$`\dot{z_1}`$ $`=`$ $`{\displaystyle \frac{ϵ}{\mathrm{\Omega }}}[{\displaystyle \frac{d}{\lambda }}(x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )^2\mathrm{cos}\mathrm{\Omega }\tau `$
$`+\gamma \left(z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)\mathrm{cos}\mathrm{\Omega }\tau `$
$`{\displaystyle \frac{\mathrm{\Omega }}{\stackrel{~}{e}}}(x_1+{\displaystyle \frac{\mathrm{\Omega }\lambda ^2}{a}}(z_1\mathrm{cos}\mathrm{\Omega }\tau z_2\mathrm{sin}\mathrm{\Omega }\tau ))(x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )\mathrm{sin}\mathrm{\Omega }\tau ]`$
$`\dot{z_2}`$ $`=`$ $`{\displaystyle \frac{ϵ}{\mathrm{\Omega }}}[{\displaystyle \frac{d}{\lambda }}(x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )^2\mathrm{sin}\mathrm{\Omega }\tau `$
$`\gamma \left(z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)\mathrm{sin}\mathrm{\Omega }\tau `$
$`{\displaystyle \frac{\mathrm{\Omega }}{\stackrel{~}{e}}}(x_1+{\displaystyle \frac{\mathrm{\Omega }\lambda ^2}{a}}(z_1\mathrm{cos}\mathrm{\Omega }\tau z_2\mathrm{sin}\mathrm{\Omega }\tau ))(x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )\mathrm{cos}\mathrm{\Omega }\tau ]`$
The proof of the existence of weak MAWs close to $`P=0`$ relies on a series of theorems from J. Hale’s monograph . We reproduce the relevant theorems in appendix B, and refer to them as the need arises.

Note that the transformation $`(\tau ,x_1,x_2,z_1,z_2)(\tau ,x_1,x_2,z_1,z_2)`$ leaves the system (3.1) invariant. So, by definition B.1 of appendix B the system has the property $`E`$ with respect to $`Q`$, with

$$Q=\mathrm{diag}(1,1,1,1).$$

As we are interested only in the solutions with definite parity, we may start the iteration with the vector

$$\stackrel{}{X_0}=(0,a_2,0,a_4).$$

According to Theorem B.4, our solution $`z(\tau ,\stackrel{}{X_0},ϵ)`$ has the property

$$Qz(\tau ,\stackrel{}{X_0},ϵ)=z(\tau ,\stackrel{}{X_0},ϵ),$$

which means that our solutions are either symmetric or antisymmetric. According to Theorem B.5, the second and the fourth determining equations are always zero for this starting vector. For the first and the third determining equations, the zeroth order solution of $`\stackrel{}{\stackrel{~}{X}}`$ , i.e. $`\stackrel{}{X_0}`$, may be substituted, and we get

$`k_2a_2+{\displaystyle \frac{\lambda d}{a}}(a_2^2+{\displaystyle \frac{1}{2}}a_4^2)`$ $`=`$ $`0`$ (30)
$`{\displaystyle \frac{\gamma }{2\mathrm{\Omega }}}a_4+{\displaystyle \frac{\lambda ^2\mathrm{\Omega }a_2a_4}{2a\stackrel{~}{e}}}{\displaystyle \frac{d}{\lambda \mathrm{\Omega }}}a_2a_4`$ $`=`$ $`0.`$ (31)
From (31), we have two possibilities: either $`a_4=0`$ or

$$\gamma +a_2\left(\frac{\lambda ^2\mathrm{\Omega }^2}{a\stackrel{~}{e}}\frac{2d}{\lambda }\right)=0.$$
(32)
When $`a_4=0`$, using $`\stackrel{}{X_0}=(0,a_2,0,0)`$ in (3.1) leads to a trivial constant solution. In the following, we consider only the second case (32). We can solve (30) and (32) for $`\gamma `$ and $`a_4`$ and prove that the system (3.1) has periodic solutions. Note that we have three free parameters $`ϵ,a_2,k_2`$. But as we will see further, $`ϵ`$ and $`a_2`$ are always combined as $`ϵa_2`$ in the first approximation controlling the amplitude and the period of the solution, and the combination will therefore be regarded here as one single free parameter. For general periodic solutions, $`a_2`$ can be interpreted as a phase control parameter, i.e., a parameter giving the initial location on the periodic orbit at $`\tau =\tau _0`$. Here, because we only consider symmetric solutions, the translational symmetry of the autonomous system is broken, and that is the reason why $`ϵ`$ and $`a_2`$ combine into a single parameter. The remaining parameter $`k_2`$ can be chosen freely, for example as to satisfy the consistency condition (25), which, when the zeroth order solution is substituted, becomes at order $`(ϵ^2)`$:

$$\frac{d^2}{4}\mathrm{\Omega }^2a_4^2+\mu \left(da_2+\frac{k_2a}{\lambda }\right)^2=0.$$
(33)
At zeroth order, $`\mathrm{\Omega }^2=\omega _1^2=4\mu `$ and $`\stackrel{~}{e}=4\mu \lambda `$. Solving the system of equations (30),(32) and (33), we get

$$\{\begin{array}{ccc}\hfill k_2& =& \frac{3\lambda }{a}da_2\hfill \\ \hfill \gamma & =& \frac{c}{a}a_2\hfill \\ \hfill a_4& =& \pm 2a_2\hfill \end{array}.$$

We can write out the Jacobian for those three equations explicitly:

$$J=\left(\begin{array}{ccc}a_2& 0& \frac{\lambda d}{a}a_4\\ 0& 1& 0\\ \frac{2\mu a}{\lambda }(da_2+\frac{k_2a}{\lambda })& 0& \frac{d^2}{2}\mathrm{\Omega }^2a_4\end{array}\right).$$

The determinant of this Jacobian is

$$detJ=\frac{1}{2}d^2\mathrm{\Omega }^2a_2a_40a_20.$$

We now invoke theorem B.2, reproduced in appendix B, and conclude our proof that system (5) and (6) has periodic solutions near $`P=0`$. We shall give approximate solutions in section 4, and show that in this case they contain defects.

### 3.2 Case II

Eigenvalue 0 has only one eigenvector. In this case, we assume that $`\frac{\lambda }{a}\stackrel{~}{e}+k_10`$ to the zeroth order in $`ϵ`$, so without loss of generality we can choose $`k_2=0`$. Then $`\frac{\lambda }{a}e=\frac{\lambda }{a}\stackrel{~}{e}+k_1`$. Implementing the transformation $`\stackrel{}{X}=D\stackrel{}{\stackrel{~}{X}}`$ with

$$D=\left(\begin{array}{cccc}0& \stackrel{~}{e}& 0& 0\\ 1& 0& 0& k_1\omega _1^2\\ 0& 1& 1& 0\\ 0& 0& 0& 1\end{array}\right)$$

we have

$$\frac{d}{d\tau }\stackrel{}{\stackrel{~}{X}}=M(\omega _1)\stackrel{}{\stackrel{~}{X}}+ϵ\left(\begin{array}{c}\frac{d}{\stackrel{~}{e}}k_1(\stackrel{~}{y}+\stackrel{~}{z})^2\\ \frac{1}{\stackrel{~}{e}}(\stackrel{~}{x}k_1\frac{\lambda }{\stackrel{~}{e}}\stackrel{~}{w})(\stackrel{~}{y}+\stackrel{~}{z})\\ \frac{1}{\stackrel{~}{e}}(\stackrel{~}{x}k_1\frac{\lambda }{\stackrel{~}{e}}\stackrel{~}{w})(\stackrel{~}{y}+\stackrel{~}{z})\\ \frac{d}{\lambda }(\stackrel{~}{y}+\stackrel{~}{z})^2\end{array}\right),$$

where

$$M(\omega _1)=D^1AD=\left(\begin{array}{cccc}\hfill 0& \frac{\lambda e}{a}& 0& 0\hfill \\ \hfill 0& 0& 0& 0\hfill \\ \hfill 0& 0& 0& 1\hfill \\ \hfill 0& 0& \omega _1^2& 0\hfill \end{array}\right).$$

As in case I, let $`\mathrm{\Omega }^2=\omega _1^2+ϵ\gamma `$ and perform the same transformation (28) into variables $`x_1,x_2,z_1,z_2`$. We then obtain a 4-$`D`$ system similar to (3.1). However, in the equation for $`\dot{x_1}`$, there is an $`ϵ`$-free term. In order to use the successive approximation method, further transformations are required. Let $`\rho `$ such that $`\rho ^2=ϵ`$. With the transformation $`x_2\rho x_2,ϵ\rho ^2`$ we recover the standard form

$`\dot{x_1}`$ $`=`$ $`{\displaystyle \frac{\rho \lambda e}{a}}x_2\rho ^2k_1{\displaystyle \frac{d}{\stackrel{~}{e}}}\left(\rho x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)^2`$
$`\dot{x_2}`$ $`=`$ $`{\displaystyle \frac{\rho }{\stackrel{~}{e}}}\left[x_1k_1{\displaystyle \frac{\mathrm{\Omega }\lambda }{\stackrel{~}{e}}}\left(z_1\mathrm{cos}\mathrm{\Omega }\tau z_2\mathrm{sin}\mathrm{\Omega }\tau \right)\right]\left(\rho x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)`$
$`\dot{z_1}`$ $`=`$ $`{\displaystyle \frac{\rho ^2}{\mathrm{\Omega }}}[{\displaystyle \frac{d}{\lambda }}(\rho x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )^2\mathrm{cos}\mathrm{\Omega }\tau `$
$`+\gamma \left(z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)\mathrm{cos}\mathrm{\Omega }\tau `$
$`{\displaystyle \frac{\mathrm{\Omega }}{\stackrel{~}{e}}}(x_1{\displaystyle \frac{\mathrm{\Omega }\lambda }{\stackrel{~}{e}}}k_1(z_1\mathrm{cos}\mathrm{\Omega }\tau z_2\mathrm{sin}\mathrm{\Omega }\tau ))(\rho x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )\mathrm{sin}\mathrm{\Omega }\tau ]`$
$`\dot{z_2}`$ $`=`$ $`{\displaystyle \frac{\rho ^2}{\mathrm{\Omega }}}[{\displaystyle \frac{d}{\lambda }}(\rho x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )^2\mathrm{sin}\mathrm{\Omega }\tau `$
$`\gamma \left(z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau \right)\mathrm{sin}\mathrm{\Omega }\tau `$
$`{\displaystyle \frac{\mathrm{\Omega }}{\stackrel{~}{e}}}(x_1{\displaystyle \frac{\mathrm{\Omega }\lambda }{\stackrel{~}{e}}}k_1(z_1\mathrm{cos}\mathrm{\Omega }\tau z_2\mathrm{sin}\mathrm{\Omega }\tau ))(\rho x_2+z_1\mathrm{sin}\mathrm{\Omega }\tau +z_2\mathrm{cos}\mathrm{\Omega }\tau )\mathrm{cos}\mathrm{\Omega }\tau ],`$
The system (3.2) has the same symmetry as identified in the case I. If we are only interested in solutions with definite parity, we may again start the iteration with $`\stackrel{}{X_0}=(0,a_2,0,a_4)`$. To the second order $`(\rho ^2)`$, the determining equations are:

$`a_2{\displaystyle \frac{\lambda e}{a}}\rho {\displaystyle \frac{da_4^2k_1}{2\stackrel{~}{e}}}+0(\rho ^3)`$ $`=`$ $`0`$ (37)
$`\rho {\displaystyle \frac{\gamma a_4}{2\mathrm{\Omega }}}\rho ^2\left({\displaystyle \frac{da_2a_4}{\lambda \mathrm{\Omega }}}+{\displaystyle \frac{\lambda \mathrm{\Omega }a_2a_4k_1}{2\stackrel{~}{e}^2}}\right)+0(\rho ^3)`$ $`=`$ $`0.`$ (38)
From the second equation we obtain either $`a_4=0`$ (trivial for our purposes, as discussed above) or

$$\gamma \rho a_2\left(\frac{2d}{\lambda }+\frac{\lambda \mathrm{\Omega }^2k_1}{\stackrel{~}{e}^2}\right)+0(\rho ^2)=0.$$
(39)
If we backtrack the transformations made, it is clear that the consistency condition requires that we keep terms up to the fourth order $`(\rho ^4)`$. We found that with the substitution

$$e=\frac{\alpha \omega +\mu }{1+\alpha \beta }+\rho ^2(\rho ^2\omega _3\rho da_2),$$

where $`\omega _3`$ is a new parameter, only the fourth or higher order terms are left in the consistency condition. From the definition $`e=R_0^2=(\omega \mu \alpha )/(\beta \alpha )`$ and the above equation, we get $`\omega \mu \beta `$ and then $`e\mu `$ to the zeroth order. So $`R_0\sqrt{\mu }`$, $`q0`$, which means that this solution bifurcates from the HOS $`A=\sqrt{\mu }\mathrm{exp}(i\omega t)`$. To the leading order ($`\rho ^4`$), we are allowed to use the following substitutions in the consistency condition (25):

$`a_20\omega \mu \beta \mathrm{\Omega }\sqrt{{\displaystyle \frac{2\mu (1+\alpha \beta )}{1+\alpha ^2}}},`$
$`k_1{\displaystyle \frac{\mu \lambda }{a}}\left(1+{\displaystyle \frac{2\lambda (1+\alpha \beta )}{1+\alpha ^2}}\right)\stackrel{~}{e}{\displaystyle \frac{2\mu \lambda (1+\alpha \beta )}{1+\alpha ^2}}.`$ (40)
The resulting equation is of a relatively simple form:

$$a_4^2(\lambda +d^2(1+\alpha \beta )(1+\alpha ^2+\lambda +\lambda \alpha \beta ))+4(1+\alpha \beta )^2\lambda \mu \omega _3=0.$$
(41)
From (37) it follows that $`a_2`$ is of order $`\rho `$, and from (39) that $`\gamma 0(\rho ^2)`$. After a change of variable $`a_2=\rho a_{22}`$ and keeping only the highest order for the equations, we can rewrite (37) and (39) as

$`a_{22}{\displaystyle \frac{\lambda e}{a}}{\displaystyle \frac{k_1da_4^2}{2\stackrel{~}{e}}}`$ $`=`$ $`0`$ (42)
$`\gamma `$ $`=`$ $`0.`$ (43)
For $`e,\stackrel{~}{e},k_1`$ we use the values in (40). From (42), (43) and (41), we can solve for $`a_{22},\gamma ,\omega _3`$. The Jacobian of those equations is

$$J=\left(\begin{array}{ccc}\frac{\lambda e}{a}& 0& 0\\ 0& 1& 0\\ 0& 0& 4(1+\alpha \beta )^2\lambda \mu \end{array}\right),$$

So, $`detJ=4(1+\alpha \beta )^2\lambda ^2\mu e/a0`$ for $`1+\alpha \beta 0`$. According to Theorem B.2, we have proved that equations (5) and (6) possess periodic solutions.

## 4 Analytic form of periodic solutions, stability analysis and numerical tests

We have proved in the preceding section the existence of symmetric periodic solutions in case I and II. In both cases, a small parameter $`ϵ`$ or $`\rho `$ ensures the convergence of successive approximations. However, we did not give a bound on the highest value of this parameter, nor did we show that the solutions which we obtain are the ones observed in numerical simulations. In this section we give the approximate analytical form of periodic solutions. We compare them with direct numerical integration of the CGLe in case II.

The solutions are shown to be independent of $`\lambda `$ to order $`ϵ`$ in case I and to order $`ϵ^2`$ in case II. In addition, these solutions should also satisfy the 3-$`D`$ ODE mentioned in section 1.3 which do not contain $`\lambda `$, so they can be matched with the solutions of the 3-$`D`$ system in a unique way, independent of the value of $`\lambda `$. Hence, we conclude that to all orders the physical solutions are identical for the two values of $`\lambda `$.

The two cases are taken separately. In this section, we reinstate $`x`$ as the spatial variable, $`R=R(x)`$.

### 4.1 Case I

Using (10), (26) and the case I calculations of the preceding section, we have after some algebra:

$`R^2`$ $`=`$ $`2ϵda_2(1\pm \mathrm{cos}\mathrm{\Omega }x)`$ (44)
$`\varphi _x`$ $`=`$ $`{\displaystyle \frac{ϵa_2}{2(1+\alpha ^2)\mathrm{\Omega }}}{\displaystyle \frac{\mathrm{sin}2\mathrm{\Omega }x\pm 2\mathrm{sin}\mathrm{\Omega }x}{1\pm \mathrm{cos}\mathrm{\Omega }x}}.`$
To the first order of $`ϵ`$, $`R`$ and $`\varphi _x`$ are independent of $`\lambda `$. The $`\pm `$ sign selects two solutions which transform into each other by translating by a half period. This is reminiscent of the spatial translational invariance in the symmetric solution space. From the definitions of $`e,\mathrm{\Omega }`$ and from (27), (3), we get to the first order:

$`\omega `$ $`=`$ $`\mu \alpha 3ϵa_2`$
$`\omega _1^2`$ $`=`$ $`4\mu +{\displaystyle \frac{6ϵda_2}{1+\alpha ^2}}(\alpha \beta +2\alpha ^2+3)`$
$`\mathrm{\Omega }`$ $`=`$ $`\omega _1+{\displaystyle \frac{ϵ\gamma }{2\omega _1}}.`$ (45)
We see that $`\omega `$ and $`\mathrm{\Omega }`$ are independent of $`\lambda `$. On the other hand, for periodic boundary conditions, we can use Fourier modes directly to transform the PDE (1) to a finite set of approximate ODE’s by Galerkin truncation. Then the stationary solution can be obtained by solving a set of nonlinear algebraic equations.

#### Numerical comparison

If we take as an example the following parameter values (previously used in ) for which defect chaos is expected:

$$\alpha =1.5,\beta =1.2$$

and fix the size of the domain to $`L=24`$, then at $`\mu =0.072644`$, $`\omega =0.097879`$, a periodic solution of period $`L/2`$ is found. This solution has $`R_{max}0.0750`$. On the other hand, if we use the same $`\alpha ,\beta ,\mu `$ and search for $`R_{max}0.075`$ by adjusting $`ϵ`$ (we always keep $`a_2=1`$), we find that

$$ϵ0.00380,\omega =0.097566,\text{ period }\frac{2\pi }{\mathrm{\Omega }}=12.0102.$$

The approximate analytic solution and the numerical solution of the exact CGLe agree very well. The profile of $`R`$ from our successive approximation is shown in Fig. 3.

#### Structure near the defect

It is easy to see from (44) that only $`ϵa_2>0`$ is the physically interesting combination. However, we may wonder whether it is really true that $`R^2=dP+e`$ remains non-negative everywhere while touching zero at some points. Fig. 3 and the first equation of (44) suggest a positive answer to this question. But since we have only an approximate solution, further justification is needed. Suppose at some instant $`x_0`$, we have $`dP+e=0`$ on the periodic orbit. From the consistency condition (25), at this transition point

$$\frac{d^2(1+\alpha ^2)}{4}Q^2+\left(\frac{a}{\lambda }\stackrel{~}{N}\lambda Q\right)^2=0,$$

so, $`Q=\stackrel{~}{N}=0`$. According to (18), $`\dot{\stackrel{~}{M}}=0`$ and $`\dot{P}=0`$. Assume that $`\dot{Q}=0`$, then $`\dot{\stackrel{~}{N}}0`$ since the point is not an equilibrium. At next instant $`x_0+\delta x`$, the consistency condition can not be satisfied as the two sides of (25) have different orders of $`\delta x`$. So we conclude that $`\dot{Q}0`$ at the point $`x_0`$, which means that $`Q(x_0+\delta x)`$ has negative sign to that of $`Q(x_0\delta x)`$. Thus, after touching the zero value plane, $`dP+e`$ returns to the positive half space again. The turning happens exactly on the $`dP+e=0`$ plane. We claim that $`dP+e0`$ always holds and the equality holds periodically. From (44), in the neighborhood of $`R=0`$ at $`x=x_0`$ on the periodic orbit, $`R`$ behaves like

$$R\left(\frac{d\dot{Q}}{2}\right)^{1/2}|xx_0|,$$

and is manifestly a non-analytic function of $`x`$.

We do not discuss the stability of the solutions in case I, as this has already been accomplished by Takáč  who has proven that these solutions are unstable.

### 4.2 Case II

To the first order of $`ϵ`$, the solutions are

$$\{\begin{array}{ccc}\hfill x_1& =& \frac{ϵk_1a_4^2}{8\stackrel{~}{e}^2\mathrm{\Omega }}(2d\stackrel{~}{e}+\lambda k_3)\mathrm{sin}2\mathrm{\Omega }x\hfill \\ \hfill x_2& =& ϵ(a_2^2\frac{\lambda k_1a_4^2}{4\stackrel{~}{e}^2}\mathrm{cos}2\mathrm{\Omega }x)\hfill \\ \hfill z_1& =& \frac{ϵa_4^2}{12\stackrel{~}{e}^2\lambda \mathrm{\Omega }^2}(3(3d\stackrel{~}{e}^2+\lambda ^2\mathrm{\Omega }^2k_1)\mathrm{sin}\mathrm{\Omega }x+(d\stackrel{~}{e}^2\lambda ^2\mathrm{\Omega }^2k_1)\mathrm{sin}3\mathrm{\Omega }x)\hfill \\ \hfill z_2& =& a_4+\frac{ϵa_4^2}{12\stackrel{~}{e}^2\lambda \mathrm{\Omega }^2}(3(\lambda ^2\mathrm{\Omega }^2k_1d\stackrel{~}{e}^2)\mathrm{cos}\mathrm{\Omega }x+(\lambda ^2\mathrm{\Omega }^2k_1d\stackrel{~}{e}^2)\mathrm{cos}3\mathrm{\Omega }x),\hfill \end{array}$$

where $`ϵ=\rho ^2>0`$, and $`a_4`$ is a free parameter. In the following, we will see that $`ϵ`$ and $`a_4`$ always emerge in the combination $`ϵa_4`$. To the second order, $`\omega `$ is

$$\omega =\mu \beta +\frac{ϵ^2a_4^2}{4\mu (1+\alpha ^2)}\left(\frac{1+\alpha \beta }{\beta \alpha }+\frac{\beta \alpha }{1+\alpha \beta }\right).$$

It is independent of $`\lambda `$, and therefore $`e,b,\mathrm{\Omega }`$ are also independent of $`\lambda `$. $`R`$ and $`\varphi _x`$ can also be calculated to the second order:

$`R^2`$ $`=`$ $`{\displaystyle \frac{d^2}{2\mu }}(ϵa_4)^2+dϵa_4\mathrm{cos}\mathrm{\Omega }x+{\displaystyle \frac{d(ϵa_4)^2}{12\mathrm{\Omega }^2}}\left({\displaystyle \frac{c}{a}}+{\displaystyle \frac{e}{b}}\right)\mathrm{cos}2\mathrm{\Omega }x+e`$ (46)
$`\varphi _x`$ $`=`$ $`{\displaystyle \frac{ϵa_4}{\mu \mathrm{\Omega }(1+\alpha ^2)}}\left[e\mathrm{sin}\mathrm{\Omega }x{\displaystyle \frac{ϵa_4}{24\mathrm{\Omega }^2}}\left(6d\mathrm{\Omega }^2+{\displaystyle \frac{7e^2}{a\mathrm{\Omega }^2}}+{\displaystyle \frac{7ce}{a}}\right)\mathrm{sin}2\mathrm{\Omega }x\right]`$
So clearly $`R`$ and $`\varphi _x`$ are independent of $`\lambda `$. Similarly, the different signs of $`a_4`$ will give the same solution up to a half-period translation. This solution is the one observed in the numerics when passing the Eckhaus instability for underlying wavevector $`q=0`$. Linear stability analysis reveals  that the $`q=0`$ state, the most stable state under the long wavelength perturbations, becomes unstable when the size of the system is such that the smallest possible nonzero wavenumber $`k`$ satisfies

$$k^2<\frac{2\mu (1+\alpha \beta )}{1+\alpha ^2}\kappa ^2.$$

It is easy to see that $`\kappa ^2=\omega _1^2`$ up to order $`(\rho ^4)`$.

For our parameter choices $`\mu =1,\alpha =1.5,\beta =1.2`$, the bifurcation size of the system is $`L_0=\frac{2\pi }{\kappa }=8.95492`$. In the following, we will first prove the stability of our solutions near the bifurcation point. Then we will compare them with the stable solutions observed in numerics.

#### Stability analysis: presentation

Assume that $`A=R\mathrm{exp}(i\varphi )`$ where $`R,\varphi `$ is an exact solution of (1). The perturbed solution is assumed to be $`\overline{A}=(R+r)\mathrm{exp}(\varphi +\theta )`$, where $`r,\theta `$ is the perturbation on the amplitude and phase, separately. Substitute it into (1), keeping only the linear terms in $`r`$ and $`\theta `$. We have

$`r_t`$ $`=`$ $`(\mu \varphi _x^2\alpha \varphi _{xx}3R^2)r+r_{xx}2\alpha \varphi _xr_x`$ (47)
$`(2R\varphi _x+2\alpha R_x)\theta _x\alpha R\theta _{xx}`$
$`R\theta _t`$ $`=`$ $`(\omega \alpha \varphi _x^2+\varphi _{xx}3\beta R^2)r+\alpha r_{xx}+2\varphi _xr_x`$ (48)
$`+(2R_x2\alpha R\varphi _x)\theta _x+R\theta _{xx},`$
where in (48) we have used $`\varphi _t=\omega `$. To study the stability of the starting solution $`A`$, we treat these equations as an eigenvalue problem for a two components vector, i.e., we let $`r_t=\sigma r`$, $`\theta _t=\sigma \theta `$ and we investigate the spectra $`\sigma `$ of the linear operator resulting from (47) and (48) in the $`C^1`$ continuous periodic function space. As the CGLe has global phase invariance, the eigenvalue equations always have solution $`(r,\theta )=(0,\theta _0)`$ with eigenvalue $`\sigma =0`$. At the same time, spatial translational invariance implies that another eigenmode has $`\sigma =0`$. As a result, saying that the solution is stable means that it is stable up to a phase and a spatial translation, and that all other eigenmodes have eigenvalues with negative real parts.

Invoking the expression for $`R,\varphi _x`$ to the second order of $`ϵ`$, the coefficients of various terms of $`r,\theta `$ and their derivatives in (47) and (48) become explicit functions of $`x`$. The resulting linear operator on $`(r,\theta )`$ has even parity due to the symmetry of our solution, and we can consider the even and odd solutions of $`r,\theta `$ separately. If we set $`ϵ=0`$, i.e., the starting state $`A`$ is a plane wave state, then $`\mathrm{cos}(n\mathrm{\Omega }x)`$ and $`\mathrm{sin}(n\mathrm{\Omega }x)`$ are the eigenfunctions of the unperturbed linear operator. They give the stability spectrum of the plane waves. Now, let us move a little (to the order of $`ϵ`$) beyond the bifurcation point. The eigenfunctions are still $`\mathrm{cos}(n\mathrm{\Omega }x)`$ and $`\mathrm{sin}(n\mathrm{\Omega }x)`$ up to $`ϵ`$ corrections. For example, if the even solutions are considered first, we assume that to the first order the eigenfunctions are (the time dependence for $`r,\theta `$ has been suppressed):

$`r`$ $`=`$ $`m_1\mathrm{cos}(n\mathrm{\Omega }x)+ϵ(m_0\mathrm{cos}((n1)\mathrm{\Omega }x)+m_2\mathrm{cos}((n+1)\mathrm{\Omega }x))`$ (49)
$`\theta `$ $`=`$ $`n_1\mathrm{cos}(n\mathrm{\Omega }x)+ϵ(n_0\mathrm{cos}((n1)\mathrm{\Omega }x))+n_2\mathrm{cos}((n+1)\mathrm{\Omega }x)),`$ (50)
where $`n`$ is a non-negative integer. Note that we do not include the terms such as $`ϵ^2\mathrm{cos}((n\pm 2)\mathrm{\Omega }x)`$ in the above expressions because they induce corrections of order $`ϵ^3`$ or higher in the eigenvalues. Now if we substitute (49) and (50) into the eigenvalue equations and identify the coefficients of $`\mathrm{cos}(n\mathrm{\Omega }x),\mathrm{cos}((n1)\mathrm{\Omega }x)\text{ and }\mathrm{cos}((n+1)\mathrm{\Omega }x)`$, a set of six homogeneous linear equations for $`m_0,m_1,m_2,n_0,n_1,n_2`$ can be derived. The determinant of the coefficient matrix will give an eigenvalue equation for $`\sigma `$. The resulting expression is too complicated to merit being displayed here.

Before bifurcation, the HOS is stable. The first instability occurs for $`n=1`$ mode, one eigenvalue of which is very close to $`0`$ near the bifurcation point, being negative before and positive after. Meanwhile, for $`n>1`$ modes, the corresponding eigenvalues have negative real parts bounded away from zero. As the bifurcating solution emerges continuously from the HOS, near the bifurcation point $`(ϵ1)`$ the perturbed linear operator has all the eigenvalues with negative real parts away from $`0`$ for $`n>1`$ and one eigenvalue close to $`0`$ for $`n=1`$. So, we only need to check the stability of our solutions for $`n=1`$.

For convenience, we can fix parameters $`\alpha `$ and $`\beta `$ to any values allowed by (24) and perform the above stability analysis of the solution.

#### Stability analysis: numerical checks

The numerical values we used are $`\mu =1.0,\alpha =1.5,\beta =1.2,a_4=1`$. The eigenvalue equation is then

$`7.9860ϵ^2\sigma +(56.42363.394ϵ^2)\sigma ^2+(82.56475.135ϵ^2)\sigma ^3`$
$`+(45.02228.859ϵ^2)\sigma ^4+(10.9234.3059ϵ^2)\sigma ^5`$
$`+(1.00.18864ϵ^2)\sigma ^6`$ $`=`$ $`0.`$
$`\sigma =0`$ corresponds to the neutral mode associated with the global phase invariance. All others solutions have negative real parts. The $`\sigma _{}=0.14154ϵ^2`$ solution is the interesting one. If we use the same parameter values to calculate the stability of the HOS, the eigenvalue equation for $`n=1`$ is

$$ϵ^2(0.211220.26402\sigma )+2.98462\sigma +\sigma ^2=0.$$

To the second order in $`ϵ`$, we have $`\sigma =2.984620.19325ϵ^2`$ or $`\sigma _+=0.07077ϵ^2`$. The later positive eigenvalue indicates that the plane wave solution is not stable. We note that $`2\sigma _+=\sigma _{}`$ to order $`ϵ^2`$ which indicates a supercritical pitchfork bifurcation. We have proved that this equality holds exactly at the bifurcation point for any values of $`\alpha `$ and $`\beta `$, and this justifies the above numerical checks. Under perturbation the HOS will evolve to the modulated amplitude solution given above. When the instability is saturated, the corresponding eigenvalue for the MAW is negative. If we change the sign of $`a_4`$ or use the other value of $`\lambda `$, the eigenvalue does not change, as expected.

If we alternatively consider the odd-parity function space $`\{\mathrm{sin}(n\mathrm{\Omega }x)\}_n`$, we obtain the following eigenvalue equation:

$`(28.211533.6568ϵ^2)\sigma +(27.176321.2703ϵ^2)\sigma ^2+(8.92308`$
$`4.04125ϵ^2)\sigma ^3+(10.19287ϵ^2)\sigma ^4=0.`$
This equation is quartic because for $`n=1`$ only two modes $`\mathrm{sin}\mathrm{\Omega }x`$ and $`\mathrm{sin}2\mathrm{\Omega }x`$ are used. Now $`\sigma =0`$ corresponds to the neutral mode associated with the spatial translation of the CGLe. Other eigenvalues of the equation have negative real parts bounded away from $`0`$.

To summarize, our solution is stable in the whole phase space of the CGLe, up to a phase and a spatial translation.

In ref , B. Janiaud et al. have investigated the stability of traveling waves near the Eckhaus instability in Benjamin-Feir stable regime. They derived a necessary condition for the bifurcation to be supercritical and located the corresponding regions as two strips in the $`\alpha ,\beta `$ parameter space. We have studied the stationary MAWs in the Benjamin-Feir unstable regime and found that the bifurcation from the HOS to MAWs is always supercritical, even when parameter values lay outside of the region given in ref .

In ref , application of the perturbation method to the zeroth order $`(ϵ^0)`$ equation gave nonzero eigenvalue $`\lambda _0=2/\beta `$. This can not be correct since the zeroth order equation just gives the stability of the unstable HOS. Furthermore, in the Galerkin projection, somewhat surprisingly the $`N=1`$ truncation was found to give a better result than the $`N=2`$ truncation. In our case, if we use only the first order expressions for $`R,\varphi _x`$ in (47) and (48), we cannot get the correct eigenvalues even near the bifurcation point, not to mention that it would not be possible to extend the result to the next bifurcation.

#### Comparaison with numerical integration of the CGLe

In our numerical simulations we employed a pseudo-spectral method to evolve equation (1) using $`128`$ modes. For system size $`L<L_0`$, we always recover the HOS ($`q=0`$). For $`L`$ slightly larger than $`L_0`$, however, the solution relaxes to the modulated amplitude solution given irrespective of the initial condition. Figure 4 depicts the stable steady solutions given by the two methods.

## 5 Conclusion

We have reformulated the stationary one-dimensional CGLe in a finite box with periodic boundary conditions as a fourth-order ODE for a variable $`P`$ that can be interpreted as the modulation of the amplitude squared of a plane wave solution. This reformulation enabled us to prove the existence of stationary MAW solutions in the two limit cases corresponding to the bifurcation of the trivial solution $`A=0`$ (case I), or to the bifurcation of the plane wave solution of zero wavenumber (case II), when those solutions are within the Benjamin-Feir-Newell regime, or more generally in a region of instability defined by (24). That region coincides with the Eckhaus domain if $`|\alpha |=|\beta |`$, but it is different otherwise. We proved the stability of MAW solutions for the full CGLe in a finite box with periodic boundary conditions in case II, where a homogeneous plane wave becomes unstable. We tested our analytical results by comparison of numerical integrations of the full CGLe with our approximate analytical solutions.

In case I, unstable periodic hole solutions were shown to exist. This could not be inferred from any phase equation: around the defect point $`A=0`$, the amplitude behaves non-analytically, namely piecewise affinely, and the phase is not defined. In case II we found the symmetric stable solutions observed in the numerical integrations of the CGLe just beyond the bifurcation point, using the box size $`L`$ as the bifurcation parameter. This bifurcation was shown to be always supercritical in the Benjamin-Feir-Newell unstable regime. The MAWs continue to exist when the size $`L`$ is increased.

The analysis of MAWs bifurcating from a plane wave with wavenumber $`0<q<1`$ should be similar to the study of case II. It would be interesting to study the higher order instabilities of MAWs when the system size is increased beyond the region in which our analysis takes place. It has been observed that stationary symmetrical MAWs bifurcate into uniformly-propagating asymmetrical ones via a drift-pitchfork bifurcation. This happens when $`L`$ is increased as a consequence of the growth of the amplitude of the modulation, and the increase of the spectral richness of the MAW solution. Moreover, MAWs are expected to be the building blocks of phase turbulence, and the analytical analysis of their global stability may lead to a characterization of the suspected transition between phase and defect chaos in the CGLe .

Acknowledgments The authors thank Georgia Tech. Center for Nonlinear Science and G.P. Robinson for support. Conversations with J. Lega are gratefully acknowledged.

## Appendix A Derivation of the governing equation

We use (15) and (10) to rewrite (12) using only $`P`$ and its spatial derivatives:

$$P(aP_{xxx}+bP_x+cPP_x)_xP_x(aP_{xxx}+bP_x+cPP_x)=(dP+e)P^3$$
(51)
Note that this equation contains even numbers of derivatives of $`P`$ in each term in parenthesis, and also that the powers of $`P`$ increase while the derivatives decrease. We now rewrite the equation in a form which take advantage of this structure. For example, the following equation is equivalent to (51) for any real $`\lambda `$:

$`P(aP_{xxx}+bP_x+(c+\lambda )PP_x)_xP_x(aP_{xxx}+bP_x+(c+\lambda )PP_x)`$
$`=P^2(\lambda P_{xx}+dP^2+eP),`$
or, put in another form and introducing another real parameter $`k`$:

$`\left({\displaystyle \frac{(aP_{xx}+bP+\frac{c+\lambda }{2}P^2)_x}{P}}\right)_x=\lambda P_{xx}+dP^2+\stackrel{~}{e}P+{\displaystyle \frac{a}{\lambda }}kP,`$
where we have written $`\stackrel{~}{e}+\frac{a}{\lambda }k=e`$.

In this equation, we have three free parameters: besides $`\omega `$, introduced by the ansatz (4) as the carrier frequency of the solution, we have introduced free parameters $`\lambda `$ and $`k`$. We now fix $`\lambda `$ by imposing the condition

$`{\displaystyle \frac{a}{\lambda }}={\displaystyle \frac{b}{\stackrel{~}{e}}}={\displaystyle \frac{c+\lambda }{2d}},`$ (52)
which allows us to write the equation in a more suggestive form:

$$\left(\frac{(\lambda P_{xx}+dP^2+\stackrel{~}{e}P)_x}{P}\right)_x=\frac{\lambda }{a}(\lambda P_{xx}+dP^2+\stackrel{~}{e}P)+kP,$$

the equation (17) that leads to the 4-$`D`$ ODE of section 3.

$`\lambda `$ is determined by (52):

$$\lambda ^2+c\lambda 2ad=0.$$
(53)
The discriminant of (53) is

$`\mathrm{\Delta }`$ $`=`$ $`c^2+8ad`$
$`=`$ $`\left({\displaystyle \frac{1+\alpha ^2}{2(\beta \alpha )}}\right)^2\left({\displaystyle \frac{9(1+\alpha \beta )^2}{(\beta \alpha )^2}}+8\right).`$
So $`\mathrm{\Delta }>0`$ for any real values of $`\alpha `$ and $`\beta `$, and the quadratic equation (53) always has two real roots

$`\lambda ={\displaystyle \frac{3(1+\alpha \beta )(1+\alpha ^2)}{4(\beta \alpha )^2}}\pm {\displaystyle \frac{1+\alpha ^2}{4(\beta \alpha )^2}}\sqrt{9(1+\alpha \beta )^2+8(\beta \alpha )^2}`$ (54)
Note that $`\lambda `$ is a function of $`\alpha `$ and $`\beta `$ only. In some applications , the two values of $`\lambda `$ correspond to two distinct solutions of the CGLe. In our case, $`\lambda `$ is an intermediate variable used in the derivation and the proofs, but our solutions to the CGLe do not distinguish the two values of $`\lambda `$.

## Appendix B Theorems used in the proofs

We use successive approximation method to prove the existence of modulated amplitude waves. Below are listed several theorems from the theory of nonlinear oscillations taken from Hale’s monograph .

Consider the system of equations

$$\dot{z}=Az+ϵZ(\tau ,z,ϵ)$$
(55)
where $`A`$ is a constant matrix, $`ϵ,\tau `$, and $`z,Z^n`$. $`Z`$ is a continuous function of $`\tau ,z,ϵ`$, periodic in $`\tau `$ of period $`T`$. In the following, we only consider the case that $`Z`$ is a smooth function. Without loss of generality, $`A`$ can always be assumed to have the standard form

$$A=\mathrm{diag}(0_p,B),$$

Where $`0_p`$ is a $`p\times p`$ zero matrix and $`B`$ is a constant matrix with the property that the equation $`\dot{y}=By`$ has no nontrivial periodic solution of period $`T`$. Under these settings, if the successive approximation is applied to (55), we have

###### Theorem B.1

Given $`d>b>0`$, there is an $`ϵ_1>0`$ such that for any given constant $`p`$ vector $`a,a<b`$ and real $`ϵ,|ϵ|<ϵ_1`$, there is a unique function

$$z^{}(\tau )=z(\tau ,a,ϵ),\text{with }\underset{𝜏}{sup}z^{}(\tau )<d$$

which has continuous first derivative with respect to $`\tau `$ and satisfies

$$\dot{z^{}}=Az^{}+ϵZ(\tau ,z^{},ϵ)ϵP_0Z(\tau ,z^{},ϵ).$$

Furthermore, $`z(\tau ,a,0)=a^{}`$, $`a^{}=\mathrm{col}(a,0)`$, $`P_0(z^{})=a^{}`$, and $`z(\tau ,a,ϵ)`$ has continuous first derivatives with respect to $`a,ϵ`$.

$`P_0`$ is defined as a projection operator on the Banach space $`S`$ of continuous periodic functions of period T. If $`fS`$, write $`f=\mathrm{col}(g,h)`$ where $`g`$ is a $`p`$ vector and $`h`$ is a $`np`$ vector, then

$$P_0(f)=\mathrm{col}(T^1_0^Tg(t)𝑑t,0)$$

So, $`P_0`$ brings an element $`f`$ in $`S`$ to a constant vector which has the average values of $`g`$ over one period as the first $`p`$ components and zeros as the rest components. The equation satisfied by $`z^{}`$ is different from (55) by a constant vector. By a proper choice of the starting vector $`a`$, we may make this constant vector zero to obtain a solution for the system (55). The mathematical statement is give by the following theorem.

###### Theorem B.2

Let $`z(\tau ,a,ϵ)`$ be the function given by the Theorem B.1 for all $`ab<d,|ϵ|ϵ_1`$. If there exist an $`ϵ_2ϵ_1`$ and a continuous function $`a(ϵ)`$ such that

$$P_0Z(\tau ,z(\tau ,a(ϵ),ϵ),ϵ)=0,\text{with }a(ϵ)b\text{ for }|ϵ|ϵ_2$$
(56)
then $`z(\tau ,a(ϵ),ϵ)`$ is a periodic solution of system (55) for $`ϵϵ_2`$. Conversely, if system (55) has a periodic solution $`\overline{z}(\tau ,ϵ)`$, of period $`T`$, $`\overline{z}(\tau ,ϵ)d,|ϵ|ϵ_2`$, then $`\overline{z}(\tau ,ϵ)=z(\tau ,a(ϵ),ϵ)`$.

Therefore, the existence of a continuous function $`a(ϵ)`$ satisfying (56) is a necessary and sufficient condition for the existence of a periodic solution of system (55) of period $`T`$. As we do not know the exact functional form of the periodic solution, the condition (56) could not be solved explicitly. But by using implicit function theorem, we can show that the substitution into (56) of a proper approximate function of $`z(\tau ,a,ϵ)`$ leads to the existence condition for periodic solutions.

###### Theorem B.3

In the system (55), let

$$Z=\mathrm{col}(X,Y),z=\mathrm{col}(x,y)$$

where $`X,x`$ are $`p`$ vectors and define

$$X_0(x,y,ϵ)=\frac{1}{T}_0^TX(\tau ,x,y,ϵ)𝑑\tau .$$

If there is a $`p`$ vector $`a_0,a_0<d`$, such that

$$X_0(a_0,0,0)=0,det\left[\frac{X_0(a_0,0,0)}{x}\right]0$$
(57)
then there exists an $`ϵ_1>0`$ and a periodic function $`z(\tau ,ϵ),|ϵ|ϵ_1`$, of system (55) of period $`T`$ with $`z(\tau ,0)=\mathrm{col}(a_0,0)`$.

If we need to determine other parameters as a function of $`ϵ`$ in practical applications, similar theorems could be derived. Specifically, in the main text we consider the period $`T`$ as a function of $`ϵ`$. It is clear that theorem B.3 applies if we suppose $`T(ϵ)`$ is continuous in $`ϵ`$ and bounded for $`|ϵ|ϵ_1`$. Furthermore, despite the use of the zeroth approximation in the above theorem, the $`n`$th approximation could be used instead. If simple (non-vanishing determinant) solutions to the determining equations can be found for $`ϵ`$ in the neighborhood of 0 then system (55) has a periodic solution.

If the system which we are studying possesses certain symmetries, we can prove the existence of particular symmetric solutions by a simplified version of determining equations. Let us define the symmetry first.

###### Definition B.1

Let $`\dot{z}=f(\tau ,z)`$, where $`z,f^n`$, be a system of differential equations. It is said to have the property $`E`$ with respect to $`Q`$ if there exists a nonsingular matrix $`Q`$ such that

$$Q^2=IQf(\tau ,Qz)=f(\tau ,z)QP_0=P_0Q$$

where $`P_0`$ is the projection operator defined before.

Under this symmetry assumption the following theorems apply:

###### Theorem B.4

Suppose $`Q=\mathrm{diag}(Q_1,Q_2)`$ where $`Q_1`$ is a $`p\times p`$ matrix. If system (55) has property $`E`$ with respect to this $`Q`$ for all $`ϵ`$. If $`a,ab,`$ is a $`p`$ vector and $`a^{}=\mathrm{col}(a,0)`$ is a $`n`$ vector, chosen in such a way that $`Qa^{}=a^{}`$, then the solution $`z(\tau ,a,ϵ)`$ satisfies the relation

$$Qz(\tau ,a,ϵ)=z(\tau ,a,ϵ)$$

and consequently,

$$Z(\tau ,z(\tau ,a,ϵ),ϵ)=Qz(\tau ,z(\tau ,a,ϵ),ϵ)$$

###### Theorem B.5

If the $`j`$-th element of the diagonal of the matrix $`Q_1`$ in Theorem B.4 is $`+1`$, then the $`j`$-th equation in the determining equations is equal to zero for every vector $`a^{}`$ in Theorem B.4.

The system (18) derived here from the 1-$`D`$ CGLe has this symmetry, so the number of determining equations can be reduced using these two theorems.


format: text

text:
# Gluino Pair Production at Linear 𝒆⁺𝒆⁻ Colliders

## I Introduction

Supersymmetry (SUSY) is generally considered to be one of the most promising extensions of the Standard Model (SM) of particle physics. Its attractive features include the cancellation of quadratic divergences in the Higgs sector, which implies that the soft SUSY breaking masses of the (yet unobserved) superpartners of the SM particles can not be much greater than the electroweak scale. If SUSY is indeed responsible for the stabilization of this scale against the Planck scale, supersymmetric particles should therefore be discovered either at Run II of the Fermilab Tevatron Carena:2000yx ; Abel:2000vs ; Culbertson:2000am ; Ambrosanio:1998zf ; Allanach:1999bf  or at the CERN LHC :1999fr ; Abdullin:1998pm . In particular, the strongly coupling squarks and gluinos should be copiously produced at hadron colliders and lead to first measurements of their masses and production cross sections Beenakker:1996ch . Precision measurements of masses, mixings, quantum numbers, and couplings must, however, be performed in the clean environment of a future linear $`e^+e^{}`$ collider because of the large hadronic SM background and theoretical scale and parton density uncertainties at the Fermilab Tevatron and CERN LHC. For example, in $`e^+e^{}`$ annihilation the center-of-mass energy of the collision is exactly known, and threshold energy scans allow for a precise mass determination of pair-produced SUSY particles. It will then be possible to establish whether the masses and couplings of the electroweak gauginos and of the gluino are indeed related, as expected. A global analysis should ultimately lead to a reconstruction of the SUSY breaking model and its parameters. Along these lines, detailed studies have recently been performed for squarks, sleptons, charginos, and neutralinos Aguilar-Saavedra:2001rg , but not for gluinos, the reason being that gluino pairs are produced at $`e^+e^{}`$ colliders only at the one-loop level, while all other sparticles are produced at tree level. At tree level, gluinos can be produced in pairs only in association with two quarks Campbell:1981hf , or they are produced singly in association with a quark and a squark Schiller:kq ; Brandenburg:2002ff . Both processes result in multi-jet final states, where phase space is limited and gluinos may be hard to isolate.

In the Minimal Supersymmetric Standard Model (MSSM) Nilles:1983ge ; Haber:1984rc , the exclusive production of gluino pairs in $`e^+e^{}`$ annihilation is mediated by $`s`$-channel photons and $`Z^0`$-bosons, which couple to the gluinos via triangular quark and squark loops. In earlier studies of this process only the very low center-of-mass energy region ($`\sqrt{s}=20`$ GeV) with pure photon exchange and no squark mixing Nelson:1982cu  or $`Z^0`$-boson decays into light ($`m_{\stackrel{~}{g}}<m_Z/2`$) Kane:xp ; Campbell:kf  and very light ($`m_{\stackrel{~}{g}}=3\mathrm{}\mathrm{\hspace{0.17em}5}`$ GeV) Kileng:1994kw ; Djouadi:1994js  gluinos have been considered. Some authors have presented results only for $`m_{\stackrel{~}{g}}=0`$ GeV and outdated top quark masses of $`20\mathrm{}\mathrm{\hspace{0.17em}50}`$ GeV Kane:xp ; Djouadi:1994js , while others have neglected the mixing of left- and right-handed squark interaction eigenstates into light and heavy mass eigenstates Nelson:1982cu ; Kane:xp ; Campbell:kf , which turns out to control the production cross section to a large extent.

It is the aim of this Article to study the potential of high-energy linear $`e^+e^{}`$ colliders for the production of gluino pairs within the MSSM. Taking into account realistic beam polarization effects, photon and $`Z^0`$-boson exchange, and current mass exclusion limits, we scan the MSSM parameter space for various $`e^+e^{}`$ center-of-mass energies to determine the regions, where gluino production should be visible. Furthermore we clarify the theoretical questions of the relative sign between the two contributing triangular Feynman diagrams, of the possible presence of an axial vector anomaly, and the conditions for vanishing cross sections – three related issues, which have so far been under debate in the literature. The remainder of this Paper is organized as follows: In Sec. II, we present our analytical results and compare them with existing results in the literature. Various numerical cross-sections for gluino pair production at future high-energy linear $`e^+e^{}`$ colliders are computed and discussed in Sec. III, and Sec. IV contains our Conclusions. Our conventions for squark mixing are defined in App. A, and a summary of all relevant Feynman rules is given in App. B.

## II Analytical results

The scattering process

$$e^{}(p_1,\lambda _1)e^+(p_2,\lambda _2)\stackrel{~}{g}(k_1)\stackrel{~}{g}(k_2)$$
(1)
with incoming electron/positron momenta $`p_{1,2}`$ and helicities $`\lambda _{1,2}`$ and outgoing gluino momenta $`k_{1,2}`$ proceeds through the two Feynman diagrams A and B in Fig. 1

with $`s`$-channel photon and $`Z^0`$-boson exchange and triangular quark and squark loops. Higgs boson exchange is not considered due to the negligibly small electron Yukawa coupling, but it could well be relevant at muon colliders. The process occurs only at the one-loop level, since the gluino as the superpartner of the gauge boson of the strong interaction couples neither directly to leptons nor to electroweak gauge bosons. Taking into account chiral squark mixing (see App. A) and using the Feynman rules in App. B, we decompose the corresponding scattering amplitude

$$=\underset{V=\gamma ,Z^0}{}L_\mu ^ViD_V^{\mu \nu }G_\nu ^V$$
(2)
into the lepton current

$$L_\mu ^V=\overline{v}(p_2,\lambda _2)\left[ie\gamma _\mu (v_e^Va_e^V\gamma _5)\right]u(p_1,\lambda _1),$$
(3)
the photon and $`Z^0`$-boson propagators

$$iD_V^{\mu \nu }=\frac{ig^{\mu \nu }}{sm_V^2+i\eta },V=\gamma ,Z^0,$$
(4)
which depend on the squared center-of-mass energy $`s=(p_1+p_2)^2`$, and the gluino current

$`G_\nu ^V`$ $`=`$ $`e\overline{u}(k_2){\displaystyle \underset{q}{}}[{\displaystyle \underset{i}{}}(i\mathrm{\Gamma }_{i,1}^aA_\nu ^{i,V}i\mathrm{\Gamma }_{i,2}^b+i\mathrm{\Gamma }^{}{}_{i,2}{}^{a}\stackrel{~}{A}_{\nu }^{i,V}i\mathrm{\Gamma }^{}{}_{i,1}{}^{b})`$
$`+{\displaystyle \underset{i,j}{}}(i\mathrm{\Gamma }_{i,1}^a\mathrm{\Gamma }_q^{ij,V}B_\nu ^{ij,V}i\mathrm{\Gamma }_{j,2}^b+i\mathrm{\Gamma }^{}{}_{i,2}{}^{a}\mathrm{\Gamma }_{q}^{ji,V}\stackrel{~}{B}_\nu ^{ij,V}i\mathrm{\Gamma }_{j,1}^b)]v(k_1)`$
with squark mass eigenstates $`i,j=\{1,2\}`$,

$$A_\nu ^{i,V}=\frac{\mathrm{d}^Dq}{(2\pi )^D}\mu ^{4D}\frac{(q/+k/_2+m_q)\gamma _\nu (v_q^Va_q^V\gamma _5)(q/k/_1+m_q)}{(q^2m_i^2+i\eta )[(qk_1)^2m_q^2+i\eta ][(q+k_2)^2m_q^2+i\eta ]}$$
(6)
and

$$B_\nu ^{ij,V}=\frac{\mathrm{d}^Dq}{(2\pi )^D}\mu ^{4D}\frac{(q/m_q)(2qk_1+k_2)_\nu }{(q^2m_q^2+i\eta )[(qk_1)^2m_j^2+i\eta ][(q+k_2)^2m_i^2+i\eta ]}$$
(7)
and the quark flavor $`q`$ flowing both ways in the corresponding diagrams A and B, so that $`\stackrel{~}{A}_\nu ^{i,V}=A_\nu ^{i,V}(v_q^Vv_q^V),\stackrel{~}{B}_\nu ^{ij,V}=B_\nu ^{ij,V},`$ and $`\mathrm{\Gamma }^{}=C\mathrm{\Gamma }^\mathrm{T}C^1=\mathrm{\Gamma }`$ Denner:kt .

Eq. (II) can be simplified with the Dirac equation and the anti-commutation relations for Dirac matrices, and the tensor loop integrals in Eqs. (6) and (7) can be expressed through the standard coefficient functions $`C_{k(l)}`$ of the metric tensor $`g_{\mu \nu }`$ and tensors constructed from the outgoing gluino momenta $`k_{1,\mu }`$ and $`k_{2,\mu }`$ Passarino:1978jh ; Denner:kt . The gluino current then reduces to

$$G_\nu ^V=ie\frac{\alpha _s}{2\pi }\frac{\delta ^{ab}}{2}\overline{u}(k_2)\gamma _\nu \gamma _5v(k_1)\underset{q}{}(A_q^V+B_q^V)$$
(8)
with

$`A_q^V`$ $`=`$ $`{\displaystyle \underset{i}{}}[C_0^{qi}(m_q^2a_{qiV}^{}m_{\stackrel{~}{g}}^2a_{qiV}^++2m_qm_{\stackrel{~}{g}}\widehat{a}_{qiV})+C_1^{qi}4m_{\stackrel{~}{g}}(m_q\widehat{a}_{qiV}m_{\stackrel{~}{g}}a_{qiV}^+)`$
$`+C_{00}^{qi}(2D)a_{qiV}^+C_{11}^{qi}2m_{\stackrel{~}{g}}^2a_{qiV}^++C_{12}^{qi}(s2m_{\stackrel{~}{g}}^2)a_{qiV}^+],`$
$`B_q^V`$ $`=`$ $`{\displaystyle \underset{i,j}{}}C_{00}^{qij}2b_{qijV},`$ (10)
where $`C_{k(l)}^{qi}=C_{k(l)}(m_{\stackrel{~}{g}}^2,s,m_{\stackrel{~}{g}}^2,m_{\stackrel{~}{q}_i}^2,m_q^2,m_q^2)`$ and $`C_{00}^{qij}=C_{00}(m_{\stackrel{~}{g}}^2,s,m_{\stackrel{~}{g}}^2,m_q^2,m_{\stackrel{~}{q}_j}^2,m_{\stackrel{~}{q}_i}^2)`$ are massive (infrared-finite) three-point functions and $`C_1^{qi}=C_2^{qi}`$ and $`C_{11}^{qi}=C_{22}^{qi}`$ in diagram A.

$`a_{qiV}^\pm `$ $`=`$ $`v_q^V(S_{i1}^qS_{i1}^qS_{i2}^qS_{i2}^q)\pm a_q^V,`$
$`\widehat{a}_{qiV}`$ $`=`$ $`a_q^V(S_{i1}^qS_{i2}^q+S_{i2}^qS_{i1}^q),\mathrm{and}`$ (11)
$`b_{qijV}`$ $`=`$ $`S_{i1}^qS_{j1}^q\mathrm{\Gamma }_q^{ij,V}S_{i2}^qS_{j2}^q\mathrm{\Gamma }_q^{ji,V}`$
are combinations of vector ($`v_q^V`$), axial vector ($`a_q^V`$), and derivative couplings ($`\mathrm{\Gamma }_q^{ij,V}`$), and elements of the squark mixing matrix $`S`$. Pairs of identical (Majorana) gluinos are therefore produced by a parity violating axial vector coupling induced by mass differences between the chiral squarks and the axial vector coupling $`a_q^Z`$ of the $`Z^0`$-boson. The (mass-independent) ultraviolet singularities contained in the $`C_{00}`$-functions cancel among $`A_q^V`$ and $`B_q^V`$ in $`D=42\epsilon `$ dimensions. As we have checked explicitly (even for complex squark mixing matrices), adding the two amplitudes induces not only a cancellation of the ultraviolet singularities and of the logarithmic dependence on the scale parameter $`\mu `$ introduced in Eqs. (6) and (7) in order to preserve the mass dimension of the loop integrals, but also a destructive interference of the finite remainders. This happens separately for each weak isospin partner, as is to be expected for triangular loop diagrams involving one axial vector and two scalar (not vector) couplings and no closed fermion loop.

The (finite) total cross section for incoming electrons/positrons with helicities $`\lambda _{1,2}=\pm 1/2`$ is then

$$\sigma _{\lambda _1\lambda _2}(s)=\frac{\alpha _e^2\alpha _s^2(N_C^21)\beta ^3s}{24\pi }\underset{V_1,V_2}{}\left[\frac{Q_{\lambda _1\lambda _2}^{V_1V_2}}{(sm_{V_1}^2)(sm_{V_2}^2)}\underset{q}{}(A_q^{V_1}+B_q^{V_1})(A_q^{V_2}+B_q^{V_2})^{}\right]$$
(12)
with

$$Q_{\lambda _1\lambda _2}^{V_1V_2}=(v_e^{V_1}v_e^{V_2}+a_e^{V_1}a_e^{V2})(14\lambda _1\lambda _2)(v_e^{V_1}a_e^{V2}+v_e^{V_2}a_e^{V_1})(2\lambda _12\lambda _2),$$
(13)
color factor $`N_C=3`$, and gluino velocity $`\beta =\sqrt{14m_{\stackrel{~}{g}}^2/s}`$, which contains the expected factors of $`\beta ^3`$ and $`s`$ for $`P`$-wave production of two spin-1/2 Majorana fermions. The distribution in the center-of-mass scattering angle $`\theta `$,

$$\frac{\mathrm{d}\sigma _{\lambda _1\lambda _2}}{\mathrm{d}\mathrm{\Omega }}(s)=\frac{3}{8\pi }(1+\mathrm{cos}^2\theta )\sigma _{\lambda _1\lambda _2}(s),$$
(14)
is independent of the gluino mass and has to be integrated over just one hemisphere, since the two final state particles are identical Kileng:1994kw ; Petcov:1984nf . As a consequence, the forward-backward asymmetry vanishes for Majorana fermions, but not for Dirac fermions.

Our result for diagram A agrees with the unpolarized result

$$\sigma (s)=\frac{1}{4}\underset{\lambda _{1,2}=\pm 1/2}{}\sigma _{\lambda _1\lambda _2}(s)$$
(15)
in Eq. (4.5) of Kileng and Osland, if we identify Kileng:1994kw ; Kileng:2002

$`C_0^{qi}`$ $`=`$ $`F_{qqi}^{00},`$
$`C_1^{qi}`$ $`=`$ $`+F_{qqi}^{01},`$
$`C_{00}^{qi}`$ $`=`$ $`G_{qqi}/2,`$ (16)
$`C_{11}^{qi}`$ $`=`$ $`F_{qqi}^{02},`$
$`C_{12}^{qi}`$ $`=`$ $`F_{qqi}^{11},`$
and reverse the sign of $`\widehat{b}_q`$ to account for opposite conventions of squark mass eigenstates. However, our result for diagram B disagrees in sign with Eq. (4.5) of Kileng and Osland, if we identify Kileng:1994kw ; Kileng:2002

$$C_{00}^{qij}=G_{ijq}/2.$$
(17)
If the sign of diagram B is reversed, the ultraviolet singularities cancel only after adding the contributions from the two weak isospin partners with opposite values of $`T_q^3`$. We trace this sign discrepancy to the Feynman rules employed in Ref. Kileng:1994kw , which exhibit a relative minus sign to those in Ref. Haber:1984rc  for the $`Z^0`$-boson coupling to quarks, but not to squarks, whereas our Feynman rules (see App. B) agree with those in Ref. Haber:1984rc  in the limit of no squark mixing. Except for the relative sign of diagrams A and B and in the limit of vanishing gluino mass, we also find agreement with Djouadi and Drees Djouadi:1994js . We confirm, however, the relative sign for diagrams A and B of Campbell, Scott, and Sundaresan Campbell:kf , who found (for non-mixing chiral squarks of different mass) that the ultraviolet singularities cancel separately for each weak isospin partner, and that there is no anomaly. This had also been claimed previously by Kane and Rolnick Kane:xp  for chiral squarks of equal mass. In their limit, the cross section depends only on the weak isospin and not on the charge of the (s)quarks Kane:xp , and the contribution of the photon vanishes Nelson:1982cu . For the contribution of the $`Z^0`$-boson to vanish, we must have Kileng:1994kw

1. mass degeneracy in each quark isospin doublet, $`m_d=m_u`$ etc.,
2. mass degeneracy in each squark isospin doublet, $`m_{\stackrel{~}{d}_1}=m_{\stackrel{~}{d}_2}=m_{\stackrel{~}{u}_1}=m_{\stackrel{~}{u}_2}`$ etc.,

which contradicts the condition $`m_q=m_{\stackrel{~}{q}}`$ found by Kane and Rolnick Kane:xp . Condition (1) is violated most strongly for the third generation, as is condition (2) for most SUSY breaking models.

## III Numerical Results

The analytical results presented in the previous Section have been obtained in two independent analytical calculations. They have been implemented in compact Fortran computer codes, which depend on the LoopTools/FF library Hahn:1998yk ; vanOldenborgh:1990yc  for the evaluation of the massive tensor three-point functions. As a third independent cross-check, we have recalculated the production of gluino pairs in $`e^+e^{}`$ annihilation with the computer algebra program FeynArts/FormCalc Hahn:2001rv  and found numerical agreement up to 15 digits.

Our calculations involve various masses and couplings of SM particles, for which we use the most up-to-date values from the 2002 Review of the Particle Data Group Hagiwara:pw . In particular, we evaluate the electromagnetic fine structure constant $`\alpha (m_Z)=1/127.934`$ at the mass of the $`Z^0`$-boson, $`m_Z=91.1876`$ GeV, and calculate the weak mixing angle $`\theta _W`$ from the tree-level expression $`\mathrm{sin}^2\theta _W=1m_W^2/m_Z^2`$ with $`m_W=80.423`$ GeV. Among the fermion masses, only the one of the top quark, $`m_t=174.3`$ GeV, plays a significant role due to its large splitting from the bottom quark mass, $`m_b=4.7`$ GeV, while the latter and the charm quark mass, $`m_c=1.5`$ GeV, could have been neglected like those of the three light quarks and of the electron/positron. The strong coupling constant is evaluated at the gluino mass scale from the one-loop expression with five active flavors and $`\mathrm{\Lambda }_{\mathrm{LO}}^{n_f=5}=83.76`$ MeV, corresponding to $`\alpha _s(m_Z)=0.1172`$. A variation of the renormalization scale by a factor of four about the gluino mass results in a cross section uncertainty of about $`\pm 25`$ %. Like the heavy top quark, all SUSY particles have been decoupled from the running of the strong coupling constant.

We work in the framework of the MSSM with conserved $`R`$\- (matter-) parity, which represents the simplest phenomenologically viable model, but which is still sufficiently general to not depend on a specific SUSY breaking mechanism. Models with broken $`R`$-parity are severely restricted by the non-observation of proton decay, which would violate both baryon and lepton number conservation. We do not consider light gluino mass windows, on which the literature has focused so far and which may or may not be excluded from searches at fixed target and collider experiments Hagiwara:pw . Instead, we adopt the current mass limit $`m_{\stackrel{~}{g}}200`$ GeV from the CDF Affolder:2001tc  and D0 Abachi:1995ng  searches in the jets with missing energy channel, relevant for non-mixing squark masses of $`m_{\stackrel{~}{q}}325`$ GeV and $`\mathrm{tan}\beta =3`$. Values for the ratio of the Higgs vacuum expectation values, $`\mathrm{tan}\beta `$, below 2.4 are already excluded by the CERN LEP experiments lhwg:2001xx . If not stated otherwise, we will present unpolarized cross sections for a $`\sqrt{s}=500`$ GeV linear $`e^+e^{}`$ collider like DESY TESLA, gluino masses of $`m_{\stackrel{~}{g}}=200`$ GeV, and squark masses $`m_{\stackrel{~}{q}}m_{\stackrel{~}{Q}}=m_{\stackrel{~}{D}}=m_{\stackrel{~}{U}}=m_{\stackrel{~}{S}}=m_{\stackrel{~}{C}}=m_{\stackrel{~}{B}}=m_{\stackrel{~}{T}}=m_{\mathrm{SUSY}}=325`$ GeV. We will consider two cases of large squark mass splittings: I.) On the one hand, the masses of the superpartners of left- and right-handed quarks need not be equal to each other. In this scenario we will vary the right-handed up-type squark mass parameters $`m_{\stackrel{~}{U},\stackrel{~}{C},\stackrel{~}{T}}`$ between 200 and 1500 GeV. II.) On the other hand, the superpartners of the heavy quarks can mix into light and heavy mass eigenstates (see App. A). This alternative is restricted by the CERN LEP limits on the light top and bottom squark masses, $`m_{\stackrel{~}{t}_1}100`$ GeV and $`m_{\stackrel{~}{b}_1}99`$ GeV lswg:2002xx , and on SUSY one-loop contributions Barbieri:1983wy ; Drees:1990dx ; Chankowski:1993eu  to the $`\rho `$-parameter, $`\rho _{\mathrm{SUSY}}<0.0012`$ Hagiwara:pw . In this case we assume the maximally allowed top squark mixing with $`\theta _{\stackrel{~}{t}}=45.2^{}`$, $`m_{\stackrel{~}{t}_1}=110`$ GeV, and $`m_{\stackrel{~}{t}_2}=506`$ GeV, which can be generated by choosing appropriate values for the Higgs mass parameter, $`\mu =500`$ GeV, and the trilinear top squark coupling, $`A_t=534`$ GeV. For small values of $`\mathrm{tan}\beta `$, mixing in the bottom squark sector remains small, and we take $`\theta _{\stackrel{~}{b}}=0^{}`$. Although the absolute magnitude of the cross section depends strongly on the gluino mass and collider energy, the relative importance of the different contributions is very similar also for higher gluino masses and collider energies.

First we examine the conditions found in Section II for vanishing of the photon and $`Z^0`$-boson contributions, restricting ourselves to the third generation. Since we expect the photon contribution to cancel for equal left- and right-handed squark masses, we vary the right-handed top squark mass parameter, $`m_{\stackrel{~}{T}}m_{\stackrel{~}{t}_R}`$, between 200 and 1500 GeV, but keep $`m_{\stackrel{~}{t}_L}m_{\stackrel{~}{Q}}=m_{\stackrel{~}{B}}=m_{\mathrm{SUSY}}=325`$ GeV fixed (case I), since top and bottom squarks generally interfere destructively due to their opposite charge and weak isospin quantum numbers. As can be seen from Fig. 2, the photon contribution cancels indeed for

$`m_{\stackrel{~}{T}}m_{\stackrel{~}{t}_R}=m_{\stackrel{~}{t}_L}m_{\mathrm{SUSY}}`$. This is due to the fact that for photons $`\widehat{a}_{qi\gamma }=0`$ and $`b_{q12\gamma }=b_{q21\gamma }=0`$ in Eq. (11), while unitarity of the squark mixing matrix leads to $`a_{q1\gamma }^\pm =a_{q2\gamma }^\pm `$ and $`b_{q11\gamma }=b_{q22\gamma }`$. Therefore, the photon contributions cancels for all flavors $`q`$ with equal squark masses. Due to their charge, top (s)quarks contribute four times as much as bottom (s)quarks, whose contribution is even more suppressed by the condition $`m_{\stackrel{~}{b}_L}m_{\stackrel{~}{b}_R}`$. The $`Z^0`$-boson contribution can never cancel, since $`m_tm_b`$, and therefore it depends only weakly on $`m_{\stackrel{~}{T}}`$, but it can become minimal for $`m_{\stackrel{~}{T}}m_{\stackrel{~}{t}_R}=m_{\stackrel{~}{t}_L}=m_{\stackrel{~}{b}_R}=m_{\stackrel{~}{b}_L}m_{\mathrm{SUSY}}`$. As $`m_{\stackrel{~}{T}}`$ gets significantly larger (or smaller) than $`m_{\mathrm{SUSY}}`$, the photon contribution starts to dominate over the $`Z^0`$-boson contribution.

If only $`m_{\stackrel{~}{T}}`$ differs from $`m_{\mathrm{SUSY}}`$, the third generation contributes almost 100% to the total cross section. However, if $`m_{\stackrel{~}{U}}=m_{\stackrel{~}{C}}=m_{\stackrel{~}{T}}`$ are varied simultaneously, all three generations contribute to the total cross section, which can therefore become significantly larger. This is shown in Fig. 3, where (s)quark loop contributions from all three generations

have been taken into account.

When $`m_{\stackrel{~}{U}}=m_{\stackrel{~}{C}}=m_{\stackrel{~}{T}}=m_{\mathrm{SUSY}}`$ and large mass splittings are generated only by mixing in the top squark sector (case II), photon contributions are suppressed by more than two orders of magnitude. Fig. 4 shows that the

$`Z^0`$-boson contributions from top and bottom squarks interfere destructively due to opposite values of their weak isospin quantum numbers, except for $`\theta _{\stackrel{~}{t}}45.2^{}`$, where the imaginary parts of the amplitudes interfere constructively. It is therefore advantageous to keep the bottom squark mass splitting small. As is also evident from Fig. 4, mixing in the bottom squark sector is of little importance. Note that the central region with maximal top/bottom squark mixing is excluded by the CERN LEP limits on $`m_{\stackrel{~}{t}_1}`$, $`m_{\stackrel{~}{b}_1}`$, and the $`\rho `$-parameter.

When $`m_{\mathrm{SUSY}}`$ and the diagonal elements of the squark mixing matrix (see App. A) become much larger than the quark masses and the off-diagonal elements of the matrix, the role of squark mixing is expected to be reduced. This is confirmed numerically in Fig. 5, where the dependence of

the gluino production cross section on $`m_{\mathrm{SUSY}}`$ is shown for the cases of maximal and vanishing top squark mixing. Squarks from the first two generations contribute at most 10% at low $`m_{\mathrm{SUSY}}`$ and are otherwise strongly suppressed.

At future linear $`e^+e^{}`$ colliders it will be possible to obtain relatively high degrees of polarization, i.e. about 80% for electrons and 60% for positrons Aguilar-Saavedra:2001rg . In Fig. 6 we therefore

investigate the effect of choosing different electron/positron polarizations on the gluino pair production process, including contributions from all (s)quarks. Since the $`++`$ and $``$ helicity amplitudes vanish for both photons and $`Z^0`$-bosons, we only show the squares of the remaining $`+`$ and $`+`$ amplitudes, which coincide for photons, but not for $`Z^0`$-bosons. The unpolarized cross section falls short of the polarized ones, so that a high degree of polarization is clearly desirable.

With the realistic degrees of polarization mentioned above, we show in Fig. 7 a scan in the center-of-mass energy of a future $`e^+e^{}`$

collider for various gluino masses and maximal top squark mixing (case II). The cross section rises rather slowly due to the factor $`\beta ^3`$ in Eq. (12) for $`P`$-wave production of the gluino pairs. For $`m_{\stackrel{~}{g}}=200`$ GeV we observe an interesting second maximum, which arises from the intermediate squark pair resonance at $`\sqrt{s}=2m_{\mathrm{SUSY}}=650`$ GeV. At threshold, the cross section depends strongly on the gluino mass and is largest for $`m_{\stackrel{~}{g}}=200`$ GeV, which we consider to be the lowest experimentally allowed value. It drops fast with $`m_{\stackrel{~}{g}}`$, so that for $`m_{\stackrel{~}{g}}>500`$ GeV no events at colliders with luminosities of 1000 fb<sup>-1</sup> per year can be expected, irrespective of their energy. Smaller squark mixing (cf. Fig. 4) or larger values of $`m_{\mathrm{SUSY}}`$ (cf. Fig. 5) will reduce the cross section even further. Far above threshold, it drops off like $`1/s`$ and becomes independent of the gluino mass.

The slow rise of the cross section can be observed even better in Fig. 8, where the sensitivity of a $`\sqrt{s}=500`$ GeV collider like

DESY TESLA to gluino masses around 200 GeV has been plotted. For the CERN LHC experiments, a precision of $`\pm 30\mathrm{}\mathrm{\hspace{0.17em}60}`$ ($`12\mathrm{}\mathrm{\hspace{0.17em}25}`$) GeV is expected for gluino masses of 540 (1004) GeV :1999fr ; Abdullin:1998pm . If the masses and mixing angle(s) of the top (and bottom) squarks are known, a precision of $`\pm 5\mathrm{}\mathrm{\hspace{0.17em}10}`$ GeV can be achieved at DESY TESLA for $`m_{\stackrel{~}{g}}=200`$ GeV and maximal top squark mixing with an integrated luminosity of 100 fb<sup>-1</sup> per center-of-mass energy point.

A center-of-mass energy scan for the scenario with no squark mixing, but large left-/right-handed squark splitting (case I) is shown in Fig. 9 for

light (full and dashed curves) and heavy (dotted and dot-dashed curves) gluino masses. Since the photon contributes now significantly to the cross section, it proves to be advantageous to choose the lepton polarization such that the $`Z^0`$-boson interferes constructively with the photon, even though it is by itself slightly smaller than for the opposite choice. Since all three generations add now to the cross section, it can become almost an order of magnitude larger than in the mixing scenario (case II), and even gluino masses of 1 TeV may be observable at a multi-TeV collider like CERN CLIC. However, also here the cross section drops sharply when the squark mass splitting is reduced from 1500 TeV (full and dotted curves) to values close to $`m_{\mathrm{SUSY}}`$ (dashed and dot-dashed curves).

## IV Conclusions

The purpose of this Paper has been two-fold: First, we have resolved a long-standing discrepancy in the literature about the relative sign of the quark and squark loop contributions to the production of gluino pairs in $`e^+e^{}`$ annihilation. We confirm the result of two older papers that the divergence cancels for each squark flavor separately and not between weak isospin partners Kane:xp ; Campbell:kf  and trace the sign problem in one case to the Feynman rules employed in the corresponding calculation Kileng:1994kw . Our results rely on two completely independent analytical calculations and one computer algebra calculation.

Second, we have investigated the prospects for precision measurements of gluino properties, such as its mass or its Majorana fermion nature, at future linear $`e^+e^{}`$ colliders. We have taken into account realistic beam polarization effects, photon and $`Z^0`$-boson exchange, and current mass exclusion limits. Previously, only light gluinos at center-of-mass energies up to the $`Z^0`$-boson mass had been investigated. Within the general framework of the MSSM, we have concentrated on two scenarios of large left-/right-handed up-type squark mass splitting and large top squark mixing, which produce promisingly large cross sections for gluino masses up to 500 GeV or even 1 TeV. Gluino masses of 200 GeV can then be measured with a precision of about 5 GeV in center-of-mass energy scans with luminosities of 100 fb<sup>-1</sup>/point. However, when both the left-/right-handed squark mass splitting and the squark mixing remain small, gluino pair production in $`e^+e^{}`$ annihilation will be hard to observe, even with luminosities of 1000 fb<sup>-1</sup>/year.

###### Acknowledgements.

We thank B. Kileng, G. Kramer, and P. Osland for useful discussions and A. Brandenburg and B.A. Kniehl for a careful reading of the manuscript. This work has been supported by Deutsche Forschungsgemeinschaft through Grant No. KL 1266/1-2 and through Graduiertenkolleg Zukünftige Entwicklungen in der Teilchenphysik.

## Appendix A Squark Mixing

The (generally complex) soft SUSY-breaking terms $`A_q`$ of the trilinear Higgs-squark-squark interaction and the (also generally complex) off-diagonal Higgs mass parameter $`\mu `$ in the MSSM Lagrangian induce mixings of the left- and right-handed squark eigenstates $`\stackrel{~}{q}_{L,R}`$ of the electroweak interaction into mass eigenstates $`\stackrel{~}{q}_{1,2}`$. The squark mass matrix Haber:1984rc ; Gunion:1984yn

$$^2=\left(\begin{array}{cc}m_{LL}^2+m_q^2& m_qm_{LR}^{}\\ m_qm_{LR}& m_{RR}^2+m_q^2\end{array}\right)$$
(18)
with

$`m_{LL}^2`$ $`=`$ $`(T_q^3e_q\mathrm{sin}^2\theta _W)m_Z^2\mathrm{cos}2\beta +m_{\stackrel{~}{Q}}^2,`$ (19)
$`m_{RR}^2`$ $`=`$ $`e_q\mathrm{sin}^2\theta _Wm_Z^2\mathrm{cos}2\beta +\{\begin{array}{c}m_{\stackrel{~}{U}}^2\mathrm{for}\mathrm{up}\mathrm{type}\mathrm{squarks},\hfill \\ m_{\stackrel{~}{D}}^2\mathrm{for}\mathrm{down}\mathrm{type}\mathrm{squarks},\hfill \end{array}`$ (22)
$`m_{LR}`$ $`=`$ $`A_q\mu ^{}\{\begin{array}{c}\mathrm{cot}\beta \mathrm{for}\mathrm{up}\mathrm{type}\mathrm{squarks}\hfill \\ \mathrm{tan}\beta \mathrm{for}\mathrm{down}\mathrm{type}\mathrm{squarks}\hfill \end{array}`$ (25)
is diagonalized by a unitary matrix $`S`$, $`S^2S^{}=\mathrm{diag}(m_1^2,m_2^2)`$, and has the squared mass eigenvalues

$$m_{1,2}^2=m_q^2+\frac{1}{2}\left(m_{LL}^2+m_{RR}^2\sqrt{(m_{LL}^2m_{RR}^2)^2+4m_q^2|m_{LR}|^2}\right).$$
(26)
For real values of $`m_{LR}`$, the squark mixing angle $`\theta _{\stackrel{~}{q}}`$, $`0\theta _{\stackrel{~}{q}}\pi /2`$, in

$$S=\left(\begin{array}{cc}\mathrm{cos}\theta _{\stackrel{~}{q}}& \mathrm{sin}\theta _{\stackrel{~}{q}}\\ \mathrm{sin}\theta _{\stackrel{~}{q}}& \mathrm{cos}\theta _{\stackrel{~}{q}}\end{array}\right)\mathrm{with}\left(\begin{array}{c}\stackrel{~}{q}_1\\ \stackrel{~}{q}_2\end{array}\right)=S\left(\begin{array}{c}\stackrel{~}{q}_L\\ \stackrel{~}{q}_R\end{array}\right)$$
(27)
can be obtained from

$$\mathrm{tan}2\theta _{\stackrel{~}{q}}=\frac{2m_qm_{LR}}{m_{LL}^2m_{RR}^2}.$$
(28)
If $`m_{LR}`$ is complex, one may first choose a suitable phase rotation $`\stackrel{~}{q}_R^{}=e^{i\varphi }\stackrel{~}{q}_R`$ to make the mass matrix real and then diagonalize it for $`\stackrel{~}{q}_L`$ and $`\stackrel{~}{q}_R^{}`$. $`\mathrm{tan}\beta =v_u/v_d`$ is the (real) ratio of the vacuum expectation values of the two Higgs fields, which couple to the up- and down-type (s)quarks. The weak isospin quantum numbers for left-handed up- and down-type (s)quarks with hypercharge $`Y_q=1/3`$ are $`T_q^3=\{+1/2,1/2\}`$, whereas $`Y_q=\{4/3,2/3\}`$ and $`T_q^3=0`$ for right-handed (s)quarks, and their fractional electromagnetic charges are $`e_q=T_q^3+Y_q/2`$. The soft SUSY-breaking mass terms for left- and right-handed squarks are $`m_{\stackrel{~}{Q}}`$ and $`m_{\stackrel{~}{U}}`$, $`m_{\stackrel{~}{D}}`$, respectively, and $`m_Z`$ is the mass of the neutral electroweak gauge boson $`Z^0`$.

## Appendix B Feynman Rules

Denoting squark mass eigenstates by $`i,j,\mathrm{}`$, Lorentz indices by $`\mu ,\nu ,\mathrm{}`$, and color indices of the fundamental (adjoint) representation of the color symmetry group SU(3) by $`l,m,\mathrm{}`$ $`(a,b,\mathrm{})`$, we obtain the following propagators in Feynman-gauge:

$`{\displaystyle \frac{ig^{\mu \nu }}{p^2m_V^2+i\eta }},m_V^2=\{0,m_Z^2\}`$ (29)
$`{\displaystyle \frac{i(p/+m_q)\delta _{lm}}{p^2m_q^2+i\eta }}`$ (30)
$`{\displaystyle \frac{i\delta _{ij}\delta _{lm}}{p^2m_{\stackrel{~}{q}}^2+i\eta }}`$ (31)
Dirac fermions carry an arrow, which indicates the fermion number flow, whereas Majorana fermions, such as gluinos, do not carry arrows. An additional arrow is depicted next to all fermion lines in order to obtain a unique orientation of the fermion flow, which is evaluated according to the rules in Denner:1992vz . The interaction vertices are given by Haber:1984rc

$`ie\gamma _\mu (v_f^Va_f^V\gamma _5)\delta _{lm}\{\begin{array}{cc}v_f^\gamma =e_f\hfill & a_f^\gamma =0\hfill \\ v_f^Z=\frac{T_f^32e_f\mathrm{sin}^2\theta _W}{2\mathrm{sin}\theta _W\mathrm{cos}\theta _W}\hfill & a_f^Z=\frac{T_f^3}{2\mathrm{sin}\theta _W\mathrm{cos}\theta _W}\hfill \end{array}`$ (34)
$`ie(p+p^{})_\mu \mathrm{\Gamma }_q^{ij,V}\delta _{lm}\{\begin{array}{ccc}\mathrm{\Gamma }_q^{ij,\gamma }\hfill & =& e_q\delta _{ij}\hfill \\ \mathrm{\Gamma }_q^{ij,Z^0}\hfill & =& \frac{(T_q^3e_q\mathrm{sin}^2\theta _W)S_{j1}S_{i1}^{}e_q\mathrm{sin}^2\theta _WS_{j2}S_{i2}^{}}{\mathrm{sin}\theta _W\mathrm{cos}\theta _W}\hfill \end{array}`$ (37)
$`i\mathrm{\Gamma }_{i,1}^a=i\sqrt{2}\widehat{g}_sT_{lm}^a(S_{i1}P_LS_{i2}P_R)`$ (38)
$`i\mathrm{\Gamma }_{i,2}^a=i\sqrt{2}\widehat{g}_sT_{ml}^a(S_{i1}^{}P_RS_{i2}^{}P_L)`$ (39)
The Feynman rules in Eq. (34) apply to photon and $`Z^0`$-boson interactions with quarks and charged leptons. The latter carry electromagnetic charge $`e_{\mathrm{}}=1`$ and weak isospin $`T_{\mathrm{}}^3=1/2`$ (left-handed) and 0 (right-handed), but no color ($`\delta _{lm}1`$). The gauge couplings $`g`$ and $`g^{}`$ of the weak isospin and hypercharge symmetries SU(2)<sub>L</sub> and U(1)<sub>Y</sub> have been expressed in terms of the electromagnetic coupling $`e=g\mathrm{sin}\theta _W`$ and the sine and cosine of the weak mixing angle $`\mathrm{sin}\theta _W/\mathrm{cos}\theta _W=g^{}/g`$. The gluino-quark-squark vertices depend on the generators of the SU(3) color symmetry group, $`T_{lm}^a`$, and on the Yukawa coupling $`\widehat{g}_s`$, which is identical to the strong gauge coupling $`g_s`$ in leading order, and on the squark mixing matrix $`S`$, but not on the orientation of the fermion flow.