Passive scalars, moving boundaries, and Newton's law of cooling

We study the evolution of passive scalars in both rigid and moving slab-like domains, in both horizontally periodic and infinite contexts. The scalar is required to satisfy Robin-type boundary conditions corresponding to Newton's law of cooling, which lead to nontrivial equilibrium configurations. We study the equilibration rate of the passive scalar in terms of the parameters in the boundary condition and the equilibration rates of the background velocity field and moving domain.


Introduction and motivation
In this paper we study the passive scalar equation where θ = θ(t, x) ∈ R measures some scalar quantity at time t ≥ 0 and position x ∈ Ω(t), u(t, x) ∈ R 3 is a given divergence-free velocity field, and κ > 0 is the diffusivity.
Here Ω(t) is a three-dimensional open domain that may depend on time. A full specification of Ω(t) will be provided later. The scalar θ is understood to be passive in the sense that u affects the dynamics of θ through (1.1), but θ plays no role in the dynamics of u: it is passive. We shall think of u as the velocity field of an incompressible fluid evolving in the (possibly moving) domain Ω(t) and θ as the temperature of the fluid. Passive scalars are a popular model of turbulent diffusion (see the exhaustive survey of Majda-Kramer [9] and references therein), mixing phenomena (see for instance Lin-Thiffeault-Doering [8] or Thiffeault [13]), as well as pollution and combustion (see Warhaft [14] and references therein). If we consider (1.1) on T 3 , then, as we will show below, the solution θ will converge exponentially as t → ∞ to its average, which is an equilibrium solution of (1.1). An interesting result of Constantin-Kiselev-Ryzhik-Zlatoš [3] shows that certain flows u can actually enhance this convergence to equilibrium. This was extended to time-periodic flows by Kiselev-Shternberg-Zlatoš [4] and to flows on R 2 by Zlatoš [15] (in which case the equilibrium is 0).
We aim to study the temporal decay properties of the scalar θ in terms of the decay properties of the given velocity u, the asymptotic behavior of Ω(t), and the the boundary conditions imposed on θ and u. We investigate the decay rate of θ for long times and also show that the decay leads to global-in-time solutions for arbitrary data.
1.1. Periodic boundary conditions. To motivate our analysis let's first consider a relatively trivial example, when Ω is taken to be a periodic box T 3 = R 3 /Z 3 and u(t, x) is a divergence-free, periodic, smooth vector field. In this case, by using the periodicity and the incompressibility condition (i.e. div u = 0), it is easy to see that a solution θ to (1.1) obeys the L 2 -energy law Similarly, we may integrate (1.1) to find that d dt T 3 θdx = 0. (1.3) is a constant determined by the initial average of the temperature. Since θ − θ ave also solves (1.1) with periodic boundary conditions, we also know that 1 2 We now recall the (sharp) Poincaré inequality in T 3 : Hence, by combining (1.5) and (1.6), we have 1 2 which in turn implies that θ(t, ·) − θ ave L 2 ≤ e −4π 2 κt θ(0, ·) − θ ave L 2 . (1.8) Therefore, we deduce that θ converges exponentially fast to θ ave . Since the inequality (1.6) is sharp, this decay result is optimal. We would like to highlight three important features of this decay analysis. First, the periodic boundary conditions on θ play a role in the decay through the derivation of (1.5). Second, the divergence-free condition and the periodic boundary conditions on u also play a role in the derivation of (1.5); in this particular setting u actually disappears entirely from (1.5) and has no ultimate impact on the decay rate in (1.8). Third, the fact that Ω(t) = T 3 is static in time gives rise to a static constant in the Poincaré inequality (1.6), which in turn plays a serious role in determining the decay rate in (1.8). Indeed, the Poincaré constant and the diffusion constant κ both appear directly in the decay rate.
Throughout the rest of this paper we will investigate how changes in these three features affect the decay of θ. We will consider both static (and rigid) domains as well as moving domains with a range of boundary conditions for θ and u. In the rest of the section we will give concrete examples of couplings between u and θ in order to provide motivation for our main results. However, in the main results of Section 2, we will consider a more general and weaker coupling between u and θ.

Rigid boundary.
We consider a fluid confined between two rigid horizontal plates located at x 3 = 0 and x 3 = −d for d > 0. We allow the horizontal cross-section to be either infinite, in which case we define Γ = R 2 , or else periodic, in which case we set Γ = (L 1 T) × (L 2 T), where L 1 , L 2 > 0 are the periodicity lengths and LT = R/(LZ). Then the rigid domain is (1.9) We shall write Σ + = {x 3 = 0} and Σ − = {x 3 = −d} (1.10) to denote the top and bottom boundaries of Ω.
We assume that the fluid obeys the gravity-driven incompressible Navier-Stokes equations in Ω u = 0 on Σ + and Σ − . (1.11) Here µ > 0 is the fluid viscosity, g > 0 is the gravitational strength, and e 3 = (0, 0, 1) ∈ R 3 . Actually, the gravitational force plays no real role in (1.11) since we may absorb this conservative force into the pressure by swapping p → p + gx 3 . We have included gravity here simply to maintain a similarity to the moving boundary case, where gravity does play a role. We assume that the temperature of the fluid satisfies the passive scalar equation (1.1) in Ω. Notice that θ does not appear in (1.11): it is passive in the dynamics of u. However, the velocity u does play a role in (1.1). It remains to specify boundary conditions for θ.
Since we are thinking of θ as the temperature of the fluid, there is a natural continuum of choices for boundary conditions at each plate. We will assume that the space adjacent to each plate is maintained at a constant temperature, with θ ext,− and θ ext,+ denoting the temperature below the bottom plate and the temperature above the top plate, respectively. One possibility is to suppose that the plate conducts heat perfectly into the fluid, which requires the temperature to match the temperature external to the plate. This is modeled through a Dirichlet condition θ = θ ext . A more general condition supposes that the plate is an imperfect heat conductor. In this case it is common (see for example Slattery's book [11] for a discussion of these conditions for heat conduction and mass transfer) to model the plate with Newton's law of cooling, which requires the Robin-type boundary condition where β ∈ [0, ∞] is the heat transfer coefficient and ν is the outward-pointing unit normal on the plate. Notice when β = ∞ we recover the Dirichlet condition θ = θ ext ; when β = 0 we recover what is known as an insulating condition, corresponding to a homogeneous Neumann boundary condition. We may reduce the number of parameters to be considered by studying the evolution of θ − θ ext,− instead of θ. The equation (1.1) is still satisfied, but this allows us to replace the temperature below the lower plate by 0 and the temperature above the upper plate bȳ (1.13) Since the materials above and below the fluid may in principle be different, we will allow for different heat transfer coefficients in Newton's law of cooling at the top and bottom. We summarize the equations for θ as follows: (1.14) Here we have used the fact that ν = ±e 3 on Σ ± . Henceforth we will write β = (β + , β − ) ∈ [0, ∞] 2 for the pair of heat conductivity coefficients.
1.3. Moving boundary. We shall also consider a fluid evolving in a moving domain Ω(t).
Since we wish to compare the decay properties of θ in moving domains and in rigid domains, we will suppose that Ω(t) is of the same general "slab-like" form as the domain Ω considered above. We again assume Γ = R 2 or Γ = (L 1 T) × (L 2 T) to allow for infinite or periodic horizontal cross sections. We then assume that where η : R + × Γ → R is the "free surface function." Then is the free surface of the fluid, while is the rigid bottom of the fluid domain.
The motion of the fluid is governed by the gravity-driven incompressible Navier-Stokes system in Ω(t): Here µ > 0 is the fluid viscosity, ν is the outward-pointing unit normal on Σ(t), I is the 3 × 3 identity matrix, (Du) ij = ∂ i u j + ∂ j u i is the symmetric gradient of u, g > 0 is the gravitational constant, e 3 = (0, 0, 1), σ ≥ 0 is the surface tension coefficient, and H is twice the mean curvature of the surface Σ(t).
The dynamics of the temperature of the fluid are described by the convection-diffusion equation (1.1) with boundary conditions of the form considered above: As in the rigid boundary case, we will consider β = (β + , β − ) ∈ [0, ∞] 2 . Notice that θ is passive in the dynamics of u: it does not appear in (1.18). On the other hand, both the velocity field u and the changing boundary (i.e. η) affect the dynamics of θ.
1.3.1. Reformulation in flattened coordinates. In moving boundary problems, it is convenient to flatten the free surface by using a suitable coordinate transformation. We will utilize a flattening coordinate transformation introduced by Beale [1]. To this end, we consider the fixed domain Ω given by (1.9), for which we will write the coordinates as x ∈ Ω. We shall again write Σ + = {x 3 = 0} for the upper boundary and Σ − = {x 3 = −d} for the lower boundary. We view η as a function on R + × Γ. We then definē Note that Φ(t, Σ + ) = {y 3 = η(t, y 1 , y 2 )} = Σ(t) and Φ(t, ·)| Σ − = Id Σ − , i.e. Φ maps Σ + to the free surface Σ(t) and keeps the lower surface fixed. We write (1.23) Here J = det ∇Φ is the Jacobian of the coordinate transformation. If η is sufficiently small (in an appropriate Sobolev space), then the mapping Φ is a C 1 diffeomorphism. This allows us to transform the problem to one on the fixed spatial domain Ω for t ≥ 0. In the new coordinates, the PDE (1.18) becomes with given initial data u(0, x) = u 0 (x), η(0, x ′ ) = η 0 (x ′ ) where x ′ = (x 1 , x 2 ). We will employ the Einstein convention of summing over repeated indices for vector and tensor operations. This allows us to define the differential operators used in (1.24) as follows. The actions of ∇ A , div A , and ∆ A are given by ( We have also written for the non-unit normal to Σ. Similarly, in the new coordinates, the PDE (1.19) reads with given initial data θ(0, x) = θ 0 (x). Here ν = N |N | . Notice that all the differential operators in (1.24) and (1.26) are connected to η and thus to the geometry of the free surface.

1.4.
Equilibria. An equilibrium solution corresponds to u = 0 and Ω(t) = Ω given by (1.9) (and hence η = 0 in the case of a moving boundary). In this case (1.14) and (1.26) reduce to the same equations, so the equilibria coincide for the rigid and moving boundaries. To denote the dependence on β = (β + , β − ) let us denote the equilibrium temperature by θ β eq . We need to solve ∆θ β eq = 0 with the given boundary conditions. The general equilibrium solution is then Notice in particular that the equilibrium with either β + = ∞ or β − = ∞, but not both, is recovered from (1.27) by taking a limit: (1.28) The equilibrium with β ± = ∞ (pure Dirichlet boundary conditions) is recovered through a further limit: (1.29) We note that unlike in the case of a periodic box, the equilibrium temperature is in general not a constant because of the boundary conditions. We will see that solutions of (1.14) and (1.26) converge asymptotically to these equilibrium solutions when β = (0, 0). When β = (0, 0) the equilibrium is not unique; this is related to the fact that solutions converge (at least in the periodic case) to a constant determined by the initial data.

Main results
2.1. Well-posedness. We begin with a discussion of the notion of weak solutions to (1.14) and (1.26), given as perturbations of the equilibrium state θ β eq . Since θ β eq is not in L 2 (Ω) when Ω is horizontally infinite, we must formulate the weak problem in terms of the perturbation. We will also discuss the construction of local-in-time solutions.
First we define the appropriate functional setting for solutions. Let H 1 β (Ω) denote the space (2.1) Now we give the notion of weak solutions for the rigid problem.
We say that a a function θ : [0, T ] × Ω → R is a weak solution to (1.14) if the following three conditions are satisfied.
(i) We have the inclusions Here ·, · * denotes the pairing between H 1 β (Ω) and its dual, and The justification of the weak formulation (2.4) is standard: if θ is a smooth solution to (1.14) then (2.4) may be derived through integration by parts and use of the boundary conditions satisfied by θ β eq . In this case we have (2.6) Remark 2.3. Standard arguments and the inclusions of (i) imply that This is the sense in which (ii) is required to hold.
follows from the standard trace theory.
Remark 2.5. We require the inclusion (2.2) in order to make the integrals Notice, though, that our definition of weak solution does not require that div u = 0. We will invoke this condition only in a global setting.
We may readily prove the following local well-posedness result for (1.14). The proof, which we omit for the sake of brevity, is a standard application of a Galerkin method based on the same energy estimates that we will use for our global existence and decay results. We refer to Section 3 for the energy estimates.
as well as the estimate We say that a a function θ : [0, T ] × Ω → R is a weak solution to (1.26) if the following three conditions are satisfied.
(i) We have the inclusions Here ·, · * denotes the pairing between H 1 β (Ω) and its dual, C β t :

13)
and F β t : H 1 β (Ω) → R is the t−dependent force term given by Some remarks are in order.
Remark 2.8. The weak formulation (2.12) is justified in the same way as (2.4). The only difference is that the computations in the integration by parts are somewhat more involved. For the sake of readability we have moved the computation to the appendix: see Lemma B.5.
Remark 2.9. The influence of the moving boundary is manifest in (2.12) through the appearance of the terms J, A, etc, all of which depend on η, the free surface function. Remark 2.11. In the weak solution definition we make no assumptions about u and η satisfying any of the equations of (1.24). We will need some of these only in the global theory.
We may again readily deduce the existence of local weak solutions to (1.26) by employing a standard Galerkin method based on the energy estimates that we will use for our global existence and decay results, and Lemmas B.1, B.2, and B.3. Notice in particular that Lemma B.2 provides for the ellipticity condition on A and for the fact that the boundary integrals are controlled via the usual H 1 (Ω) trace theory. We again refer to Section 3 for the energy estimates, and again omit the proof of local existence result. Proposition 2.12. Fix T ∈ (0, ∞). Let u, η satisfy (2.9) and (2.10), and let θ 0 ∈ L 2 (Ω). Then there exists a unique weak solution θ : [0, T ] × Ω → R to (1.26).

Global solutions and equilibration.
We now state our main results on the existence of global solutions and their decay to equilibrium. We begin with the case of a rigid domain.
(i) When β ± = 0 (pure Neumann boundary conditions) and Γ = (L 1 T) × (L 2 T), then we define Also, regardless of the decay of u, θ − θ β eq − θ avg decays to 0 exponentially fast with the following estimate (2.20) (ii) When β ± = 0 and Γ = R 2 , we have the estimate we further assume the L 2 decay of u: where g(t) → 0 as t → ∞. Then θ − θ β eq decays to 0 with the following estimate Some remarks are in order.
Remark 2.14. When β = (0, 0) and Γ = R 2 we fail to actually show that the solutions decay to equilibrium. This is due to the lack of coercivity induced by the infinite cross-section. Instead we have bounds on the deviation of θ from the equilibrium. When β = (0, 0) and Γ = (L 1 T) × (L 2 T) we find that θ converges exponentially fast to the constant θ β eq + θ avg . Recall that θ β eq is an arbitrary constant when β = (0, 0). This means that we can view the equilibrium as being determined by the data in the sense that (2.19) implies Remark 2.16. In the case of insulating boundary conditions, β = (0, 0), the decay rate is independent of the velocity field u. The same is true for any β ∈ [0, ∞] 2 ifθ = 0 since in this case θ β eq = 0.
Remark 2.17. If g(t) decays faster than e −µ β t , then the decay rate of θ is dictated by e −µ β t .
for some constant C. If g(t) decays slower than e −µ β t , the decay is dominated by g(t).
In the case of moving boundary, we have the following result.

Remark 2.21.
In the moving boundary case, the decay rate is not as explicit as the rigid case. We need the decay of not only velocity u but also moving boundary η in order to derive the decay rate of θ for β > 0.
Remark 2.22. Again ifθ = 0 then the decay is not influenced by the asymptotics of u. The asymptotics of η, and hence Ω(t), do play a role in the bounds on J, which contribute the c 0 terms to the decay rate.
Remark 2.23. The reason we include ∇ * η H 1/2 (Γ) in (2.36) but η H 3/2 (Γ) in (2.38) is that in case (iii) we also need η H 1 (Γ) to decay, which in particular requires the L 2 decay of η. This is due to the appearance of the term 1 − K |N | in F β t . This term also causes the difference in the structure of the estimates (2.39) and (2.34), (2.37). The former estimates the decay of the square of the L 2 norm while the latter two estimate the decay of the L 2 norm.
We conclude our discussion of the moving boundary problem with a couple concrete of examples. In particular, we consider cases (iii) and (iv) under the extra assumption that u and η satisfy the system (1.24). Initially suppose that σ > 0 in (1.24), which corresponds to surface tension on the moving interface. When Γ = R 2 the work of Beale-Nishida [2] showed that with the slower rate of decay for h determined by the slow decay of η(t, ·) L 2 . When Γ = (L 1 T) × (L 2 T) the work of Nishida-Teramoto-Yoshihara [10] showed that for some γ > 0. Thus we see that the faster decay in the periodic case also leads to faster decay of the passive scalar. Now consider the case σ = 0. When Γ = R 2 the work of Guo-Tice [6] showed that where λ ∈ (0, 1) is a measure of the negative Sobolev regularity of the initial data. When Γ = (L 1 T) × (L 2 T), Guo-Tice showed in [5] that for any m ≥ 4 there exists a smallness condition on the data for (1.24) that guarantee that (2.43) Again we find that the periodicity leads to faster decay rates.

L 2 -energy estimates
In this section we present the basic L 2 energy estimates.
3.1. Rigid boundary. We now derive the L 2 energy estimate in the case of a static rigid domain. In order to state these, we must first define the dissipation functional in terms of β.
Recall the spaces H 1 β (Ω) defined in (2.1). We define D β : Lemma 3.1. Assume that u satisfies the assumptions of Theorem 2.13, and suppose θ is a weak solution to (1.14). Then we have 1 2 where D θ is given by (3.1).
Proof. We may use θ − θ β eq as a test function in (2.4) to deduce that where and A standard computation shows that We may use the fact that div u = 0 in Ω and u 3 = 0 on ∂Σ ± to compute (3.9) Finally, we compute The equality (3.2) then follows by combining these.
We remark that the right-hand-side of (3.2) does not have a definite sign in general and it includes not only the temperature fluctuation but also the vertical velocity u 3 .
In the case of the Neumann boundary condition, when β = (0, 0), the equilibrium is given by θ β eq = C. In this case the right-hand-side of (3.2) is zero, and hence we obtain the equality 1 2 In order to deduce decay information we need a slightly different version of this equality.

Moving boundary.
We now seek to derive the L 2 energy estimate in the case of a moving boundary. In this case we define the time-dependent dissipation functional M t β : Here we have written M t β to emphasize the dependence on time t: J, K, A, and N are all understood to be evaluated at time t in (3.15) We can derive the L 2 energy estimate. We will employ Lemma B.4 for many of the calculations.
where M t β is given by (3.15) and F β t is given by (2.14).
Proof. We use θ − θ β eq ∈ H 1 β (Ω) as a test function in (2.12) to derive the equality where (3.20) Employing (B.16) and the identity JK = 1, we rewrite I as (3.21) Identity (B.19) shows that u j JA jk e 3 · e k = u 3 = 0 on Σ d , and we know thatd = 0 on Σ d , so we may integrate by parts to see that The incompressibility condition (2.27) is equivalent to A jk ∂ k u j = 0, and an easy computation shows that ∂ t J = ∂ 3 (∂ tη d). Since ∂ t η = u · N on Σ, we may use (B.17) and (B.18) to deduce that I 2 = 0. Hence, We simply rewrite II in (3.19) as It remains to handle III in (3.20). Since ∂ k θ β eq = δ 3k ∂ k θ β eq and ∂ 3 θ β eq is constant, we can rewrite where M t β is given by (3.15). Also, if we write Proof. To prove (3.27) we use θ(t, ·) − θ β eq − C as a test function and argue as in the proof of Lemma 3.3. To prove (3.29) we argue similarly, using 1 as a test function to derive the equality We also compute, using Lemma B.4, Then (3.29) follows by integrating these identities and using the definition (3.28).

Coercivity of the dissipation
In order for the L 2 energy estimates of the previous section to give rise to decay results, we need for the dissipation terms D β , defined by (3.1) for rigid boundaries, and M t β , defined by (3.15) for moving boundaries, to be L 2 -coercive. In this section we pursue the proof of this coercivity estimate, which is essentially a Poincaré inequality. 4.1. One-dimensional analysis. When β = (0, 0), our coercivity estimates will be based on a corresponding coercivity estimate in the vertical direction. Our aim now is to prove this estimate.

(4.2)
We now prove a variational principle.
Then there exists ζ β ∈ C such that Moreover, and in fact this inequality is an equality when β ∈ {0, ∞} 2 .
Proof. The existence of ζ β , a minimizer of D β over the constraint set C, follows from a standard application of the direct methods in the calculus of variations, so we shall omit the details. It's clear that µ β := D β [ζ β ] ≥ 0. The minimizer must then satisfy the equations with the obvious reinterpretation when β ± = ∞. When β ∈ (0, ∞) 2 , the problem (4.6) corresponds to finding eigenfunctions of the Laplacian with Robin boundary conditions in one dimension. Standard analysis (see for example, Chapter 4.3 of Strauss's book [12]) shows that µ β must be the smallest positive solution to Although this equation cannot be solved analytically in general, it is a simple matter to see that When β + = 0 and β − ∈ (0, ∞), the problem (4.6) requires that µ β be the smallest positive solution to which must lie in the interval (0, π/2). By employing the estimate 2 for x ∈ (0, π/2), (4.10) we may show that µ β satisfies the bound 1 4d 2 In the case β − = 0 and β + ∈ (0, ∞) a similar analysis allows us to deduce that 1 4d 2 The case β + = ∞ and β − ∈ [0, ∞) for problem (4.6) leads to the equation Similarly, the case β − = ∞ and β + ∈ [0, ∞) leads to the equation (4.14) In either case we find that with equality achieved if the finite β term actually vanishes. The case β + = β − = ∞ in (4.6) leads to the equation The smallest positive solution is then The remaining cases β + = β − = 0 in (4.6) leads to the equation (4.16) as well, but the solution must only be non-negative, and hence µ β = 0.
As a consequence of the minimality of µ β we deduce a corresponding one-dimensional coercivity estimate. Hence, for any ζ ∈ H 1 β ((−d, 0)) we may estimate which is the desired inequality.

Rigid boundary.
With the one-dimensional coercivity in hand, we can now derive the general coercivity estimate in the case of a rigid boundary. We recall the definition of H 1 β (Ω) from (2.1) and D β from (3.1). We also define the vertical part of D β according to When Γ = (L 1 T) × (L 2 T) this infimum is actually a minimum.
We know from Lemma 4.1 that µ β > 0 when β = (0, 0). In this case the conclusion of the proposition is non-trivial. When β = (0, 0) the result is actually trivial since µ 0 = 0. Our next goal is to derive a non-trivial estimate when β = (0, 0), under a stronger assumption on the functions in question. We will only be able to do so in the case of a periodic horizontal cross section. Moreover, the inequality is sharp in the sense that Proof. We define the constraint set The direct method in the calculus of variations allows us to produce a minimizer of D β over C.
That is, we can find ψ 0 ∈ C such that Standard arguments reveal that ψ 0 must satisfy for n 1 , n 2 , m ∈ N. However, the condition Ω ψ 0 eliminates the zero eigenvalue as a possibility, and so λ must be the second smallest eigenvalue, which is easily computed: (4.41) We might try something similar when Γ = R 2 . The problems then are two-fold. First we don't know for sure that ϕ ∈ L 1 (Ω), and so the condition Ω ϕ = 0 need not make sense. Second, and more fundamental, is that when β = (0, 0) The upshot of this analysis is that the dissipation functional D β is simply not L 2 -coercive when Γ = R 2 and β = (0, 0).

4.3.
Moving boundary. Now we derive a coercivity estimate in the case of a moving boundary. Let M t β be given by (3.15).
Again this estimate is only useful when β = (0, 0). When β = (0, 0) we must use Proposition 4.4 as the basis of our estimate. Let M t β be given by (3.15). and that the minimizer is constant function ϕ 0 . We may use the direct methods in the calculus of variations to find ϕ 1 ∈ H 1 (Ω)\{0} such that ψ ∈ H 1 (Ω)\{0} and ψ, ϕ 0 0,t = 0 =: λ 1 (t). (4.56) For any w ∈ H 1 (Ω)\{0} we will write for the co-dimension one subspace that is perpendicular to w with respect to the inner-product ·, · 0,t . We claim that ψ 1 satisfies the maximin principle To prove the claim we first set V = span{ϕ 0 , ϕ 1 }. The definitions of ϕ 0 , ϕ 1 easily imply that and hence that From this and the definition of λ 1 (t) as a minimizer over X t (ϕ 0 ), we deduce that (4.58) holds, proving the claim. Now we return to the proof of Proposition 4.4; we may use the matching homogeneities of D β [·] and · 2 L 2 (Ω) to rewritē Since JK = 1 we have that . (4.64) The inequalities (4.50) and (4.51) allow us to estimate   Since ψ 0 is constant, this implies (4.52) when ψ = 0. The case ψ = 0 is trivial.

Proof of Main results
We start with the proof of Theorem 2.13.
Proof of Theorem 2.13. We will only prove the estimates; with them in hand the existence of global solutions then follows from Proposition 2.6 and a standard continuation argument. The global solutions must then satisfy the decay estimates as well. The argument for (iii) is somewhat simpler, so we begin with it. From Lemma 3.1 and Proposition 4.3, we first obtain 1 2 Let f (t) := θ(t, ·) − θ β eq L 2 (Ω) . Since u(t) L 2 (Ω) ≤ g(t), we deduce from (5.1) that f satisfies df dt which leads to d dt e µ β t f ≤ ∂ 3 θ β eq e µ β t g(t).
We next present the proof of Theorem 2.18.
Proof of Theorem 2.18. As in the proof of Theorem 2.18, we will only prove the estimates. The global existence claim follows from these and a continuation argument. We start with (i). When β = (0, 0), we have the energy identity (3.27) and the zero-average condition (3.29) from Lemma 3.4. Together with the coercivity estimate (4.51) of Proposition 4.6, this implies that which easily implies (2.34).
To prove (ii) we use the identity (3.26), which followed from Lemma 3.3, and simply integrate in time.
For (iii), we first note that from (3.16), (4.45), and the Cauchy-Schwarz inequality, ≤ Cg(t). (5.9) Combining these, we obtain for a universal constant C > 0. Let f := (θ − θ ∞ eq ) √ J L 2 (Ω) . Then f satisfies By integrating in t, we deduce (2.37). In order to establish (iv), we first use the energy identity (3.16) to obtain =: I + II. For I we argue as in case (iii) to bound . Suppose that f belongs to eitherḢ q (Γ) or elseḢ q−1/2 (Γ) for some q ∈ N (hereḢ s is the usual homogeneous Sobolev space of order s).
Also, for q ∈ N we have the estimate Proof. These estimates are proved in [7]: Lemmas A.7 and A.8 handle the case Γ = R 2 and Lemmas A.9 and A.10 handle the case Γ = (L 1 T) × (L 2 T).
Appendix B. Some tools for the moving boundary problem In this section we record a number of results that are essential the analysis of the moving boundary problem. We begin with some applications of Lemma A.1 to η and ∂ t η.
Lemma B.1. We have the estimates for all ϕ ∈ H 1 (Ω). where ∇ * η = (∂ 1 η, ∂ 2 η) is the horizontal gradient; also  Proof. The first, third, fourth identities may be found in Lemma 2.1 of [6]. The second is in the proof of Lemma 2.2 of the same paper.
Finally, we record an energy identity for smooth solutions of (1.26) that motivates the definition of weak solution in Definition 2.7.