RELAXATION APPROXIMATIONS TO FRONT PROPAGATION by Shi Jin^{1} and Markos A. Katsoulakis^{2} 1 School of Mathematics Georgia Institute of Technology Atlanta, GA 30332 2 Department of Mathematics and Statistics University of Massachusetts Amherst, MA 01003

Abstract. We introduce a relaxation model for front propagation problems. Our proposed relaxation approximation is a semilinear hyperbolic system without singularities. It yields a direction-depedent normal velocity at the leading term and captures, in the Chapman-Enskog expansion, the higher order curvature dependent corrections, including possible anisotropies.

1. Introduction.

The relaxation approximation to systems of conservation laws, first introduced by Jin and Xin [JX1], has been an active area of research in recent years [JX1, KT, Nat]. The idea there is to introduce a relaxation approximation to a system of conservation laws that can generate entropy solutions in the zero relaxation limit. It also gives rise to a class of relaxation schemes that are total-variation-diminishing and free of Riemann solvers. For a general Hamilton-Jacobi equation a class of relaxation approximation was introduced by Jin and Xin [JX2] using the conservation law formulation of the Hamilton-Jacobi equations.

In this paper we construct a related relaxation approximation for front propagation problems. Our proposed relaxation approximation not only yields the direction-dependent normal velocity of the front at the leading term, but also captures, in the Chapman-Enskog expansion, the higher order curvature dependent corrections, including possible anisotropies. This new relaxation model is a first order semilinear hyperbolic system without singularities. First, it may constitute a new regularization for front propagation problems, including illposed problems such as evolutions corresponding to a nonconvex interfacial energy. As in

Typeset by AMS-TEX

The research of S. Jin is partially supported by the National Science Foundation No. DMS-9404157. The research of M. A. Katsoulakis is partially supported by the National Science Foundation No. DMS-9500717.

the case of conservation laws [JX1], it also provides a relaxation scheme for the computation of propagating fronts. In addition to the numerical advantages of relaxation schemes for conservation laws, our approximation does not have singularities that arise in the level set formulation of front propagation problems.

2. Relaxation approximation for general Hamilton-Jacobi equations. Hamilton-Jacobi equations arise in a variety of applications, among others in control

theory, geometric optics, image analysis and phase transitions. It is well-known that such equations in general do not have classical solutions. The correct class of suitable weak solutions, namely the viscosity solutions, was introduced by Crandall and Lions [CL].

The Cauchy problem for a Hamilton-Jacobi equation in R^{n},

(2.1)  \begin{matrix}∂_{t}u + H(∇u) = 0 , \\ u(x, 0) = u_{0}(x) , x ∈ R^{n} \end{matrix}

has a corresponding system of conservation laws obtained formally by differentiating (2.1),

(2.2)  \begin{matrix}∂_{t}p + ∇H(p) = 0 , \\ p(x, 0) = p_{0}(x) = ∇u_{0}(x) , \end{matrix}

where

(2.3)  p(x, t) = ∇u(x, t) .

The equivalence of the weak solutions obtained in the vanishing viscosity limit, of (2.1) and (2.2) for a convex Hamiltonian H, is a classical result in one space dimension [Kru]. In [JX2] this result was also proved for any space dimension. Moreover, a relaxation approximation to (2.2) and subsequently to (2.1), was introduced by Jin and Xin in [JX2]. Their relaxation approximation is

(2.4a)  ∂_{t}p + ∇w = 0 ,
(2.4b)  ∂_{t}w + a∇ · p = −\frac{1}{ϵ}(w − H(p)) ,
(2.4c)  ∂_{t}u + w = 0 ,

where ϵ > 0 is the relaxation time, a is a positive constant satisfying the stability condition

(2.5)  a > |∇_{p}H|^{2} .

In order to understand the behavior of the solution of (2.4) for small ϵ, the Chapman-Enskog expansion can be applied to (2.4) [CLL, JX2]. By ignoring the O(ϵ^{2}) terms, this expansion yields

(2.6)  ∂_{t}u + H(∇u) = ϵ [aΔu − (∇_{p}H(p))^{T} ∇p∇_{p}H(p)] .

The stability condition (2.5) guarantees that the O(ϵ) correction to the Hamilton-Jacobi equation is dissipative.

The relaxation approximation (2.4) has several advantages, such as finite speed of propagation, linear convection, and greater regularity. It also allows the construction of efficient relaxation schemes that are Riemann solver free. Moreover, it might provide an alternative way to construct the weak solution of (2.1). In the case of multidimensional scalar conservation laws, entropy solutions were rigoriously obtained as the zero relaxation limit of suitable relaxation approximations [KT, Nat]. The zero relaxation limit and the entropy properties of general nonlinear hyperbolic systems with relaxation, were studied in [CLL].

3. A relaxation approximation for fronts propagating with curvature dependent velocity.

A front Γ_{t} propagating with normal velocity V = −1 can be described by the zero level set of an auxiliary function u solving the Hamilton-Jacobi equation (2.1), where H(∇u) = |∇u|. By choosing a = 1, (2.6) becomes

(3.1)  ∂_{t}u + |∇u| = ϵ tr [(I − \frac{∇u ⊗ ∇u}{|∇u|^{2}}) ∇^{2}u] .

This equation is the level set formulation of a front Γ_{t} propagating with normal velocity

(3.2)  V = −1 − ϵκ ,

where κ is the mean curvature of the hypersurface Γ_{t}, u is positive inside Γ_{t}. The front is recovered from (3.1) as the level set Γ_{t} = {x : u(x, t) = 0}; in this case the mean curvature is given by the formula

κ = ∇ · n = −\frac{1}{|∇u|}tr [(I − \frac{∇u ⊗ ∇u}{|∇u|^{2}}) ∇^{2}u] ,

where n = −\frac{∇u}{|∇u|} is the outnormal of the front and the velocity is V = \frac{u_{t}}{|∇u|}. Of course, to generate the motion V = 1 − ϵκ, one just needs to take H(∇u) = −|∇u| and a = 1. Similar to Hamilton-Jacobi equations, such equations in general have only weak–viscosity–solutions [CGG, ESp], so the level set formulation allows for the interpretation of (3.2) in a weak sense past all possible singularities and changes of topological type of the front. For an overview to the theory of viscosity solutions for first and second order equations, we refer to [CIL].

An important feature of the level set formulation for (3.2), or any such other local law, is the geometric property of the corresponding PDE (3.1): the zero level sets of the solution u remain invariant under the transformation u ↦ ϕ(u) for an increasing function ϕ with ϕ(0) = 0. The geometric property guarantees that the evolution of the zero level set Γ_{t} depends only on the zero level set Γ_{0} of the initial datum and its sign inside and outside Γ_{0}, and not on the particular choice of the initial condition of (3.1). In addition to its analytical advantages, the level set approach was successfully used in the numerical simulation of curvature dependent front propogation [OS].

Going back to (2.4), we underline the observation that the Chapman-Enskog expansion leading to (3.1) (for the “degenerate” choice a = 1), indicates that the relaxation approximation yields a front propagating with normal velocity V = −1 − ϵκ, including the mean curvature as the next order correction to the front velocity V = −1. Curvature corrections of this type arise in a variety of physical phenomena, one of the most striking being in phase transition problems, for instance the Allen-Cahn theory [AC].

The local law V = −1 or equivalently the PDE u_{t} + |∇u| = 0, describe the isotropic evolution of an interface with normal velocity one. In an analogous way, anisotropic front propagation depending on the normal direction of the front is described by the local law

(3.3)  V = −H(n) ,

where n is the outnormal of the evolving front Γ_{t}. The PDE governing such a propagation is given by the Hamilton-Jacobi equation

(3.4)  u_{t} + H(∇u) = 0 ,

where H is a positively homogeneous function of degree one, i.e.,

H(λp) = λH(p) for λ > 0 ,

and the front is recovered as the zero level set of the viscosity solution of (3.4). Similar to the isotropic case, the higher order correction to (3.4) is given by terms depending on the curvature tensor. This fact has been demonstrated formally and rigorously in phase transition problems, both in microscopic models [Spo, KS] and phase field theories [MWBCS,ES].

If we attempt to approximate (3.4) with the relaxation model (2.4), we readily see that (2.6) is not geometric, therefore it does not constitute a legitimate level set formulation for an underlying evolving front. Here we introduce a relaxation model for (3.4) that will yield a geometric evolution in the Chapman-Enskog expansion, capturing curvature dependent terms as a higher order effect.

We propose the following hyperbolic system with relaxation

(3.5a)  ∂_{t}p + ∇w = 0 ,
(3.5b)  ∂_{t}w + ∇ · q = −\frac{1}{ϵ}(w − H(p)) ,
(3.5c)  ∂_{t}q + ∇w = −\frac{1}{ϵ}(q − J(p)) ,
(3.5d)  ∂_{t}u + w = 0 ,

where

(3.6)  J(p) = \frac{1}{2}∇_{p}H^{2}(p) = H(p)∇_{p}H(p) .

If H(p) = |p| then J(p) = p, and (3.5) recovers the relaxation approximation (2.4) with the choice of q = J(p). Notice also that the relaxation approximation (3.5), at least when q = J(p), corresponds to a perturbation of the Hamilton-Jacobi equation (3.4) by a nonlinear wave operator:

u_{t} + H(∇u) = ϵ(∇ · J(∇u) − u_{tt}) .

Before we study the behavior of (3.5) for small ϵ, we establish the relation between p and ∇u, following [JX2]. By (3.5d),

(3.7)  ∂_{t}∇u + ∇w = 0 .

Combining (3.7) with (3.5a) gives

(3.8)  ∂_{t}(∇u − p) = 0 .

By choosing the initial data to be in local equilibrium i.e.,

(3.9)  p(x, 0) = p_{0}(x) ≡ ∇u_{0}(x) ,

(3.5) has no initial layer. Moreover, (3.8) and (3.9) give

(3.10)  ∇u(x, t) = p(x, t)  for t ≥ 0 .

We now apply the Chapman-Enskog expansion to (3.5),

(3.11)

w(x, t) = H(p) − ϵ (∂_{t}w + ∇ · J(p)) + O(ϵ^{2})
= H(p) − ϵ (∂_{t}H(p) + ∇ · J(p)) + O(ϵ^{2})
= H(p) − ϵ (∇_{p}H · ∂_{t}p + ∇ · J(p)) + O(ϵ^{2})
= H(p) − ϵ (−∇_{p}H · ∇w + ∇ · J(p)) + O(ϵ^{2})
= H(p) − ϵ (−(∇_{p}H)^{T} ∇p∇_{p}H + ∇ · J(p)) + O(ϵ^{2})
= H(p) − ϵ tr [(−∇_{p}H ⊗ ∇_{p}H + ∇_{p}J(p))∇p] + O(ϵ^{2}) .

Recall that the homogeneity of H implies the following properties:

(3.12)

p∇_{p}H(λp) = H(p) , ∇_{p}H(λp) = ∇_{p}H(p) ,
p∇_{p}^{2}H(λp) = 0 , ∇_{p}^{2}H(λp) = \frac{1}{λ}∇_{p}^{2}H(p) .

Consequently

(3.13)  ∇_{p}J(p) = ∇_{p}(H(p)∇_{p}H(p)) = ∇_{p}H(p) ⊗ ∇_{p}H(p) + H(p)∇_{p}^{2}H(p) .

Substituting in (3.11) we get

(3.14)  w = H(p) − ϵ tr [H(p)∇_{p}^{2}H∇p] + O(ϵ^{2}) .

We now replace w in (3.5d) by the previous expression and set p = ∇u. After ignoring the O(ϵ^{2}) term, we get the equation

(3.15)  u_{t} + H(∇u) = ϵ H(∇u) tr [∇_{p}^{2}H(∇u)∇^{2}u] .

It easily follows from (3.12) that this equation satisfies the geometric property and the level set Γ_{t} = {x : u(x, t) = 0} evolves with normal velocity

(3.16)  V = −H(n) + ϵ H(n) tr [∇_{p}^{2}H(n)∇n] .

By a rescaling c^{2}ϵ ↦ ϵ, cH ↦ H in (3.5), one can also obtain the slightly more general anisotropic mean-curvature flow

(3.17)  V = −cH(n) + ϵ H(n) tr[∇_{p}^{2}H(n)∇n] .

Interface laws similar to (3.16) or (3.17) are obtained as singular limits of anisotropic phase field equations. In such models one first postulates the free energy [MWBCS, TC]

(3.18)  E(v) = ∫_{R^{n}}[ϵ^{2}A(∇v) + W (v)]dx ,

where A(p) = \frac{1}{2}H^{2}(p) and W is a double well potential. Then an analogue of the Allen-Cahn equation is derived as the L^{2}–gradient flow of (3.18), which after rescaling is:

(3.19)  v_{t}^{ϵ} = ϵ∇ · ∇_{p}A(∇v^{ϵ}) − \frac{1}{ϵ}W′(v^{ϵ}) .

Formal asymptotic expansions [MWBCS], as well as rigorous analysis [ES], demonstrate that (3.17) governs the asymptotic behavior of the solution of (3.19) (c in (3.17) is related to the relative minima of the double well potential W), as ϵ → 0.

4. Discussion of the relaxation model (3.5).

It is evident that the equation (3.15) is well-posed only if H is positive and convex. However, there are physical examples where H is not convex [CH], thus equation (3.15) has backward parabolic regions for suitable gradients ∇u. In this case microstructure develops [AG], and (3.15) is considered after an ad hoc convexification of H. This approach, at least in two dimensions, gives rise to the so-called crystalline motion by curvature [Tay], which describes the evolution of faceted curves. However the relaxation approximation (3.5) is welldefined in any dimension for even a non-convex H and as a result, it may provide a suitable approximation of the underlying front, both analytically and numerically. In this sense, our model is close in spirit to Suliciu’s relaxation model for van der Waals type equations describing dynamical phase changes [Sul].

Finally, note that (3.15) is a second order fully nonlinear equation with a singularity at ∇u = 0, while the relaxation model (3.5) is a first order semilinear hyperbolic system without singularities. This fact indicates that (3.5) may provide numerical advantages in the computation of propagating fronts governed by (3.15). All these issues, as well as some of their analytical aspects, will be explored in our forthcoming work [JK].

References.

[AC] Allen S.M. and Cahn J.W., 1979. A macroscopic theory for antiphase boundary motion and its application to antiphase domain coarsening, Acta Metal. 27, pp. 1085–1095.

[AG] Angenent S.B. and Gurtin M.E., 1989. Multiphase thermomechanics with interfacial structure, 2. Evolution of an isothermal interface, Arch. Rat. Mech. Anal. 108, pp. 323-391.

[CH] Cahn J.W., and Hoffman D.W., 1974. A vector thermodynamics for anisotropic surfaces -2. Curved and faceted surfaces, Acta Metal. 22, pp.1205-1214.

[CLL] Chen G.-Q,, Levermore C.D. and Liu T.-P., 1994. Hyperbolic conservation laws with stiff relaxation terms and entropy, Comm. Pure Appl. Math., 47, pp. 789-830.

[CGG] Chen Y.-G., Giga Y., and Goto S., 1991. Uniqueness and existence of viscosity solutions of generalized mean curvature flow equations, J. Diff. Geom. 33, pp. 749–786.

[CIL] Crandall M.G., Ishii I. and Lions P.L., 1992. User’s guide to viscosity solutions of second order partial differential equations, Bul. AMS, 27, pp. 1-67.

[CL] Crandall M.G. and Lions P.L., 1983. Viscosity solutions of Hamilton-Jacobi equations, Trans. AMS, 277, pp. 1-42.

[ES] Elliott C.M., and Schatzle R., 1996. The limit of the fully anisotropic double-obstacle Allen-Cahn equation in the non-smooth case, SIAM J. Math. Anal., (to appear)

[ESp] Evans L.C., and Spruck J., 1991. Motion of level sets by mean curvature I, J. Diff. Geom. 33, pp. 635–681.

[JK] Jin S. and Katsoulakis M.A. (in preparation).

[JX1] Jin S. and Xin Z.P., 1995. The relaxation schemes for systems of conservation laws in arbitrary space dimensions, Comm. Pure Appl. Math., 48, pp. 235-277.

[JX2] Jin S. and Xin Z.P., 1996. Numerical Passages from Systems of Conservation Laws to Hamilton-Jacobi Equations, SIAM J. Num. Anal., submitted.

[KT] Katsoulakis M.A. and Tzavaras A.E., 1997. Contractive relaxation systems and the scalar multidimensional conservation law, Comm. PDE, 22, pp. 195-233.

[Kru] Kruzkov S.N., 1975. Generalized solutions to the Hamilton-Jacobi equations of Eikonal type I,

Math. USSR-Sb., 27, pp. 406-466.

[MWBCS] McFadden G.B., Wheeler A.A., Braun R.J., Coriell S.R., Sekerka R.F., 1993. Phase field

models for anisotropic interfaces, Phys. Rev. E. 48, pp2016-2024.

[Nat] Natalini R., 1996. A discrete kinetic approximation of entropy solutions to multidimensional scalar conservation laws, preprint.

[OS] Osher S. and Sethian J.A., 1988. Fronts propagating with curvature dependent speed: algorithms based on Hamilton-Jacobi formulations, J. Comp. Phys., 78, pp. 12-49.

[Spo] Spohn H., 1993. Interface motion in models with stochastic dynamics, J. Stat. Phys. 71, pp. 1081–1132.

[Sul] Suliciu I., 1990. Modeling phase transitions by rate-type constitutive equations: shock wave structures, Int. J. Eng. Sci. 28, pp. 829-841.

[Tay] Taylor J.E., 1993. Motion of curves by crystalline curvature including triple junctions and boundary points, Proc. Symp. Pure Math. 54, pp. 417-438.

[TC] Taylor J.E., Cahn J.W., 1993. Linking anisotropic sharp and diffuse motion laws via gradient flows, J. Stat. Phys. 77, pp. 183-197.