Falling Drop in an Unbounded Liquid Reservoir: Steady-state Solutions

The equations governing the motion of a three-dimensional liquid drop moving freely in an unbounded liquid reservoir under the influence of a gravitational force are investigated. Provided the (constant) densities in the two liquids are sufficiently close, existence of a steady-state solution is shown. The proof is based on a suitable linearization of the equations. A setting of function spaces is introduced in which the corresponding linear operator acts as a homeomorphism.


Introduction
Consider a drop of liquid with density ρ 1 submerged into an unbounded reservoir of liquid with density ρ 2 .Assume the liquids are immiscible.We investigate the motion of the drop under the influence of a constant gravitational force and surface tension on the interface.Specifically, we shall show existence of a steady-state solution to the governing equations of motion, provided the difference ρ 1 − ρ 2 of the densities is sufficiently small.
The dynamics of a falling (or rising) drop in a quiescent fluid has attracted a lot of attention in the field of fluid mechanics.Such flows have been studied extensively both experimentally and numerically with truly fascinating outcomes (see [3] for a comprehensive overview and further references), but it remains an intriguing task to analytically validate the observations.The observed dynamics can be characterized as a series of bifurcations with respect to the Reynolds number as parameter.Broadly speaking, steady-state solutions are observed for small Reynolds numbers, with bifurcations into oscillating motions as the Reynolds number increases.Bifurcations into more complex solutions can be observed as the Reynolds number increases even further.
In the following, we shall investigate the steady-state solutions corresponding to small Reynolds numbers.A small Reynolds number is equivalent to a small density difference ρ 1 − ρ 2 .One of our aims is to develop a functional analytic framework that can be used not only to study steady states but also as foundation for further investigations into the dynamics described above.In particular, the framework should facilitate a stability analysis of the steady states, and an investigation of the Hopf-type bifurcations (into oscillating motions) observed in experiments.For this purpose, it should satisfy certain properties.First and foremost, it should be possible to identify function spaces within the framework such that the differential operator of the linearized equations of motion acts as a homeomorphism.Second, the framework should have a natural extension to a suitable time-periodic framework (recall that a steady-state solution is trivially also time-periodic).Third, the framework should adequately facilitate a spectral analysis of the operators obtained by linearizing the equations of motion around a steady state.To meet these criteria, we propose a framework of Sobolev spaces.Although a setting of Sobolev spaces seems natural, and by far the most convenient to work with, it is by no means trivial to identify one that conforms to the problem of a freely falling (or rising) drop.Indeed, one of the novelties of this article is the introduction of such a Sobolevspace setting that meets at least the first and most important criteria, and possibly also the other two, mentioned above, and in which existence of steady-state solutions can be shown effortlessly for small data.The investigation of steady-state solutions is not new, though.It was initiated by Bemelmans [5] and advanced by Solonnikov [12,13].However, the analysis carried out by Bemelmans and Solonnikov do not lead to a framework of Sobolev spaces.Indeed, for reasons that will be explained in detail below, the approaches of both Bemelmans and Solonnikov cannot be adapted to a Sobolev-space setting with the desired properties.
We shall consider the most commonly used model for two-phase flows with surface tension on the interface.It is assumed both fluids are Navier-Stokes liquids, that is, incompressible, viscous, and Newtonian.It is further assumed that the fluids are immis-cible with surface tension on their interface in normal direction proportional to the mean curvature.Moreover, we consider a system in which the drop is a ball B R 0 of radius R 0 when no external forces act on the system, that is, in its stress free configuration.If we choose a coordinate system attached to the falling drop, these assumptions lead to the following equations of motion for a steady state (see Section 2 for details on the derivation): (1.1) Here, Γ η denotes the interface between the two liquids, which we may assume to be a closed manifold parameterized by a "height" function η∶ ∂B R 0 → R describing the displacement of the drop's boundary points in normal direction.The domain Ω η ⊂ R 3 bounded by Γ η describes the domain occupied by the drop, and the exterior domain the region of the liquid reservoir.The drop velocity −λe 3 , λ ∈ R, is assumed to be directed along the axis of the (constant) gravitational force ge 3 .The first two equations in (1.1) are the Navier-Stokes equations written in a moving frame of reference, where u∶ R 3 ∖ Γ η → R 3 denotes the Eulerian velocity field of the liquids, p∶ R 3 ∖ Γ η → R the scalar pressure field, and T(u, p) denotes the corresponding Cauchy stress tensor.The density function ρ∶ R 3 ∖ Γ η → R is constant in both components of R 3 ∖ Γ η .The third equation states that the surface tension in normal direction on the interface Γ η is proportional to the mean curvature H, with σ > 0 a constant.The notation ⟦ ⋅ ⟧ is used to denote the jump of a quantity across Γ η .Immiscibility of the two liquids under a no-slip assumption at the interface is expressed via the fourth and fifth equation.Observe that the normal velocity on the interface then coincides with that of the moving frame, which moves with the same velocity −λe 3 as the falling drop.The equations are augmented with a volume condition for the drop and the requirement that the liquid in the reservoir is at rest at spatial infinity in the sixth and seventh equation, respectively.
A key part of our investigation is directed towards finding an appropriate linearization of (1.1) with respect to the unknowns u, p, λ and η.The canonical linearization, i.e., around the trivial state (0, 0, 0, 0), leads to the Navier-Stokes equations (1.1) 1-2 being replaced with the Stokes system (1.2) An analysis based on this linearization would have to be carried out in a setting of function spaces conforming to the properties of the Stokes problem.Such a setting, however, is not suitable for an investigation of the exterior domain Navier-Stokes equations in a moving frame.Since the falling drop, and thus the frame of reference, moves with a nonzero velocity −λe 3 , the appropriate linearization of the Navier-Stokes equations in the exterior domain is an Oseen system.At least in a setting of classical Sobolev spaces, the steady-state exterior-domain Navier-Stokes equations in a moving frame can only be solved in a framework of Sobolev spaces conforming to the Oseen linearization.
To resolve this issue, we propose to rewrite the system (1.1) as a perturbation around a state (u 0 , p 0 , λ 0 , η 0 ) with λ 0 ≠ 0. A subsequent linearization of (1.1) then yields the Oseen problem 3) The main challenge, and indeed novelty of this article, is to determine a suitable state (u 0 , p 0 , λ 0 , η 0 ) that renders the problem well posed in a framework of classical Sobolev spaces.
The starting point of our investigation were the articles [12,13] by Solonnikov, which contain a number of truly outstanding ideas on how to analyze (1.1).However, Solonnikov overlooks the necessity of an Oseen linearization as described above.Instead, he employs a Stokes linearization and consequently a setting of function spaces in which the nonlinear term λ∂ 3 u cannot be correctly treated on the right-hand side.Our approach resolves this issue.
We derive the steady-state equations of motion for the falling drop and state the main theorem in the following Section 2. The aforementioned framework of Sobolev spaces is then introduced in Section 4. Fundamental L r estimates are established in Section 5, and a reformulation of (1.1) in a fixed reference configuration in Section 6.The linearization around a non-trivial state is carried out in Section 7. In Section 8 we show in Theorem 8.1 that the operator corresponding to this linearization is a homeomorphism in our framework of Sobolev spaces, which finally enables us to establish a proof of the main theorem, namely the existence of a steady-state solution for ρ 1 − ρ 2 sufficiently small.

Equations of motion and statement of the main theorem
We derive the system of equations governing the motion of a freely falling drop in a liquid under the influence of a constant gravitational force.We shall express these equations in a frame of reference with origin in the barycenter of the drop.More specifically, we denote by ξ(t) the barycenter of the falling drop with respect to an inertial frame, whose coordinates we denote by y, and express the equations of motion in barycentric coordinates x(t, y) ∶= y − ξ(t).In these coordinates, the domain Ω (1) t ⊂ R 3 occupied by the drop at time t satisfies Ω (1) t x dx = 0. (2.1) We let Ω (2) denote the domain of the surrounding liquid reservoir, and put t .The surface Γ t ∶= ∂Ω (1) t describes the interface between the two liquids.Moreover, we let µ 1 , µ 2 and ρ 1 , ρ 2 denote the constant viscosities and densities of the drop and the liquid reservoir, respectively.The functions then describe the viscosity and density of the liquid occupying the point x at a given time t.Expressed in a frame of reference attached to the barycenter ξ, the conservation of momentum and mass of both liquids is described by the Navier-Stokes system where v denotes the Eulerian velocity field in the liquids, p the pressure, the Cauchy stress tensor, and b ∈ R 3 a constant gravitational acceleration.One can decompose the velocity field and pressure term into t → R describing the liquid flow in the drop, and another part describing the flow in the reservoir.We employ the notation Γ to denote the jump in a quantity on the interface between the two liquids.Concerning the physical nature of the interface, we make the basic assumption that slippage between the two liquids cannot occur, i.e., a no-slip boundary condition, and that liquid cannot be absorbed in the interface.Consequently, there is no jump in the velocity field neither in tangential nor in normal direction: Since the liquids are immiscible, the normal component of the liquid velocity at the interface coincides with the velocity of the interface itself.If Φ Γ denotes a Lagrangian description of the interface in barycentric coordinates, the immiscibility condition therefore takes the form (2.4) In the classical two-phase flow model, surface tension on the interface, i.e., the difference in normal stresses of the two liquids, is proportional to the mean curvature in normal direction and in balance in tangential direction: Since we consider the motion of a drop in a quiescent liquid, the velocity in the reservoir vanishes at spatial infinity lim Due to the incompressibility of the liquid drop, its volume is constant.Since we consider a drop that takes the shape of the ball B R 0 in its stress free configuration, this volume is prescribed by In conclusion, the system obtained by combining (2.1)-(2.8)governs the motion of a liquid drop falling freely in a liquid reservoir under the influence of a constant gravitational force.
In this article, we seek to establish existence of a steady-state solution, that is, a time-independent solution to (2.1)- (2.8).Such a solution is of course only steady with respect to the chosen frame of reference; in our case the frame attached to the barycenter.Other types of steady states can be investigated by analyzing time-independent solutions in other frames.For example, it is conceivable that falling drops can perform steady rotating motions, which should be investigated by considering the equations of motion in a rotating frame of reference.
The unknowns in (2.1)-(2.8)are the functions v, p, ξ, Φ Γ .The mean curvature H can be computed from Φ Γ .The viscosities µ 1 , µ 2 > 0, surface tension σ > 0 and the prescribed volume 4π 3 R 3 0 of the drop are constants, which may be chosen arbitrarily.Also the gravitational force b ∈ R 3 is an arbitrary constant, but upon a re-orientation of the coordinates we may assume without loss of generality that it is directed along the negative e 3 axis, i.e., b = −ge 3 with g > 0. The constant densities ρ 1 , ρ 2 > 0 shall be restricted to pairs whose difference ρ 1 − ρ 2 is sufficiently small.In this sense, we treat ρ 1 − ρ 2 as the data of the system.Since the geometry Ω (1) t , Γ t of the problem is determined by the unknown description Φ Γ of the interface, (2.1)-(2.8) is a free boundary problem.
At the outset, it is clear that (2.9) can have multiple solutions.This is best illustrated by considering ρ 1 = ρ 2 , in which case the trivial solution with v = 0, λ = 0 and constant pressures p (1) , p (2) is a steady-state solution if σH equals the constant hydrostatic pressure difference p (1) − p (2) between the drop and the reservoir.Since a constant mean curvature H is realized whenever Ω (1) is a multiple of disjoint balls, we obtain for each Φ Γ describing one or more spheres a trivial solution by adjusting the hydrostatic pressure difference accordingly (depending on the fixed volume Ω (1) ).In the case (2.9) above, the fixed volume of Ω (1) coincides with the volume of the ball B R 0 .With constant pressures satisfying p (1) − p (2) = 2 R 0 , the ball B R 0 therefore becomes an admissible steady-state drop configuration when ρ 1 = ρ 2 .We shall single out this configuration for further investigation in the sense that we investigate non-trivial steady-states with a configuration close to the ball B R 0 for ρ 1 ≠ ρ 2 with ρ 1 − ρ 2 sufficiently small.
From a physical perspective, a smallness condition is only meaningful when expressed in a non-dimensional form.In order to obtain a dimensionless formulation of (2.9), we choose R 0 as characteristic length scale, V 0 ∶= √ gR 0 as the characteristic velocity, ρ 1 + ρ 2 as characteristic density, (ρ 1 + ρ 2 )R 0 V 0 as the characteristic viscosity, and (ρ 1 + ρ 2 )R 0 V 2 0 as the characteristic surface tension.Investigating the resulting non-dimensional equations of motion, we will establish existence of a non-trivial steady-state solution with drop configuration close to the unit ball B 1 .For this purpose, it is convenient to introduce (in the non-dimensionalized coordinates) the normalized pressures We then obtain the following system of non-dimensional equations: x dx = 0. (2.10) Observe that the mean curvature now appears in the form (H + 2) that vanishes if Γ is the unit sphere, which means that (v, p, λ) = (0, 0, 0) is a trivial solution when ρ 1 −ρ 2 = 0. We shall employ a parameterization of Γ over the unit sphere S 2 ⊂ R 3 and subsequently linearize (2.10).The linearization of the operator σ(H + 2), however, has a non-trivial kernel.To circumvent an introduction of the corresponding compatibility conditions, we employ an idea from [13] and replace the two equations x dx = 0 (2.11) in (2.10) with the equations x dx (2.12) The resulting system then reads x dx (2.13) The systems (2.10) and (2.13) are equivalent.Clearly, (2.10) implies (2.13).To verify the reverse implication, observe that (2.13) 5-7 imply Ω (1)   x dx x dx x dx The matrix ∫ Γ n ⊗ n dS is symmetric positive definite and thus invertible.Consequently, the equation above implies ∫ Ω (1) x dx = 0. We conclude that (2.13) implies (2.10).
Since we investigate existence of non-trivial steady-states in a drop configuration close to the ball B 1 (in non-dimensionalized coordinates) under the restriction that the difference in densities of the two liquids is sufficiently small, it is convenient to introduce ρ ∶= ρ 1 − ρ 2 as smallness parameter.Moreover, it is convenient to parameterize the interface Γ via a height function η ∶ S 2 → R that describes the drop's displacement in normal direction with respect to its unit sphere S 2 ⊂ R 3 stress-free configuration.The geometry then becomes a function of η: η ∪ Ω (2) η .
The system of steady-state equations of motion finally takes the form with respect to unknowns (v, p, λ, η).
As the main result in the article we prove existence of a solution to the steady-state equations of motion (2.14) under a smallness condition on the density difference ρ.

Theorem 2.1 (Main Theorem).
There is an ε > 0 such that for 0 < ρ ≤ ε there is a solution to (2.14).The solution is smooth up to the interface, that is, Moreover, it possesses the integrability properties and admits the representation for all ε > 0, where Γ λ Oseen denotes the Oseen fundamental solution1 .The solution is symmetric with respect to rotations leaving e 3 invariant: and the velocity λ of the drop's barycenter is non-vanishing.
By far the most challenging part of proving Theorem 2.1 is to establish the existence of a solution.As mentioned in the introduction, via a perturbation around a non-trivial state we are able to solve the system in a setting of Sobolev spaces adapted from the 3D exterior-domain Oseen linearization of the Navier-Stokes equations.Consequently, we are led to a solution with the integrability properties (2.16).The symmetry (2.18) follows from the observation that (2.14) is invariant with respect to rotations leaving e 3 invariant.Higher-order regularity is obtained via a standard approach utilizing the ellipticity of (2.14), while the asymptotic profile (2.17) is a direct consequence of (2.16) and and a celebrated result of Babenko [4] and Galdi [6].Observe that the coefficient vector in the asymptotic expansion, which at the outset is given by η e 3 acting on the liquid drop, that is, the difference of the gravitational force and the buoyancy force.

Notation
We use capital letters to denote global constants in the proofs and theorems, and small letters for local constants appearing in the proofs.
By B R ∶= B R (0) we denote a ball in R n centered at 0 with radius R.Moreover, we let for a domain Ω ⊂ R n .Additionally, we use S 2 ∶= ∂B 1 to denote the unit sphere.By we denote the twofold half space, which is the union of the two domains We use the notation (x ′ , x 3 ) for a vector Lebesgue spaces are denoted by L q (Ω) with associated norms ⋅ q,Ω .By W k,q (Ω) we denote the corresponding Sobolev space of order k ∈ N 0 with norm ⋅ k,q,Ω , and we introduce the subspaces Moreover, W −k,q (Ω) and W −k,q 0 (Ω) denote the dual spaces of W k,q ′ (Ω) and W k,q ′ 0 (Ω), respectively, where q ′ ∶= q q−1 .We further introduce homogeneous Sobolev spaces D k,q (Ω) defined by and the corresponding seminorm In general, D k,p (Ω) is not a Banach space.However, ⋅ k,q,Ω defines a norm on C ∞ 0 (Ω), and the completion is therefore a Banach space.By Sobolev's Embedding Theorem, D k,q 0 (Ω) can be identified with a subspace of L 1 loc (Ω) if kq < 3. We denote its dual space by D −k,q ′ 0 (Ω).For a sufficiently smooth manifold Γ ⊂ R 3 and s > 0, s ∈ N, we let W s,q (Γ) denote the Sobolev-Slobodeckij space of order s with norm ⋅ s,q,Γ .

Preliminaries
In this section we introduce a bespoke framework of Sobolev spaces for the investigation of (2.14).For this purpose, we let Ω ⊂ R 3 denote a domain of the same type as in Section 2, that is, we assume Ω (1) ⊂ R 3 is a bounded domain and Ω (2) ∶= R 3 ∖ Ω (1) is a domain, Γ ∶= ∂Ω (1) , For a function u ∶ Ω → R we use the abbreviations .
The function n = n Γ denotes the unit outer normal at Γ.If u is sufficiently regular, we set where the restriction has to be understood in the trace sense.Furthermore, δ(Ω) denotes the diameter of Ω (1) .When considering a function u ∶ Ω → R, we often have to distinguish between its properties on the disjoint sub-domains Ω (1) and Ω (2) .To this end, function spaces of the type (1) , u (2) ∈ X (2)   are introduced.Equipped with the norm such a space X is isomorphic to the direct sum of the spaces X (1) and X (2) .Clearly, X is a Banach space if X (1) and X (2) are so.Let q ∈ 1, 3 2 and r ∈ (3, ∞).For λ 0 ∈ R, λ 0 ≠ 0 the space X q,r,λ 0 Oseen ∶= X q,r,λ 0 Oseen (Ω (2) equipped with the norm is the canonical solution space for solutions to the exterior domain Oseen problem for forcing terms f in L q Ω (2) ∩ L r Ω (2) ; see for example [7,Chapter VII.7].Let and The bespoke framework of Sobolev spaces we shall employ in our investigation of (2.14) is then given by In Theorem 8.1 we show that the operator corresponding to the appropriate linearization of (2.14) maps X q,r λ 0 homeomorphically onto Y q,r .The following embedding is valid:
Proof.The above estimates for the part u (1) of u defined on a bounded domain follows directly from Sobolev embedding theorems.Concerning the part u (2) defined on an exterior domain, it follows from [7, Lemma II.6.1] and the Sobolev inequality that 1,λ 0 .

Auxiliary linear problem
Let Ω be a domain of the same type as in Section 4, i.e., satisfying (4.1).We further assume that the boundary Γ is at least Lipschitz.The linear system is an integral part of the linearization of (2.14).It is a two-phase strongly coupled Oseen (λ 0 ≠ 0) or Stokes (λ 0 = 0) system.Since the coupling is strong, the questions of existence and uniqueness of solutions as well as a priori estimates hereof cannot be investigated by means of a simple decomposition into two classical Oseen/Stokes problems.In the following, we carry out an analysis of (5.1) in the framework of the Sobolev spaces introduced in the previous section.Existence and uniqueness of solutions is first shown in a setting of weak solutions, and strong a priori estimates of Agmon-Douglis-Nirenberg type are subsequently established first in the half space and then in the general case via a localization technique.The main result of the section is contained in Theorem 5.9 and Theorem 5.10.

Weak solutions
We introduce a weak formulation of (5.1) in the setting of the function spaces: Ω (1)   p dx = 0 .
In the following, we establish existence and uniqueness as well as higher-order regularity of weak solutions to (5.1) in this framework.We start with the definition of a weak solution: as well as div u = g in Ω and u ⋅ n = h 1 on Γ.
Existence of a weak solution u can be shown by standard techniques; we sketch a proof below.
Theorem 5.2.Assume that the boundary Γ is Lipschitz.For every 3 satisfying where Proof.We sketch a proof of existence following [7, Proof of Theorem VII.2.1] based on a Galerkin approximation.To this end, a Schauder basis constructed.This function space is clearly separable, whence such a basis can be constructed via a Gram-Schmidt procedure.We consider first the case (g, h 1 ) = (0, 0).Existence of an approximate solution of order m ∈ N, that is, a vector field u m ∶= ∑ m l=1 ξ l ϕ l satisfying the equation in (5.2) for all test functions in span{ϕ 1 , . . ., ϕ m }, then follows directly from the fact that the matrix an approximate solution u m .Employing u m itself as a test function in the weak formulation, one obtains a uniform bound on S(u m ) 2 , which, since u m is solenoidal, also implies a uniform bound as in (5.4) on u m 1,2 .A weak solution to (5.1) is now obtained as the limit u ∶= lim m→∞ u m in H .The general case of non-vanishing g and h 1 follows by a lifting argument.Employing a right inverse of the . The compatibility condition (5.3) ensures that ∫ Ω (1) g − div u 1 dx = 0 so that we can find u 2 ∈ D 1,2 0 (R 3 ) with div u 2 = g − div u 1 and satisfying u 2 = 0 on Γ as well as (5.4); see for example [7,Theorem III.3.1 and III.3.6].The ansatz u = v + u 1 + u 2 now reduces the problem to the case above with respect to the unknown v.We thus conclude existence of a weak solution.
A pressure p can be associated to a weak solution u such that (u, p) becomes a solution to (5.1) in the sense of distributions: and with and consider the functional which is continuous on H M by Sobolev embedding.We further introduce the space p dx = 0 and the operator div ∶ H M → L 2 0,M .The operator is surjective, which is seen by solving for arbitrary p ∈ L 2 0,M the two equations according to [7, Theorem III.3.1].It follows that the operator and hence also its adjoint div * are both closed.Since u is a weak solution, F M vanishes on the kernel of div and consequently belongs to the image of div * .We thus obtain a function After possibly adding a constant to p M +1 , we may assume p M +1 = p M in B M .The sequence {p M } ∞ M =1 then induces a pressure p ∈ L 2 loc (R 3 ) satisfying (5.5) and ∫ Ω (1) p dx = 0.It remains to establish L 2 (R 3 ) integrability of p.If λ 0 = 0, the functional F M remains continuous if H M is replaced with H.In this case the argument above directly yields a pressure p ∈ L 2 0 (R) 3 satisfying (5.5).Subsequently choosing a function ϕ ∈ H with div ϕ = p in R 3 and ϕ 1,2 ≤ c 0 p 2 , which can be done via [7, Theorem III.3.1 and Theorem III.3.6], one obtains (5.6) by inserting ϕ into (5.5).If λ 0 ≠ 0, it suffices to observe that (u, p) is a weak solution to an Oseen problem in the exterior domain Ω (2) , whence [7, Theorem VII.7.2] yields p ∈ L 2 0 (R) 3 satisfying (5.5)-(5.6).
Provided u and p are sufficiently regular, integration by parts in (5.5) reveals that (u, p) is a classical solution to (5.1).Higher-order regularity of (u, p) can be obtained via a classical approach under appropriate regularity assumptions on the data: 3 , then a weak solution u ∈ D 1,2 0 (R 3 ) 3 to (5.1) with associated pressure p ∈ L 2 loc (R 3 ) satisfying (5.5) also satisfies Proof.The proof is a standard application of a well-known technique based on difference quotients.In fact, with only minimal modification it is similar to a proof of higher-order regularity for solutions to the Stokes system with prescribed normal velocity and tangential stress on the boundary; see [14, Proof of Theorem 2].For the sake of completeness, we sketch the proof.We include only the case h 1 = 0 and k = 0.The general case h 1 ≠ 0 and k > 0 follows by a simple lifting technique and iteration procedure, respectively.Since higher-order regularity in Ω away from the boundary Γ is well known for Stokes systems (see for example [7, Section IV.2]), we focus on regularity up to the boundary Γ.To this end, consider an arbitrary x ∈ Γ and choose a cube Q r (x) ⊂ R 3 , centered at x with side length r, such that Γ ∩ Q r (x) can be parameterized by a C 3 function ω.
Without loss of generality, we may assume that x = 0 and where Q ′ r (0) ⊂ R 2 is the two-dimensional cube centered around 0, and that ∇ω(0) = 0 as well as We introduce test functions The transformed fields (U, P) satisfy the weak formulation where F 0 contains up to first-order terms of u and zeroth-order terms of p, and F 1 contains first-order terms of U multiplied with components of ∇ω.The magnitude of the latter terms can be made small by choosing r small.Difference quotients are denoted by . Importantly, difference quotients D −h l D h l U in tangential direction l = 1, 2 are admissible as test functions in W 1,2 0,Γ 0 Q r (0) and can therefore be inserted into (5.8),which yields an estimate of S(D h l U ) 2 in terms of lower-order norms of u and p as well as D h l P 2 .A similar bound on ∇D h l U 2 follows from Korn's inequality.Choosing in (5.8) a test function l P, a bound on D h l P 2 in terms of lower-order norms of u and p is obtained.Such a test function is constructed by setting Existence of solutions to the two equations above and the estimates ) follows as a combination of ∂ 3 div U = ∂ 3 G and the regularity of U 's tangential derivatives.Finally, the distributional derivative ∂ 3 P can now be isolated in (5.8) to deduce in each half of the cube , where O(x) is a neighborhood of x and O (1) (x) ∶= O(x) ∩ Ω (1) , O (2) (x) ∶= O(x) ∩ Ω (2) .Higher-order regularity of (u, p) up to the boundary Γ is thereby established.
Finally, uniqueness of a weak solution to (5.1) can be established.In fact, uniqueness can be obtained in a much larger class of distributional solutions with even less summability at spatial infinity than u ∈ L 6 (R 3 ) satisfied by a weak solution via Sobolev embedding.The theorem below is not optimal in this respect, but suffices for the purposes of this article.

Twofold half space
The main challenge towards a priori L r estimates of solutions to (5.1), i.e., a priori estimates of Agmon-Douglis-Nirenberg type, is to obtain such estimates in the halfspace case under disregard of the lower-order terms in the equations; the general case then follows via a localization argument.We therefore first consider system where n = −e 3 .We shall implicitly identify ∂ Ṙ3 with R 2 .
In the celebrated work [2] of Agmon, Douglis and Nirenberg, a priori L r estimates for strong solutions to elliptic systems with boundary values of a certain type were established.For Stokes systems, both the Dirichlet boundary condition and the slip boundary condition that make up the boundary values in (5.9) fall into the category of so-called Agmon-Douglis-Nirenberg problems covered by [2].However, since the system (5.9) is strongly coupled, it does not itself fall into this category, nor can it be decomposed into two systems to which the estimates from [2] can be applied separately.It therefore seems unavoidable that L r estimates for solutions to (5.9) have to be established without the help of [2].We present a comprehensive proof below.We shall not employ the technique from [2] based on singular integrals, but choose instead an approach based on Fourier multipliers and real interpolation that seems particularly well suited for coupled systems such as (5.9).
The main result on a priori L r estimates of solutions to (5.9) is contained in Theorem 5.8 below.The proof is divided into two lemmas, Lemma 5.6 and Lemma 5.7, which the reader may find interesting in their own right.For technical reasons, it is convenient to decompose both the data and the solution to (5.9) into one part with lower frequency support and another part with higher frequency support in tangential directions e 1 , e 2 .We shall repeatedly employ the Fourier transform F R 2 with respect to these two directions.To this end, observe that F R 2 u(⋅, x 3 ) (ξ ′ ) is well-defined in the sense of distributions S ′ (R 3 ) when u ∈ L r (R 3 ) for some r ∈ (1, ∞), which will be the case whenever such an expression appears below. where Proof.A solution to (5.10) can be constructed explicitly.To this end, consider first a sufficiently smooth right-hand side b ∈ S (R 2 ) (5.12) Therefore p satisfies ξ ′ 2 p − ∂ 2 3 p = 0 and thus We insert p into (5.12) 1 and (5.12) 2 and solve the resulting differential equations.Taking into account the boundary conditions (5.12) 4 , we obtain Inserting the above formula for û into (5.12) 3 , we find that Consequently, a solution to (5.10) is given by where Although M b has a singularity, (u, p) as defined above is a well-defined solution, smooth on Ṙ3 even, due to the assumption that b(ξ ′ ) has support away from 0. In order to provide an estimate for the solution, we let κ 1 4 ∈ C ∞ 0 (R 2 ) with κ 1 4 = 0 on B 1 4 (0) and κ 1 4 = 1 on R 2 ∖ B 1 2 (0), and consider the truncation of the solution operator.The singularity of M ϕ makes it necessary to employ the truncation κ 1 4 to ensure that K is well-defined.We shall use real interpolation to show that K extends to a bounded operator K ∶ W 2−1 r,r (R 2 ) → W 2,r ( Ṙ3 ).To this end, we observe for m ∈ N 0 and any x 3 ∈ R that the symbol ξ ′ ↦ ( ξ ′ x 3 ) m e − ξ ′ x 3 is an L r (R 2 )-multiplier.Specifically, one may verify that sup whence it follows from the Marcinkiewicz Multiplier Theorem (see for example [9, Corollary 6.2.5]) that the Fourier-multiplier operator with symbol ξ ′ ↦ ( ξ ′ x 3 ) m e − ξ ′ x 3 is a bounded operator on L r (R 2 ) with operator norm independent of x 3 , that is, sup < ∞. (5.15) We return to (5.14) and employ (5.15) to deduce where the restriction in the norm of the left-hand side to the twofold real line Ṙ in the second estimate is required since ∂ x 3 K(ϕ) has a singularity at x 3 = 0.It follows that This estimate shall serve as an interpolation endpoint.To obtain the opposite endpoint, we again employ (5.15) to infer sup where the last estimate relies on the truncation introduced in K.It follows that (5.17) Consequently, (5.16) and (5.17) imply whence K extends to a bounded operator K ∶ W 2−1 r,r (R 2 ) → W 2,r ( Ṙ3 ).Recalling the formula (5.13) for the solution u to (5.10) and that supp F 0) via a standard mollifier procedure, we conclude the lemma in its entirety. where corresponding to the Fourier multiplier appearing in (5.20) extends to a bounded operator , and we can therefore introduce the corresponding solution (u, p) ∈ W 2,r ( Ṙ3 ) 3 × W 1,r ( Ṙ3 ) to (5.10) from Lemma 5.6.By construction, u ⋅ n = H 1 on ∂ Ṙ3 .Moreover, recalling (5.13) we compute Consequently, (u, p) is a solution to (5.18).Employing (5.11) we deduce and conclude the lemma.
Theorem 5.8.Let r ∈ (1, ∞) and Then all solutions (u, p) ∈ W 2,r ( Ṙ3 ) 3 × W 1,r ( Ṙ3 ) to (5.9) satisfy where Proof.We decompose both the solution and the data into one part with lower and another part with higher frequency support in tangential directions e 1 , e 2 .For this purpose, we introduce cut-off functions κ α ∈ C ∞ 0 (R 2 ) with κ α = 0 on B α (0) and κ α = 1 on R 2 ∖ B 2α (0), and put Similarly, we introduce p # , p and f # , g # , h 0# , h 1# , h 2# .Observe that (u # , p # ) solves (5.9) with respect to data (f # , g # , h 0# , h 1# , h 2# ).We shall construct another solution satisfying estimate (5.21), and subsequently show that it coincides with (u # , p # ).To this end, we let g + # ∈ W 1,r (R 3 ) denote an extension of g # R 3 + to W 1,r (R 3 ).Specifically employing Heesten's extension operator (see for example [1,Theorem 4.26]) one readily verifies that the extension retains the property that the Fourier transform (in tangential directions) is well defined.Similarly, we introduce an extension of we then obtain a field G # ∈ W 2,r ( Ṙ3 ) 3 with div G # = g # in Ṙ3 .Moreover, a straightforward application of Marcinkiewicz's Multiplier Theorem (see for example [8, Corollary 5.2.5]) yields Owing to the fact that supported away from 0, the expressions above are well defined and yield functions with Moreover, another straight-forward application of Marcinkiewicz's Multiplier Theorem yields Utilizing Lemma 5.6, we construct a solution Finally, by Lemma 5.7 there is a solution

A priori estimates for strong solutions
We return to the linearized two-phase-flow Navier-Stokes problem (5.1),where Ω is a domain of the same type as in Section 4, i.e., satisfying (4.1).Based on the estimates obtained in the twofold-half-space case in Theorem 5.8, we shall establish L r estimates of solutions to (5.1).The Oseen case (λ 0 ≠ 0) and Stokes case (λ 0 = 0) are treated separately in Theorem 5.9 and Theorem 5.10, respectively.
Proof.We first consider data g dx = We fix an R > δ(Ω) and observe that (u, p) ∈ W 2,r (Ω 2R ) 3 × W 1,r (Ω 2R ) by Sobolev embedding.According to the regularity assumptions, Γ can be covered by a finite number of balls Γ ⊂ ⋃ m i=1 B r i (x i ) each of which upon a rotation R i can be mapped to B r i (0) by a with ∇Φ i ∞ arbitrarily small for sufficiently small radii r i , i = 1, . . ., m.The covering can clearly be augmented with bounded open sets O 1 ⊂⊂ Ω (1) and Employing a partition of unity subordinate to such a covering, we can decompose and transform the solution (u, p) into m solutions (u i , p i ) ∈ W 2,r ( Ṙ3 ) 3 × W 1,r ( Ṙ3 ), i = 1, . . ., m, to the twofold half-space Stokes problem (5.9), two solutions (u m+1 , p m+1 ), (u m+2 , p m+2 ) ∈ W 2,r (R 3 ) 3 × W 1,r (R 3 ) to a whole-space Stokes problem, and finally one solution (w, q) ∈ D 1,2 0 (R 3 ) 3 × L 2 0 (R 3 ) to the whole-space Oseen problem (5.28) In all three cases, the data contain lower-order terms of u and p supported in B 2R .Furthermore, the data in the twofold half-space Stokes equations satisfied by (u i , p i ), i = 1, . . ., m, also contain higher-order terms of u and p supported in B 2R and multiplied with components of ∇Φ i .By Sobolev embeddings, we have w ∈ D 1,2 0 (R 3 ) ↪ L 6 (R 3 ), and it is therefore easy to verify, for example by applying the Fourier transform in (5.28), that (w, q) coincides with the solution from [7, Theorem VII.4.1] and therefore satisfies (5.29) with a constant c 1 = c 1 (q, r, λ) independent of λ 0 .A similar estimate is satisfied by the solutions (u m+1 , p m+1 ) and (u m+2 , p m+2 ) to the whole-space Stokes problems by [7, Theorem IV.2.1].Moreover, Theorem 5.8 implies that (u i , p i ), i = 1, . . ., m, also satisfies the estimate, provided a covering is chosen with ∇Φ i ∞ sufficiently small in relation to C 6 so that the higher-order terms can be absorbed on the left-hand side.We thus conclude where c 2 = c 2 (Γ, q, r, λ) > 0. It remains to show that the lower-order terms of u and p on the right-hand side can be neglected.This can be achieved by a standard contradiction argument.Assuming that ∃c > 0 ∀0 < λ 0 < λ ∀solutions (u, p) ∈ X q,r 1,λ 0 × X q,r 2 w.r.t.data (5.27 does not hold, one can utilize (5.30) to construct a sequence (λ n , u n , p n ) normalized such that u n W 1,r (Ω 2R ) + p n L r (Ω 2R ) = 1 and with λ n → λ and (u n , p n ) weakly convergent in the Banach space to a solution (u, p) to (5.1) with parameter λ ∈ [0, λ] and homogeneous right-hand side.The restriction q < 3 2 is critical in this step.Theorem 5.5 implies ).We conclude (5.31).Therefore, the lower-order terms of u and p on the right-hand side in (5.30) can be neglected, which yields (5.26).Uniqueness of the solution follows from Theorem 5.5, and the theorem is thereby established for data satisfying (5.27).However, it is easy to verify that data satisfying (5.27) are dense in the space Y q,r 1 × Y q,r 2,3 × Y q,r 4 .Consequently, the general case follows by a density argument.

Reformulation on a fixed domain
The steady-state equations of motion as expressed in (2.14) in a frame attached to the barycenter of the falling drop form a classical free boundary problem.Specifically, the boundary Γ depends on the unknown height function η.For further analysis it is necessary to refer all unknowns in this so-called current configuration to a fixed domain reference configuration.This section is devoted to such a reformulation.
As mentioned in the introduction and further elaborated on in Section 2, we investigate a falling drop whose stress-free configuration, i.e., the configuration when the density in the two liquids is the same, is the unit ball B 1 in non-dimensionalized coordinates.Our aim is to establish existence of steady-state configurations close to the stress-free configuration B 1 for small density differences.Canonically, we therefore choose as the fixed liquid reference domain.
In order to refer the equations of motion to Ω 0 , we first construct a suitable coordinate transformation Φ η based on the height function η.For technical reasons, it is important that Φ η retains any rotational symmetry possessed by η.Lemma 6.1.Let r ∈ (3, ∞).There is an extension operator The extension operator is invariant with respect to rotations, that is, for all R ∈ SO(3): If r > 3, there is a δ 0 > 0 such that for any η ∈ W 3−1 r,r (S 2 ) with η W 3−1 r,r < δ 0 the mapping is continuous and maps Since the Laplace operator is rotational invariant, also the solution H η is invariant with respect to rotations of the data η.Let χ ∈ C ∞ 0 (R) be a cut-off function with χ(s) = 1 for s ≤ 2 and χ(s) = 1 for s ≥ 3. Putting we obtain an operator with the desired properties.Observe that Moreover, by (6.1) we clearly have det ∇Φ η = det I + ∇E(η) > 0 when η W 3−1 r,r (S 2 ) is sufficiently small.In this case, Φ η is a C 2 -diffeomorphism onto its image Ω by the global inverse function theorem of Hadamard.
We shall use Φ η to change the coordinates and consequently express (2.14) in the reference configuration Ω 0 .To this end, we set In order to simplify the notation, we put and introduce the transformed stress tensor Observe that an application of the Piola identity yields div T η (w, q) = J η div T(v, p) ○ Φ η and div(A η w) = J η (div w) ○ Φ η .
The normal vector n Γ at Γ expressed in the coordinates of the reference configuration is given by , and the transformed tangential projection by With this notation, the steady-state equations of motion (2.14) take the following form in the reference configuration: with respect to unknowns (w, q, λ, η).We use the notation n = n S 2 in the following.
In the next step, we exploit an inherent symmetry in (6.9) and simplify the system by replacing (6.9) 7 with We shall a posteriori verify that a solution to the simplified system exhibits axial symmetry around e 3 and consequently satisfies Consequently, a solution to the simplified system (1 + η) 3 − 1 dS = 0, lim x →∞ w(x) = 0 (6.10) with unknowns (w, q, λ, η) is also a solution to (6.9).The analysis in the remaining part of the article is carried out on the system (6.10).

Linearization
A main challenge is to identify a suitable linearization of (6.10) such that the fully nonlinear system can be solved via a perturbation technique.Indeed, as explained in the introduction, the trivial linearization obtained by neglecting all nonlinear terms is not suitable since it leads to a Stokes-type rather than an Oseen-type problem.Instead, we shall linearize the equations around a non-trivial first-order approximation.
In order to identify the first-order approximation, we utilize an idea going back to Happel and Brenner [10] and introduce as auxiliary field a solution to the system  By adding a constant to P (1) , that is, replacing P with we may assume, by choosing the constant C appropriately, that This choice of λ 0 (ρ) combined with the fact that the symmetry (7.4) implies ⟦T(U, P)n⟧ dS = 0 (j = 1, 2) means that (λ 0 (ρ)U, λ 0 (ρ)P, λ 0 (ρ), 0) is a solution to the trivial linearization of (6.10) around the zero state, that is, to the system obtained by neglecting in (6.10) all nonlinear terms with respect to (w, q, λ, η).The state (λ 0 (ρ)U, λ 0 (ρ)P, λ 0 (ρ), 0) can therefore be seen as a first-order approximation of the solution to (6.10).
To conclude the linearization, we express the mean curvature H on Γ as a function of η.As in [11, Section 2.2.5], we obtain where ∆ S 2 and ∇ S 2 denote the Laplace-Beltrami operator and the surface gradient on the unit sphere S 2 , respectively, and Then we have containing all the nonlinear terms.
Proof.Recalling (6.7), we observe that I − A η 1 contains only terms of first and second order with respect to components of ∇E(η 1 ).Utilizing that W 2,r (R 3 ∖ S 2 ) is an algebra for r > 3, and the Sobolev embedding W The first assertion of the lemma then follows from (6.1) in Lemma 6.1.The next assertions follows in a similar manner.Concerning the estimates involving F −1 η 1 , we recall from (6.5)-(6.7)that Consequently, we obtain an estimate of I − F −1 η 1 W 1,∞ as above, provided J η 1 is bounded away from 0. To this end, we recall (6.6) and choose δ 1 so small that J η 1 > 1 2 for η 1 W 3−1 r,r ≤ δ 1 .One may now verify the rest of the assertions analogously.
The linearization (7.10) is a result of expressing the velocity field and pressure term as a perturbation (7.9) around a truncated auxiliary field (U R , P R ).The truncation is necessary to avoid right-hand side terms in (7.10) with inadmissible decay properties.Instead, compactly supported right-hand side terms appear.Suitable estimates of these terms are established in the following lemma.In particular, the magnitude of their norms are estimated in terms of the distance R of the truncation χ R from the drop domain: Lemma 8.3.Let q ∈ 1, 3  2 , r ∈ (3, ∞) and δ 1 be the constant from Lemma 8.2.For all where C 12 = C 12 (q, r, δ 1 ).Moreover, where C 13 = C 13 (q, r, δ 1 ).
Proof.Let η = η 1 .Recalling from Lemma 6.1 that A η (x) = F η (x) = I for x ≥ 4, we utilize Lemma 8.2 to estimate div T Recalling the truncation (7.8), the pointwise decay of the auxiliary fields (7.3), and that supp Since r > q, we obtain an even better estimate for div T(U R , P R ) r with respect to decay in R, and thus conclude the first assertion of the lemma.The other inequalities in (8.8) follow in a similar manner.
The most critical estimate in (8.9) is the second one.Employing Lemma 8.2 together with the integrability properties (7.2) and the pointwise decay (7.3) of the auxiliary fields, we conclude since R > 4. The remaining estimates in (8.9) are verified in a similar fashion.
We are now in a position to show existence of a solution to (7.10).
Finally, we are able to prove the main theorem of the article.
A boot-strapping argument based on coercive L r estimates in the whole and half space for the principle part of the operators L 1−4 and L 7 , furnished by Theorem 5.8 in the former case and well-know estimates for the classical Laplace operator in the latter case, yields higher-order regularity.More specifically, after smoothing out the boundary in the L 1−4 part of equation (8.10), difference quotients of (u, p) can be estimated using Theorem 5.8, which implies additional regularity of (u, p).In turn, classical L r estimates for the Laplace operator in the 2D whole space yields bounds on difference quotients for η after smoothing out the interface in the L 7 part of equation (8.10).In both cases, we choose ε and thus ρ sufficiently small in order to absorb higher-order terms from the right-hand side.Bootstrapping this procedure, we conclude regularity of arbitrary order for both (u, p) and η, and thereby deduce that the solution is smooth up to the boundary.
R 2 [b] ⊂ R 2 ∖B 1 2 (0), we clearly have u = K(b).It follows that u 2,r ≤ c 8 b 2−1 r,r .In a completely similar manner, one shows that also p 1,r ≤ c 8 b 2−1 r,r .Thus the lemma follows for this particular choice of b ∈ S (R 2 ).Since any b

( 7 . 2 )
Moreover, standard regularity theory for the Stokes problem implies that both U and P are smooth in Ω 0 , and well-known decay estimates for the 3D exterior domain Stokes problem (see for example[7,  Theorem V.3.2])yield U = O( x −1 ), ∇U = O( x −2 ) and P = O( x −2 ) as x → ∞.(7.3)Additionally, both the Stokes operator and the boundary operator on the left-hand side of (7.1) are invariant with respect to rotations.Since the data on the right-hand side is clearly invariant with respect to rotations R ∈ SO(3) leaving e 3 invariant, the solution (U, P) retains this symmetry:∀R ∈ SO(3), Re 3 = e 3 ∶ R ⊺ U (Rx) = U (x), P(Rx) = P(x).(7.4)