INTRODUCTION

In this paper we study the dynamics of a two-dimensional incompressible fluid permeating a homogeneous porous medium. This physical phenomenon is modeled by the classical Darcy law [37] 𝜇𝑢(𝑥, 𝑦, 𝑡) = -∇ 𝑥,𝑦 𝑝(𝑥, 𝑦, 𝑡) -𝜌(0, 1), ∇ 𝑥,𝑦 ⋅ 𝑢(𝑥, 𝑦, 𝑡) = 0, (𝑥, 𝑦) ∈ Ω 𝑡 ⊂ ℝ 2 , 𝑡 ∈ ℝ + .

(1.1)

Here 𝑢 is the fluid velocity, 𝑝 is the fluid pressure, and the positive constants 𝜇 and 𝜌 are respectively the dynamic viscosity and fluid density. The permeability of the media and the gravitational acceleration have been normalized for the sake of notational simplicity. This problem is known as the one-phase Muskat problem [63] and, interestingly, is mathematically equivalent to the vertical Hele-Shaw problem driven by gravity. It is a fundamental parabolic free boundary problem in fluid mechanics. We should point out that the horizontal Hele-Shaw problem driven by fluid injection or suction has a different nature and has been studied extensively.

The fluid occupies a time-dependent domain Ω 𝑡 ⊂ ℝ 2 whose boundary Σ 𝑡 moves with the fluid

where 𝑛 is the outward pointing unit normal to Σ 𝑡 . Neglecting surface tension, the pressure is continuous across the free boundary

One is interested in the geometry and regularity of the free boundary Σ 𝑡 as time evolves. Regarding the geometry, there are two distinguished cases: graph and non-graph boundary. For both kinds of geometry, the existence and uniqueness of strong solution for a short time have been established even for much more general settings, including multi-phase, with rigid boundaries, with surface tension, nonconstant permeability. See for example [20,42,45], and more recent developments in [2-7, 21, 28-30, 47, 61, 64, 65].

In this paper we are concerned with long-term dynamics of the Muskat problem. Global existence and uniqueness of solutions have been obtained either when the initial free boundary is the graph of a small function (in certain function spaces) [21,26,30,44,66] or is close to a circle [20,27,49]. There are many recent developments for the former in the direction of low regularity solutions. These include "medium data" in the Wiener algebra [24,25,50,51], "medium data" in the Lipschitz norm [12,13], small data in critical Sobolev spaces allowing for arbitrarily large slopes [4,33,52] and even infinite slope [5]. We note that these norms are scaling invariant for the setting considered.

The geometry has an interesting interplay with the regularity. It is known that for the twophase problem (two fluids separated by a free interface), some initial graph interfaces with large slopes can turn over, passing from a stable regime to an unstable one [16]. Subsequently the solutions lose regularity in finite time [14]. On the other hand, there exist solutions that shift stability regime from stable to unstable and then return back to stable [32]. Particles on the free boundary cannot collide along a smooth curve (splat singularity) for both the two-phase [31] and one-phase problems [34]. However, the one-phase problem is more singular in the sense that it can develop splash singularity [15] from some non-graph initial boundary, while the two-phase problem cannot [53]. Here, splash singularity occurs when two particles collide at a single point on the free boundary while the boundary remains regular. The remaining less well-understood scenario is the one-phase problem starting from graph initial boundaries. Then, it is known that in stark contrast to its two-phase counterpart, the free boundary of the one-phase problem cannot turn over. Two fundamental questions for this scenario are: (1) Does there exist a unique global solution? (2) If yes, what is its long-term regularity?

In this paper we affirmatively answer the first question. More precisely, we prove that if the initial free boundary is the graph of a periodic Lipschitz function, then there exists a global Lipschitz solution in the strong 𝐿 ∞ 𝑡 𝐿 2 𝑥 sense, and hence almost everywhere. Moreover, it is the unique viscosity solution [36]. In order to establish this result, we first show that the graph free boundary of the one-phase Muskat problem obeys a nonlinear integro-differential equation. Here, for 𝑓, 𝑔 ∶ 𝕋 ≡ ℝ∕2𝜋ℤ → ℝ the Dirichlet-Neumann operator 𝐺(𝑓)𝑔 is given by 𝐺(𝑓)𝑔(𝑥) = 1 4𝜋 𝑝.𝑣. ∫ 𝕋 sin(𝑥-𝑥 ′ )+sinh(𝑓(𝑥)-𝑓(𝑥 ′ ))𝜕 𝑥 𝑓(𝑥) cosh(𝑓(𝑥)-𝑓(𝑥 ′ ))-cos(𝑥-𝑥 ′ ) 𝜃(𝑥 ′ )𝑑𝑥 ′ , (1.5) where 𝜃 ∶ 𝕋 → ℝ satisfies 1 2 𝜃(𝑥) + 1 2𝜋 𝑝.𝑣. ∫ 𝕋 sinh(𝑓(𝑥)-𝑓(𝑥 ′ ))-sin(𝑥-𝑥 ′ )𝜕 𝑥 𝑓(𝑥) cosh(𝑓(𝑥)-𝑓(𝑥 ′ ))-cos(𝑥-𝑥 ′ ) 𝜃(𝑥 ′ )𝑑𝑥 ′ = 𝜕 𝑥 𝑔(𝑥).

(1.6)

We shall prove using deep results from layer potential theory for 𝐶 1 and Lipschitz domains [46,68] that 𝐺(𝑓)𝑔 is well-defined in 𝐿 2 (𝕋) with a quantitative bound, provided that 𝑓 ∈ 𝑊 1,∞ (𝕋) and 𝑔 ∈ 𝐻 1 (𝕋). Our main result is stated as follows.

Theorem 1.2. For all 𝑓 0 ∈ 𝑊 1,∞ (𝕋), there exists

such that 𝑓| 𝑡=0 = 𝑓 0 , 𝑓 satisfies (1.4) in 𝐿 ∞ 𝑡 𝐿 2 𝑥 , and ‖𝑓(𝑡)‖ 𝑊 1,∞ (𝕋) ≤ ‖𝑓 0 ‖ 𝑊 1,∞ (𝕋) a.e. 𝑡 > 0.

Moreover, 𝑓 is the unique viscosity solution of (1.4) and is stable in 𝐿 ∞ (𝕋).

We note that the stability in the theorem above is in the sense of Corollary 6.4. See also Proposition 5.2.

Remark 1.3. We shall prove in Theorem 6.3 that viscosity solutions of (1.4) obey the comparison principle. Consequently, every modulus of continuity of 𝑓 0 is preserved by 𝑓(𝑡) for all 𝑡 > 0.

It turns out that sufficiently smooth solutions of the one-phase Muskat problem obey the comparison principle. See Proposition 5.1 and [3]. As a consequence, so long as the free boundary remains to be a graph, its slope is bounded by the initial slope. Thus starting from a Lipschitz graph, any sufficiently smooth solution must be a graph.

The comparison principle has been discovered for other free boundary problems including the Stefan problem [8,58,59] and the horizontal Hele-Shaw problem [11,19,58]. Motivated by the fact that these models may develop singularities in finite time, the notion of viscosity solutions, introduced by Crandall and Lions [36], has been employed to construct global solutions past singularities. In the aforementioned works, viscosity solutions for the unknown scalars are merely continuous and so are their zero level sets-the free boundary. Consequently, the dynamics of the free boundary is not satisfied in the classical sense (pointwise). This is also the case for variational weak solutions constructed in [43,57]. Higher regularity of weak (viscosity or variational) solutions to these problems has been studied in [8,22,23]. In particular, for the horizontal onephase Hele-Shaw problem, the authors in [22] proved that if the initial Lipschitz norm of the free boundary is small, then for a short time the unique viscosity solution satisfies the equations pointwise. We also mention that by using the methodology of convex integration, non-unique weak solutions have been constructed in the purely unstable [17,48,62] and partially unstable scenarios [18].

As far as the Muskat problem is concerned, to the best of our knowledge there has not been any global well-posedness result for initial data of arbitrary size, either for weak or strong solutions. We should point out that for the two-phase problem, the authors in [39] proved the existence (without uniqueness) of global weak solutions that are monotone in ℝ. In order to obtain viscosity solutions that satisfy the equation almost everywhere, one cannot appeal to the Perron method as in the aforementioned works. Instead, we construct solutions by the vanishing viscosity approach. More precisely, for small 𝜀 > 0 we consider the approximate equation 𝜕 𝑡 𝑓 𝜀 = -𝜅𝐺(𝑓 𝜀 )𝑓 𝜀 + 𝜀𝜕 2 𝑥 𝑓 𝜀 .

(1.7)

The added viscous term 𝜀𝜕 2 𝑥 𝑓 𝜀 retains the comparison principle for smooth solutions. The first difficulty is to establish global regularity for (1.7) with large data. This is done in Proposition 5.3 via the layer potential representation (1.5)-(1.6) of the Dirichlet-Neumann operator 𝐺(𝑓)𝑔 and a careful decomposition of its singular and non-singular parts. An important ingredient in its proof is the 𝐿 2 bound ‖𝐺(𝑓)𝑔‖ 𝐿 2 (𝕋) ≤ 𝐶(1 + ‖𝑓‖ Lip(𝕋) )

2 ‖𝜕 𝑥 𝑔‖ 𝐿 2 (𝕋) ,

where the slope ‖𝑓‖ Lip(𝕋) is controlled for all time. The proof of (1.8) relies on a quantitative bound for the inverse

𝐼 -𝐾 * [𝑓]

) -1 , where 𝐾 [𝑓] is the boundary double layer potential for the fluid domain {𝑦 < 𝑓(𝑥)}. This is established in Section 4. Under the 𝐶 2,𝛼 regularity condition for the boundary, a similar bound was obtained [28]. The quantitative bound is also crucial to the uniform 𝐿 2 bound for the solution 𝜃 𝜀 of (1.6) with 𝑓 = 𝑓 𝜀 . Finally, from the uniform Lipschitz bound for 𝑓 𝜀 and the uniform 𝐿 2 bound for 𝜃 𝜀 , we are able to pass to the limit 𝜀 → 0 in Equations (1.4)-(1.5) and (1.6) using a decomposition into small and large scales inspired by [24,25]. We thus obtain a global Lipschitz solution in the strong 𝐿 ∞ 𝑡 𝐿 2 𝑥 sense. An advantage of the viscosity regularization (1.7) is that one can show, in a direct manner, that the constructed solution 𝑓 is also a viscosity solution (see Definition 6.1). The proof of this makes use of the contraction estimate for the Dirichlet-Neumann operator associated to two different domains [65]. In order to obtain the uniqueness of viscosity solutions we prove that they obey the comparison principle. By the sup and inf-convolution technique, this reduces to proving the consistency of viscosity solutions. That is, if a viscosity solution is smooth (𝐶 1,1 ) at a point (𝑥 0 , 𝑡 0 ) then it satisfies Equation (1.4) classically at the same point. It is known that the harmonic extension 𝜙 of 𝑓 to the domain {𝑦 < 𝑓(𝑥)} is 𝐶 1 at the boundary point (𝑥 0 , 𝑓(𝑥 0 )) in all nontangential directions (see, for instance, [9,Lemma 11.17]), so that 𝐺(𝑓)𝑓 is classically well-defined at 𝑥 0 . However, without min-max formulas available at hand as in [10,19,67], the 𝐶 1 regularity seems insufficient. Here, we prove that in two dimensions, the harmonic extension 𝜙 is in fact 𝐶 1,𝛼 at (𝑥 0 , 𝑓(𝑥 0 )) with quantitative estimates (see Theorem 2.12 and Corollary 2.14). This pointwise 𝐶 1,𝛼 regularity for harmonic functions in Lipschitz domains is of independent interest. It allows us to reduce to the Dirichlet-Neumann operator for an interior disk tangent to {𝑦 = 𝑓(𝑥)} at (𝑥 0 , 𝑓(𝑥 0 )), so that an explicit integral formula can be used (see Proposition B.1) to conclude the consistency of viscosity solutions.

THE DIRICHLET-NEUMANN OPERATOR

Notation 2.1. For any set  ∈ ℝ 𝑑 , we denote

For 𝛼 ∈ (0, 1] we say that 𝑔 is 𝐶 1,𝛼 at 𝑥 0 ∈  if there exists 𝑣 ∈ ℝ 𝑑 and positive numbers 𝑀 and 𝛾 such that

If 𝑥 0 is an interior point of  then (2.1) implies that ∇𝑔(𝑥) exists and equals 𝑣.

Notation 2.2. For 𝑓 ∶ 𝕋 𝑑 → ℝ we denote

For 𝜙 ∶ Ω 𝑓 → ℝ and 𝑧 0 ∈ Σ we denote by

the limit of 𝜙 when Ω 𝑓 ∋ 𝑧 → 𝑧 0 in the direction of 𝑁.

Notation 2.3. For 𝑥 0 ∈ ℝ 𝑚 , 𝑚 ≥ 1 we denote the ball of radius 𝑟 centered at 𝑥 0 by 𝐵 𝑟 (𝑥 0 ). When 𝑥 0 = 0, we shall write 𝐵 𝑟 = 𝐵 𝑟 (0). (2.7)

Definition and global properties

The one-phase Muskat problem has a compact reformulation in terms of the Dirichlet-Neumann operator.

Proposition 2.5. If the fluid domain is given by Ω 𝑓 for some 𝑓 ∶ ℝ × (0, 𝑇) → ℝ such that

then 𝑓 satisfies

where 𝜅 = 𝜌∕𝜇.

This reformulation has been exploited in [3,47,64,65]. We recall its proof for completeness. which finishes the proof of (2.8). □

Proof

In order to define 𝐺(𝑓)𝑔 we first study the well-posedness of the elliptic problem (2.6). It was proved in Proposition 3.6 [65] that (2.6) has a unique variational solution when 𝑓 ∈ Lip(ℝ 𝑑 ) and

Of course, this is only nonstandard when the domain is unbounded. Adaption to horizontally periodic domains is straightforward. To fix notation for later purposes, we prove Proposition 2.6. Let 𝑓 ∈ Lip(𝕋 𝑑 ) and 𝑔 ∈ Ḣ 1 2 (𝕋 𝑑 ). Then there exists a unique variational solution 𝜙 ∈ Ḣ1 (Ω 𝑓 ) to (2.6), where

is endowed with the norm ‖𝑢‖ Ḣ1 (Ω 𝑓 ) = ‖∇𝑢‖ 𝐿 2 (Ω 𝑓 ) . Moreover, 𝜙 satisfies

.

(2.10)

We shall refer to the solution 𝜙 of (2.6) as the harmonic extension of 𝑔 to Ω 𝑓 .

for all 𝑓 𝑗 , 𝑔 ∈ 𝐻 𝑠 0 (𝕋 𝑑 ).

2.2

Pointwise properties in two dimensions

2.2.1

Pointwise 𝐶 1,𝛼 estimate for harmonic functions Suppose that 𝑈 is a Lipschitz domain in ℝ 2 . For (𝑥 0 , 𝑦 0 ) ∈ ℝ 2 and 𝑟 > 0, we denote 𝑈 𝑟 (𝑥 0 , 𝑦 0 ) = 𝐵 𝑟 (𝑥 0 , 𝑦 0 ) ∩ 𝑈 and 𝑈 𝑟 = 𝑈 𝑟 (0). We also define the half ball as

We assume that 0 ∈ 𝜕𝑈. Suppose that there exists some 𝑟 0 > 0 such that in a coordinate system, 𝜕𝑈 ∩ 𝐵 2𝑟 0 can be represented by a Lipschitz graph with Lipschitz constant 𝐿 > 0.

Let 𝑢 be a harmonic function in 𝑈, which vanishes on 𝜕𝑈.

Lemma 2.10. Under the conditions above, there exist 𝜀 0 = 𝜀 0 (𝐿) > 0 and 𝑀 1 = 𝑀 1 (𝐿) > 0 such that 𝑢 ∈ 𝐶

), 𝛾 = 𝜋∕(2𝜋 -𝛽) ∈ ( In view of (2.20), it is also easily seen that 𝑣 is a harmonic function in 𝑈 1∕2 (𝑥 0 , 𝑦 0 ). Moreover, 𝑣 ≥ 0 on 𝜕𝑈 ∩ 𝐵 1∕2 (𝑥 0 , 𝑦 0 ) and 𝑣 ≥ ‖𝑢‖ 𝐿 ∞ (𝑈 3∕2 ) on 𝑈 ∩ 𝜕𝐵 1∕2 (𝑥 0 , 𝑦 0 ). By the comparison principle, we have |𝑢| ≤ 𝑣 in 𝑈 1∕2 (𝑥 0 , 𝑦 0 ), which implies that

in 𝑈 1∕2 (𝑥 0 , 𝑦 0 ). This together with (2.19) and the interior regularity of harmonic functions complete the proof with 𝜀 0 = 𝛾 - .

□

We remark that alternatively the above Hölder estimate can be obtained by using the 𝑊 1,4+𝜀 estimate by Jerison-Kenig [56] and the Morrey embedding.

Assumption 2.11.

There exist constants 𝑀 0 , 𝑟 0 > 0 and function 𝜓 in (-𝑟 0 , 𝑟 0 ) such that in a coordinate system 𝜓(0) = 𝜓 ′ (0) = 0, 𝑈 𝑟 0 = {(𝑥, 𝑦) ∈ 𝐵 𝑟 0 ∶ 𝑦 > 𝜓(𝑥)}, and for any 𝑟 ∈ (0, 𝑟 0 ), {(𝑥, 𝑦) ∈ 𝐵 𝑟 ∶ 𝑦 > 𝑀 0 𝑟 2 } ⊂ 𝑈 𝑟 .

(2.21)

Moreover, 𝑈 satisfies the exterior ball condition at 0 with radius 1∕𝑀 0 .

Note that if 𝜕𝑈 is 𝐶 1,1 at 0 ∈ 𝜕𝑈, then Assumption 2.11 is satisfied.

Theorem 2.12. Under the conditions above, 𝑢 is 𝐶 1,𝛼 at 0, that is, for any (𝑥, 𝑦) ∈ 𝑈 such that √ 𝑥 2 + 𝑦 2 < 𝑟 0 , we have

where 𝑀 > 0 is a constant depending only on 𝑀 0 𝑟 0 , and 𝐿, and 𝛼 ∈ (0, 1) is a small constant depending only on 𝐿.

First, by scaling we may assume that 𝑟 0 = 1, so that the new 𝑀 0 becomes 𝑀 0 𝑟 0 , and 𝜀 0 and 𝐿 remain the same. Moreover, by dividing 𝑢 by a constant we may also assume that ‖𝑢‖ 𝐿 2 (𝑈 2 ) = 1.

Since 𝜕𝑈 satisfies an exterior ball condition, by using a barrier argument, we know that 𝑢 is Lipschitz at 0 and

(2.23)

By the Caccioppoli inequality, it is easily seen that for 𝑟 ∈ (0, 1),

where we used (2.23) Now we take a smooth domain 𝐸 such that 𝐵 + 2∕3 ⊂ 𝐸 ⊂ 𝐵 + 3∕4 . For any (𝑥 0 , 𝑦 0 ) ∈ ℝ 2 and 𝑟 > 0, denote 𝐸 𝑟 (𝑥 0 , 𝑦 0 ) = {(𝑥, 𝑦) ∈ ℝ 2 ∶ 𝑟 -1 (𝑥 -𝑥 0 , 𝑦 -𝑦 0 ) ∈ 𝐸}, Γ 𝑟 (𝑥 0 , 𝑦 0 ) = {(𝑥, 𝑦) ∈ 𝜕𝐸 𝑟 (𝑥 0 , 𝑦 0 ) ∶ 𝑦 = 𝑦 0 }.

Clearly, for 𝑟 sufficiently small, by (2.21) we have 𝐸 𝑟 (0, 𝑀 0 𝑟 2 ) ⊂ 𝑈 𝑟 .

Take a smooth function 𝜂 = 𝜂(𝑠) on ℝ such that 𝜂(𝑠) = 0 in (-∞, 1) and 𝜂(𝑠) = 1 in (2, ∞). Denote 𝜂 𝑟 (𝑠) = 𝜂(𝑠∕(𝑀 0 𝑟 2 )). A simple calculation reveals that 𝑢(𝑥, 𝑦)𝜂 𝑟 (𝑦) satisfies

and 𝑢(𝑥, 𝑦)𝜂 𝑟 (𝑦) = 0 on Γ 𝑟 (0, 𝑀 0 𝑟 2 ). Note that the right-hand side is supported in a narrow strip {(𝑥, 𝑦) ∈ 𝑈 𝑟 ∶ 𝑀 0 𝑟 2 < 𝑦 < 2𝑀 0 𝑟 2 }.

We decompose 𝑢𝜂 𝑟 in 𝐸 𝑟 (0, 𝑀 0 𝑟 2 ) as follows. Let 𝑤 = 𝑤 𝑟 be a weak solution to

with the zero Dirichlet boundary condition on 𝜕𝐸 𝑟 (0, 𝑀 0 𝑟 2 ). Then 𝑣 = 𝑣 𝑟 = 𝑢𝜂 𝑟 -𝑤 satisfies

and 𝑣 = 0 on Γ 𝑟 (0, 𝑀 0 𝑟 2 ). Estimates of 𝑤. Since 𝐸 𝑟 (0, 𝑀 0 𝑟 2 ) is smooth, by the 𝑊 1,𝑝 estimate, we know that for any 𝑝 < ∞, 𝑤 ∈ 𝑊 1,𝑝 (𝐸 𝑟 (0, 𝑀 0 𝑟 2 )). Notice that 𝜂 ′ 𝑟 (𝑦) ≤ 𝑀|𝑦 -𝜓(𝑥)| -1 . By using Hardy's inequality and a duality argument (see the proof of [41, Theorem 3.5]), we have ‖∇ 𝑥,𝑦 𝑤‖ 𝐿 𝑝 (𝐸 𝑟 (0,𝑀 0 𝑟 2 )) ≤ 𝑀‖∇ 𝑥,𝑦 𝑢‖ 𝐿 𝑝 (𝑈 𝑟 ∩{𝑦<2𝑀 0 𝑟 2 }) .

(2.26)

Now we fix 𝑝 = (2+𝑝 0 ) 2 and let 𝑞 > 1 be such that

(2.28)

Estimates of 𝑣. Since 𝐵 + 𝑟∕2 (0, 𝑀 0 𝑟 2 ) ⊂ 𝐸 𝑟 (0, 𝑀 0 𝑟 2 ), we know that 𝑣 is harmonic in 𝐵 + 𝑟∕2 (0, 𝑀 0 𝑟 2 ) and vanishes on the flat boundary. By the boundary estimate for harmonic functions,

, Thus by the mean value theorem and because 𝑣(0, 𝑀 0 𝑟 2 ) = 𝜕 𝑥 𝑣(0, 𝑀 0 𝑟 2 ) = 0, for any 𝜅 ∈ (0, 1∕4),

Then there exists 𝛼 ∈ (0, 1) depending only on ‖𝑓‖ 𝑊 1,∞ (𝕋) , and there exists 𝑀 > 0 depending only on ‖𝑓‖ 𝑊 1,∞ (𝕋) and 𝑀 0 𝛾, such that

Proof. Assume without loss of generality that (𝑥 0 , 𝑓(𝑥 0 )) = (0, 0). Let 𝑝 = 𝜙 -𝑦 so that 𝑝 is harmonic in Ω 𝑓 and vanishes on 𝜕Ω 𝑓 . Since 𝜕Ω 𝑓 is a Lipschitz graph globally, Theorem 2.12 implies that where 𝛼 and 𝑀 are as in the statement. Since 𝑓 ∈ 𝑊 1,∞ (𝕋), in (2.13) we can take 𝑔(𝑥, 𝑦) = 𝑓(𝑥)𝜒(𝑦 -𝑓(𝑥)), where the cutoff 𝜒 is supported on [-1, 0] and identically 1 on [-

, 0]. Then 𝜙 = 𝑢 + 𝑔, where 𝑢 ∈ Ḣ1 0 (Ω 𝑓 ) and ‖𝑔‖ 𝐻 1 (Ω 𝑓 ) ≤ 𝐶(‖𝑓‖ 𝑊 1,∞ ). In the bounded strip Ω 𝑓,𝛾 = {𝑓(𝑥) -2𝛾 < 𝑦 < 𝑓(𝑥)}, Poincaré's inequality yields

upon recalling (2.14). In addition, ‖𝑔‖ 𝐿 2 (Ω 𝑓 ∩𝐵 2𝛾 ) ≤ 𝐶 0 𝛾‖𝑓‖ 𝐿 ∞ (𝕋) , where 𝐶 0 is an absolute constant. It follows that ‖𝜙‖ 𝐿 2 (Ω 𝑓 ∩𝐵 2𝛾 ) ≤ 𝛾𝐶(‖𝑓‖ 𝑊 1,∞ ). Finally, inserting 𝑝 = 𝜙 -𝑦 into (2.36) we conclude the proof. □

Pointwise comparison principle for Dirichlet-Neumann

for all 𝑥 ∈ 𝕋 and 𝑓 1 (𝑥 0 ) = 𝑓 2 (𝑥 0 ), then 𝐺(𝑓 𝑗 )𝑓 𝑗 are classically well-defined at 𝑥 0 and

Proof. Let 𝜙 𝑗 be the harmonic extension of 𝑓 𝑗 to Ω 𝑗 ≡ Ω 𝑓 𝑗 as given by Proposition 2.6. Then 𝑝 𝑗 ∶= 𝜙 𝑗 -𝑦 are harmonic in Ω 𝑗 and identically vanish on 𝜕Ω 𝑗 = {(𝑥, 𝑦) ∶ 𝑦 = 𝑓 𝑗 (𝑥)}. Since 𝑥 0 is a 𝐶 1,1 boundary point, Theorem 2.12 implies that the normal derivatives 𝜕 𝑁 𝑗 𝑝 𝑗 , where 𝑁 𝑗 (𝑥) = (-𝜕 𝑥 𝑓 𝑗 (𝑥), 1), are classically well-defined at 𝑥 0 . Consequently, in view of the definition (2.5), we have

are classically well-defined at (𝑥 0 , 𝑓(𝑥 0 )). We first claim that

If 𝜙 𝑗 are smooth then so do 𝑝 𝑗 and (2.39) is a consequence of the maximum principle for harmonic functions. Here, following Proposition 4.3 [65] we prove (2.39) for variational solutions. Indeed, recall from (2.15) that

This implies ∫ Ω 𝑗 |∇ 𝑥,𝑦 𝑝 - 𝑗 | 2 = 0 and thus 𝑝 𝑗 ≥ 0 a.e. Ω 𝑗 . We obtain the claim (2.39). Since 𝑓 1 ≤ 𝑓 2 , we have Ω 1 ⊂ Ω 2 and 𝑝 ∶= 𝑝 1 -𝑝 2 satisfies 𝑝(𝑥, 𝑓 1 (𝑥)) = -𝑝 2 (𝑥, 𝑓 1 (𝑥)) ≤ 0 by virtue of (2.39). In addition, noticing 𝑝 = 𝜙 1 -𝜙 2 we find from (2.40) that

Then as for (2.42), upon choosing 𝜑 = max{𝑝 -sup 𝑥∈𝕋 𝑑 (𝑝(𝑥, 𝑓 1 (𝑥))), 0} we obtain

Therefore, 𝑝 attains the zero maximum at the boundary point (𝑥 0 , 𝑓 1 (𝑥 0 )). In addition, we have 𝑁 1 (𝑥 0 ) = 𝑁 2 (𝑥 0 ) since 𝜕 𝑥 𝑓 1 (𝑥 0 ) = 𝜕 𝑥 𝑓 2 (𝑥 0 ) at the maximum point 𝑥 0 of 𝑓 1 -𝑓 2 . We deduce that

Then (2.37) follows from this and (2.38). □ Remark 2.16. Proposition 2.15 is valid in all dimensions. Indeed, Theorem 2.12 was only used to conclude that the normal derivatives 𝜕 𝑁 𝑗 𝑝 𝑗 are classically well-defined at (𝑥 0 , 𝑓 𝑗 (𝑥 0 )). However, the same conclusion can be deduced from [9, Lemma 11.17] which is valid in all dimensions. On the other hand, the stronger 𝐶 1,𝛼 regularity obtained in Theorem 2.12 will be crucial in the proof of the comparison principle for viscosity solutions in Section 6.

CONTOUR DYNAMICS

In this section, we use potential theory to express the Dirichlet-Neumann operator, thereby obtaining a contour dynamics reformulation of the Muskat problem.

The Newtonian kernel for 𝕋 × ℝ is

We shall identify functions defined on 𝕋 to functions defined on the graph Σ = {(𝑥, 𝑓(𝑥)) ∶ 𝑥 ∈ 𝕋}.

Precisely, ℎ ∶ 𝕋 → ℝ is identified to h ∶ Σ → ℝ defined by h(𝑥, 𝑓(𝑥)) = ℎ(𝑥). The double layer potential associated to Ω 𝑓 for a function ℎ ∶ 𝕋 → ℝ is

where 𝑧 = (𝑥, 𝑦) ∈ (𝕋 × ℝ) ⧵ Σ and

A direct calculation using (hyperbolic) trigonometric identities gives

where in the above and throughout the paper, we assume that the range of arctan is [0, 𝜋), instead of the usual range (-𝜋∕2, 𝜋∕2). Hence,

3)

The single layer potential is

Note that the preceding single layer potential differs from the one in [46,68] by the arc length element

where 𝑧 = (𝑥, 𝑓(𝑥)) and 𝑧 ′ = (𝑥 ′ , 𝑓(𝑥 ′ )). The adjoint of 𝐾[𝑓] is

If 𝑓 ∈ Lip(𝕋) and ℎ ∈ 𝐿 𝑝 (𝕋) for some 𝑝 ∈ (1, ∞) then the nontangential limit

)

holds a.e. 𝑧 ∈ 𝕋 (see Theorem 1.10 [68]). Consequently, the unique solution 𝜙 of the Dirichlet problem (2.6) is given by

It was proved in Section 4 [68] that for Lipschitz domains, the mappings

where

) -1

(𝜕 𝑥 𝑔).

𝐼 + 𝐾[𝑓]

) -1 𝑔. Note that Θ ∈ 𝐻 1 (𝕋) in view of (3.10). Setting 𝜃 = 𝜕 𝑥 Θ ∈ 𝐿 2 (𝕋) we deduce from

) 𝜃(𝑥 ′ )𝑑𝑥 ′ , (𝑥,𝑦) ∈ 𝕋 × ℝ.

It follows from (3.5) and (3.6) that

) -1 𝜕 𝑥 𝑔.

For (𝑥, 𝑦) ∈ (𝕋 × ℝ) ⧵ Σ, we compute

Set 𝑇(𝑥) = (1, 𝜕 𝑥 𝑓(𝑥)) to be the tangent to Σ. For Ω 𝑓 ∋ (x, 𝑦) → 𝑁 (𝑥, 𝑓(𝑥)), by the definition (2.5) and the mean-value theorem, we have

Since 𝑓 ∈ Lip(𝕋) and 𝜃 ∈ 𝐿 2 (𝕋), Theorem 1.6 [68]

QUANTITATIVE BOUND FOR INVERSE OF DOUBLE LAYER POTENTIAL

According to (3.11),

Proof. In order to get this result, we will show that there exists a constant

We claim that (4.2) yields

This in turn implies the desired bound (4.1). Indeed, by the triangle inequality, we have

where in the last inequality we used the fact that 𝑀 ≥ 1. This finishes the proof of (4.3).

Recall that the single layer potential given by

With Ω = {(𝑥, 𝑦) ∶ 𝑥 ∈ 𝕋, 𝑦 < 𝑓(𝑥)}, 𝐻 satisfies (see [68])

where 𝜕 𝑁 is the normal derivative from inside the domain Ω and 𝜕 𝑐 𝑁 is to the normal derivative from inside the complementary set Ω 𝑐 = ℝ 2 ⧵ Ω. By Theorem 1.6 [68], the tangential derivative 𝜕 𝜏 𝐻 is continuous across the 𝜕Ω. Thus we shall not distinguish the tangential derivative of 𝐻 when approaching 𝜕Ω from either side. Since ∫ 𝕋 𝜂𝑑𝑥 = 0 we can write

and deduce that as |𝑦| → ∞, 𝐻(𝑥, 𝑦) is bounded and ∇𝐻 decay uniformly in 𝑥 ∈ 𝕋.

We consider a compactly supported vector field 𝑉(𝑥, 𝑦) given by

where the function 𝐺(𝑥, 𝑦) is given by

and Γ 𝛿 (𝑥, 𝑦) is a smooth approximation of the identity with compact support. The parameter 𝛿 > 0 is taken small enough so that

One could check that

where 𝐶 is independent of 𝑓. Let us recall the Rellich identity for 𝑉 and 𝐻:

where ∇𝑉 is a 2 × 2 matrix acting on ∇𝐻. Using this together with the normal and tangential decomposition for ∇𝐻 = (𝜕 𝑛 𝐻)𝑛 + (𝜕 𝜏 𝐻)𝜏 gives

Condition (4.9), Cauchy-Schwarz and Young's inequalities yield

Plugging in the bounds for 𝑉 in (4.10) we find that

Using the boundedness of 𝐻 and the decay of ∇𝐻 as |𝑦| → ∞, we can integrate by parts to obtain

By the fundamental theorem of calculus and Cauchy-Schwarz's inequality,

It follows from this and Cauchy-Schwarz's inequality for (4.13) that Plugging above inequality into (4.12) and using Young's inequality we get

Next, we employ the Rellich identity (4.11) for 𝑉 and 𝐻 in Ω 𝑐 = ℝ 2 ⧵ Ω to get

As before, using Young's inequality and (4.14) for Ω 𝑐 we obtain

.

The obvious bound |𝜕 𝜏 𝐻| 2 ≤ |∇𝐻| 2 combined with Young's inequality yields

The above argument when applied to Ω gives

Combining (4.15) and (4.16) allows us to relate both normal derivatives due to the continuity of the tangential derivative:

).

Thus we obtain (4.4) with 𝑀 = 2(12𝐶𝜋) 2 (1 + ‖𝑓‖ Lip(𝕋) ) 5∕2 . The opposite inequality

follows by the same argument, completing the proof. □ Remark 4.2. The invertibility of (

𝐼 -𝐾 * ) on 𝐿 2 0 was proved in [68] using the inequalities

The inequalities in (4.2) dispense with the second term on the right-hand of (4.18) and hence allowed us to obtain the quantitative bound (4.1).

As an interesting corollary of the proof of Proposition 4.1, we deduce a quantitative bound for the tangential derivative of the single layer potential in terms of the double layer potential. This in turn yields the continuity of 𝐺(𝑔) ∶ Ḣ1 → 𝐿 2 assuming only that the boundary is Lipschitz.

Corollary 4.3. Let 𝑓 ∈ Lip(𝕋).

There exists an absolute constant 𝐶 > 0 such that for all 𝜂 ∈ 𝐿 2 0 (𝕋), we have

Consequently,

for all 𝑔 ∈ Ḣ1 (𝕋).

Proof. Denote 𝐻(𝑥, 𝑦) = 𝑆[𝑓](𝜂)(𝑥, 𝑦) the single layer potential as given by (4.5). Now (4.17

where we recall that 𝜏 and 𝑛 are respectively the unit tangent and normal to the surface {𝑦 = 𝑓(𝑥)}.

On the other hand, it follows from the second equation in (4.6) that

VISCOSITY REGULARIZATION

We consider the viscosity regularization of the Muskat problem:

𝜕 𝑡 𝑓 𝜀 = -𝜅𝐺(𝑓 𝜀 )𝑓 𝜀 + 𝜀𝜕 2 𝑥 𝑓 𝜀 , 𝜀 ≥ 0.

(5.1)

Comparison principle for smooth solutions

Proposition 5.1 (Comparison principle). Assume that 𝑓 𝑗 ∈ 𝐶([0, 𝑇]; 𝐶 2 (𝕋)), 𝑗 = 1, 2 are smooth solutions of (5.1) with 𝜀 ≥ 0.

Proof. Assume by contradiction that 𝑀 0 = max 𝕋×[0,𝑇] (𝑓 1 -𝑓 2 ) = (𝑓 1 -𝑓 2 )(𝑥 0 , 𝑡 0 ) > 0. We have 𝑡 0 > 0 since 𝑓 1 (⋅, 0) ≤ 𝑓 2 (⋅, 0). Choose 𝜂 > 0 sufficiently small so that (𝑓 1 -𝑓 2 -𝜂𝑡)(𝑥 0 , 𝑡 0 ) = 𝑀 0 -𝜂𝑡 0 > 0, and thus 𝑀 * = max 𝕋×[0,𝑇] (𝑓 Owing to the translation invariance of (5.1), comparison principle implies maximum principles for the amplitude and the slope of solutions. (5.3)

In particular, for 𝛾(𝑧) = 𝑧‖𝑓(⋅, 0)‖ Lip(𝕋) we have

(5.4)

If 𝑓 1 and 𝑓 2 are smooth solutions in 𝐶([0, 𝑇]; 𝐶 2 (𝕋)) of ( 5.1) with 𝜀 ≥ 0, then We turn to prove (5.5). We have 𝑓 1 (⋅, 0) ≤ 𝑓 2 (⋅, 0) + ‖𝑓 1 (⋅, 0) -𝑓 2 (⋅, 0)‖ 𝐿 ∞ (𝕋) and 𝑓 2 (⋅, 𝑡) + ‖𝑓 1 (⋅, 0) -𝑓 2 (⋅, 0)‖ 𝐿 ∞ (𝕋) is a smooth solution of (5.1). By the comparison principle in Proposition 5.1,

Thus we can interchange the role of 𝑓 1 and 𝑓 2 to obtain (5.5). □

Global well-posedness of smooth solutions

By virtue of Proposition 3.1, in the class of Lipschitz solutions, (5.1) is equivalent to the system

sinh(𝑓 𝜀 (𝑥) -𝑓 𝜀 (𝑥 -𝑥 ′ )) -sin 𝑥 ′ 𝜕 𝑥 𝑓 𝜀 (𝑥) cosh(𝑓 𝜀 (𝑥) -𝑓 𝜀 (𝑥 -𝑥 ′ )) -cos 𝑥 ′ 𝜃 𝜀 (𝑥 -𝑥 ′ )𝑑𝑥 ′ = -𝜅𝜕 𝑥 𝑓 𝜀 (𝑥).

(5.6)

We prove that (5.6) is globally well-posed for any smooth initial data in Sobolev spaces.

Proposition 5.3. Let 𝜀 > 0 and 𝑠 ≥ 2. For each 𝑓 𝜀 0 ∈ 𝐻 𝑠 (𝕋), there exists a unique solution 𝑓 𝜀 ∈ 𝐶([0, 𝑇]; 𝐻 𝑠 (𝕋)) ∩ 𝐿 2 ([0, 𝑇]; 𝐻 𝑠+1 (𝕋)) a solution to (5.1) for any time 𝑇 > 0 with initial data 𝑓 𝜀 0 . In particular, 𝑓 𝜀 obeys the maximum principles

Proof. We first note that the dissipative term 𝜀𝜕 2 𝑥 𝑓 𝜀 makes the regularized Muskat problem (5.1) semilinear. Then by the contraction mapping method, it was proved in Section 4.1.5 [65] that for any 𝑠 >

where 𝑇 𝜀 is the maximal time. Moreover, if 𝑇 𝜀 < ∞ then lim sup

Therefore, to conclude that 𝑓 𝜀 is global in time, it suffices to prove that the 𝐻 𝑠 norm of 𝑓 𝜀 is bounded on [0, 𝑇 𝜀 ). To this end, we shall prove that for any 𝑇 < 𝑇 𝜀 , the 𝐶([0, 𝑇]; 𝐻 𝑠 ) norm of 𝑓 𝜀 can be controlled by ‖𝑓 𝜀 0 ‖ 𝐻 𝑠 and the 𝐿 ∞ ([0, 𝑇]; 𝑊 1,∞ ) norm of 𝑓 𝜀 , where the latter is uniformly bounded in time owing to the maximum principles (5.4) and (5.5).

For notational simplicity we shall write 𝑓 𝜀 = 𝑓. In what follows,  ∶ ℝ + → ℝ + denotes continuous nondecreasing functions that may change from line to line but only depend on 𝜅.

𝐿 2 estimate. Using (𝐺(𝑓)𝑓, 𝑓) ≥ 0 and (𝜕 2 𝑥 𝑓, 𝑓) ≤ 0, we find at once that ‖𝑓(𝑡)‖ 𝐿 2 ≤ ‖𝑓 0 ‖ 𝐿 2 ∀𝑡 > 0.

(5.8)

Ḣ1 estimate. Multiplying (5.1) by -𝜕 2 𝑥 𝑓, then integrating by parts and using the estimate (4.20), we obtain

By a Grönwall argument we find that

(5.9)

Ḣ2 estimate. We differentiate (5.1) in 𝑥, then multiply by -𝜕 We claim that for any 𝜈 ∈ (0, 1),

(5.10) Taking (5.10) for granted, we choose 𝜈 = 𝜀 4𝜅 so that the above estimates yield

(5.11)

Then applying Grönwall's lemma and (5.9) to bound

we arrive at the Ḣ2 estimate

(5.12)

The proof of (5.10) requires a careful decomposition of the singular integral representation (3.12)-(3.14) of 𝐺(𝑓)𝑓 and is postponed to Lemma 5.4 below. Higher regularity. In this step we prove that 𝐻 𝑠 regularity with 𝑠 > 2 can be propagated. To this end, we use Sobolev estimates for the Dirichlet-Neumann operator with the control norm in

+ 𝑥 provided by (5.12). Indeed, fixing 𝑠 0 ∈ (

, 2) and applying (2.17) with 𝑔 = 𝑓 and 𝜎 = 𝑠, we have ‖𝐺(𝑓)𝑓‖ 𝐻 𝑠-1 ≤ (‖𝑓‖ 𝐻 𝑠 0 )‖𝑓‖ 𝐻 𝑠 .

(5.13)

Then an 𝐻 𝑠 energy estimate for (5.1) gives

where in the last inequality we used Young's inequality and the maximum principle for the 𝐿 2 norm of 𝑓. By Grönwall's lemma, we obtain ‖𝑓(𝑡)‖ Then using the Ḣ2 estimate (5.12) together with (5.8), (5.5) and (5.4), we arrive at ‖𝑓‖ 𝐿 ∞ ([0,𝑇];𝐻 𝑠 ) ≤ 𝐶(‖𝑓 0 ‖ 𝐻 𝑠 , 𝑇, 𝜀) ∀𝑇 < 𝑇 𝜀 .

(5.16)

This concludes the proof of global regularity for (5.1). □

We now prove the claim (5.10).

Lemma 5.4.

There exists  ∶ ℝ + → ℝ + nondecreasing such that for any 𝜈 ∈ (0, 1),

Here and in what follows, 𝐶 denotes a universal constant that may change from line to line. Then the integral kernel in 𝐺 1,1 is bounded by 𝐶‖𝜕 𝑥 𝑓‖ 𝐿 ∞ , giving

Using Young's inequality for convolutions, we bound

.

Next in view of the inequalities

we deduce that To control 𝐺 1,4 , we use Taylor's expansion to have

which in conjunction with Gagliardo-Nirenberg's inequality implies

The integral in 𝐺 1,5 has a Calderón-Zygmund type kernel, whence

Gathering the above estimates yields

)

where we bounded the ‖𝜕 𝑥 𝑓‖ 𝐿 2 term in the estimates for 𝐺 1,1 and 𝐺 1,4 by

For any 𝜈 ∈ (0, 1), we choose

with  sufficiently large so that The integral kernel in 𝐺 2,1 is bounded by an absolute constant so that

The identity

allows us to get

where we applied Cauchy-Schwarz's inequality for the integral in 𝑥 ′ . Using the inequality

we can bound 𝐺 2,3 as

Finally, 𝐺 2,4 obeys the same bound as 𝐺 1,5 :

Gathering all the above estimates leads to with 𝐴(𝑥, 𝑥 ′ ) =

) , and

The integral kernel in 𝐺 3,1 is smooth enough so that

By brutal force we have

.

For the next two terms, Gagliardo-Nirenberg's inequality gives

The term 𝐺 3,5 is bounded by

As for 𝐺 1,4 , we have

Finally, the kernels in 𝐺 3,7 are of Calderón-Zygmund type in such a way that

Gathering the above estimates we find that

)

.

For any 𝜈 ∈ (0, 1), using Young's inequality for the first term and choosing

,

To proceed, we differentiate the second equation in (5.6) with respect to 𝑥 to have where Since the kernel integral in 𝐺 1 4,1 is bounded by 𝐶‖𝜕 𝑥 𝑓‖ 𝐿 ∞ , we have

Similarly, the kernel integral in 𝐺 2 4,2 is bounded by an absolute constant, yielding

Using convolution properties, we bound

.

Similarly, 𝐺 4 4,1 can be controlled as

We have proved that

)

.

For any 𝜈 ∈ (0, 1), choosing

, we obtain

) .

(5.28)

and

Since the integral kernel in 𝐺 1 4,2 is bounded by  (‖𝑓‖ 𝑊 1,∞ ), it follows that

Next we estimate 𝐺 3 4,2 and 𝐺 4 4,2 as

and

We thus obtain

For any 𝜈 ∈ (0, 1), choosing

gives

Then plugging (5.28) and (5.29) into (5.27) we obtain

(5.30)

Finally, combining (5.21), (5.23), (5.24) and (5.30) we conclude that

which finishes the proof. □

VISCOSITY SOLUTIONS AND COMPARISON PRINCIPLE

By virtue of the comparison principle for the Dirichlet-Neumann operator given by Proposition 2.15 and the fact that 𝐺(𝑓 + 𝐶)(𝑓 + 𝐶) = 𝐺(𝑓)𝑓 for any constant 𝐶, we have the following natural definition of viscosity solutions for the Muskat problem (2.8).

Definition 6.1 (Viscosity solutions).

A function 𝑓 ∶ 𝕋 × [0, 𝑇] is called a viscosity subsolution (resp. supersolution) of (2.8) on (0, 𝑇) provided that (i) 𝑓 is upper semicontinuous (resp. lower semicontinuous) on 𝕋 × [0, 𝑇], and (ii) for every 𝜓 ∶ 𝕋 × (0, 𝑇) → ℝ with 𝜕 𝑡 𝜓 ∈ 𝐶(𝕋 × (0, 𝑇)) and 𝜓 ∈ 𝐶((0, 𝑇); 𝐶 1,1 (𝕋)),

if 𝑓 -𝜓 attains a global maximum (resp. minimum) over 𝕋 × [𝑡 0 -𝑟, 𝑡 0 ] at (𝑥 0 , 𝑡 0 ) ∈ 𝕋 × (0, 𝑇) for some 𝑟 > 0, then 𝜕 𝑡 𝜓(𝑥 0 , 𝑡 0 ) ≤ -𝜅(𝐺(𝜓)𝜓)(𝑥 0 , 𝑡 0 ) (resp. ≥). (

A viscosity solution is both a viscosity subsolution and viscosity supersolution.

The next proposition shows that if a viscosity solution is regular (𝐶 1,1 ) at a point (𝑥 0 , 𝑡 0 ) then it satisfies the equation classically at the same point. Proposition 6.2 (Consistency). Let 𝑓 be a viscosity subsolution of (2.8) on (0, 𝑇). Assume that 𝑓 ∈ 𝑊 1,∞ (𝕋 × (0, 𝑇)) and 𝑓 is 𝐶 1,1 at (𝑥 0 , 𝑡 0 ) ∈ 𝕋 × (0, 𝑇). Then 𝐺(𝑓)𝑓 is classically well-defined at (𝑥 0 , 𝑡 0 ) and 𝜕 𝑡 𝑓(𝑥 0 , 𝑡 0 ) ≤ -𝜅(𝐺(𝑓)𝑓)(𝑥 0 , 𝑡 0 ).

(6.

2)

The corresponding statement for viscosity supersolutions holds true.

Proof. According to Proposition 2.15, 𝐺(𝑓)𝑓 is classically well-defined at (𝑥 0 , 𝑡 0 ). We shall write 𝑋 = (𝑥, 𝑡) and 𝑋 0 = (𝑥 0 , 𝑡 0 ). Since 𝑓 is 𝐶 1,1 at 𝑋 0 , for some 𝑟 0 ∈ (0, we have

Let 𝜓 be the tangent parabola with double opening:

We can find a family of functions 𝜓 𝑟 ∈ 𝐶 ∞ (𝕋 × ℝ), 𝑟 ∈ (0, 𝑟 0 ) satisfying

The construction of 𝜓 𝑟 is postponed to Appendix C. The properties in (6.5) imply that 𝜓 𝑟 (𝑡 0 )'s are uniformly 𝐶 1,1 at 𝑥 0 , that is, there exists 𝑀 0 > 0 independent of 𝑟 ∈ (0, 𝑟 0 ) such that

Since each 𝜓 𝑟 is a valid test function for the viscosity subsolution 𝑓, we have

At the maximum point 𝑋 0 of 𝑓 -𝜓 𝑟 , we have 𝜕 𝑡 𝑓 = 𝜕 𝑡 𝜓 𝑟 . Thus (6.2) will be a consequence of (6.7) and lim 𝑟→0 (𝐺(𝜓 𝑟 )𝜓 𝑟 )(𝑋 0 ) = (𝐺(𝑓)𝑓)(𝑋 0 ). (

We shall skip the time variable in the remainder of the proof because 𝑡 = 𝑡 0 is fixed in (6.8). The proof of (6.8) proceeds in two steps.

Step 1. Let 𝜙 𝑟 (resp. 𝜙) be the harmonic extension of 𝜓 𝑟 to Ω 𝑟 ≡ Ω 𝜓 𝑟 (resp. Ω 𝑓 ). We claim that lim 𝑟→0 𝜙 𝑟 (𝑥, 𝑦) = 𝜙(𝑥, 𝑦) ∀(𝑥, 𝑦) ∈ Ω 𝑓 , (

where we note that Ω 𝑓 ⊂ Ω 𝑟 for all 𝑟. From (2.41), (2.42) and the uniform boundedness of ‖𝜓 𝑟 ‖ 𝑊 1,∞ (𝕋) , we have

Fix a sequence of bounded subsets Ω 𝑗 ⋐ Ω 𝑓 satisfying Ω 𝑗 ⊂ Ω 𝑗+1 and ∪ 𝑗≥1 Ω 𝑗 = Ω 𝑓 . Let 𝑟 𝑛 be an arbitrary sequence of positive numbers converging to 0. We relabel 𝜙 𝑛 = 𝜙 𝑟 𝑛 and Ω 𝑛 ≡ Ω 𝑟 𝑛 . The uniform bound (6.10) implies that every subsequence of 𝜙 𝑛 has a subsequence converging weakly-* in 𝐿 ∞ (Ω 𝑓 ). For notational simplicity, we write 𝜙 𝑛 → 𝜙 ∞ weakly-* in 𝐿 ∞ (Ω 𝑓 ). We have the gradient estimate for the harmonic functions 𝜙 𝑗

By virtue of the Arzelà-Ascoli theorem, for every 𝑗 there exists a subsequence of 𝜙 𝑛 converging to 𝜙 ∞ in 𝐶(Ω 𝑗 ). Using the Cantor diagonal argument, we can find a subsequence 𝜙 𝑛 𝑘 → 𝜙 ∞ in 𝐶(Ω 𝑗 ) for all 𝑗. Thus (6.9) will follow once we can prove that

This can be achieved by using the variational characterization of 𝜙 𝑛 𝑘 and 𝜙. Indeed, we recall from (2.13) and (2.14) that 𝜙 𝑛 = 𝑢 𝑛 + 𝜓 𝑛 , where 𝜓 𝑛 ∈ Ḣ1 (Ω 𝑛 ) with 𝜓 𝑛 (𝑥, 𝜓 𝑛 (𝑥)) = 𝜓 𝑛 (𝑥) (in the trace sense) and 𝑢 𝑛 ∈ Ḣ1 0 (Ω 𝑛 ) satisfies

In particular, the uniform bounds

hold. Since Ω 𝑓 ⊂ Ω 𝑛 , upon passing to a subsequence of 𝑛 𝑘 , we have

As for the convergence of 𝜓 𝑛 , we recall from Theorem A.2 that we can choose 𝜓 𝑛 (𝑥, 𝑦) = 𝜓 𝑛 ♯ (𝑥, 𝑦 -𝜓 𝑛 (𝑥)), where 𝜓 𝑛 ♯ (𝑥, 𝑧) = 𝑒 𝑧|𝐷 𝑥 | 𝜓 𝑛 (𝑥) for (𝑥, 𝑧) ∈ 𝕋 × ℝ -. Since 𝜓 𝑛 → 𝑓 in 𝐶 𝛼 (𝕋) for all 𝛼 ∈ (0, 1), 𝜓 𝑛 ♯ → 𝑓 ♯ ∶= 𝑒 𝑧|𝐷 𝑥 | 𝑓(𝑥) in 𝐿 ∞ (𝕋 × ℝ -). Consequently, 𝜓 𝑛 (𝑥, 𝑦) → 𝑓(𝑥, 𝑦) ∶= 𝑓 ♯ (𝑥, 𝑦 -𝑓(𝑥)) in 𝐿 ∞ (Ω 𝑓 ). Then for all 𝜃 ∈ 𝐶 ∞ 𝑐 (Ω 𝑓 ), using integration by parts and the dominated convergence theorem gives

where 𝑁 = (-𝜕 𝑥 𝑓, 1). As 𝐶 ∞ 𝑐 (Ω 𝑓 ) is dense in Ḣ1 (Ω 𝑓 ), it follows that

Recalling that 𝜙 𝑛 𝑘 = 𝑢 𝑛 𝑘 + 𝜓 𝑛 𝑘 , we deduce from (6.13) and (6.14) that

On the other hand, since Ḣ1 0 (Ω 𝑓 ) ⊂ Ḣ1 0 (Ω 𝑛 ), (6.12) implies

By virtue of (6.13) and (6.14), letting 𝑛 = 𝑛 𝑘 → ∞ in the preceding equality we obtain

This together with (6.15) proves that 𝜙 ∞ is the harmonic extension of 𝑓 to Ω 𝑓 in the variational sense. By the uniqueness of variational solutions we conclude that 𝜙 ∞ = 𝜙, proving (6.11).

Step 2. Recall that 𝑓 ∈ 𝑊 1,∞ (𝕋) is 𝐶 1,1 at 𝑥 0 . In addition, 𝜓 𝑟 ∈ 𝐶 ∞ (𝕋) is tangent to 𝑓 from above at 𝑥 0 , so that 𝑁 𝑟 (𝑥 0 ) = 𝑁(𝑥 0 ) and there exists an interior disk 𝐵 tangent to Ω 𝑓 at 𝑧 0 = (𝑥 0 , 𝑓(𝑥 0 )). Clearly, 𝐵 is also tangent to all Ω 𝑟 . Assume without loss of generality that 𝐵 = 𝐵 1 (0). In addition, in view of the uniform 𝐶 1,1 property (6.6), Corollary 2.14 implies that 𝜙 𝑟 's are uniformly 𝐶 1,𝛼 at (𝑥 0 , 𝑓(𝑥 0 )). In particular, there exists 𝑀 > 0 such that for all 𝑟 < 𝑟 0 , the trace 𝑔 𝑟 = 𝜙 𝑟 | 𝐵 satisfies The convergence (6.9) implies that 𝑔 𝑟 (𝑥) → 𝑔(𝑥) for all 𝑥 ∈ 𝐵. The uniform bounds (6.10) and (6.16) then allow us to apply the dominated convergence theorem to conclude the claim (6.8). □ Using Proposition 6.2, we prove that viscosity solutions obey the comparison principle. An immediate consequence is the uniqueness of viscosity solutions. Proof. For classical sub and supersolutions, one can follow the proof of Proposition 5.1. Here the source difficulty comes from the low regularity of 𝑓 and 𝑔. We employ the regularization using sup-and inf-convolutions. For small 𝛿 > 0, let 𝑓 𝛿 and 𝑔 𝜀 be respectively the sup-convolution and inf-convolution of 𝑓 and 𝑔: 𝑓 𝛿 The boundedness of 𝑓 and 𝑔 is only used in (i).

Next we prove that each 𝑓 𝛿 (resp. 𝑔 𝛿 ) is a viscosity subsolution (resp. supersolution) of (2.8). This follows from the translation invariance of (2.8). Indeed, let 𝜓 be a test function as in Definition 6.1 such that 𝑓 𝛿 -𝜓 has a global maximum over 𝕋 × [𝑡 0 -𝑟, 𝑡 0 ] at (𝑥 0 , 𝑡 0 ) for some 𝑟 > 0. Let (𝑦 0 , 𝑠 0 ) ∈ 𝕋 × [0, 𝑇] be the point where the maximum in the definition of 𝑓 𝛿 (𝑥 0 , 𝑡 0 ) is attained. We have that 𝑓 𝛿 (𝑥 0 , 𝑡 0 ) = 𝑓(𝑦 0 , 𝑠 0 ) -1 2𝛿

𝑓 𝛿 (𝑥, 𝑡) ≥ 𝑓(𝑥 -𝑥 0 + 𝑦 0 , 𝑡 -𝑡 0 + 𝑠 0 ) -1 2𝛿

Therefore, the function ψ(𝑥, 𝑡) = 𝜓(𝑥 + 𝑥 0 -𝑦 0 , 𝑡 + 𝑡 0 -𝑠 0 ) + 1 2𝛿

is a valid test function for which 𝑓 -ψ has a global maximum over 𝕋 × [𝑠 0 -𝑟, 𝑠 0 ] at (𝑦 0 , 𝑠 0 ). Therefore, 𝜕 𝑡 ψ(𝑦 0 , 𝑠 0 ) ≤ -𝜅(𝐺( ψ) ψ)(𝑦 0 , 𝑠 0 ), and hence 𝜕 𝑡 𝜓(𝑥 0 , 𝑡 0 ) ≤ -𝜅(𝐺(𝜓)𝜓)(𝑥 0 , 𝑡 0 ). Thus each 𝑓 𝛿 is viscosity subsolution of (2.8). By an analogous argument, each 𝑔 𝛿 is viscosity supersolution of (2.8). We claim that for very 𝜀 > 0, there exists 𝛿(𝜀) > 0 such that for all 𝛿 ≤ 𝛿(𝜀), 𝑓 𝛿 (𝑥, 𝑡) ≤ 𝑔 𝛿 (𝑥, 𝑡) + 𝜀 ∀(𝑥, 𝑡) ∈ 𝕋 × [0, 𝑇]. (6.18) Taking this for granted, the half-relaxed limits in (6.17) yield 𝑓(𝑥, 𝑡) ≤ 𝑔(𝑥, 𝑡) for all (𝑥, 𝑡) ∈ 𝕋 × [0, 𝑇], proving the comparison principle. To prove (6.18), we assume by contradiction that for some 𝜀 0 > 0, there exists a sequence 𝛿 𝑛 → 0 such that for all 𝑛, 𝑀 𝑛 ∶= max

PROOF OF Theorem 1.2

Let 𝑓 0 be an arbitrary initial surface in 𝑊 1,∞ (𝕋). We construct a global solution to the contour dynamics formulation (1.4)-(1.6) by regularizing initial data and the Muskat problem. Precisely, for 𝜀 ∈ (0, 1), let 𝑓 𝜀 0 = Γ 𝜀 * 𝑓 0 (𝑥) be the mollification of 𝑓 0 where Γ 𝜀 is an approximation of the identity. Since 𝑓 𝜀 0 ∈ 𝐻 𝑠 (𝕋) for any 𝑠 > 0, Proposition 5.3 yields the unique global solution 𝑓 𝜀 to the regularized Muskat Equation (5.1) with initial data 𝑓 𝜀 0 . Moreover, 𝑓 𝜀 is smooth and in view of (5.7),

Setting 𝜃 𝜀 = ( Let 𝜀 𝑛 → 0 and relabel 𝑓 𝑛 = 𝑓 𝜀 𝑛 and 𝜃 𝑛 = 𝜃 𝜀 𝑛 . From the uniform bounds (7.1) and (7.2), we obtain the weak* convergences (upon extracting subsequences)

In particular, (7.1) implies that

We now prove that 𝑓 and 𝜃 satisfy (1.4)-(1.6) for all 𝑡 > 0. Fix an arbitrary time 𝑇 > 0.

Step It then follows from (5.1) that 𝜕 𝑡 𝑓 𝜀 is uniformly bounded in 𝐿 ∞ ([0, 𝑇]; 𝐻 -1 (𝕋)). We have the continuous embedding 𝑊 1,∞ (𝕋) ⊂ 𝐶(𝕋) ⊂ 𝐻 -1 (𝕋) where the first one is compact by the Azelà-Ascoli theorem. Thus the Aubin-Lions Lemma [60] implies that 𝑓 𝑛 → 𝑓 in 𝐶(𝕋 × [0, 𝑇]). (7.5)

In particular, lim 𝑡→0 + 𝑓(⋅, 𝑡) = 𝑓 0 (⋅).

Step 2: 𝑓 is the unique viscosity solution. It suffices to prove that 𝑓 is a viscosity subsolution since its uniqueness then follows at once from Theorem 6.3. Indeed, assume that for 𝜓 ∶ 𝕋 × (0, 𝑇) → ℝ with 𝜕 𝑡 𝜓 ∈ 𝐶(𝕋 × (0, 𝑇)) and 𝜓 ∈ 𝐶((0, 𝑇); 𝐶 1,1 (𝕋)), 𝑓 -𝜓 attains a global maximum over 𝕋 × [𝑡 0 -𝑟, 𝑡 0 ] at (𝑥 0 , 𝑡 0 ) ∈ 𝕋 × (0, 𝑇) for some 𝑟 > 0. Setting 𝜓(𝑥, 𝑡) = 𝜓(𝑥, 𝑡) + |𝑡 -𝑡 0 | 2 we have that 𝑓 -𝜓 attains a strict global maximum over 𝕋 × [𝑡 0 -𝑟, 𝑡 0 at (𝑥 0 , 𝑡 0 ). By the uniform convergence (7.5) and 𝜕 𝑡 𝜓 ∈ 𝐶(𝕋 × (0, 𝑇)), letting 𝑛 → ∞ in the preceding inequality yields

where we used the fact that 𝐺(𝜓)𝜓 = 𝐺(𝜓)𝜓. We claim that lim 𝑛→∞ (𝐺(𝜓)𝜓)(𝑥 𝑛 , 𝑡 𝑛 ) = (𝐺(𝜓)𝜓)(𝑥 0 , 𝑡 0 ).

Taking this for granted and noticing that 𝜕 𝑡 𝜓(𝑥 0 , 𝑡 0 ) = 𝜕 𝑡 𝜓(𝑥 0 , 𝑡 0 ) we obtain

whence 𝑓 is a viscosity subsolution. Analogously, 𝑓 is also a viscosity supersolution. To prove the claim (7. For 𝜓(𝑡 0 ) ∈ 𝐶 1,1 (𝕋), 𝐺(𝜓(𝑡 0 ))𝜓(𝑡 0 ) ∈ 𝐶 𝛼 (𝕋) for all 𝛼 ∈ (0, 1) and thus 𝑁 𝑛 2 → 0. As for 𝑁 𝑛 1 we apply the continuity and contraction estimates in Propositions 2.8 and 2.9 to have

since 𝜓 ∈ 𝐶((0, 𝑇); 𝐶 1,1 (𝕋)). This finishes the proof of (7.8). We shall prove in the following steps that 𝑓 and 𝜃 solve (1.4)-(1.6) by passing to the limit in the equivalent formulation (5.6) of (5.1).

Step 3: Equation (1.6)

where

The first and the last integrals in (7.9) converge to their corresponding limit in view of the weak* convergences in (7.3). We show now that as 𝑛 → ∞,

We split 𝐼 𝑛 1 -𝐼 1 = 𝐷 𝑛 1 + 𝐷 𝑛 2 , where

)

From the uniform bound (7.

)) 𝑑𝑥𝑑𝑥 ′ 𝑑𝑡.

Proof. By the density of 𝐶 ∞ 𝑐 (Ω 𝑓 ) in Ḣ1 (Ω 𝑓 ), it suffices to prove the estimate (A.1) for 𝑢 ∈ 𝐶 ∞ 𝑐 (Ω 𝑓 ).

Step 1: 𝑓 = 0, that is, Ω 𝑓 = Ω 0 = 𝕋 𝑑 × ℝ -. We extend 𝑢 to 𝕋 𝑑 × ℝ by even reflection. That is, 𝑢(𝑥, 𝑦) = 𝑢(𝑥, -𝑦) for 𝑦 > 0. Clearly, ‖𝑢‖ Ḣ1 (𝕋 𝑑 ×ℝ) ≤ 2‖𝑢‖ Ḣ1 (𝕋 𝑑 ×ℝ -) . We have Tr(𝑢)(𝑥) = 𝑢(𝑥, 0) and

By the Cauchy-Schwartz inequality, for 𝑘 ∈ ℤ 𝑑 ⧵ {0},

.

Step 2: 𝑓 ∈ Lip(𝕋 𝑑 ). Setting 𝑣(𝑥, 𝑧) = 𝑢(𝑥, 𝑧 + 𝑓(𝑥)) for (𝑥, 𝑧) ∈ 𝕋 𝑑 × ℝ -, we have Tr(𝑢)(𝑥) = 𝑣(𝑥, 0) and 𝑣 ∈ Ḣ1 (𝕋 𝑑 × ℝ -) with ‖𝑣‖ Ḣ1 (𝕋 𝑑 ×ℝ -) ≤ (1 + ‖𝑓‖ Lip(𝕋 𝑑 ) )‖𝑢‖ Ḣ1 (Ω 𝑓 ) .

Using In the integral defining 𝐽 2 , it suffices to consider 𝑥 ′ near 0 because away from 0 the integrand is uniformly bounded as 𝑟 → 1. For sufficiently small 𝑎 ∈ (0, 𝛾), let 𝐽 0 2 be the contribution of |𝑥 ′ | < 𝑎. In view of (B.2), we have On the other hand, 𝐿 2 can be treated as in (B.8). The proof is complete. □ APPENDIX C: CONSTRUCTION OF 𝝍 𝒓 Let 𝑓 ∈ 𝑊 1,∞ (𝕋 × (0, 𝑇)) be 𝐶 1,1 at 𝑋 0 = (𝑥 0 , 𝑡 0 ) ∈ 𝕋 × (0, 𝑇). Recall the definitions (6.3) and (6.4) of the parabolas 𝜓 and 𝜓 tangent to the graph {𝑦 = 𝑓(𝑥, 𝑡)} from above at (𝑥 0 , 𝑓(𝑥 0 , 𝑡 0 )). We construct a family of functions 𝜓 𝑟 ∈ 𝐶 ∞ (𝕋 × ℝ), 𝑟 ∈ (0, 𝑟 0 ) satisfying the properties in (6.5). Assume without loss of generality that 𝑓(𝑋 0 ) = 0 so that Let Γ 𝑟 be an approximation of the identity in ℝ 2 as 𝑟 → 0 and denote ℎ (𝑟) = ℎ * Γ 𝑟 for any ℎ ∶ ℝ 2 → ℝ. For any 𝑈 ⊂ ℝ 2 , ‖ℎ (𝑟) -ℎ‖ 𝐿 ∞ (𝑈) ≤ 𝑟‖ℎ‖ Lip(𝑈+𝐵 𝑟 (0)) .

For 𝑟 ∈ (0, 𝑟 0 ), let 𝑔 𝑟 = 𝑓 (𝑟 2 ) + (𝑀 + 2)𝑟 2 , so that 𝑔 𝑟 -𝑓 ≥ 2𝑟 2 . By increasing 𝐶 if necessary, we can assume 𝐶 > 4(𝑀 + 2). Then for |𝑋 -𝑋 0 | ∈ [𝑟, 𝑟 0 ), we have Let 𝐹(𝑎, 𝑏) = min{𝑎, 𝑏}. Note that ‖𝐹‖ Lip(ℝ 2 ) = 1 and 𝐹 (𝛿) (𝑎, 𝑏) = 𝑎 if 𝑎 < 𝑏 -𝛿. In addition, 𝐹 (𝛿) (𝑎, 𝑏) ≤ 𝑎 (𝛿) = 𝑎 and likewise 𝐹 (𝛿) (𝑎, 𝑏) ≤ 𝑏.

Define 𝜓 𝑟 = 𝐹 (𝑟 2 ) (𝑔 𝑟 , 𝜓) in 𝐵 𝑟 (𝑋 0 ) and 𝜓 𝑟 = 𝑔 𝑟 in ℝ 2 ⧵ 𝐵 𝑟 (𝑋 0 ). In fact, (C.1) shows that 𝜓(𝑋) -𝑔 𝑟 (𝑋) ≥ 𝑟 2 in an open neighborhood of 𝜕𝐵 𝑟 (𝑋 0 ), so that 𝜓 𝑟 = 𝑔 𝑟 in the same neighborhood. It follows that 𝜓 𝑟 is a smooth function on ℝ 2 . We have 𝜓 𝑟 ≤ 𝜓 in 𝐵 𝑟 (𝑋 0 ) and in view of (C.1), 𝜓 𝑟 (𝑋) ≤ 𝜓(𝑋) for |𝑋 -𝑋 0 | ∈ [𝑟, 𝑟 0 ). Consequently, 𝜓 𝑟 ≤ 𝜓 in 𝐵 𝑟 0 (𝑋 0 ). In 𝐵 𝑟 (𝑋 0 ), when 𝜓 𝑟 (𝑋) < 𝜓(𝑋), we must have 𝜓(𝑋) ≥ 𝑔 𝑟 (𝑋) -𝑟 2 , whence 𝜓 𝑟 (𝑋) ≥ min{𝑔 𝑟 (𝑋), 𝜓(𝑋)} -𝑟 2 ≥ 𝑔 𝑟 (𝑋) -2𝑟 2 ≥ 𝑓(𝑋).

Therefore, 𝜓 𝑟 ≥ 𝑓 in 𝐵 𝑟 (𝑋 0 ).

Let 𝛿 ∈ (0, 𝑟 where 𝐼 = [𝑥 0 -𝜋, 𝑥 0 + 𝜋] and 𝑋 0 = (𝑥 0 , 𝑡 0 ). We have proved that 𝜓 𝑟 satisfies all the properties in (6.5) but with 𝕋 replaced by 𝐼. To finish the proof, we note that 𝑥 0 is the midpoint of 𝐼 and 𝜓 𝑟 = 𝑔 𝑟 for |𝑋 -𝑋 0 | > 𝑟. Therefore, since 𝑔 𝑟 ∈ 𝐶 ∞ (𝕋 × ℝ) and 𝐼 has length 2𝜋, we can periodize 𝜓 𝑟 as 𝜓 𝑟 (𝑥 + 𝑘2𝜋, 𝑡) = 𝜓 𝑟 (𝑥, 𝑡) ∀(𝑥, 𝑡) ∈ 𝐼 × ℝ, ∀𝑘 ∈ ℤ and obtain 𝜓 𝑟 ∈ 𝐶 ∞ (𝕋 × ℝ).