Solitary solutions to the steady Euler equations with piecewise constant vorticity in a channel

We consider a two-dimensional, two-layer, incompressible, steady flow, with vorticity which is constant in each layer, in an infinite channel with rigid walls. The velocity is continuous across the interface, there is no surface tension or difference in density between the two layers, and the flow is inviscid. Unlike in previous studies, we consider solutions which are localised perturbations rather than periodic or quasi-periodic perturbations of a background shear flow. We rigorously construct a curve of exact solutions and give the leading order terms in an asymptotic expansion. We also give a thorough qualitative description of the fluid particle paths, which can include stagnation points, critical layers, and streamlines which meet the boundary.


Introduction
The incompressible Euler equations have been studied for over 250 years, but even when considering time-independent, two-dimensional flows, there is still much being discovered.See [5], or [48] for a tighter focus on the steady case.The time-independent case is of interest both for its own sake, and to better understand possible end states of the time-dependent problem.For instance, steady states give counterexamples to the phenomenon of invicid damping, where solutions of the time-dependent Euler equations decay in some sense, despite energy being conserved [4,42].
The study of solutions to the Euler equations with piecewise constant vorticity is classical.There are broadly speaking two cases: vortex patches, where the vorticity has compact support, and vorticity fronts, with unbounded regions of vorticity.Studies of these phenomena date back to the 19th century, for instance by Rayleigh [47], Kelvin [51], and Kirchhoff [39] with a more recent wave of interest in vortex patches including work such as [9,31,11,6,18].Time-dependent results have also been shown for vorticity fronts [33,7,32,34].
A common setting for study of the Euler equations is in a two-dimensional channel.It is shown in [26] that in a channel, any solution to the steady Euler equations is either a shear flow, or has a stagnation point.This contrasts with water waves with a free upper boundary, which exhibit many interesting solutions without stagnation [30].Indeed, many of the standard techniques in the classical study of water waves assume the stream function is strictly monotone in the vertical variable, so cannot allow for stagnation points.However, there has been great progress in understanding non-monotone stream functions in recent years [53,43,54,56,13], often using the weaker assumption that the vorticity can be expressed as a vorticity function of the stream function.
We consider a two-dimensional flow in a channel with piecewise constant vorticity, and find explicit expansions of exact solutions, while thoroughly describing their qualitative behaviour.This behaviour includes stagnation points, which complements [26], and since the boundary at which the vorticity changes is free, this problem also fits in with the literature on water waves with stagnation.However, for some values of our parameters, there are solutions for which even the weaker assumption that a vorticity function exists does not hold.We focus on localised solutions, which have many challenges not present in the periodic case.Reformulating our problem as an evolution equation, with the horizontal spatial direction playing the role of time, allows us to use powerful tools from spatial dynamics.These are particularly useful when working outside a periodic regime, especially a centre manifold theorem of Mielke [44], which, for sufficiently small localised solutions, allows us to reduce our PDE to an ODE.
1.1.The problem.We consider a two-dimensional flow which is incompressible and inviscid.The flow has two layers, each with the same constant density, but with different constant vorticities.There is no surface tension between the two layers.The fluid domain is infinite in width, but bounded above and below by rigid walls of constant height.It is also steady, that is, the wave profile moves with constant speed.In other words, if X is the horizontal coordinate, and t is time, a steady solution depends on (X, t) only through the combination X − ct, for some constant wave speed c.Then moving to a co-moving frame of reference gives a time-independent flow.We use dimensionless units for distance and time, so that the channel has height 1, and the difference in vorticity between the two layers is 1.
More precisely, given real constants ω 0 , ω 1 with ω 0 − ω 1 = 1, and h ∈ (0, 1), we seek η ∈ C 2 (R) which in turn defines regions Ω 0 = {(X, Y ) | 0 < Y < h + η(X)} and Ω 1 = {(X, Y ) | h + η(X) < Y < 1}.We also seek a velocity field U, V ∈ C 1 (Ω 0 ) ∩ C 1 (Ω 1 ) ∩ C 0 (Ω 0 ∪ Ω 1 ).These must satisfy where (1.1a) is the incompressibility condition, and (1.1b) enforces the piecewise constant vorticity.The vorticity of each fluid element is constant, so fluid cannot enter or leave either of Ω 0 and Ω 1 .This gives the kinematic boundary conditions (1.1e) We call {(X, Y ) | Y = h+η(X)} the interface.Points at which U = V = 0 are often of interest, and we refer to these as stagnation points, or points at which the flow stagnates.We use critical layers to mean curves along which the horizontal velocity U vanishes [53].A solution is shear, or a shear flow if it has no X dependence.We say a solution is localised, or homoclinic, or forms a solitary wave if lim X→∞ η(X) = lim X→−∞ η(X) = 0.A solitary wave with η(X) > 0 for all X is called a wave of elevation, and one with η(X) < 0 for all X, a wave of depression.
1.2.Main results.We construct localised solutions to (1.1) by bifurcating from a shear flow.Localised solutions present many challenges that periodic solutions do not.For instance, since periodic solutions can be considered to have a compact domain, one has access to compact embeddings between Hölder and Sobolev spaces.Furthermore, smallamplitude periodic solutions are a linear phenomenon, whereas small-amplitude solitary waves are weakly-nonlinear; see (4.26) and the surrounding discussion.
We now, slightly informally, state our main result.Let the trivial solution, which we will bifurcate from, be given by where the constant c is the speed of the flow at the interface; see Figure 1.
In the terminology of Lin and Zeng [42], there are steady structures which are an arbitrarily small, but non-zero, localised perturbation from a trivial shear flow.
Remark 1.2.Each of the solutions described in Theorem 1.1 in fact gives rise to an infinite family of solutions through phase shifts, i.e., through transformations of the form X → X + a.
on the other hand, then there does exist a vorticity function.See Proposition 2.4 for more details.
We can also give a relatively complete picture of what the flow looks like, and where and how it stagnates.All the solutions constructed in Theorem 1.1 have a unique interior stagnation point, and some have two further stagnation points on the boundary.Much of the qualitative behaviour can be determined by the sign of θ, and whether ω 0 is larger or smaller than the critical value 1 − h.Theorem 1.4.Every non-shear solution constructed in Theorem 1.1 has either a unique stagnation point in Ω 0 ∪ Ω 1 , or one stagnation point in Ω 0 ∪ Ω 1 and two saddle points on the boundary closest to the interior stagnation point.The uniqueness, nature and location of the interior stagnation point, and whether η ε gives a wave of elevation or depression, depends only on θ and ω 0 , as laid out in Table 1.Moreover η ε is strictly monotone on X > 0 and on X < 0. The solution has a critical layer, which is unbounded if ω 0 ̸ = 1 − h, and bounded if ω 0 = 1 − h.If ω 0 = 1 − h, there is a streamline which connects the saddle points on the boundary.See Figure 2.
Table 1.Behaviour of solutions in different regions of parameter space; see Figure 2. The stagnation point being referred to is the interior stagnation point.
Figure 2. Streamlines of solutions in the six different regions of parameter space from Table 1.The dashed streamlines correspond to the interface Y = h + η(X).The dotted line in parameter space corresponds to θ = 0, we have not constructed solutions here.
1.3.Related work.One body of work that this problem is related to is that of travelling water waves with critical layers, as we are considering the Euler equations with a free boundary.See [30] and the references within for a general overview of the steady water waves literature, including in particular, problems with solutions exhibiting critical layers.Wahlén [53] gave the first exact construction of a rotational flow with critical layers and a free surface when considering the one-layer, constant vorticity case.Since then, new results have taken this line of inquiry in different directions.For example, Ehrnström, Escher, Villari, and Wahlén [17,16] found solutions with an arbitrary number of critical layers, with stagnation points exhibiting "cat's eye" structures.Matioc then moved away from these constant or affine vorticity functions, and considered a twolayer problem, where the two layers had distinct, constant vorticities and densities [43].Walsh, Bühler, and Shatah considered gravity-vorticity waves, [54], where they took the lower layer to be of finite depth and zero vorticity, and the upper layer to be either of infinite height and constant vorticity, or finite height and a more general vorticity.The study of critical layers is not limited to local perturbations, for instance, Constantin, Strauss, and Vȃrvȃrucȃ [13] considered gravity waves of constant vorticity, and used global bifurcation to construct large amplitude, overhanging waves.Wheeler considered a flow with an arbitrary number of layers, with waves of vorticity [56].Whereas [53,17,16,13] considered the case when a vorticity function exists, we can work outside this scope.More generally, water wave problems usually have a discontinuity in velocity across free surfaces, due to surface tension, or a change in density [43,54,56].We do not have such a "jump" condition.As already mentioned, our work here is also connected to invicid damping, where small smooth perturbations from a shear flow decay to a shear flow when evolved under the Euler equations.Bedrossian and Masmoudi were the first to demonstrat this phenomenon at the non-linear level [4].Recently [49], Sinambela and Zhao constructed less regular flows for which linear invicid damping occurs.For a general overview see the survey [35].Invicid damping is an analogue of Landau damping; see [46].
Results on invicid damping are often complemented by results on the flexibility or rigidity of a flow.If a given steady solution of the Euler equations, often a shear flow, has nearby steady solutions, this is said to be flexible, if not it is rigid.If one can show flexibility near a shear flow, that is, that steady non-shear flows exist arbitrarily close to a shear flow, then invicid damping cannot hold.We go on to show that such a family of solutions exists.Lin and Zeng [42] showed that any periodic travelling steady solution with vorticity sufficiently close to a constant must be a shear flow.They also showed a flexibility result, that for arbitrary horizontal period, there exists steady flows in any neighbourhood of the Couette flow in H 5 2 .Constantin, Drivas, and Ginsburg [14] showed that given modest regularity assumptions, shear, steady, non-stagnating flows in a channel can be perturbed to give non-shear, steady, non stagnating flows in a perturbed channel.Hamel and Nadirashvili [26] went on to show that in fact, in a channel, any steady flow which does not stagnate is a shear flow, and in the half plane, any flow which does not stagnate and has bounded velocity is a shear flow.Zelati, Elgindi, and Widmayer consider waves of vorticity, [15], and construct steady analytic flows arbitrarily close to a shear flow.In [49], Sinambela and Zhao considered a large family of shear flows, and constructed steady flows arbitrarily close in H s for arbitrarily large s.The recent work by Franzoi, Masmoudi, and Montalto [21] constructs spatially quasiperiodic steady flows near the Couette flow.
Most of the work mentioned so far considered the periodic case, apart from [21] which considers the quasiperiodic case.As is often done in water wave problems, in order to consider localised solutions instead, we turn to spatial dynamics.This asks us to reformulate our PDE as an evolution equation, with X playing the role of time, and turns out to yield many tools which are especially good at finding non-periodic solutions.One such spatial dynamical tool is a centre manifold reduction used by Kirchgässner [37].However this can only be applied to semi-linear problems.Ours is a quasilinear problem, and for that, we need a stronger centre manifold reduction of Mielke [44].This allows us to reduce to an ODE.Note, this is an exact reduction, not an approximation; the evolution in X of all sufficiently small bounded solutions to our infinite dimensional problem is confined to a two-dimensional manifold.Although we only encounter twodimensional centre manifolds in this paper, with more than two layers we would expect to have higher dimensional manifolds such as those in [8].There are great many papers applying spatial dynamics to water wave problems; see the surveys [23,24] and the references therein for those focusing on gravity-capillary waves.In one of the more classical papers, Kirchgässner used spatial dynamics to construct capillary-gravity waves [38].Groves, Toland, and Buffoni [8] considered the free surface problem and found solitary solutions lying on a four-dimensional manifold, from this centre manifold reduction of Mielke.Groves and Wahlén constructed Stokes waves with vorticity using spatial dynamics in [25].More recently, Wang [55] also uses this reduction, when considering a problem similar to ours, except with free upper boundary and fluids of different densities.Kozlov, Kuznetsov, and Lokharu [40] found solitary gravity waves with vorticity and critical layers.They too considered the rotational case, and like us, they went on to find solutions with streamlines attached to the boundary.In addition to classical spatial dynamics, recently-developed centre manifold techniques 'without a phase space' [20,10] have been applied to problems related to water waves [52,36,12]; also see [1].There are also many classical constructions of solitary water waves which employ fixed-point methods rather than spatial dynamics [41,22,2,3,50].
One body of work that also considers solutions to the two-dimensional Euler equations with continuous velocity but a discontinuity in vorticity is that of vortex patches.These are solutions of the Euler equations with compactly supported vorticity.See the work by Burbea [9], by Hmidi, Mateu and Verdera [31], by Hassainia and Wheeler [29], and by Castro, Córdoba and Gómez-Serrano [11].Hassainia, Masmoudi and Wheeler [28] also find vortex patches with cat's eye structures, which also appear in many papers on water waves with critical layers.
A vorticity front can be thought of as a generalisation of vortex patches, where the vorticity has unbounded support.Hunter, Moreno-Vasquez, Shu and Zhang consider a time dependent problem similar to ours, but without rigid walls [34]: they seek vorticity fronts, as do we, and they also consider perturbations from the same trivial flow, but extended to their infinite-depth domain.It differs in that where we ultimately reduce to a perturbation of the Korteweg-De Vries equation, which is local in x, they reduce to equations involving a Hilbert transform, which is non-local in x.Furthermore, they consider the case where both layers are semi-infinite, and they consider the time-dependent case.They also find that their problem is non-dispersive.Subsequent work [32,33,7] considers an equation that approximates the evolution of the interface between vorticity fronts, and show that despite quadratic nonlinearity in the equation, the approximation holds for cubically non-linear time scales.1.4.Outline of the paper.Firstly, in Section 2 we prove various preliminary results about the background shear flow, conserved quantities, and the stream function.In Section 3 we change coordinates to reformulate the problem into a system of PDEs and boundary conditions on known domains.We then write our PDEs in the form of an evolution equation: ∂ ∂x (u, v, η) = F(u, v, η; c), where the only derivatives to appear in F are with respect to y.Finally, we write F as L + R, where L is the Fréchet derivative of F, and R is the non-linear remainder term.In Section 4 we prove Theorem 1.1.First, in Section 4.1, we verify some hypotheses about the spectrum of L. This allows us to apply the centre manifold theorem in Section 4.2, which lets us reduce the problem to solving an ODE rather than a PDE.Roughly speaking, the centre manifold theorem works by finding two things.The first is a reduction function, ψ : E 0 × R → W, where W is some space of functions in y, and E 0 is the generalised kernel of L. Two functions, ξ 0 (y) and ξ 1 (y), span E 0 , so we often identify E 0 with R 2 , and consider ψ a function on R 3 .The second is an ODE of the form (a ′ (x), b ′ (x)) = F (a(x), b(x), ε), where ε is a small parameter related to c.These are such that if (a, b, ε) solves the ODE, then (x, y) → a(x)ξ 0 (y) + b(x)ξ 1 (y) + ψ(a(x), b(x), ε)(y) is a solution to our problem.The reduction function ψ is not related to, and not to be confused with, the stream function Ψ.In Section 4.2 we find solutions to leading order, and prove Theorem 1.1.
In Section 5 we prove Theorem 1.4.We show in Section 5.1 that η is monotone on x > 0 and x < 0, and that the sign of V can be similarly characterised.In Sections 5.2 and 5.3 we show, given h, ω 0 , ω 1 , how and where the flow stagnates, and whether η is deflected towards or away from this stagnation.In Section 5.3, we deal with certain critical values, where we have that U has the same sign in the entire flow, except for a bounded region of space which contains a stagnation point.
In the appendix we have the details of an argument of Kirchgässner which we allude to in Section 4.2.

Preliminaries
We first introduce some notation and terminology.Sometimes we write a function that depends on, say, (X, Y, ω).This unsubscripted ω should be taken to be equal to For an open set U, then as usual C n (U) denotes the set of n-times continuously differentiable functions with domain U.By C n (U), we mean the Banach space of functions which have domain U, and whose derivatives up to order n exist, can be continuously extended to U, and are bounded.To specify a codomain V, we use the notations C n (U, V) and C n (U, V) respectively.
As usual, incompressibility (1.1a) guarantees the existence of a stream function Ψ which is a first integral of the fluid particle motion.
which clearly vanishes on Y = h + η(X).Using the regularity of U and η, it is straightforward to check that Ψ ∈ C 0 (Ω 0 ∪ Ω 1 ) with Ψ Y = U .Differentiating (2.1) with respect to X and then integrating by parts using (1.1a), we discover that where the boundary terms at Y = h + η(X) cancel thanks to (1.1e).The regularity One important symmetry this problem has is that it is reversible in X.More precisely, if (U, V, η) solves (1.1), then the functions ( Ǔ , V , η) defined by also solve (1.1).Physically, this corresponds to reflecting the particle motion across the Y axis and reversing time.Note that the solutions in Theorem 1.1 are invariant under this symmetry.
We next introduce several well-known invariants for (1.1).Since X plays a time-like role, these are called 'conserved quantities' in what follows.
Lemma 2.2 (Conserved quantities).Let (U, V, η) be a solution to (1.1).Then the quantities We call Q 0 the mass flux of the lower layer and Q 1 the mass flux of the upper layer.The flow force S does not play a role in our later arguments, but it is included here for completeness.
Proof.First consider the mass flux of the lower layer, Q 0 , and find that The second equality is due to (1.1a), and the last is due to (1.1e).An almost identical proof shows the conservation of Q 1 .Now consider the flow force.We see that where the last line uses (1.1a) and (1.1b).Performing integration by parts on the last term in the integrand, and using (1.1e) to eliminate the boundary terms, gives us Shear flows, that is, solutions of (1.1) which do not depend on X, play a key role in the analysis.As we are viewing X as a time-like variable, we call these solutions equilibria.
, where c is an arbitrary constant equal to the horizontal velocity at y = h + η. □ If, as in (1.2), an equilibrium solution has η = 0, we call it a trivial solution.However, notice that for any general equilibrium of the form in (2.3), h can be redefined to give an equilibrium of the form in (1.2).
Many studies on steady solutions of the two-dimensional Euler equations assume the existence of a so-called 'vorticity function' γ such that the ω = U Y − V X = γ(Ψ).This then implies that the stream function Ψ satisfies the semilinear equation ∆Ψ = γ(Ψ).Interestingly, as mentioned in Remark 1.3, many of the solutions that we construct in Theorem 1.1 do not possess a vorticity function.
Proposition 2.4.Let (U * , 0, 0) be the shear solution given by U for some sufficiently small δ depending only on ω 0 and h.If , then there does exist such a function.

Proof of Proposition 2.4.
In what follows all stream functions are defined using (2.1) so that they vanish at their respective interfaces.Thus the stream function corresponding to (U * , 0, 0) is First consider the case ω 0 < 1 − 2h.Using (2.5) and ω 1 = ω 0 − 1, we find are distinct, and hence there cannot exist a global vorticity γ such that ω = γ(Ψ).The ω 0 > 2h − 2 follows by an analogous argument.
It remains to consider the case where 1 − 2h < ω 0 < 2 − 2h.We claim that Ψ < 0 in Ω 0 and Ψ > 0 in Ω 1 .This will imply that ω = γ(Ψ) where γ(t) = ω 0 − H(t) and H is the Heaviside step function.For |Y − h| ≤ √ δ, these strict signs for Ψ follow from the fact that Ψ(X, η(X)) = 0 together the uniform lower bound which holds for δ sufficiently small.For |Y − h| > √ δ, suppose first that ω 0 ≥ 1 − h.Then Ψ * has at most one critical point in (0, h), which is a local minimum, and is strictly increasing on (h, 1).Thus, for δ > 0 sufficiently small, we have max The arguments for ω 0 ≤ 1 − h are similar but with the role of the two layers reversed, and the claim is proved.□

3.1.
Pointwise flattening with a Möbius map.One of main challenges of (1.1) is that it is a free-boundary problem, i.e., the domains Ω 0 and Ω 1 are unknowns.To overcome this, we introduce new coordinates to map the problem onto a known domain, at the cost of making the PDE much more nonlinear.
In particular, we seek a coordinate transformation which flattens interface to a straight line of constant height h, pointwise in X, and fixes the upper and lower boundaries.Thus the new independent variables x = x(X, Y ) and y = y(X, Y ) should satisfy The obvious choice is to let Y depend on y in a piecewise linear fashion 1 .For our purposes, however, it is more convenient to have a globally smooth transformation, and so instead we use Möbius functions to smoothly map Y to y, defining For independent variables we take where c is a parameter.Notice that u = v = η = 0 always gives a solution, in particular the trivial solution from (1.2) with c = c.We now differentiate (3.1b) to get ) and rearrange and differentiate (3.1c) and (3.1d) to find 1 Thus the new independent variables would be for instance with ũ(x, y) = U (X, Y ) and ṽ(x, y) = V (X, Y ) as dependent variables.
for y ∈ (0, h) ∪ (h, 1), while the boundary conditions ) where η ∈ C 2 (R), and u, v 3.2.Properties of this coordinate change.One benefit of this change of coordinates is that there are conserved quantities related to the mass fluxes Q 0 , Q 1 in Lemma 2.2 which have a particularly simple form.We let and call q 0 the pseudoflux of the lower layer, and q 1 the pseudoflux of the upper layer.Notice that for solutions, these are conserved quantities, since and similarly for q 1 , where in the last step we have used Lemma 2.2.We will ultimately restrict our attention to solutions with vanishing pseudofluxes, i.e. to solutions satisfying the constraints ).Furthermore, since the pseudofluxes are conserved quantities, the same is true for any flows which tend to one of these shear flows as x → −∞ or x → +∞, in particular the solitary waves constructed in Theorem 1.1.For a general solution (U, V, η) in the original variables, without welldefined limits as x → ±∞, we can always choose the value of c in (3.1c) so that one of (3.5) and (3.6) is satisfied, but imposing both is an additional restriction.While we do not need the flow force S from Lemma 2.2 in our arguments, we nevertheless record its expression in terms of the new variables for completeness, Note while the pseudofluxes are related to -but not exactly -the mass fluxes, this is exactly the flow force, just written in the new coordinates.The reversibility in (2.2) is also preserved.That is, if (u, v, η; c) satisfies (3.4)-(3.6),and we define then (ǔ, v, η; c) also satisfy (3.4)-(3.6).Lastly, let us describe the x-independent solutions.In physical coordinates, with fixed ω and h, (1.1) has a two-parameter family of equilibria given by (2.3).However, once we enforce the constraints (3.5)-(3.6)on the psuedofluxes, we are left with are most three equilibria.
Lemma 3.1.For fixed h, ω, c, the number of equilibrium solutions of (3.4)-(3.6) is given by the number of real roots η of ) There is always one equilibrium with η = 0. Supposing as in Theorem 1.1 that θ ̸ = 0, there is one non-zero equilibrium when ε = 0, while for 0 < ε < 1 8 θ 2 there are two non-zero equilibria.
Proof of Lemma 3.1.By Lemma 2.3, any equilibrium solution of (3.4) must correspond to an equilibrium solution (2.3) of (1.1).Applying the coordinate transformation (3.1) we deduce that v = 0, η is constant, and (3.9) Enforcing (3.5)-(3.6),i.e., that the pseudofluxes are zero, yields a system of algebraic equations Solving (3.4a) for u x , and substituting (3.11a) to eliminate η x we find while similarly rearranging (3.4b) yields We then abbreviate (3.11) as Since the formula for F does not contain any x derivatives, we can think of x as being fixed, meaning u and v are functions of y only, and η is just a real number.Another way of thinking of this change of viewpoint is that F is a map between Banach spaces of functions in y.We now define these function spaces: let (3.13)More plainly, U is the subset of W × R in which none of the denominators that appear in (3.11) vanish.We see that U is open and for any c * ̸ = 0, U contains (0, 0, 0, c * ).
), and so the derivatives u y and v y are well-defined, as are the pointwise values u(h), v(0), v(1), and v(h).Comparing with (3.11), we conclude that F(u, v, η; c) is a well-defined element of L 2 ((0, 1)) × L 2 ((0, 1)) × R. To verify the integral conditions in the definition of X , let F(u, v, η; c) = (f, g, α).For the first integral condition, we see that The second integral condition follows by a similar argument.Now we show analyticity.Pointwise multiplication by a bounded function, and differentiation are both linear and bounded from W to X , therefore analytic.Similarly, the relevant trace maps are bounded and linear from U → C. As F is a composition of these analytic maps with appropriate rational functions, we conclude that it is also analytic. □ We now separate F into linear and non-linear parts around (u, v, η; c) = (0, 0, 0; c * ) with c = c * ̸ = 0, i.e., around a shear solution which does not stagnate on the interface.More precisely, we define a linear mapping L and nonlinear remainder R by For the moment the value of c * is unspecified, but we will end up focusing on the case where the coefficient function p is defined by As F is analytic by Lemma 3.3, L : W → X is a bounded linear operator while R : U → X is analytic.

The centre manifold
We now state a centre manifold theorem.This is a result due to Mielke [44]; see Theorem 3.3 of [27], and Theorem 2.1 of [45].It considers the general differential equation In our application, the variable w(x, y) corresponds to (u(x, y), v(x, y), η(x)), and ε = c − c * .
Theorem 4.1 (Centre manifold theorem).Let X , W be Hilbert spaces with W continuously embedded in X .Suppose a bounded linear operator L : W → X has spectrum σ(L), and satisfies the following three hypotheses: consists only of finitely many eigenvalues, all of which have finite algebraic multiplicity.(ii) There exist R > 0, C > 0, such that for all k ∈ R with |k| > R, we have that L − ikI : X → X is invertible, and satisfies (iii) There exists δ > 0 such that there are no z ∈ σ(L) satisfying 0 < | Re(z)| < δ.Suppose we also have a nonlinear function R and an integer N ≥ 2, such that there exists a neighbourhood U ⊂ W × R of 0 such that R ∈ C N (U, X ), and that R(0, 0) = 0, D w R(0, 0) = 0.
Define E 0 ⊂ W to be the generalised eigenspace corresponding to the purely imaginary eigenvalues of L. Let P 0 be a continuous projection onto E 0 which commutes with L, and whose kernel is called E h .
(b) M 0 (ε) contains the set of bounded solutions of (4.1a) staying in V w for all x ∈ R, i.e., if w is a solution of (4.1a) satisfying w(x) ∈ V w for all x ∈ R, then w(0) ∈ M 0 (ε).(c) Suppose w 0 satisfies the reduced equation and w 0 (x) ∈ V w ∩ E 0 for all x.Then w 0 + ψ(w 0 , ε) is a solution to the full problem, i.e., satisfies (4.1a).
4.1.Analysis of the linearised operator.We show that the linear operator L defined in (3.14) satisfies Hypotheses (i)-(iii) of Theorem 4.1.
Finding the spectrum of L is made easier by the fact that L is Fredholm.To show this, we first consider the slightly simpler operator L given by Proof.First we show L has trivial kernel.This is not a particularly complicated calculation, but variants of it will appear frequently, so this will serve as a simple example.Suppose (u, v, η) ∈ W is in the kernel of L. This means we seek continuous functions of y, namely u and v, and a real number η satisfying p(y)v(h) − v y = 0 for y ∈ (0, h) ∪ (h, 1) (4.2a) where (4.2a)-(4.2c)are just L(u, v, η) = 0 rewritten while (4.2d)-(4.2e)and the continuity of u and v are imposed by (u, v, η) ∈ W.
Notice also that L is of full range.In particular, the solution to L(u, v, η) = (f, g, α) can be shown by a direct calculation, which we leave to the reader, to be It can also be shown straightforwardly that this (u, v, η) ∈ W. Therefore, L is invertible as desired.□ Corollary 4.3.The spectrum of L consists only of eigenvalues.
Proof.We have just shown L is invertible, and hence in particular that it is Fredholm with index 0. Therefore, since L − L has finite-dimensional range, and therefore is a compact operator, L is also Fredholm of index 0. Similarly, L − zI is Fredholm of index 0 for all z ∈ C. Therefore L − zI is injective if and only if it is surjective.This means the spectrum of L is exactly those z ∈ C for which L − zI has non-trivial kernel.
□ We are now ready to verify Hypothesis (i) of Theorem 4.1.Writing z = ik, we ask for which k ∈ C the operator L − ikI has non-trivial kernel.The answer turns out to be largely captured by the dispersion relation: where However, not all eigenvalues are given by solutions to (4.3).If then L − zI has non-trivial kernel, no matter the values of c * , ω 0 , ω 1 .We are primarily concerned with (4.3) for two reasons.Firstly, (4.4) is very rarely satisfied; in particular, it requires h to be rational, which generically is not true.Secondly, even if there exists a z * satisfying (4.4), this does not contradict our hypotheses.This is because this z * must be real with |z * | > π, and the hypotheses are only concerned with eigenvalues on, or arbitrarily close to, the imaginary axis.Proof.Suppose (u, v, η) for some constants A and B. Therefore, inserting the formula for v into the first component of (L − ikI)(u, v, η) = 0 yields The conditions on v(h) and the continuity of u give the equations Proof.We seek solutions to L(u, v, η) = 0.It can be easily shown using the equations that come directly from L(u, v, η) = 0, the conditions on v, and the integral conditions on u, that if a non-trivial kernel exists, it must be spanned by a vector of the form The choice c * = h(1 − h) then guarantees p is such that u is continuous at h, so that we indeed have a one-dimensional kernel.We now show that there are no other purely imaginary eigenvalues.Lemma 4.4 means that this is equivalent to showing that (4.3) has no non-zero real roots.We differentiate d with respect to k, and see Therefore d(k) is strictly decreasing for k < 0 and strictly increasing for k > 0, so it achieves its minimum value only at k = 0. Hence, for k ∈ R \ {0}, we have Therefore (4.3) has no nonzero real roots, thus 0 is the only purely imaginary eigenvalue.We now seek second order generalised eigenvectors, i.e., solutions to L(u, v, η) = (u * , 0, 1).Without loss of generality, we can take η = 0. Therefore we need to solve which has solution (0, v * , 0) where Performing a similar process to find a third order generalised eigenvector shows that none exist, so we conclude that when c * = h(1−h), 0 is an eigenvector, and has algebraic multiplicity 2. □ Corollary 4.6.For c * = h(1 − h), L satisfies Hypothesis (i) of Theorem 4.1.
Proof.We have shown that for this value of c * , L has one purely imaginary eigenvalue, and it has finite algebraic multiplicity.Furthermore, since L is Fredholm, its spectrum consists only of eigenvalues, therefore σ 0 = {0}.Thus, Hypothesis (i) is satisfied.□ We now verify Hypothesis (ii).Since σ 0 is finite, the inverse (L − ikI) −1 , which we view as a linear operator X → X , exists for sufficiently large k.Proof.Suppose (L − ikI)(u, v, η) = (f, g, α).For the remainder of this proof, ∥ • ∥ should be taken to mean the L 2 ((0, 1)) norm, and P = 1 + ∥p∥.
Our first step will be to find a bound on |v(h)| which grows sub-linearly in k.We see that

Thinking of this as a quadratic in |v(h)|, then applying Cauchy-Schwarz gives us
which is a bound of the desired form.Now we multiply the first component of (L−ikI)(u, v, η) = (f, g, α) by ū, the complex conjugate of the second by v, and subtract to get Integrating over [0, 1] yields The boundary conditions on v mean that ū(0)v(0) = ū(1)v(1) = 0, and so using this and Cauchy-Schwarz, Then, applying Cauchy-Schwarz shows that The third component of (L − ikI)(u, v, η) = (f, g, α) tells us that v(h)/c * − ikη = α.Applying this, then (4.10) to the final term of the previous line yields Finally, rearranging, we see For sufficiently large |k|, the denominator is larger than 1 2 |k|.Therefore there exists R > 0 such that for |k| > R, |k|∥(u, v, η)∥ X ≤ C∥(f, g, α)∥ X as required, where the constant C depends only on R, h, ω 0 , ω 1 , and c * . □ We now verify Hypothesis (iii), and show a spectral gap around the imaginary axis.This could in principle be done by only studying the solutions to (4.3).However, we use an easier proof which makes use of Hypothesis (ii) as well.
Proposition 4.8.There exists δ > 0 such that there are no z in the spectrum of L satisfying 0 < | Re(z)| < δ.
Proof.We show this by contradiction.Suppose this spectral gap does not exist.Hence, there must exist a sequence (z n ) in the spectrum such that Re(z n ) → 0. Define where σ(L) is the spectrum of L. None of the z ∈ S can satisfy (4.4), so for all z ∈ S, we have that (4.3) holds, i.e., that d(−iz) = c −1 * .Since d(−iz) is a meromorphic function in z, S is an isolated set.In other words, (z n ) is a sequence taking values in an isolated set, with Re(z n ) → 0, so we conclude However, by Hypothesis (ii) there exist constants C, R such that for all z ∈ iR with |z| > R, we have is also injective.Therefore, for all n > n * , we have that z n is in the resolvent set, but we assumed it was in the spectrum, and so arrive at our contradiction.□ We have now verified Hypotheses (i)-(iii) of Theorem 4.1 in the case c * = h(1 − h), and so proceed for the rest of the paper with c * = h(1 − h).This specific value of c * is of particular interest, as it corresponds to k = 0, so we might expect to find waves of infinite period, i.e., solitary waves.
The final step in our analysis of L is to seek the projection P 0 from Theorem 4.1.Notice that because 0 is the only purely imaginary element of the spectrum, E 0 is the generalised kernel of L. Lemma 4.5 gives us a basis {ξ 0 , ξ 1 } of E 0 , where Thus we can write P 0 as where A, B : X → R are bounded linear functions.Calculating we see that L and P 0 commute if and only if A • L = B and B • L = 0.
Proof.This can be verified by brute force.However, we now show a little of the construction, in order to give some intuition about why the projection is of this form.We where G is any primitive of g.This implies Now construct B such that it satisfies B • L = 0. Motivated by (4.14), let V be a primitive of v, and define If L had an inverse, solving A • L = B for A would give a formula for A. Doing exactly this is impossible, as L is not invertible, but we can use the fact that B depends only on v to do something very similar.Notice that L(u, v, η) = (f, g, a) implies that Notice that Therefore A and B exhibit the required behaviour when composed with L. It remains only to check that P 0 is indeed a projection onto the generalised kernel of L. We see that A(u * , 0, 1) = B(0, v * , 0) = 1, and B(u * , 0, 1) = A(0, v * , 0) = 0, so we are done.□

4.2.
Applying the centre manifold theorem.We quickly verify the hypotheses of Theorem 4.1 on the non-linear operator R. First, fix an integer N ≥ 4 which will not change for the rest of this paper.Although Theorem 4.1 only needs this N to be greater than or equal to 2, we need it to be greater than or equal to 4 to make some arguments about the regularity of the solutions we find.See, for example, (4.23b) and Appendix A.
We are now ready to apply Theorem 4.
then it must be of the of the form where ξ 0 , ξ 1 are defined in (4.11) and a and b are scalar functions.
We now introduce and clarify some points of notation.Firstly, since E 0 ×R, the domain of ψ, is a 3 dimensional space, we abuse notation slightly, and sometimes consider ψ to be a function on R 3 given by ψ(a, b, ε) = ψ(aξ 0 + bξ 1 , ε).
Secondly, since ψ takes values in W, we denote its u, v and η components as e 1 • ψ, e 2 • ψ, and e 3 • ψ respectively.This means (4.16) can be written completely equivalently as Reversibility immediately gives a symmetry on ψ.Recall we defined reversibility in (2.2).Using Theorem 3.15 of [27], we see that i.e., that the first and third components of ψ are even in b, and the second component is odd in b.
In order for a solution of the form in (4.16) to exist, a and b must satisfy a differential equation, which is found using P 0 .Inserting (4.16) into (4.15)yields We now apply P 0 to find the reduced equation (4.1d).The composition P 0 • ψ = 0, therefore P 0 • Dψ = 0.By construction, P 0 commutes with L, therefore applying P 0 to (4. 19) gives We can rewrite (4.20) more compactly as Let e 1 • F and e 2 • F denote the first and second components of F respectively.Since F is a composition of C N functions, it follows that F ∈ C N (R 3 , R 2 ).Note that F inherits the reversibility of ψ shown in (4.18).Specifically, e 1 • F is odd in b, and e 2 • F is even in b.This means that if a(x) and b(x) solve (4.21) and we define then ǎ(x) and b(x) also satisfy (4.21).
We now find expansions for ψ and F .Since for all ε, the zero function u = v = η = 0 is a solution of (3.4), Theorem 4.1 tells us that ψ(0, 0, ε) = 0. Since ψ a (0, 0, 0) = ψ b (0, 0, 0) = 0, we have Therefore, 22b) Notice that we know P 0 R(aξ 0 + bξ 1 , ε) explicitly.Substituting (4.22b) into (4.20),we obtain an expansion for F .However, appealing to parity can sharpen the order of the error term.Since e 1 • F is odd in b, and e 2 • F is even in b, we deduce where recall θ was defined in (1.3).Notice that these really are the dominating terms, in that their coefficients are non-zero.The terms written here can be found by examining P 0 R(ax 0 + bξ 1 ).Higher order terms can be found explicitly by a recursive method, but for our purposes this is unnecessary.
We now solve (4.21) for a and b.For fixed ε, by Picard-Lindelöf, (4.21) with initial condition (a(0), b(0)) = (a * , 0) has a unique solution.By reversibility, this solution curve in the a, b plane is unchanged under reflection about the a axis.In other words, for fixed ε, solutions to (4.21) that cross the a axis are symmetric about the a axis.This motivates us to focus on solutions where a is even, and b is odd.
We now rescale with the aim of eventually eliminating higher order terms, and proceed having picked ε to be strictly positive.Let Notice that F satisfies

.26)
When ε = 0, (4.26) is exactly the KdV equation for travelling waves, a nonlinear equation.We expect the solution to persist as ε moves away from 0, and describe such a solution as weakly nonlinear.There are several ways to show that the solutions do indeed persist, for example, by considering the stable and unstable manifolds of F at the fixed point at the origin.Instead, we adapt an argument of Kirchgässner [38].
(a) ãε is even and bε is odd.
(d) The explicit formula for (ã 0 , b0 ), is given by ã0 , where the O is with respect to the C 1 (R) norm.
Proof.We argue as in Proposition 5.1 of [38].For more details, see Appendix A. □ Remark 4.11.Since the phase portrait of F is symmetric about the line b = 0, the only zeroes of bε are at x = 0. Indeed, if bε vanished at two distinct points, then bε would be periodic rather than homoclinic.The oddness of bε implies that bε (0) = 0, therefore this zero is unique.We explicitly find the leading order terms in C 0 (R, H 1 ((0, 1))), and show that the solution is small in C 1 (R, H 1 ((0, 1))).Proposition 4.12.For all sufficiently small ε ≥ 0, there exist functions u ε (x, y), v ε (x, y), η ε (x) such that: 2 ) (4.28b) (d) These functions satisfy the estimates 2 ), sup (e) These solutions are homoclinic to 0 in H 1 , i.e., for any fixed ε, we have Proof.Applying the inverse of the scaling from (4.24), and appealing to Proposition 4.10(f) yields where the O is with respect to the C 1 norm.We now consider (4.17).This, the parity of the components of ψ with respect to b, and the parity of a ε and b ε (x) gives us the required parity on u ε , v ε and η ε .Using (4.22a), and the fact that the oddness of e 2 • ψ with respect to b gives us that e 2 • ψ = bO(|a| + |b| + |ε|), we see that ) where the Os are with respect to the C 0 (R, H 1 ((0, 1))) norm.Differentiating (4.17) with respect to x gives We deduce (4.29a) from this, and the fact Dψ = O(|a| + |b| + |ε|).For the other two estimates, notice that we can solve (3.4a) and (3.4b) for u y and v y , as the determinant of the linear system is nonzero so long as η, η x are sufficiently small.Doing so shows that (u ε y (x, y), v ε y (x, y)) is an analytic function of ε, u ε (x, y), u ε (x, h), u ε x (x, y), v ε (x, y), v ε (x, h), v ε x (x, y), η ε (x), and η ε x (x).Calling this function f : R 9 → R 2 , we know all the arguments of f are O(ε) uniformly in x and y.Since f (0) = 0, we see we have the estimates (4.29b).
Finally, these solutions are homoclinic to 0 because (ã ε , bε ) is homoclinic to 0, thus (a ε , b ε ) is homoclinic to 0. Since ψ(0, 0, ε) = 0, we have (u ε , v ε , η ε ) is homoclinic to 0. □ We would like to convert Proposition 4.12 into statements about U ε , V ε , and η ε .Before we do this, however, it is useful to determine their regularity.Viewing (u, v, η) as a function of x, the centre manifold theorem tells us (u ε , v ε , η ε ) ∈ C N (R, X ) for ε sufficiently small.This gives us immediately that η ∈ C N (R).Showing U ε and V ε are analytic is a little more involved.Firstly we show that u ε and v ε are in H 1 .We can then use the smoothness of the coordinate change (3.1) to show that U ε and V ε are in H 1 .This is enough regularity on the functions to mean we can apply elliptic regularity results, in particular, by examining (1.1), we see U ε and V ε are harmonic, and thus real analytic, away from the interface.
Proof.We know that for each x, the function y → u ε (x, y) has a weak derivative in L 2 ((0, 1)), which we call u y (x, y).We also know from the centre manifold theorem that ∥u(x, • )∥ H 1 ((0,1)) varies continuously with x.We have that u by continuity.A straightforward application of Fubini's theorem then shows that u ε y is the weak derivative of u ε in H 1 ((−M, M ) × (0, 1)).
We now find the x derivative.The fact that u ε ∈ C 1 (R, H 1 ((0, 1)) means that for each x, there exists a function of y in H 1 ((0, 1)), which we call u ε x (x, y), such that ∥u ε x (x, • )∥ H 1 ((0,1)) varies continuously with x, and that lim h→0 Therefore, since the trace map is bounded on H 1 ((0, 1)), we have that for all y, or in other words, u ε x (x, y) is the classical partial derivative of u ε (x, y), so it is certainly the weak derivative as well.We have that u by continuity.Rewriting (4.32) with u ε x replaced by u ε shows that u ε is in L 2 ((−M, M )× (0, 1)) as well.Therefore, since first order weak derivatives of u ε exist, and they and u ε are in L 2 ((−M, M ) × (0, 1)), we have that u ε ∈ H 1 ((−M, M ) × (0, 1)) as required.
The case for v ε follows similarly.□ We now use this to show a similar result for U ε and V ε .
Lemma 4.14.For all M > 0, we have Proof.Using Lemma 4.13, it is straightforward to show that Now in order to show that U ε and V ε are in H 1 ((−M, M ) × (0, 1)), all we need to show is that composition with the change of coordinates gives a bounded linear map from It can be seen that that where the second equality uses the formula for y, (3.1b), and the estimates on η ε in (4.28c) and (4.29a).We deduce the change of coordinates is a bounded linear map on H 1 ∩ C ∞ , but this is a dense subspace of H 1 , so we have our required result.□ We are now ready to apply elliptic regularity to show U ε and V ε are analytic.
Proof.Notice that U ε and V ε are weakly harmonic in Ω 0 and in Ω 1 , as for all ϕ ∈ C ∞ c (Ω i ), we have and similarly for V ε .By standard elliptic regularity arguments, U ε and V ε are therefore smooth harmonic functions, and hence real-analytic.□ We are now ready to prove results on U ε and V ε of a similar form to those in Proposition 4.12.
(a) U ε , η ε are even in x, and 2 ), and as required, and we already know η ε (X) is homoclinic to 0.
For conclusion (c), the estimate on η ε X comes immediately from (4.29a).As for the other two, again, we argue similarly to in Lemma 4.14, and consider (3.3) to see □ We now have everything we need to show that U ε , V ε , η ε satisfy the estimates in Theorem 1.1.
Proof.Recall (3.1), and the definition of η ε from Theorem 1.1.It will be useful here to define the quantities Notice that h + = h + O(ε 2 ), and similarly for h − .We have just shown the required result for η in Conclusion (c) of Proposition 4.12.We deal with U first.In what follows, X and Y should be understood to mean X(x, y) and Y (x, y) respectively.Using (3.1c) and (4.28), we see that Notice that despite the discontinuities in ω and H, we have that ωy(1 − y) + h 2 (1 − h) 2 H is continuous at y = h, and it is in fact Lipschitz continuous.Therefore, for y between 0 and h + , replacing ω with ω 0 and H with −1/h in (4.34) introduces errors of O(ε 2 ).More precisely, for (x, y) ∈ R × [0, h + ], we have that where as usual, the error term is uniform in (x, y).Notice that if (X, Y ) satisfies Y ≤ h + η ε (X), then it also satisfies Y ≤ Y (x, h + ).This means that when we apply the coordinate transformation (3.1b), we see that for all X, Y satisfying Y ≤ h + η ε (X), we have We argue similarly to conclude that for all X, Y satisfying Y ≥ h + η ε (X), we have Showing the result for V follows a similar but more straightforward argument.We see for (X, Y ) with Y ≤ h + η ε (X), we have 2 ), and similarly, for (X, Y ) with Y ≥ h + η ε (X), we have 2 ).□ In a sense, (4.28) is a more precise approximation of the solution than U ε , V ε , η ε .This is because in (4.28), the discontinuity in the derivatives happens at y = h, which is where it occurs in the exact solution, but the discontinuity in the derivatives of U ε ,

Streamline patterns
In this section we will prove Theorem 1.4.Both the qualitative nature of the results we show, and the general method used are similar in spirit to those in, for example [53].From now onward we assume without loss of generality that ω 0 ≤ 1 − h.We lose no generality from this, as if it is not the case, then performing the vertical reflection and relabelling Ω 0 , Ω 1 , puts us in the regime where once again ω 0 ≤ 1 − h holds.After performing this reflection, we still have the vorticity of lower layer is one greater than the vorticity of the upper layer.
5.1.Signs of components of the flow.One thing we can deduce quite straightforwardly is a sign for U ε on the interface, as (5.2) We will next show that the interface is monotone, that is that η ε x (x) has one sign for x > 0 and the opposite sign for x < 0. Proposition 5.1.For all x, η ε x (x) has the same sign as θx.Proof.At x = 0 we have that a ε x (0) = η ε x (0) = 0 by evenness of a ε and η ε .We now show the x ̸ = 0 case.A corollary of Proposition 4.10 is that for all x, bε (x) shares a sign with x.Therefore, for x ̸ = 0, b ε (x) θx > 0. (5.3) In particular, by Remark 4.11, b ε (x) = 0 if and only if x = 0. We also know using (4.23a), that for x ̸ = 0, (5.4) Therefore, when x ̸ = 0, we know that a ε x (x) ̸ = 0. Consider (4.31).We are interested in η ε x so we take the third component of each vector, then divide the resulting equation by a ε x .Since ψ a and ψ b are 3-component vectors, we continue our notational convention, and denote their third components by e 3 • ψ a and e 3 • ψ b respectively.We also know from (4.18) that e 3 • ψ is even in b and C N , hence e 3 • ψ b is odd in b and C N −1 , therefore Putting all this together, we get We know from differentiability of ψ a that e 3 where we have used the fact that b ε x (x) = O(ε 2 ), which can be deduced from the form of F in (4.21), the rescaling (4.24) and the boundedness of ãε and bε .Therefore, we see that η ε x , b ε , and θx share a sign.□ Corollary 5.2.For all x, we have θη ε (x) < 0.
Proof.We know by conclusion (b) of Corollary 4.16 that η ε tends to 0 as x tends to ±∞.This and Proposition 5.1 give us the result.□ We now show that V ε also shares a sign with θX, except on the upper and lower boundaries, where we know V ε = 0. We use a maximum principle argument, one possible reference for which would be [19, Section 6.4, Theorem 3].Proposition 5.3.For all (X, Y ) ∈ R × (0, 1), we have that V ε (X, Y ) and θX share a sign.
Proof.Notice that by oddness, we have immediately that V ε (0, Y ) = 0. We now consider the case where X > 0 and θ > 0, but the other cases follow an almost identical argument.
we can apply the strong maximum principle, which says V ε attains a minimum on Ω M i , and does so only on ∂Ω M i .We know that V ε (X, 0) = 0 from (1.1c), and that V ε (0, Y ) = 0 by oddness.We also know from (1.1e) that We know from Corollary 4.16 that ∥V ε (M, • )∥ H 1 ((0,1)) tends to 0 as M tends to infinity.This means we can use an argument by contradiction to show that V ε > 0 on Ω ∞ i .Suppose there exists ( X, Ŷ ) with X > 0 such that V ε ( X, Ŷ ) < 0. Since ∥V ε (X, • )∥ H 1 ((0,1)) → 0, there exists M > X such that for all Y ∈ [0, 1], we have Hence, the minimum V ε attains on ∂Ω M i is no less than 1 2 V ε ( X, Ŷ ), and so V ε does not attain its minimum on the boundary, contradicting the maximum principle.Thus, for all M > 0, we have that V ε ≥ 0 on Ω M i .However, we know that V ε does not attain its minimum in Ω M i , therefore V ε > 0 on Ω M i .Therefore, for all X > 0, Y ∈ (0, 1), we have The cases where θ < 0 or X < 0 follow very similarly, and we conclude that for all (X, Y ) ∈ R × (0, 1), we have that V ε (X, Y ) and θX share a sign.□ We are now in a position to start describing the critical layers and stagnation points of the flow.Recall the critical layer of the solution is the set on which U ε = 0. We start by identifying the stagnation points in the trivial solutions.If ω 0 = 1 − h, then the trivial solution stagnates at y = 0 and at y = 1.Otherwise the trivial solution stagnates in the upper layer at y = h − h(1 − h)/ω 1 , since recall we are now assuming ω 0 ≤ 1 − h.
We now investigate the stagnation points of the non-trivial solutions.We separate into two cases, which occupy different regions of parameter space.It turns out that the critical layers are qualitatively very different; in the first case they are unbounded, and in the second, bounded.5.2.Unbounded critical layer.We first assume that 1 − h > ω 0 , i.e., that stagnation occurs in Ω 1 , rather than on the upper or lower boundaries.Notice that if we also have 1 − 2h > ω 0 , then we are outside the ω = γ(Ψ) regime.
We now show there is a critical layer, and investigate some of its properties.
Proof.First notice that in Ω 0 , we have so any zeros of U ε must be in Ω 1 .Notice for all X, we have and (5.2) tell us that U ε (X, h + η ε (X)) > 0. Therefore, by the intermediate value theorem, for all X, there exists Y * ∈ (h + η ε (X), 1) such that U ε (X, Y * ) = 0.
We now show this Y * is unique.By (1.1b), in Ω 1 , This implies that for each X, U ε is a strictly decreasing function of Y , so in fact, for every X, there exists a unique Y * (X) ∈ (h + η ε (X), 1) such that U ε (X, Y * ) = 0. Since U ε is an analytic function, and U ε Y (X, Y * (X)) ̸ = 0, we can apply the analytic implicit function theorem to see that Y * is an analytic function of X. □ Remark 5.5.Solving (3.1c), we see that if . By smoothness of the coordinate transformations, we see that , or in other words, the distance between the critical layer of our solution, and of the background shear flow is of order ε, uniformly in X.In particular, for ω 0 < 1 − h, the critical layer Y * does not touch the boundary Y = 1 or the interface.
We can deduce from this lemma that we have a stagnation point at (0, Y * (0)), and will now examine its behaviour.
Proof.We know that U ε (0, Y ) = 0 if and only if Y = Y * (0), and that for Y ̸ = 0, 1, the vertical velocity V ε (X, Y ) = 0 if and only if X = 0 by Propostition 5.3.Therefore we have a unique stagnation point at (0, Y * (0)).We now consider the nature of this stagnation point.Let s, s ′ be the eigenvalues of the derivative of (U ε , V ε ) at a stagnation point.We know from incompressibility, and the fact we have real solutions that s+s ′ = 0, and ss ′ ∈ R. Hence, if the determinant of the derivative of (U ε , V ε ) is positive, s and s ′ are both purely imaginary, and we have a centre.If the determinant is negative, they are both real, and we have a saddle.Notice that by evenness of U ε , we have U ε X (0, Y ) = 0 for all Y .Using this and (1.1b), we have We know that in (0, ∞) × (0, 1), we have θV ε > 0, that V ε (0, Y ) = 0, and that V ε is harmonic, so we can apply the Hopf lemma to deduce that for all Y ∈ (0, 1), we have θV ε X (0, Y ) > 0. We also see from Conclusion (c) of Corollary 4.16 that and in particular, this is true at Y = Y * (0), the stagnation point.In other words, if θ > 0, we have a centre, if θ < 0, we have a saddle.□ 5.3.Bounded critical layer.We now consider the case ω 0 = 1 − h, ω 0 − ω 1 = 1.This implies θ = 2h − 1.Notice that it is sufficient to consider cases when h ̸ = 1 2 , as θ = 0 if and only if h = 1 2 .Notice also, that we can reflect in Y , as we did in (5.1), to insist that h > 1 2 .This will imply that θ > 0. We will now show that the critical layer is bounded, and ends where it intersects with the upper boundary.
(5.7) Note, this is not true in the unbounded critical layer case.We conclude that by the intermediate value theorem, and evenness of U ε , there exists X * > 0 such that U ε (X * , 1) = U ε (−X * , 1) = 0.
We now note that for all X ∈ (−X * , X * ), we have U ε (X, 1) < 0 < U ε (X, h + η(X)), therefore there exists Y * (X) ∈ (h + η ε (X), 1) such that U ε (X, Y * (X)) = 0. We see from (1.1b) and the fact that ω 1 < 0 that in Ω 1 , ) < 0, and similarly in Ω 0 , we have U ε Y (X, Y ) > 0. Therefore, since U ε Y (X, 0) > 0, for each X this Y * (X) is unique.We can apply the analytic implicit function theorem to conclude that Y * (X) is an analytic function.Next, we seek to show that Y * (X) → 1 as X → X * and as X → −X * .Notice that Thus, there exists δ > 0 such that for all t ∈ (0, δ), we have U ε (X * − t, 1) < 0 < U ε (X * − t, 1 − t), therefore 1 − t < Y * (X * − t) < 1, so we are done.Define the bounded region R as in (5.6).We know that U ε (X, Y * (X)) = 0, and that U (X, 1) < 0 for |X| < X * .Therefore we can apply the maximum principle to show that in R, we have U ε < 0. We can also show that U ε > 0 outside R. To do this we first recall that we deduced in (5.7) that there exists M > X * such that for all X with |X| ≥ M , and all Y , we have U ε (X, Y ) > 0.
For X ∈ (−M, M ), we can apply the Hopf lemma to V ε at (X, 0), and deduce with (1.1a) that XU ε X (X, 0) > 0. We already showed at the very start of this proof that U ε (0, 0) > 0, and so we see that for all X, U ε (X, 0) > 0. Therefore we apply the maximum principle to U ε on ((−M, M ) × (0, 1)) \ R, and conclude that on this set too, U ε > 0.
Finally, we discuss the stagnation point.We deduce from what we have just discussed that there is a stagnation point at (0, Y * (0)), and that there are no others away from the boundaries, and apply the same arguments as in Theorem 5.6 to conclude that since the stagnation happens in Ω 1 , and θ > 0, we have a centre. □ We now show that in this region of parameter space, we have a streamline which is attached to the upper layer, as shown in Figure 2.

3 . 3 .
) for η ∈ (−h, 1 − h) and c ∈ R. Either equation can be uniquely solved for c as a function of η, and eliminating c yields (3.8) as desired.□The evolution equation.We now reformulate (3.4) as an evolution equation in x.The first step is to algebraically solve (3.4) for the derivatives u x , v x , and η x .Solving the boundary condition (3.4e) for η x gives us

4. 3 .
Regularity and estimates.We now undo the scaling from (4.24) to find functions a ε (x), b ε (x) satisfying(4.21).By conclusion (c) of Theorem 4.1, inserting these into the reduction function does indeed give a solution to(3.4).Then undoing the coordinate change (3.1) gives us solutions to (1.1).