Patterns formed in a thin film with spatially homogeneous and non-homogeneous Derjaguin disjoining pressure

We consider patterns formed in a two-dimensional thin film on a planar substrate with a Derjaguin disjoining pressure and periodic wettability stripes. We rigorously clarify some of the results obtained numerically by Honisch et al. and embed them in the general theory of thin-film equations. For the case of constant wettability, we elucidate the change in the global structure of branches of steady state solutions as the average film thickness and the surface tension are varied. Specifically we find, by using methods of local bifurcation theory and the continuation software package AUTO, both nucleation and metastable regimes. We discuss admissible forms of spatially non-homogeneous disjoining pressure, arguing for a form that differs from the one used by Honisch et al. , and study the dependence of the steady state solutions on the wettability contrast in that case.


Introduction
Thin liquid films on solid substrates occur in many natural situations. For example, they appear in tear films in the eye which protect the cornea [6] or in streams of lava from a volcanic eruption [19]. Moreover, thin liquid films occur in many technological applications, such as coatings [22] and lubricants (e.g. oil films which lubricate the piston in a car engine [41], drying paint layers [18], and in the manufacture of microelectronic devices [21]). For extensive reviews of thin-film flow see, for example, Oron et al. [28] and Craster and Matar [10].
As these liquid films are thin, the Navier-Stokes equation governing their flow can be reduced to a single degenerate fourth-order quasi-linear parabolic partial differential equation (PDE) usually known as a thin-film equation. In many applications a choice of a disjoining pressure, which we denote by Π, must be made. Such a term describes the action of surface forces on the film [37]. In different situations, different forms of disjoining pressure are appropriate; these may incorporate long-range van der Waals forces and/or various types of short-range interaction terms such as Born repulsion; inclusion of a particular type of interaction can have significant effects on the wettability of the surface and the evolution of the film film, sometimes leading to dewetting phenomena, i.e. the rupture of the film and the appearance of dry spots. (Here and subsequently by "wettability" of the surface we mean the ability of a solid surface to reduce the surface tension of a liquid on contact with it such that it spreads over the surface and wets it. ) Witelski and Bernoff [43] were one of the first authors to analyse mathematically the rupture of three-dimensional thin films. In particular, considering a disjoining pressure of the form Π = −1/(3h 3 ) (we use the sign convention adopted in Honisch et al. [17]), they analysed planar and axisymmetric equilibrium solutions on a finite domain. They showed that a disjoining pressure of this form leads to finite-time rupture singularities, that is, the film thickness approaches zero in finite time at a point (or a line or a ring) in the domain. In a related more recent paper, Ji and Witelski [20] considered a different choice of disjoining pressure and investigated the finite-time rupture solutions in a model of thin film of liquid with evaporative effects. They observed different types of finitetime singularities due to the non-conservative terms in the model. In particular, they showed that the inclusion of non-conservative term can prevent the disjoining pressure from causing finite-time rupture.
A pioneering theoretical study of a thin-film equation with a disjoining pressure term given by a combination of negative powers of the thin film thickness is that by Bertozzi et al. [4], who studied the formation of dewetting patterns and the rupture of thin liquid films due to long-range attractive and short-range Born repulsive forces, and determined the structure of the bifurcation diagram for steady state solutions, both with and without the repulsive term.
Aiming to quantify the temporal coarsening in a thin film, Glasner and Witelski [15] examined two coarsening mechanisms that arise in dewetting films: mass exchange that influences the breakdown of individual droplets and spatial motion that results in droplet rupture as well as merging events. They provided a simple model with a disjoining pressure which combines the effects of both short-and long-range forces acting on the film. Kitavtsev et al. [23] analysed the long-time dynamics of dewetting in a thin-film equation by using a disjoining pressure similar to that used by Bertozzi et al. [4]. They applied centre manifold theory to derive and analysed an ordinary differentual equation model for the dynamics of dewetting.
The recent article by Witelski [42] presents a review of the various stages of dewetting for a film of liquid spreading on a hydrophobic substrate. Different types of behaviour of the film are observed depending on the form of the disjoining pressure: finite-time singularities, self-similar solutions and coarsening. In particular, he divides the evolution of dewetting processes into three phases: an initial linear instability that leads to finitetime rupture (short time dynamics), which is followed by the propagation of dewetting rims and instabilities of liquid ridges (intermediate time dynamics), and the eventual formation of quasi-steady droplets (long time dynamics).
Most of the previous studies of thin liquid films focussed on films on homogeneous substrates. However, thin liquid films on non-homogeneous chemically patterned substrates are also of interest. These have many practical applications, such as in the construction of microfluidic devices and creating soft materials with a particular pattern [30]. Chem-ically patterned substrates are an efficient way to obtain microstructures of different shapes by using different types of substrate patterning [34]. Chemical modification of substrates can also be used to avoid spontaneous breakup of thin films, which is often highly undesirable, as, for example, in printing technology [5,1].
Due to their many applications briefly described above, films on non-homogeneous substrates have been the object of a number of previous theoretical studies which motivate the present work. For example, Konnur et al. [24] found that in the case of an isolated circular patch with wetting properties different from that of the rest of the substrate that chemical non-homogeneity of the substrate can greatly accelerate the growth of surface instabilities. Sharma et al. [35] studied instabilities of a liquid film on a substrate containing a single heterogeneous patch and a substrate with stripes of alternating less and more wettable regions. The main concern of that paper was to investigate how substrate patterns are reproduced in the liquid film, and to determine the best conditions for templating.
Thiele et al. [39] performed a bifurcation analysis using the continuation software package AUTO [11] to study dewetting on a chemically patterned substrate by solving a thin-film equation with a disjoining pressure, using the wettability contrast as a control parameter. The wettability contrast measures the degree of heterogeneity of the substrate; it is introduced and defined rigorously in (5.1) in Section 5. Honisch et al. [17] modelled an experiment in which organic molecules were deposited on a chemically non-homogeneous silicon oxide substrates with gold stripes and discuss the redistribution of the liquid following deposition.
In a recent paper, Liu and Witelski [27] studied thin films on chemically heterogeneous substrates. They claim that in some applications such as digital microfluidics, substrates with alternate hydrophilic and hydrophobic rectangular areas are better described by a piecewise constant function than by a sinusoidal one. Therefore, in contrast to other studies, including the present one, they study substrates with wettability characteristics described by such a function. Based on the structure of the bifurcation diagram, they divide the steady-state solutions into six distinct but connected branches and show that the only unstable branch corresponds to confined droplets, while the rest of the branches are stable.
In the present work, we build on the work of Thiele et al. [39] and Honisch et al. [17]. Part of our motivation is to clarify and explain rigorously some of the numerical results reported in these papers. In the sinusoidally striped non-homogeneous substrate case, we offer a justification for using a form of the disjoining pressure that differs from the one used in these two papers. A detailed plan of the paper is given in the last paragraph of Section 2.

Problem Statement
Denoting the thickness of the thin liquid film by z = h(x, y, t), where (x, y, z) are the usual Cartesian coordinates and t is time, Honisch et al. [17] consider the thin-film equation where Q(h) = h 3 /(3η) is the mobility coefficient with η being the dynamic viscosity. The generalized pressure P (h, x, y) is given by where γ is the coefficient of surface tension and we follow [17] in taking the Derjaguin disjoining pressure Π(h, x, y) in the spatially homogeneous case to be of the form suggested, for example, by Pismen [29]. Here A and B are positive parameters that measure the relative contributions of the long-range forces (the 1/h 3 term) and the short-range ones (the 1/h 6 term). However, we will see that both of these constants can be scaled out of the mathematical problem. Equation In what follows, we study thin films on both homogeneous and non-homogeneous substrates. In the non-homogeneous case, we will modify (2.3) by assuming that the Derjaguin pressure term Π changes periodically in the x-direction with period L. The appropriate forms of Π in the non-homogeneous case are discussed in Section 5.
Hence, in order better to understand solutions of (2.1), we study its one-dimensional version, We start by characterising steady state solutions of (2.4) subject to periodic boundary conditions at x = 0 and x = L. In other words, we seek steady state solutions h(x) of (2.4), satisfying the non-local boundary value problem where the constant h * (> 0) denotes the (scaled) average film thickness, and the periodic boundary conditions Now we non-dimensionalise. Setting , h = Hh, and x = Lx, (2.8) in (2.5) and removing the tildes, we obtain and subject to the periodic boundary conditions (2.12) and the volume constraint Note that the problem (2.9)-(2.14) is very similar to the corresponding steady state problem for the Cahn-Hilliard equation considered as a bifurcation problem in the parametersh and ǫ by Eilbeck et al. [12]. The boundary conditions considered in that work were the physically natural double Neumann conditions. The periodic boundary conditions (2.12) in the present problem slightly change the analysis, but our general approach in characterising different bifurcation regimes still follows that of Eilbeck et al. [12], though the correct interpretation of the limit as ǫ → 0 + is that now we let the surface tension γ go to zero. In particular, we perform a Liapunov-Schmidt reduction to determine the local behaviour close to bifurcation points and then use AUTO (in the present work we use the AUTO-07p version [11]) to explore the global structure of branches of steady state solutions both for the spatially homogeneous case and for the spatially non-homogeneous case in the case of an x-periodically patterned substrate.
We first investigate the homogeneous case and, having elucidated the structure of the bifurcations of non-trivial solutions from the constant solution h =h in that case in Sections 3 and 4, we study forced rotational (O(2)) symmetry breaking in the nonhomogeneous case in Section 5. In Appendix A, we present a general result about O(2) symmetry breaking in the spatially non-homogeneous case. It shows that in the spatially non-homogeneous case, only two steady state solutions remain from the orbit of solutions of (2.9)-(2.14) induced by its O(2) invariance. We concentrate on the simplest steady state solutions of (2.9)-(2.14), as by a result of Laugesen and Pugh [25,Theorem 1] only such solutions, that is, constant solutions and those having only one extremum point, are linearly stable in the homogeneous case. For information about dynamics of one-dimensional thin-film equations the reader should also consult Zhang [46].

Liapunov-Schmidt Reduction in the Spatially Homogeneous Case
We start by performing an analysis of the dependence of the global structure of branches of steady state solutions of the problem in the spatially homogeneous case given by (2.9)-(2.14) on the parametersh and ǫ.
In what follows, we do not indicate explicitly the dependence of the operators on the parametersh and ǫ, and all of the calculations are performed for a fixed value ofh and close to a bifurcation point ǫ = ǫ k for k = 1, 2, 3, . . . defined below. We has zero mean, and rewrite (2.9) as If we set The linearisation of G at v applied to w is defined by We denote dG(0) by S, so that S applied to w is given by To locate the bifurcation points, we have to find the nontrivial solutions of the equation Sw = 0 subject to The kernel of S is non-empty and two dimensional when and is spanned by cos(2kπx) and sin(2kπx). That these values of ǫ correspond to bifurcation points follows from two theorems of Vanderbauwhede [40, Theorems 2 and 3].
In a neighbourhood of a bifurcation point (ǫ k , 0) in (ǫ, v) space, solutions of G(v) = 0 on H are in one-to-one correspondence with solutions of the reduced system of equations on R 2 , for some functions g 1 and g 2 to be obtained through the Liapunov-Schmidt reduction [16].
To set up the Liapunov-Schmidt reduction, we decompose D(G) and H as follows: and Since S is self-adjoint with respect to the L 2 -inner product denoted by ·, · , we can choose and denote the above basis for M by {w 1 , w 2 } and for N by {w * 1 , w * 2 }. We also denote the projection of H onto range S by E.
Since the present problem is invariant with respect to the group O(2), the functions g 1 and g 2 must have the form for some function p(·, ·) [9], which means that in order to determine the bifurcation structure, the only terms that need to be computed are g 1,xǫ and g 1,xxx , as these immediately give g 2,yǫ and g 2,yyy and all of the other second and third partial derivatives of g 1 and g 2 are identically zero. Following Golubitsky and Schaeffer [16], we have 14) for r = 1, 2, 3, . . . , (3.16) and we choose where w 1 ∈ ker S and w * 1 ∈ (range S) ⊥ . In particular, from (3.16) we have and so , we use the definition of ǫ k given in (3.8) and solve the second order ordinary differential equation satisfied by R(x), subject to in ker S. We obtain and hence In addition, from (3.16) we have and therefore Putting all of the information in (3.23) and (3.25) into (3.15) we obtain In addition, G ǫ (v) = v xx , so that G ǫ (0) = 0 at v = 0, and hence we have Furthermore, since dG ǫ (w) = w xx , from (3.14) we obtain (3.28) Referring to (3.13) and the argument following that equation, the above analysis shows that as long as f ′ (h) > 0 at ǫ = ǫ k a circle of equilibria bifurcates from the constant solution h ≡h. The direction of bifurcation is locally determined by the sign of g 1,xxx given by (3.26). Hence, using 1/ǫ as the bifurcation parameter, the bifurcation of nontrivial equilibria is supercritical if g 1,xxx is negative and subcritical if it is positive.
By finding the values ofh where g 1,xxx given by (3.26) with f (h) given by (2.10) is zero, we finally obtain the following proposition: so that g 1,xxx < 0 ifh ∈ (1.289, 1747) giving the result. Proof Assume that such a nontrivial solution exists. Then, sinceh ≤ 1, its global minimum, achieved at some point x 0 ∈ (0, 1), must be less than 1. (We may take the point x 0 to be an interior point by translation invariance.) But then so the point x 0 cannot be a minimum.

Two-Parameter Continuation of Solutions in the Spatially Homogeneous Case
To describe the change in the global structure of branches of steady state solutions of the problem (2.9)-(2.14) ash and ǫ are varied, we use AUTO [11] and our results are summarised in Figure 1. As Figure 1 shows, a curve of saddle-node (SN) bifurcations which originates from  Figure 1. The global structure of branches of steady state solutions with a unique maximum, including both saddle-node (SN) (shown with dash-dotted curves) and pitchfork (PF) bifurcation branches (shown with solid curves). The nucleation regime 1 <h < 2 1/3 ≈ 1.259 (Regime I), the metastable regime 2 1/3 <h < 1.289 andh > 1.747 (Regime II), and the unstable regime 1.289 <h < 1.747 (Regime III) are also indicated. Figure 1 identifies three different bifurcation regimes, denoted by I, II and III, with differing bifurcation behaviour occurring in the different regimes, namely (using the terminology of [12] in the context of the Cahn-Hilliard equation): • a "nucleation" regime for 1 <h < 2 1/3 ≈ 1.259 (Regime I), • a "metastable" regime for 2 1/3 <h < 1.289 andh > 1.747 (Regime II), and • an "unstable" regime for 1.289 <h < 1.747 (Regime III).
In Regime I, the constant solution h(x) ≡h is linearly stable which follows from analysing the spectrum of the operator S for f ′ (h) < 0 in (3.6) and (3.7), but under sufficiently large perturbations the system will evolve to a non-constant steady state solution. See [25] for an extensive discussion of the stability analysis of steady state solutions to thin-film equations.
In Regime II, as ǫ is decreased the constant solution h(x) ≡h loses stability through a subcritical bifurcation.
In Regime III, as ǫ is decreased, the constant solution h(x) ≡h loses stability through a supercritical bifurcation.

The Spatially Non-Homogeneous Case
Honisch et al. [17] chose the Derjaguin disjoining pressure Π(h, x, y) to be of the form where the function G(x, y) models the non-homogeneity of the substrate and the parameter ρ, which can be either positive or negative, is called the "wettability contrast". Following Honisch et al. [17], in the remainder of the present work we consider the specific case corresponding to an x-periodically patterned (i.e. striped) substrate. There are, however, some difficulties in accepting (5.1) as a physically realistic form of the disjoining pressure for a non-homogeneous substrate. The problems arise because the two terms in (5.1) represent rather different physical effects. Specifically, since the 1/h 6 term models the short-range interaction amongst the molecules of the liquid and the 1/h 3 term models the long-range interaction, assuming that both terms reflect the patterning in the substrate in exactly the same way through their dependence on the same function G(x, y) does not seem very plausible. Moreover, there are other studies which assume that the wettability of the substrate is incorporated in either the shortrange interaction term (see, for example, Thiele and Knobloch [38] and Bertozzi et al. [4]) or the long-range interaction term (see, for example, Ajaev et al. [2] and Sharma et al. [35]), but not both simultaneously. Hence in what follows we will consider the two cases Π(h, x) = Π LR (h, x) and Π(h, x) = Π SR (h, x), where LR stands for "long range" and SR stands for "short range" where in both of which G(x) is given by (5.2) and ρ is the wettability contrast. For small wettability contrast, |ρ| ≪ 1, we do not expect there to be significant differences between the influence of Π LR or Π SR on the bifurcation diagrams, as these results depend only on the nature of the bifurcation in the homogeneous case ρ = 0 and on the symmetry groups under which the equations are invariant. To see this, consider the spatially non-homogeneous version of (2.9), that is, the boundary value problem Seeking an asymptotic solution to (5.5)-(5.7) in the form h(x) =h + ρh 1 (x) + O(ρ 2 ) in the limit ρ → 0, by substituting this anzatz into (5.5) we find that in the case of Π LR (h, x) we have while in the case of Π SR (h, x) the corresponding result is For non-zero values of ρ, in both the Π LR and Π SR cases, the changes in the bifurcation diagrams obtained in the homogeneous case (ρ = 0) are an example of forced symmetry breaking (see, for example, Chillingworth [8]), which we discuss further in Appendix A. More precisely, we show there that when ρ = 0, out of the entire O(2) orbit only two equilibria are left after symmetry breaking.   shows how for small wettability contrast |ρ| ≪ 1, the resulting spatial nonhomogeneity introduces imperfections [16] in the bifurcation diagrams of the homogeneous case ρ = 0 discussed in Section 4. It presents bifurcation diagrams in Regimes I, II and III when the disjoining pressure Π LR is given by (5.3) for a range of small values of ρ together with the corresponding diagrams in the case ρ = 0. The corresponding bifurcation diagrams when the disjoining pressure Π SR is given by (5.4) are very similar and hence are not shown here.
For large wettability contrast, specifically for |ρ| ≥ 1, significant differences between the two forms of the disjoining pressure are to be expected. When using Π LR , one expects global existence of positive solutions for all values of |ρ|; see, for example, Wu and Zheng [44]. On the other hand, when using Π SR , there is the possibility of rupture of the liquid film for |ρ| ≥ 1; see, for example, [4,44], which means in this case we do not expect positive solutions for sufficiently large values of |ρ|. In Figure 3 we plot the branches of the positive solutions of (5.5)-(5.7) with a unique maximum when the disjoining pressure is Π LR forh = 3 and 1/ǫ = 50, so that when ρ = 0, we are in Regime II above the curve of pitchfork (PF) bifurcations from the constant solution shown in Figure 1. The work of Bertozzi et al. [4] and of Wu and Zheng [44], shows that strictly positive solutions exist for all values of |ρ|, i.e. beyond the range ρ ∈ [−2, 2] for which we have performed the numerical calculations. Figure 4 shows that the situation is different when the disjoining pressure is Π SR (with the same values ofh and ǫ). At |ρ| = 1, the upper branches of solutions disappear due to rupture of the film, so that at some point x 0 ∈ [0, 1] we have h(x 0 ) = 0 and a strictly positive solution no longer exists, while the other two branches coalesce at a saddle node bifurcation at |ρ| ≈ 5.8.
Note that in Figures 3 and 4, the non-trivial "solution" on the axis ρ = 0 is, in fact, a whole O(2)-symmetric orbit of solutions predicted by the analysis leading to Figure 1. Note that when the disjoining pressure is Π SR , given by (5.4), we were unable to use AUTO to continue branches of solutions beyond the rupture of the film at |ρ| = 1.  Figure 5. Solutions h(x) when the disjoining pressure is Π SR forh = 2 and 1/ǫ = 30 for ρ = 0, 0.97, 0.98, 0.99 and 1, denoted by "1", "2", "3", "4" and "5", respectively, showing convergence of strictly positive solutions to a non-strictly positive one as ρ → 1 − . Figure 5 shows the film approaching rupture as ρ → 1 − at the point where the coefficient of the short-range term disappears when ρ = 1, i.e. when 1 + sin(2πx) = 0 and hence at x = 3/4. These numerical results are consistent with the theoretical arguments of Bertozzi et al. [4].
Investigation of the possible leading-order balance in (2.9) when the disjoining pressure is Π SR and ρ = 1 in the limit x → 3/4 reveals that h(x) = O (x − 3/4) 2/3 ; continuing the analysis to higher order yields (5.10) Figure 6 shows the detail of the excellent agreement between the solution h(x) when ρ = 1 and the two-term asymptotic solution given by (5.10) close to x = 3/4.   Each curve in Figure 7(c) is a curve of saddle-node bifurcations. Similar to Figure 2(c), Figure 9 shows the bifurcation diagram of steady state solutions withh = 2 (Regime II) for ρ = 0 and ρ = 0.005. As 1/ǫ is increased, we pass successively through saddle-node bifurcations with the multiplicity of the steady state solutions changing from 1 to 3 to 5 and then back to 3 again. In Figure 10, we plot the five steady state solutions with a unique maximum indicated in Figure 9 by (i)-(v).

Conclusions
In the present work we have investigated the steady state solutions of the thin-film evolution equation (2.1) both in the spatially homogeneous case (2.9)-(2.12) and in the spatially non-homogeneous case for two choices of spatially non-homogeneous Derjaguin disjoining pressure given by (5.3) and (5.4). We have given a physical motivation for our choices of the disjoining pressure. For reasons explained in the last paragraph of Section 2, we concentrated on branches of solutions with a unique maximum.
In the spatially homogeneous case (2.9)-(2.14), we used the Liapunov-Schmidt reduction of an equation invariant under the action of the O(2) symmetry group to obtain local bifurcation results and to determine the dependence of the direction and nature of bifurcation in bifurcation parameter 1/ǫ on the average film thicknessh; our results on the existence of three different bifurcation regimes, (namely nucleation, metastable, and unstable) are summarised in Propositions 1 and 2 and in Figure 1 obtained using AUTO.
In the spatially non-homogeneous case (5.5)-(5.7), we clarified the O(2) symmetry breaking phenomenon (see Appendix A) and presented imperfect bifurcation diagrams in Figure 2 and global bifurcation diagrams using the wettability contrast ρ as a bifurcation parameter for fixed ǫ andh in Figures 3 and 4.
Let us discuss Figure 3 in more detail; for convenience, it is reproduced in Figure 11 with labelling added to the different branches of strictly positive steady state solutions with a unique maximum. Below we explain what these different labels mean.  Figure 11. A sketch of the bifurcation diagram plotted in Figure 3 with the different solution branches labelled.
Our results can be described using the language of global compact attractors. For more information on attractors in dissipative infinite-dimensional dynamical systems please see the monograph of Hale [14] and Wu and Zheng [44] for applications in the thin-film context. In systems such as (2.4), the global compact attractor of the PDE is the union of equilibria and their unstable manifolds. In Figures 12 and 13 we sketch the global attractor of (2.4) for ǫ = 1/50 andh = 3, the values used to plot Figure 3. For these values of the parameters the attractor is two-dimensional and we sketch its projection onto a plane. When ρ = 0, for 1/ǫ = 50, the attractor is two-dimensional; the constant solution h ≡h denoted by O has a two-dimensional unstable manifold and X corresponds to a whole O(2) orbit of steady state solutions. A sketch of the attractor in this case is shown in Figure 12.
When |ρ| takes small positive values, only two steady state solutions, denoted by A and C remain from the entire O(2) orbit, as discussed in Appendix A, while the constant solution continues to B without change of stability. The resulting attractor is sketched in Figure 13.
Increasing |ρ| causes the steady state solutions B and C to coalesce in a saddle node bifurcation, so that the attractor degenerates to a single asymptotically stable steady state solution. It would be interesting to understand why this collision of steady state solution branches occurs.
We have also explored the geometry of film rupture which occurs as ρ → 1 − when the disjoining pressure is given by Π SR ; this phenomenon is shown in detail in Figures 5 and  6. Finally, in Figures 7 and 8, we showed the results of a two-parameter continuation study in the (1/ǫ, ρ) plane, showing how the multiplicity of positive steady state solutions changes as parameters are varied, and, in particular, indicating in the case of disjoining pressure Π SR shown in Figure 8 regions in parameter space where no such solutions exist. We conjecture that in these regions the solution of the unsteady problem with any positive initial condition converges to a weak solution of the thin-film equation with regions in which h(x) = 0, i.e. solutions with film rupture. For a discussion of such (weak) solutions of thin-film equations in the homogeneous case the reader is referred to the work of Laugesen and Pugh [26].
In the case of disjoining pressure Π SR , we could not use the AUTO-07p version of AUTO to continue branches of solutions beyond rupture. It would be an interesting project to develop such a capability for this powerful and versatile piece of software. Figures 8(b) and 8(c) provide numerical evidence for the existence of a curve of saddlenode bifurcations converging to the point (0, 1) in the (1/ǫ, ρ) plane; an explanation for this feature of the global bifurcation diagrams requires further study.
To summarise: our study was primarily motivated by the work of Honisch et al. [17]. While we have clarified the mathematical properties of (2.9)-(2.12) and (5.5)-(5.7), so that the structure of bifurcations in Figure 3(a) of that paper for non-zero values of ρ is now understood, many of their other numerical findings are still to be explored mathematically. For example, the stability of ridge solutions shown in their Figure 5 in the context of the full two-dimensional problem of a substrate with periodic wettability stripes. There is clearly much work still to be done on heterogeneous substrates with more complex wettability geometry.
A final remark that might be of interest to the reader is that the solutions of equations (5.5)-(5.7), the steady state solutions of (2.4), a degenerate quasi-linear fourth-order PDE, can also be thought of as the steady state solutions of a much simpler (Rubinstein-Sternberg type [31]) second-order semi-linear non-local equation, It would be interesting to compare the dynamics of (2.4) and (6.1), for example using the spectral comparison principles of Bates and Fife [3]. For other work on non-local reaction-diffusion equations such as (6.1), please see Budd et al. [7] and the review of Freitas [13].
if and only if the right-hand side of (A 5) is orthogonal in L 2 ([0, 1]) to the solutions of Ku = 0.
Next, we show that u := h ′ 0,s solves Ku = 0. Indeed, by differentiating (A 1) with ρ = 0 with respect to x, we find that u solves the equation Now recall that h 0,s (x) = h 0 (x + s), so if we write F (x + s) = f 2 (h 0 (x + s))h ′ 0 (x + s), the function F (·) is 1-periodic in x with zero mean. Hence which has two solutions s ∈ R/[0, 1], from which we conclude there is a solution h 1 (x) only for two choices of s ∈ R/[0, 1], that is, that only two solutions to (A 1) remain from the entire O(2) orbit that exists for ρ = 0.