Effect of driving on coarsening dynamics in phase-separating systems

We consider the Cahn-Hilliard (CH) equation with a Burgers-type convective term that is used as a model of coarsening dynamics in laterally driven phase-separating systems. In the absence of driving, it is known that solutions to the standard CH equation are characterized by an initial stage of phase separation into regions of one phase surrounded by the other phase (i.e., clusters or drops/holes or islands are obtained) followed by the coarsening process, where the average size of the structures grows in time and their number decreases. Moreover, two main coarsening modes have been identified in the literature, namely, coarsening due to volume transfer and due to translation. In the opposite limit of strong driving, the well-known Kuramoto-Sivashinsky (KS) equation is recovered, which may produce complicated chaotic spatio-temporal oscillations. The primary aim of the present work is to perform a detailed and systematic investigation of the transitions in the solutions of the convective CH (cCH) equation for a wide range of parameter values, and, in particular, to understand in detail how the coarsening dynamics is affected by an increase of the strength of the lateral driving force. Considering symmetric two-drop states, we find that one of the coarsening modes is stabilized at relatively weak driving, and the type of the remaining mode may change as driving increases. Furthermore, there exist intervals in the driving strength where coarsening is completely stabilized. In the intervals where the symmetric two-drop states are unstable they can evolve, for example, into one-drop states, two-drop states of broken symmetry or even time-periodic two-drop states that consist of two traveling drops that periodically exchange mass. We present detailed stability diagrams for symmetric two-drop states in various parameter planes and corroborate our findings by selected time simulations.


Introduction
In recent years, there has been a renewed interest in the convective Cahn-Hilliard (cCH) equation as a model of coarsening dynamics in driven phase-separating systems. In the present study, we consider the following one-dimensional cCH equation that contains an additional nonlinear driving term of Burgers type: Here, u(x, t) is the order parameter field, with x and t denoting the spatial coordinate and time, respectively, and D is the driving strength. This equation was derived, for example, by Golovin et al. [16,17] as a model for a kinetically controlled growing crystal surface with a strongly anisotropic surface tension. In such a context, u is the surface slope and D is the growth driving force proportional to the difference between the bulk chemical potentials of the solid and fluid phases (see also Liu and Metiu [29] for modelling of growing crystal surfaces). Equation (1) was also obtained by Watson [51] as a small-slope approximation of the crystal-growth model obtained by Di Calro et al. [9] and Gurtin [20]. Related models have also been derived, for instance, in the context of epitaxial growth (see, for example,Šmilauer et al. [39]) and liquid droplets on inclined planes (see, for example, Thiele and Knobloch [44,45], Thiele [41]).
In the absence of driving, the cCH equation reduces to the standard CH equation [1,33], that was proposed as a model to describe phase separation (or spinodal decomposition) of two-component mixtures (see, for instance, Cahn [3,4,5], Cahn and Hilliard [6,7]). Note that the standard CH equation can be written in the following general gradient-dynamics form: where δ/δu denotes the variational derivative. The free energy F [u] is given by and ϕ(u, u x ) = 1 2 u 2 x + f (u) is the energy density, with the first term being the squaregradient term that penalizes interfaces and with the double-well potential f (u) = 1 4 u 4 − 1 2 u 2 as the local free energy. The initial dynamics of the solutions of the standard CH equation from a perturbed homogeneous state is characterized by separation into regions corresponding to different components, i.e., clusters (drops/holes or islands) of one phase surrounded by the other phase, or labyrinthine patterns of the two phases. However, after this initial stage of evolution, these structures slowly grow in size and their number decreases, i.e., the structure coarsens. In the following we refer to the structures as "drops".
Two main modes of coarsening have been identified, namely, coarsening by volume transfer and by translation. In coarsening by the volume transfer mode (which is also known as Ostwald ripening [35]), the centres of the drops remain fixed in space, while the sizes of the drops change -some grow in time, while others decrease in size and, eventually, disappear. In coarsening by the translation mode, the centres of the drops do not remain fixed, and coarsening occurs due to motion and merging of the drops. The coarsening process continues until only a single large drop remains. For a more detailed discussion of coarsening for the CH and related equations see, for example, Onuki [34], Desai [8], Thiele et al. [42], and Pototsky et al. [37]).
In the limit of strong driving, the cCH equation reduces to the well-known Kuramoto-Sivashinsky (KS) equation [25,38]. Indeed, substituting u =ũ/D into (1) and taking the limit D → ∞, one obtains the KS equation forũ (see, for example, Golovin et al. [18]). In contrast to the solutions of the CH equation, the long-time dynamics of the solutions of the KS equation is characterized by complicated chaotic spatio-temporal oscillations [21,22,40]. Thus, as the driving force is increased from zero to large values, there must appear transitions leading from the coarsening dynamics typical of the standard CH equation to complicated chaotic oscillations typical of the KS equation. We note that coarsening dynamics for the cCH equation has been studied in the limit of a weak driving force numerically by Emmott and Bray [13] and Golovin et al. [18] and analytically by Watson et al. [52], and for moderately large driving force by Podolny et al. [36], and scaling laws for the average separation between the successive phases as a function of time have been obtained. Zaks et al. [53] reported that driving can be used to stop coarsening for certain parameter values. Some stationary solutions of the cCH equation have been analysed by Korzec et al. [24]. We also note that Eden and Kalantarov [12] demonstrated the existence of a finitedimensional inertial manifold for the cCH equation. The main aim of the present work is to perform a detailed and systematic investigation of the transitions in the solutions of the cCH equation for a wide range of parameter values as the driving force is increased and to construct detailed stability diagrams in the parameter planes. Finally, note that similar transitions with increasing lateral driving strength have been investigated for various thin-film equations [45,47]. The place of the cCH and thin-film equations in a classification of one-field equations based on mass conservation and variational character is discussed in the introduction of [14].
The rest of the present work is organized as follows. In Sect. 2, we discuss basic background on the cCH equation. In Sect. 3, we discuss some theory behind singleinterface (i.e., front) and double-interface (i.e., drop) solutions. We present the results of numerical continuation of periodic drop solutions in Sect. 4. First, we discuss the results of numerical continuation with respect to the domain size for different values of the mean concentration, and then we analyze how the driving force affects inhomogeneous solutions of the CH equation. In Sect. 5, we present a systematic study of the linear stability properties of various spatially periodic traveling solutions of the cCH equation, and analyze the effect of driving on the coarsening modes of symmetric two-drop states. We produce detailed bifurcation diagrams additionally including two-drop states of broken symmetry and time-periodic two-drop states that consist of two drops that periodically exchange mass. We present detailed stability diagrams for symmetric two-drop states in various parameter planes. In addition, we support the numerical continuation results by selected time simulations. Finally, in Sect. 6 we present our conclusions.

The convective Cahn-Hilliard equation
As we focus on analyzing solutions that are stationary or time-periodic in a moving frame, it is convenient to rewrite equation (1) in a frame moving with velocity v, i.e., We are primarily interested in analyzing solutions on a spatially periodic domain, say x ∈ [0, L], and we note that u(x, t) is a conserved quantity, i.e., the mean valuē u = 1 2L L −L u dx is constant. Note that due to the symmetry (D, u) → (−D, −u), it is sufficient to only consider nonnegative values of D. In addition, the symmetry (x, u) → (−x, −u) implies that it is sufficient to only consider nonnegative mean values u. For the rest of the manuscript, we therefore assume that D ≥ 0 andū ≥ 0.
To analyze the linear stability of a spatially uniform solutionū, we consider a small perturbation of the form ∝ exp(ikx + βt) and, after linearizing equation (4), obtain the following dispersion relation: Thus, the growth rate w(k) = Re β(k) of a small-amplitude sinusoidal wave of wavenumber k is as for the standard CH equation, and the phase speed is −Imβ(k)/k = Dū − v. By solving equation w(k c ) = 0, we find the cutoff wavenumber k c : This solution exists only when 1 − 3ū 2 > 0, i.e., when |ū| < 1/3. In this case, there is a band of unstable wavenumbers, k ∈ (0, k c ). Otherwise, if |ū| ≥ 1/3, we find that w(k) < 0 for all k > 0, and we obtain the linearly stable case. Note that these uniform states may still be nonlinearly unstable.

Front and one-drop solutions
In this section, we discuss single-interface solutions (i.e., kinks and anti-kinks, or fronts) and double-interface solutions (i.e., one-drop solutions) of the standard and convective CH equations. For this purpose, we consider the cCH equation on an infinite domain. A front solution is a solution that approaches two different constants as x → ±∞. Let us denote these constants by u a and u b for x → −∞ and x → +∞, respectively. If u a < u b we obtain a so-called kink solution. If u a > u b , we obtain an anti-kink solution.
Here we call both "front". A double-interface (or one-drop) solution, is a solution that approaches the same constant (say u b ) as x → ±∞, but has a region where it approaches a different constant (say u a ), so that this region is macroscopic, i.e., sufficiently long compared to the lengths of the regions where the solution first transitions from u b to u a and then from u a to u b . Such a solution may be considered as a superposition (with small correction) of well-separated kink and anti-kink solutions. If u a > u b , we obtain a solution in the form of a drop, otherwise, we obtain a solution in the form of a hole. We note that our discussion of single-and double-interface solutions (i.e., front and one-drop states) partly follows the discussions of Emmott and Bray [13], Golovin et al. [18], Korzec et al. [24], Zaks et al. [53]. For the standard CH equation, it is well-known that front solutions have zero speed and are of the form (see Novick-Cohen and Segel [32]) There also exist periodic drop and hole solutions of drops/holes of arbitrarily large size.
In the course of our work we consider domain sizes where one or two periods of a periodic solution fit. Note that in the latter case the solution has a discrete translation symmetry with respect to a shift of half the domain size. We refer to the respective solutions as "one-period" and "two-period" states. Alternatively we refer to them as "one-drop" and "symmetric two-drop" states. For the cCH equation, a solution u 0 that is stationary in a frame moving at speed v satisfies the equation which, when integrated once, becomes where C 0 is a constant of integration that corresponds to the flux in the moving frame. Equation (10) can be rewritten as a three-dimensional dynamical system by introducing the functions y 1 = u 0 , y 2 = u 0x and y 3 = u 0xx : We note that this dynamical system preserves phase space volume, since the divergence of the corresponding vector field (or, equivalently, the trace of the Jacobian matrix) is identically zero. The fixed points of (11)-(13) satisfy y 2 = y 3 = 0 and Assuming that there exists a front solution that connects uniform solutions u a and u b we obtain that A front solution then corresponds to a heteroclinic orbit connecting the fixed point (u a , 0, 0) along the unstable manifold of u a , denoted by W u (u a ), to the fixed point (u b , 0, 0) along the stable manifold of u b , denoted by W s (u b ). In fact, it is known that equation (1) has exact kink and anti-kink solutions which have v = 0 and which are given by (see Golovin et al. [18]) for ±, respectively. Thus, for these solutions, u a = −u + and u b = u + for the case of the kink, and u a = −u − and u b = u − for the case of the anti-kink. Note that these solutions reduce to the front solutions of the standard CH equation when D = 0. Note also that kink solutions exist only for D <D ≡ √ 2. The eigenvalues for the fixed points (u a,b , 0, 0) satisfy Figures 1(a) and 1(b) show the dependence on D of the real and imaginary parts, respectively, of the eigenvalues for u + , which can be found analytically (see Zaks et al. [53]): It can be seen that λ 1 is real and negative for all D ∈ (0,D). The other two eigenvalues, λ 2,3 , have positive real parts and are real for D ∈ (0, D) and complex conjugate for D ∈ ( D,D), where D = √ 2/3. This was first pointed out by Podolny et al. [36]. Note that as D →D, u + → 0, and λ 1 → 0, λ 2,3 → ±i. The eigenvalues for −u + are −λ 1,2,3 . We conclude that dim(W u (u + )) = 2, dim(W s (u + )) = 1, dim(W u (−u + )) = 1, dim(W s (−u + )) = 2. Therefore, there is a neighbourhood of the point ( − u + , u + ) in the (u a , u b )-plane in which the kink solution exists only for u a = −u + and u b = u + , and this kink solution is u + 0 (x), given by (16). Note that there may exist other isolates points in the (u a , u b )-plane which correspond to kink solutions, and some of these solutions were computed by Zaks et al. [53]. Regarding anti-kink solutions, we conclude that there exists a one-parameter family of such solutions corresponding to a curve in some neighbourhood of the point (u + , −u + ) in the (u a , u b )-plane for each D ≥ 0. We also note that although kink solutions exist for D ∈ [0,D) and anti-kink solutions exist for any D ≥ 0, the flat parts of such solutions become linearly unstable (in the sense of temporal linear stability analysis) on a sufficiently long spatial domain when D >Ď = 2 √ 2/3. Double-interface (and, in fact, many-interface) solutions can be analysed, for instance, by using the Shilnikov-type approach, see, e.g., Glendinning and Sparrow [15], Guckenheimer and Holmes [19], Knobloch and Wagenknecht [23], Kuznetsov [26], Tseluiko et al. [50]. Indeed, let us consider, for example, the case D ∈ ( D,D). For the points (−u + , 0, 0) and (u + , 0, 0), there exists a heteroclinic orbit connecting the first point to the second (corresponding to the kink solution) and heteroclinic orbits connecting the second point to the first (corresponding to the anti-kink solutions). Then, we expect that there exists an infinite but countable number of the values of u + k , k ∈ N, in the neighbourhood of u + for which there exist homoclinic orbits for the fixed points (−u + k , 0, 0) that pass near (u + , 0, 0). Such orbits then correspond to drop solutions, and such drop solutions differ by their lengths. Then, since Re λ 2,3 < −λ 1 (note that, given that the dynamical system (11)-(13) preserves phase space volume, this inequality is automatically satisfied), Shilnikov's theory implies the existence of an infinite but countable number of subsidiary homoclinic orbits in the vicinity of the primary orbit that pass near (u + , 0, 0) several times before achieving homoclinicity. Such subsidiary homoclinic orbits correspond to multi-drop solutions. In addition, Shilnikov's theory implies the existence of an infinite number of periodic orbits in the vicinity of the primary homoclinic orbits. Such periodic orbits correspond to periodic arrays of drops. In a similar way, we can analyze hole solutions and can obtain finite or periodic arrays of hole solutions (of course, periodic arrays of hole solutions are equivalent to periodic arrays of drop solutions). We note, however, that for D >D, kink solutions do not exist, and, therefore, the double-interface or multi-interface solutions that are typical of the standard CH equation do not exist for such values of D. Nevertheless, there may still exist homoclinic orbits corresponding to pulse or anti-pulse solutions (also referred to as hump or hollow solutions, respectively). Shilnikov's theory then implies the existence of bound states or (a)periodic arrays of such pulses or anti-pulses. These solutions may still be characterized as localized drops or holes, but the nature of these solutions is different from that for the standard CH equation.

One-drop solutions for the standard CH equation
In this section, we review the structure of one-period solutions for the standard CH equation, when D = 0, for different values ofū, and, in particular, we compute solutions forū = 0.4, 0.55 and 0.6. Note that much more exhaustive results are available in the literature for the standard CH equation, e.g., [30,31,32,33,46]. We characterize the solutions by their norms δu 0 = (1/L) Forū = 0.55, the primary bifurcation at L c = 20.66 is subcritical. The branch of nonuniform solutions initially follows to decreasing values of the domain size L and is unstable up to the saddle-node bifurcation at L = L s ≈ 13.818. After this point, the branch turns back and becomes stable. The exact value ofū at which the bifurcation switches from supercritical to subcritical can be obtained by the weakly nonlinear analysis given in the Appendix, see equation (35). It turns out that this value isū * = 1/ √ 5 ≈ 0.45. Note that forū = 0.55 the energy of the nonuniform solution first increases monotonically, up to the saddle-node bifurcation, and then decreases monotonically. It remains positive up to a certain value of the domain size, L m ≈ 14.30 between L s and L c , and then becomes negative. The point L = L m is the so-called Maxwell point. At this point, both linearly stable solutions, i.e., the uniform solution and the nonuniform solution with the larger value of the norm, have the same value of the energy. For L ∈ (L s , L m ), the uniform solution has lower free energy, whereas for L > L m , the nonuniform solution has lower free energy. In Fig. 2(f), we can see that, as L increases, the solution profiles forū = 0.55 behave as in the case ofū = 0.4, except that now the width of the drop grows as 0.775 L and the width of the hole grows to 0.225 L as L increases, so that the mean value remains equal toū = 0.55.
As mentioned above, sinceū = 0.6 > 1/ √ 3, the flat solutionū = 0.6 is linearly  stable for any L, i.e., there is no primary bifurcation on the uniform solution. To produce the branch of nonuniform solutions, we can first compute the branch of nonuniform solutions for, e.g.,ū = 0, and then select a solution on this branch at a sufficiently large value of L (e.g., L = 100). We then keep L fixed and perform a continuation inū, until we reach the valueū = 0.6. This produces the nonuniform solution forū = 0.6 at L = 100. After that, we again keepū fixed and perform a continuation in L, going in both directions, which produces the whole branch of nonuniform solutions. We can observe that the branch of nonuniform solutions has a turning point at L = L s ≈ 16.327. For each L > L s , there are two nonuniform solutions, one is unstable and is of smaller norm while the other one is stable and is of larger norm. The energy of the linearly unstable nonuniform solution monotonically decreases from some positive value to zero as L increases from L s . Whereas the energy of the linearly stable nonuniform solution decreases monotonically from a positive value to negative values crossing zero at the Maxwell point, L m ≈ 17.466. In Fig. 2(i), we can see that, as L increases, the solution profiles of the upper branch of nonuniform solutions forū = 0.6 behave as in the previous cases, except that now the width of the drop grows as 0.8 L and the width of the hole grows to 0.2 L as L increases, so that the mean value remains equal toū = 0.6. The behaviour of the solutions of the lower branch of nonuniform solutions is, however, different. As L increases, their amplitude decreases approaching a constant value, and the width in the rescaled variable x/L also decreases approaching a constant value in the original variable x, so that the solution tends to an anti-pulse shape. Note that a recent study investigates how the Maxwell construction at phase coexistence emerges from bifurcation diagrams like the ones in Fig. 2 for finite-size systems when approaching the thermodynamic limit [49].

One-drop solutions for the cCH equation
We now consider how the driving force affects the one-period steady traveling-wave solutions of the CH equation. We use both the driving force, D, and the domain size, L, as the control parameters and consider three cases,ū = 0.4, 0.55 and 0.6, as we did for the standard CH equation. Note that bifurcations of periodic one-drop solutions of the cCH equation were previously analysed in detail by Zaks et al. [53] but only for u = 0. We first notice that the changeover from supercritical to subcritical primary bifurcation that we have discussed in the previous section, is affected by D. Using the weakly nonlinear analysis presented in the Appendix (see equation (35)) we can show that the line separating the regions in the (ū, D)-plane where the primary bifurcation is supercritical or subcritical is given by the equation see Fig. 3. To be more precise, for a fixed D > 0, the primary bifurcation (if it exists) is subcritical ifū ∈ (ū * ,ū * * ) and supercritical otherwise, wherē  We remind here that we consider only nonnegative values ofū and D. Equivalently, for a fixedū, the bifurcation is subcritical if D < D c and supercritical otherwise, where Note that the expression under the square root in (21) is positive only when 1/ √ 5 <ū < 1/ √ 3 (considering nonnegative values ofū), i.e., the driving force can switch the type of the bifurcation only when 1/ √ 5 <ū < 1/ √ 3. If 0 ≤ū < 1/ √ 5, the primary bifurcation is supercritical for any value of the driving force. We also remind that ifū > 1/  Fig. 4(a) shows that for all the considered values of L, the norm δu 0 is a monotonically decreasing function of D. In Fig. 4(b), we can observe that for u = 0.4 all the branches of spatially nonuniform solutions (when L is used as the control parameter) bifurcate supercritically from the homogeneous branch at L = L c , consistent with the weakly nonlinear analysis discussed above. We can also observe that for small values of D, the norm increases monotonically and tends to a constant as L increases. As D increases, the norm becomes a nonmonotonic function of L but still tends to a constant as L increases (see, for example, the line for D = 0.8). For even larger values of D this behaviour changes -the norm first monotonically increases, then it may undergo a few oscillations before monotonically decreasing. This is consistent with the fact that the one-drop or multi-drop solutions that are typical of the standard CH equation do not exist for D >D ≡ √ 2, as discussed at the end of Sect. 3. Instead, we obtain solutions of a different nature, namely, localized traveling-wave solutions, whose width remains almost unaffected by the increasing domain size, and whose norm, therefore, tends to zero according to the law 1/ √ L as L increases (this has been verified numerically). In Fig. 4(c), we can see that for smaller values of D, the solution profile has a drop shape. As D increases, the solution becomes flatter and the drop is deformed, namely, a ridge develops at the right-hand side of the drop. For larger values of D, the ridge first becomes more pronounced and then decreases in amplitude. Further, there appear additional visible oscillations in the profile that decay upstream. The appearance of such oscillations can be understood through the spatial linear stability analysis. Also, it can be observed that for any value of D, the width of the drop in the rescaled coordinate x/L increases as D increases and the cavity narrows down. In fact, as discussed above, proper drop solutions exist only for D < √ 2, and the solution profiles for D > √ 2 should rather be classified as localized anti-pulse or hollow (or as pulse or hump solutions for negative values ofū) than as drop solutions. Forū = 0.55, using D as the control parameter, we can see in Fig. 4(d) that for L < L c the branches start at D = 0, then have saddle-node bifurcations at some positive values of D, and then return to D = 0. As L increases, the saddle-node bifurcation shifts to the left. For L = L c , the branch starts at D = 0, then has one saddle-node bifurcation at a positive value of D. However, it does not go back to D = 0. Instead, the branch terminates at the horizontal axis, where δu 0 = 0, at some positive value of the driving force, D = D c ≈ 1.3064. For L > L c , the branches start at D = 0, but are characterized by two saddle-node bifurcations. After the second saddle-node bifurcation, the branch continues to infinity. For sufficiently large L, both saddle-node bifurcations annihilate each other, as is below discussed in more detail. In fact, the value D c is precisely the value at which the primary bifurcation changes from subcritical to supercritical when the domain size L is used as the control parameter, as given by equation (21). In Fig. 4(e), when L is used as the control parameter, we can observe that forū = 0.55 the primary bifurcation is indeed subcritical for D < D c while it is supercritical otherwise, in agreement with the weakly nonlinear analysis. When D < D c , there is only one saddle-node bifurcation. On the other hand, when D > D c , there are two saddle-node bifurcation -the branch bifurcates supercritically from the uniform solution, then turns back at the first saddle-node bifurcation, and then turns again at the second saddlenode bifurcation and goes off to infinity. This is consistent with the results presented in Fig. 4(d), which show that for moderately large values of L > L c there exist three different solutions for a certain range of the driving force D. In Fig. 4(f), we can see that for L = 25 andū = 0.55 there exist three different solutions at the same values of D between the two saddle-node bifurcations that occur at D s1 ≈ 1.96 and D s2 ≈ 1.62. For D = 1.75 the solutions with larger and smaller amplitudes belong to the respective upper and the lower parts of the branch shown in Fig. 4(e) and are stable, whereas the solution with the intermediate value of the amplitude belongs to the middle part of the branch and is unstable. As D increases further, we can see that the solution becomes flatter, and the ridge that was pronounced for smaller values of D decreases in amplitude. We also remind here that the solution profiles that we observe for D > √ 2 and sufficiently large L should be classified rather as anti-pulse or hollow solutions than drop solutions. Forū = 0.6, using D as the control parameter, we can see in Fig. 4(g) that for all the considered values of L, the branches start at D = 0 then have one saddlenode bifurcation at some positive values of D and return to D = 0. In Fig. 4(h), when L is used as the control parameter, we can observe that forū = 0.6 there are no primary bifurcations for all the values of D, and we always find a saddle-node bifurcation. For smaller values of D, the upper parts of the branches monotonically increase as L increases, whereas for larger value of D, the upper parts of the branches first monotonically increase and then monotonically decrease. In Fig. 4(i), we can see that whenū = 0.6 and L = 25 there are two different solutions for D < D s ≈ 1.31. In particular, for D = 0.1 and 1 the solutions with larger amplitudes belong to the upper part of the branch for L = 25 shown in Fig. 4(g) (these solutions are stable), whereas solutions with smaller amplitudes belong to the lower part of this branch (these solutions are unstable).
From Fig. 4(d), it is difficult to infer where exactly the saddle-nodes appear. To understand this process better, we follow in Fig. 5(a) the loci of saddle-node bifurcations forū = 0.55 in the (D, L)-plane. The horizontal dotted line indicates the cutoff period L c = 2π/k c for the linear stability of the uniform solutionū = 0.55. We see that for L < L c there is only one saddle-node bifurcation. On the other hand, for L > L c , there are two saddle-node bifurcations. For sufficiently large L, the two saddlenode bifurcations annihilate each other. Figure 5(b) shows the loci of the saddle-node bifurcations forū = 0.6 in the (D, L) plane. We see that for all the values of L ≥ L sn , where L sn is the locus of the saddle-node bifurcation at D = 0 (cf. Fig. 2(g)-(i)), there is one saddle-node bifurcation.

Linear stability, coarsening and time-periodic behaviour of two-drop solutions
In this section, we construct detailed bifurcation diagrams of one-and two-drop solutions of the standard CH and cCH equations and study in detail linear stability properties and coarsening behavior of such solutions. We note that formerly coarsening dynamics of the cCH equation was analysed by Watson et al. [52] (for D 1) and Podolny et al. [36] (for any D < √ 2/3) who derived a nearest-neighbour interaction theory for phase boundaries (kinks and anti-kinks) and revealed an important role of kink triplets in the coarsening process. Namely, they showed that due to mass conservation binary coalescence of phase boundaries is not possible. However, when an anti-kink is located between two kinks, it attracts them leading to simultaneous annihilation of the triplet and formation of a single kink. Note that Watson et al. [52] and Podolny et al. [36] considered the cCH equation in the form where the sign in front of the convective term is flipped. Thus, due to the symmetry (D, u) → (−D, −u), this implies for our case annihilation of a triplet where a kink is located between two anti-kinks resulting in a single anti-kink. In our study, we consider a periodic systems, i.e., in the simplest coarsening process a two-period or symmetric two-drop solution transforms into a one-period or one-drop solution. We take a computational approach with the aim to construct detailed stability diagrams in the parameter planes and to analyse transitions in the behaviour of the solutions not only for D < √ 2/3 but also for larger values of D. Assuming that u 0 is a steady solution of (4) (i.e., a steady traveling-wave solution of (1)) and thatũ is a small perturbation, we obtain the following linearized problem forũ:ũ where L is the following linear differential operator with nonconstant coefficients: The stability of u 0 then depends on the spectrum of L, which typically consists of isolated eigenvalues of finite multiplicity, if L is defined on a finite periodic domain. Numerically, the eigenvalues can be computed directly using, e.g., a Fourier spectral method, or via numerical continuation, e.g., utilizing the continuation and bifurcation software Auto07p [10]. In addition to analysing steady traveling-wave solutions, we also construct branches of solutions that are time-periodic in a moving frame (modulated traveling waves, here referred to as time-periodic branches). Such solutions are also known as relative period orbits, see, e.g., [11]. We construct such branches using the procedure described in [28]. This allows us to obtain a more complete understanding of the various transitions in the solutions.

The case of the standard CH equation
First, we note that branches of two-period solutions for the standard CH equation can be obtained from the branches of one-period solutions discussed in Sect. 4.1 by considering the identical periodic solution in a domain twice as large as the period. Our calculations show that for the standard CH equation the resulting two-drop branch has no side branches. Therefore, there is actually no need to recompute the primary branches. However, one still needs to individually analyze the linear stability as it may change when going from the one-drop to the two-drop states. First, we note that zero is always an eigenvalue of the linearized problem with the eigenfunction given by For the translation mode, the drops move towards each other, and for the volume mode the volume of one of the drops decreases while the volume of the other one increases accordingly. The eigenvalues for these modes correspond to the perturbed zero eigenvalue. The larger the separation distances between the fronts are, the closer to zero these eigenvalues become. It is also interesting to note that the translation be shown that the linearized operator does not change. For the mean value −ū, we, therefore, again obtain two coarsening modes (which are exactly the same as for the mean valueū). However, when the steady solutions are superimposed with the eigenfunctions, the roles of the coarsening modes are interchanged, namely, the dominant coarsening mode is now the volume one and the other one is now the translation mode.

The case of the cCH equation
5.2.1. Symmetry breaking First, we employ continuation to compute branches of twoperiod solutions in dependence of the driving force D for several fixed values of L and u. As for the standard CH equation, branches of two-period solutions can in fact be obtained from the branches of one-period solutions (that were discussed in Sect. 4.2) by considering domain sizes that are twice the solution period. We call the resulting solution branches two-drop primary branches. The symmetric two-drop states on such branches have the discrete internal translation symmetry. Solution branches bifurcating from these primary branches in secondary bifurcations we call secondary branches. Secondary pitchfork bifurcations break the discrete translation symmetry and, therefore, result in solutions with a larger spatial period. Hence, if such solutions are stable, the corresponding secondary bifurcations are associated with coarsening of the pattern. However, we emphasize here that at least for D < √ 2 for a two-drop solution given on a domain of certain length there exists a one-drop solution of the period equal to that domain length, and true coarsening would correspond to evolution towards such a one-drop solution. For completeness of the bifurcation diagrams, we also include the branches of one-drop states.   Fig. 8(a) that for L = 25 there are two bifurcation points on the two-drop primary branch, and the secondary branches that start at these bifurcation points continue towards large values of D. There is also one Hopf bifurcation on the two-drop primary branch, and the time-periodic branch starting at this point also extends to large values of D. Figure 8(b) shows that for L = 35 there are five bifurcation points on the two-drop primary branch. Some of the secondary branches that start at these points reach large values of D and may continue to infinity, whereas secondary branches starting at other bifurcation points reconnect to the same primary branch. In particular, the secondary branches starting at bifurcation points 1, 2 and 3 continue to infinity, while bifurcation points 4 and 5 are connected to each other by a secondary branch.
Regarding the one-drop branches, we find that for L = 25, there exists one Hopf bifurcation, and the time-periodic branch emanating at this Hopf bifurcation extends to Regarding the two-drop branches, we observe in Fig. 11(a) that for L = 35 there are four bifurcation points and one saddle-node bifurcations on the two-drop primary branch. The secondary branches that start at these bifurcation points reconnect to the two-drop primary branch. Also, we denote the upper and the lower parts of the primary branch by letters α and β, respectively. We can observe that points 1 and 2 on the upper part are connected to points 4 and 3, respectively on the lower part. On the one-drop branch we find two saddle-nodes, but there are no other bifurcation points. Figure 11(b) shows that for L = 50 there are five bifurcation points and two saddle-node bifurcations on the two-drop primary branch. Some of the secondary branches that start at these points, reach large values of D and may continue to infinity, whereas secondary branches starting at other bifurcation points reconnect to the primary branch. We call the upper part of the primary branch (up to the first saddle node) part α, the part connecting the two saddle nodes part β, and the lower part (starting from the second saddle node) part γ. We find that the secondary branch starting at bifurcation point 1 on part α continues to infinity, while bifurcation point 2 on part α is connected to point 5 on part β, and bifurcation point 3 on part α is connected to point 4 on part β. For L = 50, we additionally find that there are two Hopf bifurcations on the two-drop primary branch, denoted by symbols I and II. It is interesting to note that these bifurcation points are not connected to each other by a time-periodic branch, and the time-periodic branches that emerge from these points do not extend to large values of D. Instead, these time-     Fig. 12(a).
periodic branches are connected to side branches (the dashed blue and red branches, respectively). This is confirmed in Figs. 12(a) and (b) for the time-periodic branches starting at points I and II, respectively. Moreover, the inset in Fig. 12(b) indicates a possible exponential snaking behaviour of the time-periodic branch -one saddlenode is clearly visible, and one more can be obtained by another zoom. We conjecture that the time-periodic branch starting at point I results from a Takens-Bogdanov-type codimension-2 bifurcation at the pitchfork bifurcation point 3 (we note that for the  usual Takens-Bogdanov bifurcation the time-periodic branch emerges from a saddlenode bifurcation, not from a pitchfork bifurcation, see, for example, Kuznetsov [26]). Similarly, the time-periodic branch starting at point II results from such a codimension-2 bifurcation, but at a pitchfork bifurcation that has, at the shown value of L, moved to larger values of D (or to infinity). The time evolutions over one period of solutions corresponding to points 1 and 2 shown by red diamonds in Fig. 12(a) are shown in Figs. 13(a) and(b), respectively. In both cases, the solution behaves as a superposition of two drops periodically exchanging mass. Panel (b) confirms that as the homoclinic bifurcation is approached, the temporal period increases, and now the mass-exchange events happen burst-like over relatively short time intervals while for most of the time the solution is a quasi-steady superposition of two drops of different sizes. Figures 14 and 15 show the real parts of the dominant eigenvalues along the one-period primary branches presented in Figs We generally observe that sufficiently strong driving may destabilize one-drop solutions if the domain size is sufficiently large. Figures 16(a) and (b) show the real parts of the dominant eigenvalues along the two-drop primary branches presented in Figs. 8(a) and (b). Figure 16(a) shows that for L = 25 there are two pitchfork bifurcation points to side branches, one Hopf bifurcation to a branch of time-periodic solutions, and there is a stable interval between the second bifurcation point to a side branch and the Hopf bifurcation point, i.e., between D ≈ 1.41 and D ≈ 2.21. Interestingly, this means that sufficiently strong driving D can prevent coarsening, resulting in a stable two-drop traveling-wave solution. However, increasing D further may again destabilize such a solution resulting in two drops periodically interacting with each other (note that coarsening is still prevented). These observations are corroborated by the time-dependent simulations shown in Fig. 17  although it should be pointed out that for D = 0 this functional is not anymore a Lyapunov functional and should not necessarily be minimized in the time evolution.) It can be observed that for D = 0.3, the solution initially evolves into a two-drop solution, but around t = 1500 the drops coarsen and a one-drop solution is obtained (a one-drop solution is linearly stable for this value of D, see Fig. 14(b)). In contrast, for D = 3, a two-drop solution remains stable during the course of evolution, which agrees with the theoretical prediction (a one-drop solution is also linearly stable for this value of D, see Fig. 14(b), so the long-time evolution of solutions depends on initial conditions). For D = 5, the solution again initially tends to evolve into a two-drop solution. But as is evident from the energy and norm plots, around t = 150, the drops start to oscillate, and the solution eventually evolves into a time-periodic state resembling two drops periodically exchanging mass. We note that a one-drop solution is also linearly stable for this value of D, see Fig. 14(b). So we expect that different initial conditions can lead to time-periodic solutions or one-drop traveling-wave solutions. Figure 18 shows the most unstable eigenmode u 1 superimposed with the primary two-drop solution u 0 forū = 0.4 and L = 25. The arrows indicate the directions in which the fronts are shifted (in the same way as in Fig. 7 for the standard CH equation). Panels (a) and (b) correspond to D = 0.005 and 0.1. An interesting observation is that for the smaller value of D the most unstable mode appears to be translational (in agreement with the D = 0 case), whereas for the larger value of D the mode seems to change into a volume mode. Thus, the driving force can change the type of coarsening.

Linear stability of two-drop primary branches and coarsening
In Fig. 16(b), we can see that for L = 35 there are five pitchfork bifurcation points to side branches and no Hopf bifurcations. We also see that there are two stable intervals in D, namely, 0.82 ≤ D ≤ 1.23 and 2.32 ≤ D ≤ 8.28. Figure 19 shows the most unstable eigenmode u 1 superimposed with the primary two-drop solution u 0 forū = 0.4 and L = 35. Panels (a) and (b) correspond to D = 0.1 and 9. We observe that both modes are apparently volume modes. Figures 20(a) and (b) correspond to L = 35 and 50, respectively, atū = 0.55 (cf. Figs. 11 (a) and (b)). In Fig. 20(a), the solid and dashed lines correspond to the real and complex (having nonzero imaginary parts) eigenvalues, respectively, for part a (the upper part) of the bifurcation curve shown in Fig. 11(a). However, we additionally introduce the dot-dashed lines that correspond to the real eigenvalues for part b (the lower part) of the bifurcation curve shown in Fig. 11(a). Note that for part b, the eigenvalues are real in the shown range. Note that the green dot-dashed line corresponds to the unstable eigenvalue that is inherited from the one-drop primary branch (that is unstable). We can see that for L = 35 there are no stable intervals for the driving force D, and, therefore, in this case coarsening cannot be stabilized by sufficiently strong driving. Figure 20(b) shows that for L = 50 there are two stable intervals on part a (the  upper part) of the bifurcation diagram shown in Fig. 11(b), namely, 0.72 D 0.90 and 1.21 D 1.76. Part b (the middle part of the bifurcation diagram) is unstable, and there is a stable interval on part c (the lower part) of the bifurcation diagram, namely, D 2.15.

Linear stability of secondary branches
In this section, we analyze the linear stability of the secondary branches. Figures 21(a) and (b) correspond to L = 25 at u = 0.4 (cf. Fig. 8(a)). Panels (a) and (b) correspond to the first and second secondary branches shown by the red and green dashed lines, respectively, in Fig. 8(a). For the first secondary branch, there are two saddle-node bifurcations, while for the second secondary branch there are no saddle-node bifurcations. We can observe that for both secondary branches there is at least one eigenvalue with a positive real part for all the values of D. Therefore, both branches are unstable for all D values. So, in a time evolution, a solution does not evolve into a solution on the secondary branch. Instead, it can evolve into a two-drop solution (if D belongs to the stable interval), or a one-drop solution, or a time-periodic solution -such time evolutions are shown in Fig. 17. Figure 22 corresponds to L = 35 atū = 0.4 (cf. Fig. 8(b)). Panels (a), (b), (c) It is interesting to note that there may exist other stable solutions, and, in particular, for the initial condition chosen for panel (c), we observe that the solution evolves into a three-drop solution (that appears to be stable, at least in the time interval presented in Fig. 23(c)). In this work, we do not investigate in detail branches of n-drop solutions with n > 2. Forū = 0.55 and L = 35, we have verified that both secondary branches (shown by the red and green dashed lines in Fig. 11(a)) are unstable for all the values of D. For u = 0.55 and L = 50, we have verified that the only secondary branch that has a stable interval in D is the one connecting points 2 and 5 in Fig. 11(b). The dominant eigenvalue for this branch are shown in Fig. 24. This branch has one Hopf bifurcation and there is a stable interval between D ≈ 0.90 and D ≈ 1.68. Taking into account the fact that for the two-drop primary branch the stable intervals are 0.72 D 0.90 and D 2.12, we can conclude that for D ∈ (0.72, 0.90) a two-drop solution is stable, for D ∈ (0.9, 1.68) a symmetry-broken solution is stable, for D ∈ (1.68, 2.12) both a two-drop solution and a symmetry-broken solution are stable, for D 2.12 a two-drop solution is stable. For other values of D, neither a two-drop solution nor a symmetry-broken solution are stable. Then, as also discussed above for other cases, a time-dependent solution can, for example, evolve into a one-drop solution (that is stable for D 7.13), a time-periodic or multi-drop or quasi-periodic or chaotic solution. These observations can be corroborated by time-dependent simulations, however, we decided not to present such calculations here, as the results agree with the expectations and are generally qualitatively similar to the already presented time-dependent simulations.    figure). This point can be obtained using the weakly nonlinear analysis, see the Appendix. Indeed, for a given domain size L for a two-drop solution, using (7) we find that the value ofū at which the spatially-uniform solution changes its stability and a nonuniform solution emerges isū c = (1 − k 2 )/3, where k = 4π/L (the wavenumber is equal to 4π/L but not to 2π/L, since the value of L that we consider corresponds to a two-drop solution). For this value ofū, we can then find the value D c of D using (21) at which the nature of the primary bifurcation changes (between subcritical and supercritical). Thus, we expect (and, in fact, observe in our numerical results, that we decided not to show here) that when L is fixed and D is used as the principal continuation parameter, forū slightly greater thanū c the primary branch has a single saddle-node bifurcation and returns to D = 0, forū =ū c the primary branch has a single saddle-node bifurcation but it does not return to D = 0 and instead hits the D-axis at D = D c , and forū slightly smaller thanū c there appears one more saddle-node bifurcation out of (D c ,ū c ), and the branch extends to large values of D. For L = 35, we find that k ≈ 0.3590,ū c ≈ 0.5389 and D c ≈ 1.4480. This is in agreement with the results presented in Fig. 25(c) (see the inset showing point (1.4480, 0.5389) by a black circle -the branch showing the locations of saddle-node bifurcations appears exactly from this point). Figure 26 shows the loci of the bifurcation points on the two-drop primary branch in the (D, L)-plane forū = 0.55. We have split this figure into several parts. Panels (a) and (b) correspond to L < L c ≈ 41.32. For these values of L, the primary branch has one saddle-node bifurcation (when D is used as the principle continuation parameter, see Fig. 11(a)), and the branch returns to D = 0. Thus, the branch consists of an upper and a lower part, in Fig. 11  In panels (a) and (b) L < L c = 41.32 (so that the primary branch has a single saddle-node bifurcation, see Fig. 11(a)), and these panels correspond to parts (a) and (b), respectively, of the primary branch shown in Fig. 11(a). In (c), (d) and (e) L > L c = 41.32 (so that the primary branch has two saddle-node bifurcations, see Fig. 11(b)), and these branches correspond to parts (a), (b) and (c), respectively, of the primary branch shown in Fig. 11(b).
belong to L > L c ≈ 41.32. For these values of L, the primary branch has a pair of saddle-node bifurcations (when D is used as the principle continuation parameter, see Fig. 11(b)), and consists of three parts, the upper one denoted by "α", the middle one (connecting the two saddle-nodes) denoted by "β", and the lower one (starting from the second saddle-node and extending to infinity) denoted by "γ". Panels (c), (d) and (e) of Fig. 26 correspond to these parts α, β and γ, respectively. We note that for a more complete picture, it would be beneficial to more precisely indicate which solutions (e.g., one-drop, symmetry-broken or time-periodic solutions) are stable in the various regions where two-drop solutions are unstable. However, we do not present such a detailed 'morphological phase diagram' here and leave this as a topic for future investigation.
Finally, we would like to point out that linear stability of periodic one-drop solutions of the cCH equation was previously analysed by Zaks et al. [53] but only forū = 0 and on the infinite domain, using a Floquet-Bloch-type analysis. We study stability also for nonzero values ofū focusing on the analysis of coarsening modes, i.e., we consider stability on finite domains with lengths equal to twice the period of the onedrop solutions. This implies that forū = 0 the stability regions computed by Zaks et al. must be subsets of those computed here. Therefore, direct comparison is only appropriate forū = 0. In particular, we find good agreement with the stability results of Zaks et al. [53] given in their Table 1 on p. 715, where stability intervals (in terms of solution wavenumbers) for D = 0.5, 0.8, 1, 2, 5 and D → ∞ are presented. Consider, for example, our stability diagram Fig. 25(b) for L = 25. This value of L corresponds to the wavenumber of the one-drop solution K = 2π/(L/2) ≈ 0.503. The results of Zaks et al. [53] indicate that for D = 0.5, 0.8, 1, 2, 5 the solution with K = 0.503 must be linearly unstable. This fully agrees with the results presented in Fig. 25(b). Moreover, an interpolation of results of Zaks et al. indicates that forū = 0 and K = 0.503 there must exist a stability interval in D between D = 1 and D = 2, which agrees with the results in Fig. 25(b). Similarly, the results of Zaks et al. imply forū = 0 and L = 35 (corresponding to K ≈ 0.359) the existence of a stability interval D ∈ (D 1 , D 2 ), where D 1 ∈ (0.5, 0.8) and D 2 ∈ (0.8, 1). This agrees with the results given in Fig. 25(c). An important difference for L = 35 is that for larger values of D (say 1, 2, 5) we predict stability whereas the results of Zaks et al. imply instability. This is not a contradiction, since the stability regions of Zaks et al. must only be subsets of those computed here, as mentioned above. Finally note that our results also show good agreement with related studies for thin-film equations, in particular, when comparing the respective regimes of moderately strong driving. For instance, the sequence of instabilities and their dependence on driving strength for D 3 in our Fig. 25(a) is very similar to the corresponding behaviour in Fig. 22(b) of Ref. [45]. However, the regimes of weak driving notably differ as then the different underlying energies have a crucial influence.

Conclusions
We have analysed the effect of the driving force on the solutions of the cCH equation. Initial insight was obtained by temporal and spatial linear stability analyses of homogenous solutions and we concluded that for the driving force parameter D in the interval [0, √ 2/3) the "horizontal" parts of the fronts and drops/holes are expected to be monotonic, while for D ∈ ( √ 2/3, √ 2) spatial, decaying oscillations are expected. For D > √ 2, we do not expect to see proper drop or hole solutions. Instead, we expect to observe, for example, localized positive/negative-pulse solutions. In addition, for D ∈ (2 √ 2/3, √ 2), the horizontal parts of front-and drop/hole-solutions are linearly unstable, and thus, the solutions on large spatial domains are expected to break up into smaller structures.
Next, we presented the results of numerical continuation of single-and doubleinterface solutions (i.e., fronts and drops/holes). We first discussed the results of numerical continuation with respect to the domain size L for the standard CH equation for several values of the mean solution thicknessū and showed that for smaller values of u the primary bifurcation from the branch of homogeneous solutions is supercritical, whereas at some value ofū it changes to subcritical. The value ofū at which the type of the primary bifurcation switches can be found by the weakly nonlinear analysis. At some even larger value ofū (that, in fact, follows from the linear stability analysis), the primary bifurcation disappears, and beyond a certain value of the domain size, linearly stable homogeneous and inhomogeneous solutions and a linearly unstable inhomogeneous solution coexist. After that, we studied the effect of the driving force on inhomogeneous solutions of the CH equation. For smaller values ofū, we found that when continuation is performed in the driving force parameter D, branches of solutions extend to infinity for all sufficiently large values of the domain size. Whereas for larger values ofū the branches of solutions exhibit saddle-nodes and return to D = 0, if L is sufficiently small. For larger values of L, the branches exhibit an additional saddle-node and extend to infinity. The transition from one type of the bifurcation diagram to the other type of the bifurcation diagram happens at L = L c , where L c is the wavelength of a small-amplitude neutrally stable sinusoidal wave. For this value of L, the branch of solutions terminates at the horizontal axis at D = D c , where D c can be found by the weakly nonlinear analysis. So, for L just beyond L c , there is a range of D values for which two different stable spatially inhomogeneous solutions and one unstable inhomogeneous solution coexist. For even larger values of L, the saddle-nodes annihilate each other, and the branches extend to infinity. Also, ifū becomes sufficiently large, the branches of inhomogeneous solutions exhibit a saddle-node and return to D = 0 for all sufficiently large values of L.
Finally, we studied in detail the linear stability properties of the various possible spatially periodic traveling solutions of the cCH equation by performing numerical continuation of inhomogeneous solutions along with the dominant eigenvalues. To obtain more complete bifurcation diagrams, we also implemented a numerical procedure for continuation of time-periodic solutions. Our primary interest was in the study of the stability of symmetric two-drop solutions, and coarsening of such solutions in particular. Without driving force, the two-drop solutions have two real positive (unstable) eigenvalues that correspond to two different coarsening modes -volume and translation modes. For the volume mode, the corresponding eigenfunction tends to increase the volume of one of the drops and decrease the volume of the other one. For the translation mode, the corresponding eigenfunction tends to shift both drops in the opposite directions, so that they move towards each other. When driving is introduced, we found that one of the coarsening modes is stabilized at relatively small values of D. In addition, our results indicate that the type of a coarsening mode can change as D increases. We also found that there may be intervals in the driving force D, where there are no unstable eigenvalues, and, therefore, driving can be used to prevent coarsening. We, in addition, computed side branches of symmetry-broken solutions and analysed the stability of such solutions, and also branches of time-periodic solutions, and presented detailed stability diagrams in the (D, L)-and (D,ū)-planes. The predictions from the numerical continuation results have been confirmed by time simulations for the cCH equation. In the future, it will be of interest to undertake similar studies for related equations, such as, for example, the various variants of the Kuramoto-Sivashinsky equation and related thin-film models and to extend the study to two-dimensional and three-dimensional solutions. where