L¹-FRAMEWORK FOR CONTINUOUS DEPENDENCE AND ERROR ESTIMATES FOR QUASILINEAR ANISOTROPIC DEGENERATE PARABOLIC EQUATIONSL¹-FRAMEWORK FOR CONTINUOUS DEPENDENCE AND ERROR ESTIMATES FOR QUASILINEAR ANISOTROPIC DEGENERATE PARABOLIC EQUATIONS

Pure Mathematics

ISBN 82–553–1368–0 No. 6 ISSN 0806–2439 March 2003

L¹–FRAMEWORK FOR CONTINUOUS DEPENDENCE AND ERROR ESTIMATES FOR QUASILINEAR ANISOTROPIC

DEGENERATE PARABOLIC EQUATIONS

GUI-QIANG CHEN AND KENNETH H. KARLSEN

Abstract. We develop a general L¹–framework for deriving continuous dependence and error estimates for quasilinear anisotropic degenerate parabolic equations with the aid of the Chen-Perthame kinetic approach [9]. We apply our L¹–framework to establish an explicit estimate for continuous dependence on the nonlinearities and an optimal error estimate for the vanishing anisotropic viscosity method, without the requirement of bounded variation of the approximate solutions. Finally, as an example of a direct application of this framework to numerical methods, we focus on a linear convection-diffusion model equation and derive anL¹ error estimate for an upwind-central finite difference scheme.

Contents

1. Introduction 1

2. Entropy Solutions and Kinetic Formulation 4

3. GeneralL¹–Framework 6

4. Estimates for Continuous Dependence on the Nonlinearities 19 5. Error Estimates for the Vanishing Anisotropic Viscosity Approximation 19 6. Error Estimates for a Finite Difference Approximation 22

Acknowledgments 25

References 25

1. Introduction

We are concerned with the Cauchy problem for quasilinear anisotropic degenerate parabolic equations of second order with the form

(1.1) ∂_tu+ divf(u) =∇ ·(A(u)∇u), u(0, x) =u₀(x),

where (t, x) ∈ R+ ×R^d, div and ∇ are with respect to x ∈ R^d, u = u(t, x) is the scalar unknown function that is sought,

(1.2) u₀∈L¹(R^d)∩L^∞(R^d)

Date: March 25, 2003.

1991Mathematics Subject Classification. 35K65,35B35,35G25,35D99.

Key words and phrases. L¹–framework, degenerate parabolic equations, quasilinear, anisotropic, entropy solutions, kinetic formulation, continuous dependence, error estimates, vanishing viscosity, difference schemes.

1

is the initial function,

(1.3) f = (f1, . . . , f_d)∈(Lip_loc(R))^d is the vector–valued flux function, and

(1.4) A(u) =σ^A(u)σ^A(u)^>≥0, with σ^A∈(L^∞_loc(R))^d×K, 1≤K ≤d,

is the matrix–valued diffusion function. The symmetric d×d matrix A(u) = (a_ij(u)) has entries of the form

aij(u) =

K

X

k=1

σ_ik^A(u)σ^A_jk(u), i, j= 1, . . . , d.

On the space of symmetric matrices, we employ the usual ordering in the sense of quadratic forms. Note that the scalar hyperbolic conservation law (A≡0) is a special case of (1.1).

Nonlinear partial differential equations of type (1.1) model convection-diffusion motions in nature and occur in a variety of applications. Being very selective, we mention here only flow in porous media (see, e.g., [12] and the references therein) and sedimentation- consolidation processes [6]. It is well known that equation (1.1) possesses discontinuous solutions, and weak solutions are not uniquely determined by their initial data; hence (1.1) must be interpreted in the sense of entropy solutions [19, 30, 31]. The uniqueness of entropy solutions was proved in the one-dimensional context by Wu-Yin [32] and B´enilan-Tour´e [2]. In the multidimensional context with isotropic diffusion (that is,A(·) ≥0 is a scalar function), a general uniqueness result is much more recent and was proved by Carrillo [7] by using Kruˇzkov’s doubling of variables technique; and various extensions of his result can be found in B¨urger-Evje-Karlsen [5], Eymard-Gallou¨et-Herbin-Michel [15], Karlsen-Risebro [18], Karlsen- Ohlberger [16], Mascia-Porretta-Terracina [24], Michel-Vovelle [25], and Rouvre-Gagneux [28].

Chen-DiBenedetto [8] proved the uniqueness of unbounded entropy solutions by using the doubling of variables technique. Chen-Perthame [9] finally introduced the notion of kinetic solutions and established an L¹ well-posedness theory for the general anisotropic diffusion case by developing a kinetic approach for (1.1). Let us also mention the earlier work by Tassa [29], who proved the uniqueness for piecewise smooth weak solutions. There are also several recent studies concerned with the convergence of various numerical schemes: see [12]

for operator splitting methods, [13, 17] for monotone finite difference schemes, [15, 26, 25] for monotone finite volume schemes, and [1, 4] for BGK schemes. All these papers provide the L¹ convergence of approximate solutions without a rate of convergence (an error estimate).

As is well known,L¹ error estimates are more desirable for robust scientific computation and prediction, which have been an open problem for the general anisotropic case in numerical analysis.

In the hyperbolic context (i.e., A≡0), error estimates for the vanishing isotropic viscosity method were derived first in Kuznetsov [20] and more recently in Cockburn-Gremaud [10]

and Bouchut-Perthame [3, 27], while various estimates for continuous dependence on the nonlinearity (i.e., the flux functionf) were obtained first in Lucier [22] and later in Bouchut- Perthame [3]. Regarding degenerate parabolic problems with isotropic diffusion (that is,A(·) is a scalar function), continuous dependence estimates for semigroup solutions, and hence also error estimates for the vanishing isotropic viscosity method, were obtained by Cockburn- Gripenberg [11]; see also [18, 14] for a different approach for the case that the flux function f also depends on (t, x).

We are concerned with explicit estimates for continuous dependence on the nonlinearities and error estimates for the vanishing anisotropic viscosity method for (1.1). We mention that, even in the isotropic case, continuous dependence estimates have never been derived directly for entropy solutions. The purpose of this paper is to use the Chen-Perthame kinetic approach [9] to develop an abstract L¹–framework for continuous dependence and error estimates for (1.1) and to present several applications of this framework.

More precisely, we are interested in comparing an entropy solutionu=u(t, x) of (1.1) with an entropy solution v=v(t, x) of

(1.5) ∂tv+ divg(v) =∇ ·(B(v)∇v) + error terms, v(0, x) =v0(x), where

(1.6) v₀∈L¹(R^d)∩L^∞(R^d),

(1.7) g= (g₁, . . . , g_d)∈(Lip_loc(R))^d, and

(1.8) B=σ^B(v)σ^B(v)^>≥0, σ^B ∈(L^∞_loc(R))^d×K˜ , 1≤K ≤d.

The symmetric d×dmatrix B=B(v) = (bij(v)) has entries bij(v) =

K˜

X

k=1

σ_ik^B(v)σ^B_jk(v), i, j= 1, . . . , d.

Similar to the treatment of hyperbolic problems [3, 27], the error terms will take the form of

“partial derivatives” for applications, which will be specified later in Section 3.

The first application of our general L¹–framework is an explicit estimate for continuous dependence on the nonlinearities in (1.1). If g ≡ f (see Section 4 for the general case), u₀ ∈BV(R^d), and the error terms are zero in (1.5), we obtain that, for any t >0,

(1.9)

ku(t,·)−v(t,·)k_L1(R^d)

≤ ku₀−v0k_L1(R^d)+C√ t

s

√ A−√

B √

A−√ B

>

_∞

,

where the ∞ - norm is taken componentwise (see Section 3 for the precise definition). We must emphasize that the proof of a result like (1.9) depends in a fundamental way on using the parabolic dissipation/defect measure identified in Chen-Perthame [9], which is also the cornerstone of the uniqueness proof in [9].

The second application of ourL¹–framework is an error estimate for the vanishing anisotropic viscosity method for (1.1):

(1.10) ∂tv+ divf(v) =∇ ·

A(v)∇v

+µ∇ ·

B(v)∇v

, v(0, x) =v0(x),

where the matrix B(v)>0 is of the same type as in (1.8). If u₀ ∈BV(R^d), we prove that, for anyt >0,

(1.11) ku(t,·)−v(t,·)k_L1(R^d)≤ ku₀−v0k_L1(R^d)+C√ tµ, whereC depends only on theL^∞ norms of the matrices A andB.

Within our L¹–framework, there are two ways to obtain an L¹ estimate for u−v. A traditional way is to view the equation for the anisotropic viscous approximate solutions as

the original equation perturbed by the error terms taking the form of partial derivatives. Ifv is uniformlyBV bounded in space variables, one obtains the optimal ¹₂ rate of convergence.

However, ifv is notBV bounded, only a sub-optimal rate of convergence can be obtained in this way. The more efficient way is to derive the optimal rate of convergence from an estimate like (1.9) for continuous dependence with B properly chosen, without the requirement of bounded variation of v. Indeed, in this paper we apply the second way to establish the optimal rate of convergence for the vanishing anisotropic viscosity method for (1.1).

While the vanishing anisotropic viscosity method has received almost no attention in the literature, the vanishing isotropic viscosity method for the purely hyperbolic case (A≡0) is well-studied [3, 10, 19, 20, 27]. After our main results were finished, we noticed a preprint by Makridakis and Perthame [23], whose main result is the optimal rate of convergence for the vanishing anisotropic viscosity method for the hyperbolic problem, with the aid of the kinetic approach in Chen-Perthame [9] for (1.1) and an estimate technique via an auxiliary parabolic equation with constant diffusion. As we can see from our previous discussion, their result can also be obtained directly from our general L¹–framework (see Section 5 for the details). One motivation for studying the vanishing anisotropic viscosity method is that anisotropic viscosity approximations are closely related to finite volume numerical schemes on unstructured grids, for which uniformBV bounds are not available for finite volume schemes, and the standard error estimate theory for hyperbolic problems provides only a sub-optimal rate of convergence.

Although the significant applications of our L¹–framework are the estimate for continuous dependence on the nonlinearities and the error estimate for the vanishing anisotropic viscosity method, as an example of direct applications of this framework to numerical methods, we focus in Section 6 on a linear convection-diffusion model equation and derive an L¹ error estimate for a upwind-central difference scheme. We will present further applications of our L¹–framework to numerical methods for nonlinear degenerate parabolic-hyperbolic equations elsewhere. Also we remark that the results in this paper can be extended to more general equations with (t, x)–dependent coefficients; the details will be presented elsewhere.

This paper is organized as follows. We first establish the L¹–framework for continuous dependence and error estimates in Sections 2 and 3. Then we apply our generalL¹–framework to obtain the following results: (i) an explicit estimate for continuous dependence on the nonlinearities in Section 4; (ii) an optimal error estimate for the anisotropic vanishing viscosity method in Section 5; (iii) an error estimate for an upwind-central finite difference scheme for a linear convection-diffusion equation in Section 6.

2. Entropy Solutions and Kinetic Formulation For any entropy functionη :R→R, the corresponding entropy fluxes

q = (q₁, . . . , q_d) :R→R^d and R= (r_ij) :R→R^d×d are defined by

q⁰(u) =η⁰(u)f⁰(u), R⁰(u) =η⁰(u)A(u).

We will refer to (η, q, R) as an entropy-entropy flux triple.

For i= 1, . . . , d and k= 1, . . . , K, we let ζ_ik^A(u) =

Z u 0

σ_ik^A(w)dw

and

ζ_ik^A,ψ(u) = Z u

0

pψ(w)σ_ik^A(w)dw, for ψ∈C0(R).

According to Chen-Perthame [9], entropy solutions can now be defined as follows.

Definition 2.1 (Entropy Solutions). A function u∈L^∞(R₊;L¹(R^d))∩L^∞(R₊×R^d) is an entropy solution of the Cauchy problem (1.1)if the following conditions are satisfied:

(D.1) For any k= 1, . . . , K,

d

X

i=1

∂_xiζ_ik^A(u)∈L²(R₊×R^d).

(D.2) For any k= 1, . . . , K and ψ∈C0(R) withψ≥0,

d

X

i=1

∂_xiζ_ik^A,ψ(u) =p ψ(u)

d

X

i=1

∂_xiζ_ik^A(u)∈L²(R₊×R^d), and the parabolic dissipation measure n^u,ψ(t, x), defined by

n^u,ψ(t, x) =ψ(u(t, x))

K

X

k=1 d

X

i=1

∂xiζ_ik^A(u(t, x))

!2

, satisfies

n^u,ψ(t, x) =

K

X

k=1 d

X

i=1

∂xiζ_ik^A,ψ(u(t, x))

!2

a.e. in R+×R^d.

(D.3) There exists an entropy dissipation measure m^u,ψ(t, x) of the form m^u,ψ(t, x) =

Z

R

m^u(ξ, t, x)ψ(ξ)dξ, for anyψ∈C₀(R),

for some nonnegative entropy defect measurem^u(ξ, t, x)such that, for anyC²entropy- entropy flux triple (η, q, R) with η⁰⁰ ∈C₀(R), there holds

(2.1) ∂_tη(u) + divq(u)− ∇ · R⁰(u)∇u

=−

m^u,η00+n^u,η00

in D⁰(R₊×R^d), with initial data η(u)|_t=0=η(u0). That is, for any test functionφ∈ D(R₊×R^d),

Z

R+×R^d



η(u)∂_tφ+

d

X

i=1

q_i(u)∂_xiφ+

d

X

i,j=1

r_ij(u)∂_x²_ixjφ



dt dx +

Z

R^d

η(u0(x))φ(0, x)dx= Z

R+×R^d

m^u,η00+n^u,η00

φ dt dx.

Remark 2.1. The nonnegative parabolic defect measuren^u(ξ, t, x) can be defined as (2.2) n^u(ξ, t, x) =δ(ξ−u(t, x))

K

X

k=1 d

X

i=1

ζ_ik^A(u(t, x))

!² . Using the duality C0(R);M(R)

, the parabolic dissipation measuren^u,ψ then takes the form n^u,ψ(t, x) =

Z

R

n^u(ξ, t, x)ψ(ξ)dξ, ψ∈C0(R).

In the “diagonal case”a_ij ≡0 for alli6=j, the chain rule (D.2) is automatically satisfied.

We also follow Chen-Perthame [9] to give the equivalent kinetic formulation of entropy solutions for (1.1) which can be derived essentially from duality and the representation formula

(2.3) h(u)−h(0) =

Z

R

h⁰(ξ)χ(ξ;u)dξ, for any h∈C¹, where the indicator functionχ(ξ;u) is defined by

χ(ξ;u) =







10<ξ<u, whenu >0,

0, whenu= 0,

−1u<ξ<0, whenu <0;

see also Lions-Perthame-Tadmor [21].

For later use, we note that the following formulas are valid:

(2.4) ∂uχ(ξ;u) =δ(ξ−u), ∂_ξχ(ξ;u) =δ(ξ)−δ(ξ−u).

Definition 2.2 (Kinetic Formulation). Let u be an entropy solution of (1.1) in the sense of Definition 2.1. Then the kinetic formulation of (1.1)reads

(2.5)

∂_tχ(ξ;u) +f⁰(ξ)· ∇_xχ(ξ;u)

=

d

X

i,j=1

a_ij(ξ)∂_x²_ixjχ(ξ;u) +∂_ξ(m^u+n^u) (ξ, t, x) in D⁰_ξ,t,x, χ(ξ;u)|_t=0=χ(ξ;u0),

for some nonnegative entropy defect measurem^u, which measures “hyperbolicity” in the solu- tion, and some nonnegative parabolic defect measuren^u with the form (2.2), which measures

“parbolicity” in the solution.

3. General L¹–Framework

Let u be an entropy solution of the original problem (1.1). Let g be the flux function defined in (1.7) andB = (bij) be thed×dsymmetric matrix defined in (1.8). We then letv solve the “approximate” kinetic problem

(3.1)

∂tχ(ξ;v) +g⁰(ξ)· ∇_xχ(ξ;v)

=

d

X

i,j=1

b_ij(ξ)∂²_xixjχ(ξ;v) +∂_ξ(m^v+n^v+E) (ξ, t, x) inD_ξ,t,x⁰ , χ(ξ;v)|_t=0 =χ(ξ;v0),

for some nonnegative entropy defect measure m^v and nonnegative parabolic defect measure n^v taking the particular form

n^v(ξ, t, x) =δ(ξ−v(t, x))

K

X

k=1 d

X

i=1

ζ_ik^B(v(t, x))

!2

, with

ζ_ik^B(v) = Z v

0

σ^B_ik(w)dw.

Correspondingly, forψ∈C₀(R), define the function ζ_ik^B,ψ(v) =

Z v 0

pψ(w)σ_ik^B(w)dw.

Motivated by Bouchut-Perthame [3] and Perthame [27] in their treatment of the hyperbolic problem, we assume that the error termE(ξ, t, x) takes the form of “partial derivatives”:

E(ξ, t, x) =∂_t^J0e0(ξ, t, x) + X

J=(J1,...,Jd) 0≤|J|≤J∗

D^J_xe^J₁(ξ, t, x)

for some error termse0 and e^J₁ withJ0, J∗ ≥0 integers and J multi-indices. We assume that the error termse0 ande^J₁ satisfy

sup

ξ

|e₀(ξ,·,·)|, sup

ξ

e^J₁(ξ,·,·)

!

L¹_loc(R+;L¹(R^d))

<∞, 0≤ |J| ≤J∗, where sup is taken over all

ξ∈I(v) := [infv,supv]. Define the d×dsymmetric matrix

S(ξ) =p

A(ξ)−p

B(ξ) p

A(ξ)−p B(ξ)

>

= σ^A(ξ)−σ^B(ξ)

σ^A(ξ)−σ^B(ξ)>

. Then the entries ofS(ξ) = (sij(ξ)) take the form:

(3.2) s_ij(ξ) =

K

X

k=1

σ_ik^A(ξ)σ_jk^A(ξ)−2σ_ik^A(ξ)σ_jk^B(ξ) +σ_ik^B(ξ)σ_jk^B(ξ)

. To state the following theorem, we use the notations:

S∞:=kSk_∞= sup

ξ∈I(u₀)

kS(ξ)k_∞= sup

ξ∈I(u₀) i,j=1,...,d

|s_ij(ξ)|,

and

f⁰−g⁰

_∞:= sup

ξ∈I(u0)

f⁰(ξ)−g⁰(ξ)

_∞= sup

ξ∈I(u0) i=1,...,d

f_i⁰(ξ)−g_i⁰(ξ) .

Hereafter, C will denote positive constants, not necessarily the same at different occur- rences, which are independent of the small parameters and time variablet.

The main result of this section is the following abstractL¹–framework for error estimates.

Theorem 3.1 (General L¹–Framework). Let u ∈ C(R₊;L¹(R^d)) be an entropy solution of (1.1), and suppose v ∈ L^∞(R₊;L¹(R^d))∩L^∞(R₊×R^d)∩C(R₊;L¹(R^d)) solves the

“approximate” kinetic problem (3.1). Then, for any t >0 and any ε1,ε˜0,ε˜1>0,

(3.3)

ku(t,·)−v(t,·)k_L1(R^d)≤ ku₀−v₀k_L1(R^d)

+C E_u,t^x (ε₁) +E_v,t^t (˜ε₀) +E_v,t^x (˜ε₁) +E_t^f−g(ε₁) +E_t^A−B(ε₁) +E_v,t(˜ε₀,ε˜₁)

! , where

E_u,t^x (ε₁) = sup

|y|<ε1

τ=0,t

ku(τ,·+y)−u(τ,·)k_L1(R^d), E_v,t^t (˜ε0) = sup

0<s−τ <˜ε0

τ=0,t

kv(s,·)−v(τ,·)k_L1(R^d), E_v,t^x (˜ε₁) = sup

|y|<˜ε1

τ=0,t

kv(τ,·+y)−v(τ,·)k_L1(R^d),

E_t^f−g(ε₁) =

(ku₀k_L1(R^d)

tkf⁰−g⁰k_∞

ε1 , u0 ∈/ BV(R^d), ku₀k_BV(Rd)tkf⁰−g⁰k_∞, u₀ ∈BV(R^d), E_t^A−B(ε1) =

(ku₀k_L1(R^d) tS∞

ε²₁ , u0∈/ BV(R^d), ku₀k_BV(Rd)tS∞

ε1 , u0∈BV(R^d), and

E_v,T(˜ε0,ε˜1) = 1

˜ ε^J₀⁰

sup

ξ

|e₀(ξ,·,·)|

L¹(0,t+˜ε0;L¹(R^d))

+ X

J=(J1,...,Jd) 0≤|J|≤J∗

1

˜ ε^|J₁ ^|

sup

ξ

e^J₁(ξ,·,·)

L¹(0,t+˜ε0;L¹(R^d))

.

If g ≡ f and B ≡ A, then the terms E_u,t^x (ε₁), E_t^f−g(ε₁), and E_t^A−B(ε₁) in (3.3) can be dropped, that is, there holds

(3.4)

ku(t,·)−v(t,·)k_L1(R^d)≤ ku₀−v₀k_L1(R^d)

+C

E_v,t^t (˜ε₀) +E_v,t^x (˜ε₁) +E_v,t(˜ε₀,ε˜₁) , for anyt >0 and ε˜0,ε˜1 >0.

Proof. Some arguments in this proof follow Chen-Perthame [9] closely, for which we are very concise here and refer instead to [9] for more details.

We set ε = (ε₀, ε₁), ε₀ > 0 for the forward time regularization and ε₁ > 0 for the space regularization. We then define

ω_ε(t, x) :=ω_ε0(t)ω_ε1(x),

where

ω_ε0(t) := 1 ε0

ω₀ t

ε0

, ω_ε1(x) := 1 ε^d₁ω₁

x₁ ε1

· · ·ω₁ x_d

ε1

, andω<sub>`</sub>≥0,`= 0,1, denote the normalized regularization kernels with

Z

R

ω`(τ)dτ = 1, supp(ω0)⊂(−1,0), supp(ω1)⊂(−1,1).

We use the notations

χ:=χ(ξ, t, x) =χ(ξ;u(t, x)), χ˜:= ˜χ(ξ, t, x) =χ(ξ;v(t, x)), χ_ε:=χ_ε(ξ, t, x) =

χ ?

(t,x)ω_ε

(ξ, t, x), χ˜_ε˜:= ˜χ_ε˜(ξ, t, x) =

˜ χ ?

(t,x)ω_ε˜

(ξ, t, x), where ˜ε= (˜ε₀,ε˜₁)>0 is another pair of time-space regularization parameters. Moreover, we use the notations

m^u_ε :=m^u_ε(ξ, t, x) =

m^u ?

(t,x)ω_ε

(ξ, t, x), n^u_ε :=n^u_ε(ξ, t, x) =

n^u ?

(t,x)ω_ε

(ξ, t, x), m^v_˜ε :=m^v_ε˜(ξ, t, x) =

m^v ?

(t,x)ωε˜

(ξ, t, x), n^v_ε˜:=n^v_ε˜(ξ, t, x) =

n^v ?

(t,x)ωε˜

(ξ, t, x), andE_ε˜=E_ε˜(ξ, t, x), which is similarly defined.

We intend to study the microscopic functional

0≤Qε,ε˜(ξ, t, x) =|χ_ε|+|˜χε˜| −2χεχ˜ε˜. More precisely, we will calculate

(3.5) d

dt Z

Rξ×R^d_x

Q_ε,ε˜(ξ, t, x)dx dξ.

Note thatχ_ε(ξ, t, x) satisfies

(3.6) ∂tχε+f⁰(ξ)· ∇_xχε=

d

X

i,j=1

aij(ξ)∂_x²_ixjχε+∂ξ(m^u_ε+n^u_ε), and that ˜χ_ε˜(ξ, t, x) satisfies

(3.7) ∂tχ˜ε˜+g⁰(ξ)· ∇_xχ˜˜ε=

d

X

i,j=1

bij(ξ)∂_x²_ixjχ˜ε˜+∂_ξ(m^v_ε˜+n^v_˜ε+Eε˜).

Multiplying (3.6) by sign(ξ), using sign(ξ)χε=|χ_ε|, and then integrating in (ξ, x)∈Rξ×R^d_x yield

(3.8) d

dt Z

Rξ×R^d_x

|χ_ε|dx dξ =−2 Z

R^d_x

(m^u_ε +n^u_ε) (0, t, x)dx.

Similarly,

(3.9) d

dt Z

Rξ×R^d_x

|˜χ_ε˜|dx dξ =−2 Z

R^d_x

(m^v_ε˜+n^v_ε˜+E_ε˜) (0, t, x)dx.

We now consider the quadratic term. To this end, we need an additional regularization in the kinetic/velocity variableξ:

χ_ε,δ(ξ, t, x) :=

χε?

ξψ_δ

(ξ, t, x), χ˜_ε,δ˜ (ξ, t, x) :=

˜ χ˜ε?

ξψ_δ

(ξ, t, x),

for a standard regularization kernelψ_δ. We also need aξ-truncationT_L(ξ), which is a smooth nonnegative function with bounded support. That is, T_L(ξ) =T(ξ/L)→1 asL→ ∞ with

0≤ T(ξ)≤1, forξ ∈(−∞,∞), T(ξ) = 1, for |ξ| ≤1/2, T(ξ) = 0, for |ξ| ≥1.

The destiny of these additional parameters is thatδ↓0 first andL↑ ∞ second.

Then χ_ε,δ satisfies

(3.10)

∂tχ_ε,δ+f⁰(ξ)· ∇_xχ_ε,δ

=

d

X

i,j=1

∂_x²_ixj

(a_ijχ_ε)?

ξψ_δ +∂_ξ

(m^u_ε +n^u_ε)?

ξψ_δ

+R^u_ε,δ,

and ˜χ_˜ε,δ satisfies

(3.11)

∂tχ˜ε,δ˜ +g⁰(ξ)· ∇_xχ˜ε,δ˜

=

d

X

i,j=1

∂_x²_ixj

(b_ijχ˜_ε˜)?

ξψ_δ +∂_ξ

(m^v_ε˜+n^v_ε˜+E_˜ε)?

ξψ_δ

+R^v_ε,δ.

In (3.10) and (3.11),

R^u_ε,δ = div_x

f⁰(ξ)χ_ε,δ− f⁰χ_ε

?ξψ_δ

, R^v_ε,δ= div_x

g⁰(ξ) ˜χ_ε,δ˜ − g⁰χ˜_ε˜

?ξψ_δ .

A simple calculation reveals d

dt Z

Rξ×R^d_x

T_L(ξ)χ_ε,δχ˜_˜ε,δdx dξ (3.12)

=− Z

Rξ×R^d_x

T_L(ξ) ˜χ_ε,δ˜ f⁰(ξ)−g⁰(ξ)

· ∇_xχ_εdx dξ

+ Z

Rξ×R^d_x

T_L(ξ) ˜χ_˜ε,δ

d

X

i,j=1

∂_x²_ixj

(aijχε)?

ξψ_δ

dx dξ +

Z

R_ξ×R^d_x

T_L(ξ) ˜χ_˜ε,δ∂_ξ

(m^u_ε+n^u_ε)?

ξψ_δ

dx dξ +

Z

Rξ×R^d_x

T_L(ξ)χ_ε,δ

d

X

i,j=1

∂_x²_ixj

(bijχ˜ε˜)?

ξψ_δ

dx dξ +

Z

R_ξ×R^d_x

T_L(ξ)χε,δ∂ξ

(m^v_ε˜+n^v_ε˜+Eε˜)?

ξψδ

dx dξ +

Z

Rξ×R^d_x

T_L(ξ) ˜χ_ε,δ˜ R_ε,δ^u +χ_ε,δR_ε,δ^v dx dξ

=:

6

X

`=1

I`(t;ε,ε, δ, L).˜

As in Chen-Perthame [9], we have

(3.13) lim

δ↓0I6(t;ε, δ, L) = 0 inL^p(0, T) for any 1≤p <∞.

Writing out the convolution products explicitly, we have I1(t;ε,ε, δ, L)˜

=− Z

T_L(ξ)(f⁰(ξ)−g⁰(ξ))· ∇_xωε(t−s, x−y)ωε˜(t−s⁰, x−y⁰)

×ψδ(ξ−η)ψδ(ξ−η⁰)χ(η;u(s, y))χ(η⁰;v(s⁰, y⁰))ds dy dη ds⁰dy⁰dη⁰dx dξ.

Sending first δ ↓0 and second L↑ ∞, we get

(3.14)

L↑∞limlim

δ↓0I₁(t;ε,ε, δ, L)˜

=− Z

(f⁰(ξ)−g⁰(ξ))· ∇_xωε(t−s, x−y)ω˜ε(t−s⁰, x−y⁰)

×χ(ξ;u(s, y))χ(ξ;v(s⁰, y⁰))ds dy ds⁰dy⁰dx dξ

=:−E^f−g(t;ε,ε).˜

Integrating by parts yields I₃(t;ε,ε, δ, L) =˜ −

Z

Rξ×R^d_x

T_L⁰(ξ) ˜χ_ε,δ˜

(m^u_ε +n^u_ε)?

ξψ_δ

dx dξ (3.15)

− Z

Rξ×R^d_x

T_L(ξ)ψ_δ(ξ)

(m^u_ε+n^u_ε)?

ξψ_δ

dx dξ +

Z

Rξ×R^d_x

T_L(ξ)

δ(ξ−u) ?

(ξ,t,x)(ωεψ_δ) m^u_ε ?

ξψ_δ

dx dξ +

Z

R_ξ×R^d_x

T_L(ξ)

δ(ξ−u) ?

(ξ,t,x)(ωεψδ) n^u_ε ?

ξψδ

dx dξ

=:

4

X

`=1

I_3,`(t;ε,ε, δ, L),˜ and

I5(t;ε,ε, δ, L) =˜ − Z

Rξ×R^d_x

T_L⁰(ξ)χ_ε,δ

(m^v_ε˜+n^v_ε˜+Eε˜)?

ξψ_δ

dx dξ (3.16)

− Z

Rξ×R^d_x

T_L(ξ)ψδ(ξ)

(m^v_ε˜+n^v_˜ε+Eε˜)?

ξψδ

dx dξ +

Z

R_ξ×R^d_x

T_L(ξ)

δ(ξ−v) ?

(ξ,t,x)(ωε˜ψδ) m^v_ε˜?

ξψδ

dx dξ +

Z

Rξ×R^d_x

T_L(ξ)

δ(ξ−v) ?

(ξ,t,x)(ω_ε˜ψ_δ) n^v_ε˜?

ξψ_δ

dx dξ +

Z

Rξ×R^d_x

T_L(ξ)

δ(ξ−v) ?

(ξ,t,x)(ω_ε˜ψ_δ) E_˜ε?

ξψ_δ

dx dξ

=:

5

X

`=1

I_5,`(t;ε,ε, δ, L).˜ As in [9], we have

(3.17) lim

L↑∞lim

δ↓0I3,1(t;ε,ε, δ, L) = lim˜

L↑∞lim

δ↓0I5,1(t;ε,ε, δ, L) = 0,˜ in L^p(0, T) for any 1≤p <∞.

Moreover,

limδ↓0

I_3,2(t;ε,ε, δ, R) +˜ I_5,2(t;ε,ε, δ, R)˜ (3.18)

=− Z

R^d_x

(m^u_ε+n^u_ε+m^v_ε˜+n^v_ε˜+E˜ε) (0, t, x)dx, in L^p(0, T) for any 1≤p <∞.

Clearly,

I_3,3(t;ε,ε, δ, L), I˜ _5,3(t;ε,ε, δ, L)˜ ≥0 for any t >0.

Integrating by parts also yields I₂(ε, δ, L) =−

Z

Rξ×R^d_x

T_L(ξ)

K

X

k=1 d

X

i,j=1

∂_xi

σ_ik^Aχ˜_˜ε

?ξψ_δ

∂_xj

σ_jk^Aχ_ε

?ξψ_δ dx dξ (3.19)

+ Z

Rξ×R^d_x

T_L(ξ)

K

X

k=1 d

X

i,j=1

(

∂xi

σ_ik^Aχ˜ε˜

?

ξψδ

∂xj

σ_jk^Aχε

?

ξψδ

−∂_xiχ˜_˜ε,δ∂_xj

(a_ijχ_ε)?

ξψ_δ )

dx dξ

=:I_2,1(t;ε,ε, δ, L) +˜ I_2,2(t;ε,ε, δ, L).˜ Similarly,

(3.20) I4(t;ε,ε, δ, L) =˜ I4,1(t;ε,ε, δ, L) +˜ I4,2(t;ε,ε, δ, L),˜ where

I_4,1(t;ε,ε, δ, L)˜ (3.21)

=− Z

Rξ×R^d_x

T_L(ξ)

K

X

k=1 d

X

i,j=1

∂_xi

σ_ik^Bχ˜_ε˜

?ξψ_δ

∂_xj

σ_jk^Bχ_ε

?ξψ_δ dx dξ and

I_4,2(t;ε,ε, δ, L) =˜ Z

Rξ×R^d_x

T_L(ξ)

K

X

k=1 d

X

i,j=1

(

∂_xi

σ^B_ikχ˜_ε˜

?ξψ_δ

∂_xj

σ^B_jkχ_ε

?ξψ_δ (3.22)

−∂xi

(bijχ˜ε˜)?

ξψ_δ

∂xjχ_ε,δ )

dx dξ.

As in Chen-Perthame [9], we have

(3.23) lim

δ↓0I_2,2(t;ε,ε, δ, L) = lim˜

δ↓0I_4,2(t;ε,ε, δ, L) = 0˜ in L^p(0, T) for any 1≤p <∞.

We now study the new term E^A−B(t;ε,ε, δ, L)˜ (3.24)

:=I2,1(t;ε,ε, δ, L) +˜ I3,4(t;ε,ε, δ, L) +˜ I4,1(t;ε,ε, δ, L) +˜ I5,4(t;ε,ε, δ, L).˜ From now on in this proof, for notational simplicity, we drop writing the domains of integration. Writing out explicitly the convolution products, we have

I_2,1(t;ε,ε, δ, L)˜

=−

K

X

k=1 d

X

i,j=1

Z

T_L(ξ)∂xiωε(t−s, x−y)∂xjωε˜(t−s⁰, x−y⁰)ψδ(ξ−η)ψδ(ξ−η⁰)

×σ^A_ik(η)χ(η;u(s, y))σ_jk^A(η⁰)χ(η⁰;v(s⁰, y⁰))ds dy dη ds⁰dy⁰dη⁰dx dξ.

Similarly,

I_4,1(t;ε,ε, δ, L)˜

=−

K

X

k=1 d

X

i,j=1

Z

T_L(ξ)∂_xiω_ε(t−s, x−y)∂_xjω_ε˜(t−s⁰, x−y⁰)ψ_δ(ξ−η)ψ_δ(ξ−η⁰)

×σ^B_ik(η)χ(η;u(s, y))σ_jk^B(η⁰)χ(η⁰;v(s⁰, y⁰))ds dy dη ds⁰dy⁰dη⁰dx dξ.

Note that

I3,4(t;ε,ε, δ, L) +˜ I5,4(t;ε,ε, δ, L)˜

=

K

X

k=1

Z

T_L(ξ)ωε(t−s, x−y)ωε˜(t−s⁰, x−y⁰)ψ_δ(ξ−u(s, y))ψ_δ(ξ−v(s⁰, y⁰))

×





d

X

i=1

∂_yiζ_ik^A(u(s, y))

!2

+





d

X

j=1

∂_y⁰

jζ_jk^B(v(s⁰, y⁰))





2

ds dy ds⁰dy⁰dx dξ

≥2

K

X

k=1

Z

T_L(ξ)ω_ε(t−s, x−y)ω_˜ε(t−s⁰, x−y⁰)ψ_δ(ξ−u(s, y))ψ_δ(ξ−v(s⁰, y⁰))

×

d

X

i=1

∂yiζ_ik^A(u(s, y))

d

X

j=1

∂_y⁰

jζ_jk^B(v(s⁰, y⁰))ds dy ds⁰dy⁰dx dξ

= 2

K

X

k=1

Z

T_L(ξ)ω_ε(t−s, x−y)ω_˜ε(t−s⁰, x−y⁰)

×

d

X

i=1

∂_yiζ_ik^{A,ψδ(ξ−·)2}(u(s, y))

d

X

j=1

∂_y⁰

jζ_jk^{B,ψδ(ξ−·)2}(v(s⁰, y⁰))ds dy ds⁰dy⁰dx dξ

= 2

K

X

k=1 d

X

i,j=1

Z

T_L(ξ)∂xiωε(t−s, x−y)∂xjωε˜(t−s⁰, x−y⁰)ψδ(ξ−η)ψδ(ξ−η⁰)

×σ^A_ik(η)χ(η;u(s, y))σ_jk^B(η⁰)χ(η⁰;v(s⁰, y⁰))ds dy dη ds⁰dy⁰dη⁰dx dξ,

where we have used the chain rule (D.2), integration by parts, and (2.3). From this and the previous calculations, we find

E^A−B(t;ε,ε, δ, L)˜

≥ −

K

X

k=1 d

X

i,j=1

Z

T_L(ξ)∂xiωε(t−s, x−y)∂xjωε˜(t−s⁰, x−y⁰)ψ_δ(ξ−η)ψ_δ(ξ−η⁰)

×

σ_ik^A(η)σ^A_jk(η⁰)−2σ_ik^A(η)σ_jk^B(η⁰) +σ^B_ik(η)σ_jk^B(η⁰)

×χ(η;u(s, y))χ(η⁰;v(s⁰, y⁰))ds dy dη ds⁰dy⁰dη⁰dx dξ.

After performing the changes of variables:

z=ξ−η, z⁰ =ξ−η⁰, dz dz⁰ = dη dη⁰,

we get

E^A−B(t;ε,ε, δ, L)˜

≥ −

K

X

k=1 d

X

i,j=1

Z

T_L(ξ)∂xiωε(t−s, x−y)∂xjω˜ε(t−s⁰, x−y⁰)ψδ(z)ψδ(z⁰)

×

σ_ik^A(ξ−z)σ_jk^A(ξ−z⁰)−2σ^A_ik(ξ−z)σ_jk^B(ξ−z⁰) +σ_ik^B(ξ−z)σ_jk^B(ξ−z⁰)

×χ(ξ−z;u(s, y))χ(ξ−z⁰;v(s⁰, y⁰))ds dy dz ds⁰dy⁰dz⁰dx dξ.

From this, it easily follows that

L↑∞limlim

δ↓0E^A−B(t;ε,ε, δ, L)˜ (3.25)

≥ −

d

X

i,j=1

Z

∂_x²_ixjω_ε(t−s, x−y)ω_ε˜(t−s⁰, x−y⁰)

×sij(ξ)χ(ξ;u(s, y))χ(ξ;v(s⁰, y⁰))ds dy ds⁰dy⁰dx dξ

=:−E^A−B(t;ε,ε),˜

where we have also performed integration by parts and used the definition of s_ij(ξ) in (3.2).

Writing out the convolution products, we have I_5,5(t;ε,ε, δ, L)˜

(3.26)

= Z

T_L(ξ)ψδ(ξ−v(s, y))ψδ(ξ−η⁰)ωε˜(t−s, x−y)

× ∂_t^J0ωε˜(t−s⁰, x−y⁰)e0(η⁰, s⁰, y⁰)

+ X

J=(J1,...,Jd)≥0

|J|≤J∗

D_x^Jω_˜ε(t−s⁰, x−y⁰)e^J₁(η⁰, s⁰, y⁰)

!

ds dy dη ds⁰dy⁰dη⁰dx dξ

−→

Z

ωε˜(t−s, x−y) ∂_t^J0ωε˜(t−s⁰, x−y⁰)e0(v(s, y), s⁰, y⁰)

+ X

J=(J1,...,Jd)≥0

|J|≤J∗

D^J_xω_ε˜(t−s⁰, x−y⁰)e^J₁(v(s, y), s⁰, y⁰)

!

ds dy ds⁰dy⁰dx

=:−E_v(t; ˜ε), when L↑ ∞and δ ↓0.

Summarizing our calculations from (3.12) to (3.26), we obtain that, for any ε,ε >˜ 0, d

dt Z

R_ξ×R^d_x

χ_εχ˜_ε˜dx dξ ≥ − Z

R^d_x

(m^u_ε +n^u_ε+m^v_ε˜+n^v_ε˜+E_ε˜) (0, t, x)dx (3.27)

+E^f−g(t;ε,ε) +˜ E^A−B(t;ε,ε) +˜ E_v(t; ˜ε).

Then the estimates (3.5), (3.8), (3.9), and (3.27) yield that, for any ε,ε >˜ 0, Z

Rξ×R^d_x

Q_ε,˜ε(ξ, t, x)dx dξ

≤ Z

Rξ×R^d_x

Q_ε,˜ε(ξ,0, x)dx dξ+ Z t

0

E^f−g(τ;ε,ε)˜ dτ+ Z t

0

E^A−B(t;ε,ε)˜ dτ+ Z t

0

E_v(t;ε,ε)˜ dτ.

Similarly, we have Z

Rξ×R^d_x

Qε,˜ε(ξ,0, x)dx dξ

= Z

|u(s, y)|+

v(s⁰, y⁰)

−2 min |u(s, y)|,

v(s⁰, y)

1{sign(u(s,y)v(s⁰,y⁰))>0}

×ωε(−s, x−y)ωε˜(−s⁰, x−y⁰)ds dy ds⁰dy⁰dx

= Z

u(s, y)−v(s⁰, y⁰)

ω_ε(−s, x−y)ω_ε˜(−s⁰, x−y⁰)ds dy ds⁰dy⁰dx.

A standard calculation reveals Z

Rξ×R^d_x

Qε,˜ε(ξ,0, x)dx dξ

≤ ku₀−v₀k_L1(R^d)+ sup

0<s<ε0

ku(s,·)−u₀(·)k_L1(R^d)+ sup

0<s<˜ε0

kv(s,·)−v₀(·)k_L1(R^d)

+ sup

|y|<ε₁

ku₀(·+y)−u0(·)k_L1(R^d)+ sup

|y|<˜ε1

kv₀(·+y)−v0(·)k_L1(R^d). Similarly, we find

Z

Rξ×R^d_x

Q_ε,˜ε(ξ, t, x)dx dξ

≥ ku(t)−v(t)k_L1(R^d)− sup

0<s−t<ε₀

ku(s,·)−u(t,·)k_L1(R^d)− sup

0<s−t<ε˜0

kv(s,·)−v(t,·)k_L1(R^d)

− sup

|y|≤ε1

ku(t,·+y)−u(t,·)k_L1(R^d)− sup

|y|≤˜ε1

kv(t,·+y)−v(t,·)k_L1(R^d). Sending ε0 ↓0, we conclude

ku(t,·)−v(t,·)k_L1(R^d)

≤ ku₀−v₀k_L1(R^d)+ 2E_u,t^x (ε₁) + 2E_v,t^t (˜ε₀) + 2E_v,t^x (˜ε₁) + lim

ε0↓0

Z _t

0

E^f−g(τ;ε,ε)˜ dτ + lim

ε0↓0

Z t 0

E^A−B(τ;ε,ε)˜ dτ+ Z t

0

E_v(τ; ˜ε)dτ, 0≤t≤T.

It remains to estimate the three terms on the second line.