Pattern formation with a conservation law

Pattern formation in systems with a conserved quantity is considered by studying the appropriate amplitude equations. The conservation law leads to a large-scale neutral mode that must be included in the asymptotic analysis for pattern formation near onset. Near a stationary bifurcation, the usual Ginzburg-Landau equation for the amplitude of the pattern is then coupled to an equation for the large-scale mode. These amplitude equations show that for certain parameters all roll-type solutions are unstable. This new instability differs from the Eckhaus instability in that it is amplitude-driven and is supercritical. Beyond the stability boundary, there exist stable stationary solutions in the form of strongly modulated patterns. The envelope of these modulations is calculated in terms of Jacobi elliptic functions and, away from the onset of modulation, is closely approximated by a sech profile. Numerical simulations indicate that as the modulation becomes more pronounced, the envelope broadens. A number of applications are considered, including convection with fixed-flux boundaries and convection in a magnetic field, resulting in new instabilities for these systems.


Introduction
The mathematical theory of pattern formation has a wide range of applications. In the field of fluid mechanics, the most widely studied example is Rayleigh-Bénard convection and its variants, but there are many other examples such as Taylor-Couette flow and solidification of alloys. Pattern-forming systems can also be found in many chemical and biological systems (see the review by Cross and Hohenberg 1993).
In this paper we consider pattern formation in one space dimension in systems sufficiently large to allow a large-scale modulation of the pattern. It is assumed that the system has some reflection symmetry, that the onset of instability is stationary and that the onset of pattern formation occurs with a non-zero wavenumber. In this case, the equation governing the modulation of the small-amplitude pattern is in general the real Ginzburg-Landau equation. This well known result has recently been given rigorously by Melbourne (1999). The behaviour of the real Ginzburg-Landau equation is straightforward: a regular periodic pattern is stable provided its wavenumber lies within the Eckhaus band (Eckhaus 1965); patterns outside this band are unstable and the instability is subcritical (Fauve 1987, Tuckerman andBarkley 1990). Physically, patterns in which the wavenumber is either too large or too small adjust to a different wavenumber by adding or losing waves.
However, there is a large class of problems in which the onset of pattern formation at finite wavenumber is not governed by the Ginzburg-Landau equation. Such problems arise when

Derivation of equations from symmetry
Consider a PDE or system of PDEs with a conserved quantity and a basic zero state. If a quantity w(x, t) is conserved, it must be possible to write ∂w/∂t as the divergence of a flux. It follows that when the expression for ∂w/∂t is linearized, there can be no terms linear in w except those involving x-derivatives. If we assume, furthermore, that the system possesses the reflection symmetry x → −x w→ w (2.1) then only even x-derivatives occur. Thus w is density-like (a velocity-like quantity would satisfy x → −x, w → −w instead; see Matthews and Cox (2000), Tribelsky and Velarde (1996)). On large spatial scales, the leading linear term is the second derivative w xx , so w(x, t) obeys the diffusion equation and the growth rate λ of modes with wavenumber k is proportional to k 2 . This indicates that large-scale modes are almost neutral, evolving on a long timescale, and must be included in the analysis of pattern formation. Suppose that the system supports stationary patterns with a wavenumber k, which we take to be the critical wavenumber according to linear stability theory; we may set k = 1 with no loss of generality. Near the onset of pattern formation the appropriate ansatz is then w(x, t) ∼Ã(X, T ) e ix +Ã * (X, T ) e −ix +B(X, T ) ( 2.2) where X and T are rescaled forms of x and t, andÃ andB are small amplitudes. Note that the amplitudeÃ is complex, with A * denoting the complex conjugate of A, butB is real. The amplitude equations forÃ andB can be written down by imposing the symmetries of translation in x where χ is any real number, and reflection X → −XÃ →Ã * B →B. (2.4) In the absence of the large-scale modeB,Ã obeys the real Ginzburg-Landau equation.
Coupling terms in this equation must includeÃ because of (2.3), so the leading coupling term in theÃ equation that is consistent with the symmetries isÃB. The coefficient of this term is forced to be real by (2.4). In the absence of coupling terms, the equation forB is the diffusion equation, and (2.3) forces coupling terms to involve equal numbers ofÃ andÃ * terms. SinceB is conserved, all terms in theB equation involve derivatives. The leading coupling term in theB equation that is consistent with (2.4) is therefore (|Ã| 2 ) XX . These considerations give rise to amplitude equations of consistent asymptotic order if we introduce a small parameter such thatÃ = A,B = 2 B, X = x, T = 2 t. Furthermore, we take the system to be a distance of O( 2 ) from the marginal surface in parameter space. With these scalings, the amplitude equations truncated to leading order are then where we have assumed that, in the absence of B, the bifurcation of the pattern mode is supercritical. All coefficients in (2.5) have been set to unity by rescaling T , X, A and B. There are two coefficients, σ and µ, in (2.6) that cannot be removed by rescaling. We impose the condition σ > 0 to ensure that the long-wavelength modes are weakly damped in the absence of A, but µ may have either sign. Equations analogous to (2.5) and (2.6) have previously been derived in the context of the secondary stability of cellular patterns (Coullet and Iooss 1990).

Consider the PDE
where the terms inside the brackets are just those of the much-studied Swift-Hohenberg equation (Swift and Hohenberg 1977) with a symmetry-breaking quadratic term (Sakaguchi andBrand 1996, Matthews 1998). This equation has the required properties of reflection symmetry and a conservation law for w(x, t).
The growth rate λ of modes of wavenumber k in (2.7) is given by so for 0 < r 1 a narrow band of wavenumbers near k = k c = 1 is unstable (figure 1).

Figure 1.
The growth rate λ of the linear modes of (2.7) as a function of wavenumber k, with r = 0.1.
We now follow the standard asymptotic methods to describe the behaviour of (2.7) for small r, by writing The appropriate asymptotic expansion for w(x, t) is then w(x, t) = A(X, T ) e ix + A * (X, T ) e −ix + 2 B(X, T ) where the large-scale mode B has been introduced at order 2 . The expansion (2.10) is substituted into (2.7) and the system is then solved at successive orders of . At O( ), equation (2.7) is automatically satisfied. At O( 2 ), terms involving e 2ix appear, and from these we deduce that C = −sA 2 /9. At O( 3 ), the amplitude equation for A is deduced from the terms in e ix . The equation for B follows from the terms independent of x at O( 4 ). These equations are The sign of s is thus immaterial, and the bifurcation is supercritical if s 2 < 27 2 . If this condition is satisfied then (2.11) and (2.12) can be rescaled to give (2.5) and (2.6).
Although the system (2.5) and (2.6) is generic for pattern-forming systems with a conservation law, in fact, a conservation law alone is not sufficient to yield these equations. Only the first of the two nonlinear terms in (2.7) generates the coupling term in (2.5), so if we set s = 0 there is no coupling and the dynamics is then governed by the Ginzburg-Landau equation alone. Furthermore, if the original PDE has the symmetry w → −w then all nonlinear terms in (2.11) are cubic and again no coupling terms arise. This restriction will apply when we consider physical examples in section 6.

Analysis of the amplitude equations
In this section we analyse the steady solutions of the amplitude equations (2.5) and (2.6). We begin by considering small-amplitude, periodic solutions (rolls). We then investigate the super/subcritical nature of the instability, and, beyond the point where rolls become unstable, we calculate fully nonlinear solutions of (2.5) and (2.6).

Stability of patterns
In this section we analyse the stability of periodic solutions in the form of rolls. Rolls with wavenumber k c + q correspond to solutions of (2.5) and (2.6) of the form A = a 0 e iqX and B = 0, where the amplitude a 0 may be taken to be real and to satisfy a 2 0 = 1 − q 2 . To analyse the stability of these rolls, we set and linearize (2.5) and (2.6) in a and b. Upon Solutions of (3.2)-(3.4) proportional to e λT have a growth rate λ that satisfies the cubic equation For modes with λ real, there is instability if The most dangerous disturbance therefore arises in the limit l → 0, and so the criterion for instability is or, equivalently, (3.9) Figure 2. Region of stable rolls in (2.5) and (2.6). The full curve shows the stationary bifurcation, which is valid for any σ , and the broken line shows the Hopf bifurcation for σ = 1 4 . For µ > σ all rolls are unstable, while for µ < σ there is a band of stable rolls that widens as µ is decreased. The usual Eckhaus result of stability for q 2 < 1 3 corresponds to µ = 0.
If µ < σ then rolls of sufficiently small amplitude are unstable; the criterion for instability can be expressed in terms of q as It is clear that when the coupling between the mean and pattern modes is removed (µ = 0), or more generally when |σ/µ| 1, the usual Eckhaus criterion q 2 > 1 3 for instability is recovered. In general, the band of stable rolls is narrowed if µ > 0 and widened if µ < 0.
The stability results are summarized in figure 2. Note that the stationary bifurcation given by (3.10) depends only on the ratio µ/σ but this is not true for the Hopf bifurcation.
The mechanism of instability in (2.5) and (2.6) is rather different from that of the Eckhaus instability. Even rolls with q = 0 are unstable for µ > σ, and in this case the unstable mode has V = W so that a is real and the instability is amplitude-driven rather than phase-driven. If there are regions where |A| 2 has a minimum, then, according to (2.6), B grows in these regions if µ > 0. This leads to a positive value of B which then reduces the size of A further, through the coupling term in (2.5).

Nonlinear development of the instability
The results of section 3.1 show that for µ > σ, all periodic patterns in (2.5) and (2.6) are unstable. We now consider the nonlinear development of the instability, to determine whether the bifurcation is subcritical or supercritical. The Eckhaus instability, for example, is known to give rise to a subcritical bifurcation (Fauve 1987, Tuckerman andBarkley 1990).
For simplicity we consider the case q = 0, corresponding to rolls with wavenumber k c . These are the last rolls to become unstable as µ is increased. In this case the stationary solution is A = 1, B = 0. Writing A = 1 + a(X, T ), B = b(X, T ), cf (3.1), we find that the imaginary part of a decays so that we may take a to be real. The nonlinear equations for the real-valued quantities a and b are then, from (2.5) and (2.6), A weakly nonlinear analysis is now employed, by writing where the a i and b i are functions of X and T 2 . By considering a single mode, we find that the terms of order in (3.12) and (3.13) give b 11 + (2 + l 2 )a 11 = 0 and µ 0 = (1 + l 2 /2)σ so the instability condition (3.9) is recovered for small l. At O( 2 ) we deduce that where a 20 = (l 2 − 1)a 2 11 /4 a 22 = a 2 11 /4 and b 22 = −(1 + l 2 /2)a 2 11 . Applying the solvability condition at O( 3 ) gives the equation for the evolution on the slow timescale T 2 : (σ l 2 + 2 + l 2 ) da 11 dT 2 = 2l 2 µ 2 a 11 + 1 2 (σ l 4 (l 2 − 1))a 3 11 .
Hence the bifurcation is a supercritical pitchfork provided that l < 1. In a large periodic domain, the first modes to become unstable have small l, and so the bifurcation will in practice be supercritical.

Non-uniform steady solutions
We now extend the analysis of the previous subsection to consider fully nonlinear stationary solutions of the amplitude equations (2.5) and (2.6). Writing A = R exp(i ), we find from the imaginary part of (2.5) that (R 2 X ) X = 0 and so X = R −2 , where is constant. From (2.6) we find, for solutions periodic in X, that where · · · denotes a spatial average. In deriving (3.14), we have taken B = 0, which corresponds to an initial condition such that w = 0 (and hence w = 0 for all time). Any other initial value for w is equivalent to w = 0 under a rescaling of the parameters r and s in (2.7). From the real part of (2.5), it then follows that R satisfies This equation may be integrated once, and is then more conveniently expressed in terms of Q = R 2 : where H is a constant of integration. The solution is readily determined in terms of the Weierstrass elliptic function P (Abramowitz andStegun 1965, Whittaker andWatson 1962) to be and g 2 and g 3 are the usual parameters of P. The amplitude and period of the solution are implicitly related through the appearance of the quantity Q in (3.18)-(3.20). Although (3.17) provides an exact solution, it requires specification of µ , H , and the spatial period to permit analysis of its solutions. Rather more analytical progress is possible if we restrict attention to the invariant subspace in which A is real. We suppose that the solution has spatial period L. Equations (2.5) and (2.6) become

Integrating (3.22) twice and imposing periodicity gives the result
corresponding to (3.14). Eliminating B from (3.21) gives Before describing the periodic solutions of (3.21) and (3.22), it is instructive to consider a solution that is a good asymptotic approximation: where µ A 2 = 1 + c 2 and µ = 1 + 2c 2 /A 2 0 . Note that the solution only exists for µ > 1, which is consistent with the stability boundary (3.9) and the result from section 3.2 that the bifurcation at µ = 1 is supercritical. For large values of cL, cL A 2 ∼ 2A 2 0 . Using this approximation for A , we obtain for c the quadratic equation which has real solutions only for indicating that a saddle-node bifurcation occurs at µ = L/(L−2) and that for µ > L/(L−2) no such solution exists. At large amplitude the solution branch approaches µ = 1 and the amplitude is given by A 0 ∼ 4 √ 2 µ /L(µ − 1) 3/2 . The solution (3.25) is of particular significance because, as we shall see, it can be stable and is observed in the numerical simulations described in section 4.
Spatially periodic solutions of (3.24) can be written in terms of the Jacobi elliptic functions sn, cn and dn. In all that follows, m denotes the standard parameter of these functions (Whittaker andWatson 1962, Abramowitz andStegun 1965), and we use the notation · · · as above, so that where K is the real quarter-period, given by 3.3.1. dn. When µ > 1, there exist solutions of (3.24) in the form Given the parameter µ and the length L of the domain, the three equations in (3.28) serve to determine the three unknowns A 0 , c and m. However, the computation of solutions is more easily accomplished if L, say, is fixed, while m is varied in the range 0 m < 1 to yield corresponding values of A 0 , c and µ .
For small m, the solution takes the form  indicating that these solutions branch from the uniform state A ≡ 1; the branching takes place at µ = 1 + 2π 2 /L 2 (corresponding to m = 0). Furthermore, it is readily shown from (3.28) that for small m confirming the result of section 3.2 that the bifurcation is supercritical when L > 2π (l < 1 in the notation of section 3.2) and subcritical when L < 2π (l > 1). As m is increased from zero to unity, the corresponding values of µ increase to a maximum µ d (given approximately in (3.26) by µ d ∼ L/(L − 2)), and then decrease monotonically so that µ ∼ 1 + 4 ln(16/m 1 ) (3.31) as m → 1 − , where m 1 = 1 − m. We have used the well known result that in this limit. Thus a saddle-node bifurcation takes place at µ = µ d , and for µ > µ d no solutions of the form (3.27) exist. Figure 3 shows the bifurcation diagram for dn-type solutions of (3.21) and (3.22) subject to periodic boundary conditions with L = 5π and 10π . The saddle-node bifurcations occur at µ d = 1.146 and µ d = 1.068, respectively, in agreement with our approximation (3.26). Since the bifurcation is supercritical in each case, the lower branch of solutions is stable but the upper branch is unstable.
When m is small, the solution for A takes the form of small, almost sinusoidal perturbations to the uniform state A ≡ 1. As m is increased, the peak of A becomes narrow and the trough broadens. As m → 1 − , the peak increases greatly in amplitude and the solution takes the form where µ ∼ 1 + 2/K, and K → ∞ according to (3.32). In this limit, the solution for A thus takes the form of a single sech-type pulse as given by (3.25). The dn solutions are illustrated in figure 4. These results indicate that for µ d > µ > 1 a soliton-like solution exists to (2.5) and (2.6). For µ > µ d , numerical simulations of (2.5) and (2.6) indicate that the solution blows up, signalling a transition to a fully nonlinear solution of the original PDEs (for example, equation (2.7)). This blow-up can be avoided by adding higher-order terms to the amplitude equations, which we discuss in section 5. When m is small, A ∼ √ 2(1 − 4π 2 /L 2 ) 1/2 cos(2π X/L), and hence we require L > 2π for the existence of this solution in this limit. Further analysis indicates that for general values of m a solution exists only when

cn. When µ > 1, solutions also exist in the form
(3.36) The function g(m) is shown in figure 5. It is clear from (3.36) and figure 5 that if L > 2π the solution (3.34) exists for all m between zero and unity. If, on the other hand, L < 2π then solutions exist for m L < m < 1, where g(m L ) = L 2 . Regardless of the value of L, these solutions exist for m > m 0 ≈ 0.826 115. If L > 2π then as m is increased from zero, µ increases from unity to a maximum value µ c (given approximately by µ c = L/(L − 4)) and then decreases monotonically according to (3.31) as m → 1 − . If L < 2π then as m is increased from m L , µ decreases monotonically from infinity, again according to (3.31) as m → 1 − . As m is increased, the peaks of A become more pronounced, larger in amplitude and narrower, while the regions in which A has small amplitude become broader. In the limit m → 1 − , the solution is approximately given by two sech-type pulses, one positive and one negative; with an appropriate choice of origin in X, the solution for 0 < X < L satisfies with K and m related through (3.32). The cn solutions are illustrated in figure 6.

sn.
The final type of solution that may be expressed in the form of Jacobi elliptic functions applies when µ < 1, and is For small m, these solutions match with the cn solutions (allowing for an unimportant translation in X), with A ∼ √ 2(1 − 4π 2 /L 2 ) 1/2 sin(2π X/L). As this expression indicates, the small-m solutions require L > 2π. In fact, for general values of m, the length L of the domain must satisfy (3.39) As m increases towards m L , µ decreases monotonically to −∞, the peaks and troughs in A becoming broader. Provided L is large, the profile near X = 0 is approximately a front of the form as m → m − L , the front near X = L/2 being given by symmetry. The sn solutions are illustrated in figure 7.
A bifurcation diagram showing all three types of solution (cn, dn and sn) is given in figure 8. Note that only the dn solutions bifurcate from the periodic rolls represented by A = 1.
The solutions described above are those containing the minimal number of peaks and troughs within the periodic domain. Other solutions, containing n sets of peaks and troughs, can also be constructed by replacing K by nK in the last equation of (3.28), (3.35) or (3.38).

Numerical simulations of a model PDE
In this section we illustrate the application of the results of the previous sections through numerical simulations of the PDE (2.7), assuming periodic boundary conditions. From the analysis of section 2.2, the system is governed by the amplitude equations (2.5) and (2.6) with A pseudospectral method is employed in which the nonlinear terms are computed in physical space using fast Fourier transforms. Time advancement is carried out by approximating the nonlinear terms as a linear function of time and then obtaining exact solutions using the known growth rate of the linear terms (2.8).
In order to approach the asymptotic regime, the size of the periodic domain must be large and the driving parameter r must be small. Initially, the length of the domain is chosen to be 62π, so that 31 periods of the pattern are contained in the box, the forcing parameter is r = 0.01 and the number of grid points in physical space is 256. The effective length of the domain for the amplitude equations (2.5) and (2.6) in terms of the length scale X is then L = 62π √ r/2 = 3.1π . Note that since L > 2π , the result of section 3.2 indicates that the bifurcation from a regular pattern to a modulated pattern is supercritical. In the limit of large L the condition (3.9) for rolls to be unstable is µ > σ , which is equivalent to s 2 > 27 38 . For finite L the condition is (3.6), which yields µ > σ (1 + 2π 2 /L 2 ), or |s| > 0.9215 for L = 3.1π . In our numerical simulations of (2.7) at r = 0.01, we find good agreement with this asymptotic result, rolls being unstable for |s| > 0.927.
A sequence of numerical solutions of (2.7) with r = 0.01 for different values of s is shown in figure 9. In all cases, the solution shown is the final steady state after all transients have decayed; to reach this state the computation must be continued for several thousand time units. In figure 9(a), s = 0.9 and a stable periodic solution is found. For s = 0.93 (figure 9(b)), slightly beyond the bifurcation, a modulated pattern appears. The amplitude of the modulation is small, which is consistent with the result of section 3.2 that the bifurcation is supercritical. As s is increased, the amplitude of the modulation steadily increases (figures 9(c) and (d)), until there are regions in which the pattern is almost completely suppressed and regions in which the amplitude of the pattern becomes large. The envelope of the pattern takes the form of the dn or sech solutions described in section 3.3 (see figure 4); we have not found any of the cn or sn solutions, which suggests that they are unstable.
The increasing amplitude of these solutions means that the numerical simulations are no longer within the asymptotic regime of validity of equations (2.5) and (2.6); in particular, for this value of r we do not find the saddle-node bifurcation, which, according to (3.26), should occur when s = 0.939 in (2.7).
To find this saddle-node bifurcation it is necessary to carry out further sequences of numerical simulations at smaller values of r. In order to keep the value of L constant, the length of the domain (and the number of grid points) must be increased in proportion to r −1/2 . The numerically obtained bifurcation diagrams, with r successively reduced by a factor of four, are compared with the solution of the amplitude equations in figure 10. The function A is obtained numerically from the part of the spectrum of the numerical solution near k = 1 and scaled appropriately with r 1/2 . For a domain size of 62π with r = 0.01 there is no saddle-node bifurcation and the numerical and theoretical bifurcation diagrams are not in close agreement. With a domain size of 124π and r = 0.0025, agreement is better near the Figure 11. Numerical simulations of (2.7) showing w(x), with r = 0.01, domain size 62π (using 512 grid points) and s = 1.6. The two solutions are started from different initial conditions. The single-pulse solution appears to be a steady state, while in (a) the two pulses are slowly separating. The amplitude appears to be almost uniform in each of the humps. pitchfork bifurcation, but there is still no saddle-node bifurcation. Finally, for a domain size of 248π with r = 0.000 625 a saddle-node bifurcation appears at s = 0.944, close to the theoretical value of 0.939. Beyond this point the numerical results jump to a solution of much larger amplitude, and there is a region of hysteresis where two stable solutions exist.
For larger values of s, finite-amplitude solutions of (2.7) are found, that cannot be described in terms of the amplitude equations. In some cases there are multiple stable solutions. For example, in figure 11 we show two solutions computed at s = 1.6, from different initial states. In one plot there are two regions of pattern, and in the other only one. The amplitude of the solution appears to be approximately the same in each simulation. However, the single pulse has roughly twice the width of each of the smaller pulses, so in each simulation an equal proportion of the computational domain seems to be occupied by the pattern. Furthermore, it is clear that the envelope is rather broader than a sech profile, particularly in figure 11(b), where the envelope consists of two nearly uniform regions separated by a pair of fronts. In addition to those displayed in figure 11, we have found further steady (or quasi-steady) solutions of (2.7), including a single-pulse solution containing six large waves (rather than the seven in figure 11(b)) and a two-pulse solution where the pulses are of unequal size, one containing three and the other four large waves. Closer examination of the numerical results indicates that the multi-pulse solutions are not quite stationary; the pulses are gradually moving.

Beyond the leading-order amplitude equations
In the previous section we have seen that agreement between the amplitude equations (2.5) and (2.6) and the behaviour of the full PDE (2.7) only occurs for extremely small values of r. This agreement can be improved by introducing higher-order terms into (2.5) and (2.6). Rather than including all possible terms that might arise at the next order in (2.5) and (2.6), we focus on one term.
To motivate our choice of higher-order term, we note that when we compute the x-average of w 3 , the leading-order contribution is 4 ( 2 3 s|A| 4 + 6|A| 2 B). For moderate values of s, the term 6 4 |A| 2 B dominates. After rescaling the equations for A and B so that A satisfies (2.5), we find that B satisfies where σ and µ are given by (4.1) and For future convenience we introduce the quantity M, defined by When µ > 0, the term in δ plays a stabilizing role. For steady solutions, the governing equation for A becomes The bifurcation diagram for this equation, subject to the constraint (5.3), is indicated in figure 12. The bifurcation diagram is changed significantly even by very small values of δ. For example, whereas there is a single saddle-node bifurcation for δ = 0, when δ > 0 there is a pair of saddle-node bifurcations. This pair exists only for δ/σ less than approximately 0.0009, and for larger values of δ/σ the curve is monotonic. We now wish to determine the effect of small δ upon the envelope of the pattern, since our numerical simulations suggest that the envelope will broaden from a sech profile as the envelope becomes more developed. For small δ, equation (5.4) is approximately and so to determine the effect of δ upon the envelope, it is instructive to examine the model equation where b = 2/(1 − 16δ/3) 1/2 . When δ = 0, f = √ 2 sech x, but as δ increases the profile broadens and flattens (see figure 13), giving envelopes reminiscent of those in figure 11.

Examples
Since conservation laws arise very commonly in nature, it turns out that a wide variety of pattern formation problems are governed by the system (2.5) and (2.6) rather than the Ginzburg-Landau equation alone.
In thermal convection in a rotating layer with stress-free boundaries, the conserved quantity is horizontal momentum and the corresponding long-wave mode is a large-scale flow along the axes of the convection rolls. In the analysis of Cox and Matthews (2000), equations (2.5) and (2.6) appear as a special case, and it can be shown that rotating convection is unstable for moderate values of the Prandtl number.
In Marangoni convection (driven by surface tension), the fluid has a free surface and the equation for the layer depth is constrained by conservation of mass. Equations similar to (2.5) and (2.6), but with one extra term, were derived by Golovin et al (1994); in this case, however, the parameter µ is negative and so no instability arises.
In sections 6.1 and 6.2 below, we derive an equation similar to (2.7) for two problems in convection; in section 6.3, we directly derive amplitude equations equivalent to (2.5) and (2.6) for magnetoconvection between 'ideal' horizontal boundaries.

Application to convection with a nonlinear equation of state
A physical problem that gives rise to amplitude equations of the form considered here is convection between boundaries at which the heat flux is fixed (Roberts 1985, Matthews 1988a. The fluid has a nonlinear equation of state, relating the fluid density ρ * to the temperature T * . The coefficients α, β and T 0 are constants. This equation of state is appropriate for water at around 4 • C, for example. The crucial feature of this problem is that the total heat in the layer is conserved; the fixed-flux boundary conditions allow a large-scale, slowly evolving depth-averaged temperature mode, . We assume that the temperature profile is nonlinear, of the form (Matthews 1988a, b): where z * is the vertical coordinate. The problem is expressed in dimensionless variables by adopting the scales d = (B/A) 1/2 , Bd and d 2 /κ for length, temperature and time. We then find the equations governing two-dimensional motions to be where the xand z-components of velocity may be expressed in terms of the streamfunction ψ as (u, w) = (−ψ z , ψ x ), and the Jacobian is defined by J (f, g) = f x g z − f z g x . The perturbation to the temperature profile (6.2) is denoted by θ . The Rayleigh number is where ν and κ are, respectively, the coefficient of kinematic viscosity and the coefficient of thermal diffusivity, g is the gravitational acceleration; the Prandtl number is P = ν/κ, and = 2Bdβ/α. In deriving the PDE for , we assume that the equation of state is only mildly nonlinear, so that | | 1. We break the up-down symmetry of the problem by imposing rigid boundary conditions at z = −H but stress-free conditions at z = H . Thus The important mathematical consequence of these boundary conditions is that the equation for then contains a quadratic nonlinear term of the form ( 2 X ) XX , where X is a scaled xcoordinate. This term tends to lead to a supercritical bifurcation to rolls, unlike the quadratic term ( 2 ) XX , which arises from the nonlinear equation of state, and which tends to lead to a subcritical bifurcation.
We suppose that the heat flux through the planes z = ±H is fixed, so that the perturbation heat flux ∂θ/∂z = 0 at each boundary. The key consequence of this choice of boundary condition is that, for sufficiently small half-depth H , the critical mode of the linear stability problem has zero wavenumber (see figure 14).
We begin by considering the simpler case of a linear equation of state, with = 0. First we examine equations (6.3) and (6.4), linearized about the trivial state ψ = θ = 0. We find, making asymptotic expansions of the forms given by Roberts (1985) and Matthews (1988a), that modes with wavenumber k have linear growth rate λ, where for small |R − R 0 | and k. The coefficient a 4 is As a consistency check on these results, we note that as H → 0, R 0 ∼ 20/H 4 and a 4 ∼ 232H 2 /693. This limit corresponds to the problem with a linear equation of state and linear basic-state temperature profile, such as considered by Chapman and Proctor (1980).

Figure 14.
Marginal stability curves for convection with a nonlinear equation of state and nonlinear basic-state temperature profile. When the depth H of the fluid layer is less than some threshold H 0 , the critical wavenumber is zero and the critical Rayleigh number is R 0 , given by (6.7); when H > H 0 , the critical wavenumber k 1 is positive, and the critical Rayleigh number R 1 is less than R 0 .
After appropriately rescaling the depth of our fluid layer, we find that their results agree with ours in the limit as H → 0. When a 4 > 0, the largest wavenumbers in (6.7) are damped, and k = 0 corresponds to a minimum of the marginal stability curve. This condition is satisfied for 0 < H < H 0 , where H 0 ≈ 1.2709 is the smaller positive root of a 4 = 0 (see figure 14). The amplitude equation of interest applies when H is close to H 0 . We thus follow Roberts (1985) and Matthews (1988a) and let R = R 0 + 2 R 2 + 4 R 4 = 2 D H = H 0 + 2 H 2 (6.8) expanding the streamfunction and temperature perturbation as θ = 2 θ 2 (X, T ) + 2 θ 4 (X, z, T ) + · · · (6.9) ψ = 3 ψ 3 (X, z, T ) + 2 ψ 5 (X, z, T ) + · · · (6.10) where X = x and T = 6 t. These scalings have been adopted to ensure a balance of terms in the equation derived below. Of course, it is reasonable to question the physical relevance of such a large power of in the scaling of time, and we offer two partial justifications. The first is that, in fact, solutions to model PDEs of the form to be derived below are known to provide a good approximation to the full convection problem even when is not particularly small (and hence 6 represents a physically achievable scaling). The second justification is that in deriving amplitude equations of the form (2.5) and (2.6) it is not necessary to pass through an intermediate stage involving the time scale 6 t, this is just a convenient way to construct in an asymptotic fashion a long-wavelength model PDE.
At O( 2 ) in (6.3) and (6.4), we find that ∂ 2 θ 2 /∂z 2 = 0, subject to the boundary conditions ∂θ 2 /∂z = 0 at z = ±H . The solution for θ 2 is therefore independent of depth, as indicated in (6.9). After some algebra, we find from a solvability condition at O( 8 ) that the amplitude equation for = θ 2 (X, T ) is where the numerical values of the coefficients in (6.11) are given approximately by a 2 = 7.588 × 10 −5 (6.12) a 4 = 8.966 (6.13) a 4 = 0.1208 (6.14) a 6 = 0.3067 (6.15) a q = 5.571 − 1.086/P (6.16) To ensure that all terms in the equation for T arise at the same asymptotic order, R 2 = a 7 D, where a 7 ≈ 0.4338. After rescaling, we can write (6.11) as which is very similar in form to (2.7) considered in section 2.2. Following the asymptotic procedure outlined in section 2.2, we find that with w(x, t) as in (2.10), the amplitudes of the pattern and large-scale modes satisfy If −2 < s < 1, the coefficient of the term A|A| 2 in (6.19) is negative, and the bifurcation to pattern (in the absence of B) is supercritical. Equations (6.19) and (6.20) are then equivalent to (2.5) and (2.6), after an appropriate rescaling.
From this analysis we find that the bifurcation to rolls is supercritical if −2 < s < 1, and that all rolls are unstable if −2 < s < −

Application to thermosolutal convection
As a second physical example, we examine thermosolutal convection between boundaries at which the heat flux is fixed. We suppose that the fluid is slightly non-Boussinesq, so that the thermal diffusivity is a weak function of the local temperature: The layer is heated from below and in the basic state the saltier water lies over the fresher water, so both mechanisms are destabilizing. Both the thermal and the solutal gradients are vertical and uniform in the absence of convection. The thermal gradient is characterized by a (thermal) Rayleigh number R > 0 and the salinity gradient by a solutal Rayleigh number S < 0. When 0 > S > S c , the bifurcation to convection takes place at zero wavenumber; when S c > S, the critical wavenumber is finite. We analyse the case S = S c + 2 S 2 .
The natural scale for lengths is now the depth of the layer, and that for temperature perturbations is the temperature difference across the layer in the basic state (similarly for the salt concentration). With these scales, the dimensionless governing equations may be written in the form where is the perturbation salinity field, τ = κ s /κ and κ s is the salt diffusivity. The thermal Rayleigh number is now where α is the coefficient of cubic expansion and dT * /dz * is the (dimensional) temperature gradient in the absence of convection; the solutal Rayleigh number S is defined similarly. As in the previous example, we use free-rigid horizontal boundary conditions, with (6.5) and (6.6) applied for H = 1 2 , since the fluid layer now has unit dimensionless depth. With these boundary conditions, S c = −8352τ/41 ≈ −203.707τ .

Application to magnetoconvection
Convection in a vertical magnetic field has been a subject of a great deal of research (see, for example, Weiss 1982, Matthews 1999). The problem is motivated by the convection zone of the Sun, where the magnetic field is moved around by the fluid and regions of strong magnetic field inhibit the fluid motion. In the nonlinear regime magnetoconvection exhibits a very wide range of complicated phenomena. In magnetoconvection, the conservation law that leads to the large-scale mode is that the total flux of the magnetic field through the layer is conserved. The governing equations for two-dimensional magnetoconvection are (Proctor and Weiss 1982), in the notation of the preceding section, The magnetic field is defined through a perturbation flux function φ by (B x , B z ) = (−φ z , 1 + φ x ). The parameter Q is the Chandrasekhar number, a dimensionless measure of the magnetic field strength, while ζ is the ratio of magnetic diffusivity to thermal diffusivity. The most commonly used boundary conditions are that the upper and lower boundaries are stress-free, fixed temperature and with no horizontal magnetic field, so These boundary conditions are convenient in the sense that the linear eigenfunctions are trigonometric.
The large-scale mode in the problem arises through a rearrangement of the vertical magnetic field. There is a linear eigenfunction in which φ = φ(x, t) which obeys the diffusion equation φ t = ζ φ xx , so the growth rate λ of a mode with wavenumber k is simply λ = −ζ k 2 . Note that this mode is independent of the velocity and temperature perturbation, since if φ = φ(x, t) there is no forcing term in (6.37). Although this mode is usually ignored in the analysis of magnetoconvection, it plays a crucial role when the horizontal size of the domain is large. Note that although there are nearly marginal modes at small k and at order-one values of k, the graph of growth rate versus wavenumber is not as shown in figure 1, because the largescale mode and the order-one mode correspond to two distinct branches of eigenfunctions.
The weakly nonlinear analysis proceeds by including the large-scale mode. The general analysis of section 2 suggests that the mean mode should be included at order 2 . In fact, the mean vertical magnetic field is of order 2 and so the corresponding flux function φ is of order . The appropriate asymptotic expansions near the stationary bifurcation to convection are therefore where θ 1 and φ 1 are constants, readily determined from the linear stability problem, X = x, T = 2 t and R = R c + 2 R 2 . At order we obtain θ 1 = k/(π 2 + k 2 ), φ 1 = π/ζ(π 2 + k 2 ) and the well known result (Chandrasekhar 1961) R c = (π 2 + k 2 ) 3 + Qπ 2 (π 2 + k 2 ) k 2 .
(6.44) From (6.44) it follows that as R is increased, convection first occurs with a critical wavenumber k c which is related to Q by Henceforth it is assumed that k = k c ; note that the envelope function A(X, T ) permits convection with a wavenumber k = k c + O( ). At order 2 , a term in θ proportional to sin(2πz) and a φ term proportional to exp(2ikx) are obtained. At third order, evolution equations for A and are found by applying solvability conditions. Full details of the calculation will be given elsewhere. The resulting amplitude equations are where a 1 , . . . , a 4 are known constants, given by cumbersome functions of the parameters; a 2 and a 4 are always positive. After setting X = B, equations (6.46) and (6.47) are equivalent to (2.5) and (2.6). Equations (6.46) and (6.47) have not been derived previously, but the underlying processes of the coupling between the convection rolls A and the large-scale mode are well understood (see, for example, Proctor and Weiss 1982). Since a 4 > 0, the coupling term in (6.46) simply reflects the fact that the magnetic field suppresses convection. The coupling term in (6.47) describes the process of flux expulsion: the magnetic field strength decreases in regions where the convection is strong.
The bifurcation at the onset of convection is supercritical (i.e. a 3 > 0) provided that k 4 ζ 2 (π 2 + k 2 ) > π 2 (2k 2 − π 2 )(π 2 − k 2 ) (6.48) where Q has been eliminated in favour of k using (6.45). The result (3.9) can now be applied to determine whether a regular pattern of rolls is unstable. It is found that all rolls are unstable if k 4 ζ 2 (π 2 + k 2 ) < π 2 (2k 2 − π 2 )(k 2 + 3π 2 ) (6.49) where Q has again been eliminated using (6.45). Note that neither of the conditions (6.48) and (6.49) involves the Prandtl number P . Using k as a parameter, the region of stability of rolls can be plotted in the (Q, ζ )-plane ( figure 15). In the small region where the rolls are subcritical the weakly nonlinear analysis gives no useful information about the pattern. However, there is a large region of parameter space in which the bifurcation is supercritical but all rolls are unstable. In this region we expect to find that the stable pattern at onset consists of rolls that are modulated on a large length scale. Solutions of this type, but fully nonlinear, have recently been obtained through numerical simulations of magnetoconvection (Blanchflower 1999). Note that our analysis is confined to the onset of convection at a stationary bifurcation. For certain parameter values (in general, when Q is large and ζ is small) the onset of magnetoconvection is oscillatory; since this depends on the value of P , this condition is not shown in figure 15.

Discussion
We have derived amplitude equations (2.5) and (2.6) for the primary bifurcation to a stationary pattern in a one-dimensional system with a conservation law for a quantity that does not change sign under reflection. In this case the conservation law leads to a neutral large-scale mode which must be taken into account for a correct description of the behaviour near onset. The usual Ginzburg-Landau equation (Melbourne 1999) is coupled to an equation describing the large-scale evolution of the conserved quantity. We have shown how these equations may be obtained from symmetry considerations alone, and have also derived them directly from the governing equation(s) of the system.
In considering the small-amplitude solutions of the amplitude equations, we have shown how periodic patterns may become unstable to a new type of instability. This new instability is quite different from the usual Eckhaus instability in three ways. It is an amplitude-driven instability (whereas the Eckhaus instability is phase-driven) and it is in general supercritical. Most importantly, for certain parameter values all patterns are unstable. Our instability is also different from that which occurs in the case where the conserved quantity is velocity-like, and the reflection symmetry acts as −1 (such as in systems with Galilean invariance), where periodic patterns are always unstable and chaotic dynamics is observed at onset (Tribelsky and Velarde 1996, Tribelsky and Tsuboi 1996.
Beyond the instability, there are stable solutions in the form of a strongly amplitudemodulated periodic pattern, the envelope of which we have calculated in terms of Jacobi elliptic functions. The envelope may commonly be approximated by a sech profile. Numerical simulations reveal that these strongly modulated patterns may be stable in the full governing PDE, not just within the framework of the amplitude equations. Numerical simulations of the governing PDE also demonstrate that the amplitude equations represent a good approximation provided the amplitude is sufficiently small, but that higher-order terms play an important role as the amplitude increases; these higher-order terms can lead to significant modifications in the bifurcation diagram and a broadening of the envelope of the pattern.
Our work demonstrates that localized solutions are possible when a Ginzburg-Landau equation is coupled to an equation for a mean field, even when the coefficients of the equations are real and when the bifurcation is supercritical. Previous work has required complex coefficients (Riecke 1992a, b, Sakaguchi 1993), a subcritical bifurcation (Herrero andRiecke 1995, Sakaguchi andBrand 1996) or both. Our amplitude equations apply when the bifurcation from the basic state is stationary; the case of oscillatory bifurcation to travelling-wave solutions has been analysed in a series of papers by Riecke and co-workers (Riecke 1992a, b, 1996, Herrero and Riecke 1995 and by others (Sakaguchi 1993, Barthelet andCharru 1998), who find a variety of localized and front-type solutions. They find that the coupling to a mean field (our B) can have a significant effect on solutions and their stability, and can stabilize localized travelling-wave trains and holes.
We have applied our amplitude equations to three convection problems, leading to novel instabilities in these systems. Further details of these applications will be given elsewhere. In the first two, the heat flux is fixed on the boundaries, and the total heat in the layer is conserved, which leads to a large-scale mode in the form of a depth-averaged temperature perturbation. However, in order for our amplitude equations (2.5) and (2.6) to apply, this large-scale mode must influence the convection rolls; this is achieved by using either a nonlinear equation of state or a temperature-dependent diffusivity. In the third example, magnetoconvection, the conserved quantity is the vertical magnetic field; in this case the large-scale mode and the pattern mode correspond to two distinct branches of eigenfunctions of the linear problem. Our modulated solutions may be related to the isolated 'convectons' found recently in magnetoconvection (Blanchflower 1999). Another application of our work is in rotating convection, where the conserved quantity is momentum .
Our amplitude equations (2.5) and (2.6) are identical to those derived by Coullet and Iooss (1990), who considered the secondary stability of a one-dimensional cellular pattern. In that case, A corresponds to the amplitude of a perturbation to the basic pattern, and B to a gradient in the phase of the pattern; the amplitude equations then describe evolution near a stationary bifurcation. Coullet and Iooss (1990) found that when a pattern of wavelength L becomes unstable to a perturbation of wavelength 2L, the amplitude equations (2.5) and (2.6) apply, with A real; when the spatial period of the perturbation and the basic pattern are irrationally related, A may be complex. Our results suggest that if a pattern of wavelength L undergoes a supercritical subharmonic bifurcation then there may be stable states where the original pattern is largely undisturbed in some regions but where the perturbation is strong in others. Furthermore, provided the secondary bifurcation is supercritical, it should be possible to find stable states consisting of two domain types: in one domain the pattern is essentially the basic cellular state, in the other the pattern is strongly perturbed.
It would be of interest to investigate how this theory may be extended to problems in two or three space dimensions; in particular, we would like to know whether solutions may be obtained that are localized in space. Such localized solutions have been found experimentally by VanHook et al (1995VanHook et al ( , 1997. Localized solutions after an oscillatory bifurcation have recently been examined in two dimensions by Riecke and Granzow (1998), who found that such solutions may be obtained in a pair of Ginzburg-Landau equations coupled to a mean-field equation. Of particular interest is their result that stable localized states can be found, as here, when the primary bifurcation is supercritical and the coefficients of the amplitude equations are real.