Approximation of Hamilton-Jacobi equations with Caputo time-fractional derivative

TL;DR

This paper develops a numerical scheme for solving Hamilton-Jacobi equations involving Caputo time-fractional derivatives, demonstrating stability and convergence to the unique viscosity solution.

Contribution

It introduces an explicit time discretization method for Caputo derivatives combined with a finite difference scheme for Hamiltonians, ensuring stability and convergence.

Findings

The scheme is stable under specific discretization conditions.

The scheme converges to the unique viscosity solution.

Numerical experiments confirm theoretical results.

Abstract

In this paper, we investigate the numerical approximation of Hamilton-Jacobi equations with the Caputo time-fractional derivative. We introduce an explicit in time discretization of the Caputo derivative and a finite difference scheme for the approximation of the Hamiltonian. We show that the approximation scheme so obtained is stable under an appropriate condition on the discretization parameters and converges to the unique viscosity solution of the Hamilton-Jacobi equation.

Figures

Equations

$${\partial }_t^{\alpha }u(t,x)+H(t,x,u(t,x),Du(t,x))=0(t,x)\in Q_T:=(0,T]\times {\mathbb{T}}^d,$$

$${\partial }_t^{\alpha }u(t,x)=\frac{1}{\Gamma (1-\alpha )}{\int }_0^t\frac{{\partial }_{s}u(s,x)}{(t-s{)}^{\alpha }}ds$$

$$u(0,x)=u_0(x),x\in {\mathbb{T}}^d.$$

$$|H(t,x,r,p)-H(t,y,r,p)|\le \omega (|x-y|(1+|p|))$$

$${\partial }_t^{\alpha }f(t)=\frac{1}{\Gamma (1-\alpha )}{\int }_0^t\frac{f^{\prime }(s)}{(t-s{)}^{\alpha }}ds,$$

$${\partial }_t^{\alpha }f(t)=J[f](t)+K_{(0,t)}[f](t),$$

$$J[f](t):=\frac{f(t)-f(0)}{t^{\alpha }\Gamma (1-\alpha )},$$

$$K_{(0,t)}[f](t):=\frac{\alpha }{\Gamma (1-\alpha )}{\int }_0^t\frac{f(t)-f(t-\tau )}{{\tau }^{\alpha +1}}d\tau .$$

$$u^{\ast }(x)=\underset{r\to 0}{\lim }\sup \{u(y):y\in A\cap \overline{B(x,r)}\}$$

$$J[\phi ](\hat{t},\hat{x})+K_{(0,\hat{t})}[\phi ](\hat{t},\hat{x})+H(\hat{t},\hat{x},u^{\ast }(\hat{t},\hat{x}),D\phi (\hat{t},\hat{x}))\le 0,$$

$$\underset{[0,T]\times O}{\max }(u^{\ast }-\phi )=(u^{\ast }-\phi )(\hat{t},\hat{x}).$$

$$J[\phi ](\hat{t},\hat{x})+K_{(0,\hat{t})}[\phi ](\hat{t},\hat{x})+H(\hat{t},\hat{x},u^{\ast }(\hat{t},\hat{x}),D\phi (\hat{t},\hat{x}))\ge 0,$$

$$\underset{[0,T]\times O}{\min }(u_{\ast }-\phi )=(u_{\ast }-\phi )(\hat{t},\hat{x}).$$

$$D_{\Delta t}^{\alpha }{U}_{i,j}^{n+1}$$

$$=\frac{1}{\Gamma (1-\alpha )(1-\alpha )}\sum _{m=0}^n\frac{U_{i,j}^{m+1}-U_{i,j}^m}{\Delta t}\left(-\frac{1}{(t_{n+1}-t_{m+1}{)}^{\alpha -1}}+\frac{1}{(t_{n+1}-t_m{)}^{\alpha -1}}\right)$$

$$=\frac{1}{\Gamma (2-\alpha )}\sum _{m=0}^n\frac{(n+1-m{)}^{1-\alpha }-(n-m{)}^{1-\alpha }}{\Delta t^{\alpha }}(U_{i,j}^{m+1}-U_{i,j}^m),$$

$${\rho }_{\alpha }=\Gamma (2-\alpha ){\Delta t}^{\alpha },$$

$${\rho }_{\alpha }{D}_{\Delta t}^{\alpha }{U}_{i,j}^{n+1}=$$

$$\displaystyle=-\Big{(}(n+1)^{1-\alpha}-n^{1-\alpha}\Big{)}U_{i,j}^{0}$$

$$\displaystyle-\sum_{m=1}^{n}\Big{(}2(n+1-m)^{1-\alpha}-(n+2-m)^{1-\alpha}-(n-m)^{1-\alpha}\Big{)}U_{i,j}^{m}+U_{i,j}^{n+1}$$

$$=U_{i,j}^{n+1}-\sum _{m=0}^n{c}_m^{n+1}{U}_{i,j}^m,$$

$$c_0^{n+1}$$

$$c_m^{n+1}$$

$$D_{\Delta t}^{\alpha }{U}_{i,j}^{n+1}=\frac{1}{{\rho }_{\alpha }}\left(U_{i,j}^{n+1}-\sum _{m=0}^n{c}_m^{n+1}{U}_{i,j}^m\right).$$

$$r_{\Delta t}^{n+1}\le c_u{\Delta t}^{2-\alpha }$$

$$\frac{(n+2-m{)}^{1-\alpha }+(n-m{)}^{1-\alpha }}{2}<(n+1-m{)}^{1-\alpha }.$$

$$c_0^{n+2}-c_0^{n+1}=(n+2{)}^{1-\alpha }-2(n+1{)}^{1-\alpha }+n^{1-\alpha }=-c_1^{n+2}.$$

$$c_{m+1}^{n+2}$$

$$=2(n+1-m{)}^{1-\alpha }-(n+2-m{)}^{1-\alpha }-(n-m{)}^{1-\alpha }=c_m^{n+1}.$$

$$\sum _{m=0}^n$$

$$=(n+1{)}^{1-\alpha }-n^{1-\alpha }+\sum _{m=1}^{n}2(n+1-m{)}^{1-\alpha }-\sum _{m=0}^{n-1}(n+1-m{)}^{1-\alpha }-\sum _{m=2}^{n+1}(n+1-m{)}^{1-\alpha }$$

$$=(n+1{)}^{1-\alpha }-n^{1-\alpha }+2n^{1-\alpha }+2-(n+1{)}^{1-\alpha }-n^{1-\alpha }-1=1.$$

$$(D_1^{+}U{)}_{i,j}=\frac{U_{i+1,j}-U_{i,j}}{h}and(D_2^{+}U{)}_{i,j}=\frac{U_{i,j+1}-U_{i,j}}{h},$$

$$[D_{h}U{]}_{i,j}={\left((D_1^{+}U{)}_{i,j},(D_1^{+}U{)}_{i-1,j},(D_2^{+}U{)}_{i,j},(D_2^{+}U{)}_{i,j-1}\right)}^T.$$

$$g(x,q_1,q_1,q_2,q_2)=H(x,q),\forall x\in {\mathbb{T}}^2,\forall q=(q_1,q_2)\in {\mathbb{R}}^2.$$

Full text

AMS subject classification:

35R11, 65L12, 49L25.

Keywords:

Fractional Hamilton-Jacobi equation, Caputo time derivative, finite difference, convergence.

1 Introduction

We define a class of finite difference schemes for the time-fractional Hamilton-Jacobi equation

[TABLE]

where $\mathbb{T}^{d}$ is unit torus in $\mathbb{R}^{d}$. The symbol $\partial_{t}^{\alpha}$, for $0<\alpha\leq 1$, denotes the Caputo derivative

[TABLE]

(note that $\partial_{t}^{\alpha}$ reduces to the standard time derivative $\partial_{t}$ for $\alpha=1$) and $Du$ the gradient with respect to $x$. Equation (1.1) is completed with the initial condition

[TABLE]

In recent times, there has been an increasing interest in the study of differential equations with time-fractional derivatives. Indeed, this kind of differential operators allows us to introduce new phenomena in differential models such as memory and trapping effects [13, 14, 17, 20]. Also, the numerical approximation of differential equations with fractional time-derivative has been extensively analyzed [3, 11, 12].

Since in general smooth solutions to Hamilton-Jacobi equations are not expected to exist, for equation (1.1) a theory of weak solutions, in viscosity sense, has been introduced, in [9, 15, 21] (see also [19] for a different, but closely related, notion of weak solution). Most of the results and techniques which hold in the classical case, i.e. for $\alpha=1$, have been extended to the fractional case in order to prove the well-posedness of the Hamilton-Jacobi equation (1.1).

In the classical case, one of the most important properties of the viscosity solution theory is the stability with respect to the uniform convergence (see [2]). Starting with the seminal paper [6], this property has generated an enormous literature concerning the numerical approximation of Hamilton-Jacobi equations (see for example [7, 16, 18] and reference therein). Stability with respect to the uniform convergence is inherited by viscosity solutions of the Hamilton-Jacobi equation (1.1). Following [6], we define a general class of finite difference schemes for (1.1). We show that, under an appropriate Courant–Friedrichs–Lewy (CFL) condition of the type $\Delta t^{\alpha}=O(\Delta x)$, these schemes are monotone, stable and consistent. Moreover, relying on an adaptation of the classical Barles-Souganidis convergence Theorem [4], we prove that the numerical solutions generated by these schemes converge to the unique viscosity solution of the limit problem. In order to verify the properties of the proposed schemes, we perform several numerical tests, and to analyze the order of the approximation error, we also compute exact solutions for some time-fractional Hamilton-Jacobi equations. We have only recently become aware that a similar problem was considered in [8].

The rest of the paper is organized as follows. In Section 2, we shortly review some basic properties of the theory of viscosity solution for (1.1). Section 3 is devoted to the description of a class of finite difference schemes and their properties. In Section 4, we prove a convergence result and in Section 5 we carry out some numerical tests.

2 Viscosity solutions for Hamilton-Jacobi equation with time-fractional derivative

In this section, we briefly review definitions and some results for the continuous problem (1.1) (we refer to [9, 15] for more details).

Throughout the paper, a function $u$ on $\mathbb{T}^{d}$ will be equivalently regarded as a function defined on $\mathbb{R}^{d}$ which is $\mathbb{Z}^{d}$ periodic. We consider the following assumptions on the Hamiltonian $H$ and on the initial datum $u_{0}$.

1. (H1)

$H:\overline{Q_{T}}\times\mathbb{R}\times\mathbb{R}^{d}\to\mathbb{R}$ is continuous; 2. (H2)

there exists a modulus $\omega:[0,\infty)\to[0,\infty)$ such that

[TABLE]

for all $(t,x,r,p),\,(t,y,r,p)\in[0,T]\times\mathbb{T}^{d}\times\mathbb{R}\times\mathbb{R}^{d}$; 3. (H3)

$r\mapsto H(t,x,r,p)$ is nondecreasing for all $(t,x,p)\in{\overline{Q}_{T}}\times\mathbb{R}^{d}$; 4. (H4)

$u_{0}:\mathbb{T}^{d}\to\mathbb{R}$ is a continuous function.

For a function $f:[0,T]\to\mathbb{R}$ such that $f\in C^{1}((0,T])\cap C([0,T])$ and $f^{\prime}\in L^{1}((0,T))$, the Caputo time fractional derivative is defined by

[TABLE]

for any $t\in(0,T]$. Using integration by parts and change of variables, (2.1) can be rewritten as

[TABLE]

where

[TABLE]

By natural extension, we also define $K_{(a,b)}[f](t)$ for any $a,b$ with $0\leq a<b\leq t$. The advantage of rewriting the Caputo derivative in the form (2.2) is explained in [1, 9, 21].

For a set $A\subset\mathbb{T}^{d}$ and a function $u:A\to\mathbb{R}$, we denote by $u^{*}$ and $u_{*}$ the upper and the lower semi-continuous envelopes of $u$, i.e.

[TABLE]

and $u_{*}(x)=-(-u)^{*}$, where $B(x,r)$ is a the open ball of centre $x$ and radius $r$. We also denote by $USC(\overline{Q_{T}})$ (resp., $LSC(\overline{Q_{T}})$) the class of the upper semi-continuous (resp., lower semi-continuous) functions in $\overline{Q_{T}}$.

We give the definition of viscosity solution for (1.1) (see [15, Definiton 2.2]).

Definition 2.1**.**

• Let $O\subset\mathbb{T}^{d}$. Then
• (i)

A function $u:[0,T]\times O\to\mathbb{R}$ is said a viscosity subsolution of (1.1) in $(a,T]\times O$ if $u^{*}<+\infty$ in $(a,T]\times O$ and

[TABLE]

whenever $(\hat{t},\hat{x})\in(a,T]\times O$ and $\varphi\in C^{1,1}((a,T]\times O)\cap C([0,T]\times O)$ satisfy

[TABLE]

• (ii)

A function $u:[0,T]\times O\to\mathbb{R}$ is said a viscosity supersolution of (1.1) in $(a,T]\times O$ if $u_{*}>-\infty$ in $(a,T]\times O$ and

[TABLE]

whenever $(\hat{t},\hat{x})\in(a,T]\times O$ and $\varphi\in C^{1,1}((a,T]\times O)\cap C([0,T]\times O)$ satisfy

[TABLE]

• (iii)

A function $u:[0,T]\times O\to\mathbb{R}$ is said a viscosity solution of (1.1) in $(a,T]\times O$ if it is both a viscosity sub- and supersolution of (1.1) in $(a,T]\times O$.

In the previous definition, the notation $C^{1,1}((a,T]\times O)$ denotes the space of functions $\varphi$ such that $\varphi_{t}$, $\partial_{t}\varphi$ and $D\varphi$ are continuous in $(a,T]\times O$.

For other equivalent definitions of viscosity solutions for (1.1), we refer to [9].

The first result is a comparison principle for (1.1) (see [9, Theorem 3.1]).

Theorem 2.2**.**

Assume (H1)-(H3). Let $u\in USC(\overline{Q_{T}})$ and $v\in LSC(\overline{Q_{T}})$ be a subsolution and a supersolution of (1.1), respectively. If $u(0,x)\leq v(0,x)$ for $x\in\mathbb{T}^{d}$, then $u\leq v$ on $\overline{Q_{T}}$.

We also recall an existence result for viscosity solutions of (1.1) (see [9, Theorem 4.2]).

Theorem 2.3**.**

Assume (H1). Let $u^{-}\in USC(\overline{Q_{T}})$ and $u^{+}\in LSC(\overline{Q_{T}})$ be a subsolution and a supersolution of (1.1) such that $(u^{-})_{*}>-\infty$ and $(u^{+})^{*}<+\infty$ on $\overline{Q_{T}}$. If $u^{-}\leq u^{+}$, then there exists a solution $u$ of (1.1) that satisfies $u^{-}\leq u\leq u^{+}$ in $\overline{Q_{T}}$.

By Theorems 2.2 and 2.3, it follows an existence and uniqueness result for the solution of (1.1), (1.2) (see [9, Corollary 4.3]).

Corollary 2.4**.**

Assume (H1)-(H4). Then there exists a unique viscosity solution of (1.1) which satisfies the initial condition (1.2).

Existence and uniqueness results for the problem of (1.1), (1.2) in a bounded domain with boundary conditions in viscosity sense are discussed in [15].

3 A class of finite difference schemes

In this section we describe a finite difference scheme for the approximation of (1.1). For simplicity of notations, we assume that the Hamiltonian $H$ depends only on the state and gradient variables, i.e. $H=H(x,p)$, and that the dimension $d$ is equal to $2$. The extension for general $H$ and $d$ will be clear from this special case. Moreover, we will always identify a function $\mathbb{T}^{d}$ with its $\mathbb{Z}^{d}$-periodic extension defined in all $\mathbb{R}^{d}$.

Let $\mathbb{T}^{2}_{h}$ be a uniform grid on the torus with step $h$, (this supposes that $1/{h}$ is an integer), and denote by $x_{i,j}$ a generic point in $\mathbb{T}^{2}_{h}$ (an anisotropic mesh with steps $h_{1}$ and $h_{2}$ is possible too and we have taken $h_{1}=h_{2}$ only for simplicity). The value $U^{n}_{i,j}$ denotes the numerical approximation of the function $u$ at $(x_{i,j},t_{n})=(ih,jh,n{\Delta t})$, $0\leq i,j\leq 1/h$, $n=0,\dots,N$ (assuming that $N=T/{\Delta t}$ is an integer). We also denote by $U^{n}$ the grid function taking the value $U^{n}_{i,j}$ at $x_{i,j}\in\mathbb{T}^{2}_{h}$.

We start by describing the numerical approximation of the Caputo time-fractional derivative $\partial_{t}^{\alpha}$ introduced in [11]. The numerical derivative is obtained by approximating the time-derivative inside the fractional integral in (2.1) via finite difference and writing in compact form the expression so obtained. We approximate $\partial_{t}^{\alpha}u(x_{i,j},t_{n+1})$ by

[TABLE]

since $t_{n}-t_{m}=(n-m){\Delta t}$. Defined

[TABLE]

we obtain by (3.1)

[TABLE]

where

[TABLE]

for $1\leq m\leq n$. Thus, the approximation of the Caputo time-derivative is given by

[TABLE]

Remark 3.1**.**

Denoted by $r_{{\Delta t}}^{n+1}$ the truncation error, in [11] it is proved that

[TABLE]

where $c_{u}$ is a constant depending on the second order time-derivative of $u$. Hence the temporal accuracy of the scheme is of order $2-\alpha$.

In the following we summarize some properties of the coefficients $c_{m}$ in (3.3)

Lemma 3.2**.**

1. (i)

$c_{m}^{n+1}>0$ for $0\leq m\leq n$. 2. (ii)

$c_{0}^{n+2}-c_{0}^{n+1}=-c_{1}^{n+2}$. 3. (iii)

$c_{m+1}^{n+2}=c_{m}^{n+1}\text{ for }1\leq m\leq n$. 4. (iv)

$\sum_{m=0}^{n}c_{m}^{n+1}=1$.

Proof.

1. (i)

When $m=0$, it is clear that $c_{0}^{n+1}>0$. Consider the case where $1\leq m\leq n$. Because of the strong concavity of the function $x^{1-\alpha}$ for $x\geq 0$, by Jensen’s inequality, we have

[TABLE]

Thus, it follows that $c_{m}^{n+1}>0$. 2. (ii)

By definition,

[TABLE] 3. (iii)

By definition,

[TABLE] 4. (iv)

We have

[TABLE]

∎

For the approximation of the Hamiltonian in (1.1) we follow the approach in [6]. We introduce the finite difference operators

[TABLE]

and define

[TABLE]

In order to approximate the Hamiltonian $H$ in equation (1.1), we consider a numerical Hamiltonian $g:\mathbb{T}^{2}\times\mathbb{R}^{4}\to\mathbb{R}$, $(x,q_{1},q_{2},q_{3},q_{4})\mapsto g\left(x,q_{1},q_{2},q_{3},q_{4}\right)$ satisfying the following conditions:

• (G1)

$g$ is nonincreasing with respect to $q_{1}$ and $q_{3}$, and nondecreasing with respect to $q_{2}$ and $q_{4}$.

• (G2)

$g$ is consistent with the Hamiltonian $H$, i.e.

[TABLE]

• (G3)

$g$ is locally Lipschitz continuous.

• (G4)

There exists a constant $C$ such that

[TABLE]

Hence, recalling the approximation (3.3) of the Caputo time derivative, we consider the explicit finite difference scheme

[TABLE]

for $i,j=1,\dots,1/h$, $n=0,\dots,N-1$, where $\rho_{\alpha}$ is defined as in (3.2) and

[TABLE]

In (3.7), $[U^{n}]_{i,j}$ represents the set of the values of $U^{n}$ used to compute the scheme at $x_{i,j}$, except that the value $U^{n}_{i,j}$ itself, and $h$ is the space discretization step. The scheme is completed with the initial condition

[TABLE]

Note that $U^{n+1}$ depends on all the past history $U^{m}$, $m=0,\dots,n$ of the solution.

For $\alpha=1$, the scheme (3.6) reduces to the standard finite difference approximation

[TABLE]

of the Hamilton-Jacobi equation

[TABLE]

3.1 Stability properties of the scheme

We set $Q^{h,{\Delta t}}_{n}=\mathbb{T}^{2}_{h}\times\{0,\dots,n{\Delta t}\}$ and we denote by $\mathcal{G}$ the space of the grid functions on $\mathbb{T}^{2}_{h}$ and by $\mathcal{G}^{n}$, $n=0,\dots,N$, the set of the grid function on $Q^{h,{\Delta t}}_{n}$, i.e.

[TABLE]

Moreover, we set $\|U\|_{\infty}=\sup_{i,j}|U_{i,j}|$ for $U=\{U_{i,j}\}_{i,j=0}^{1/h}\in\mathcal{G}$, and $\|\mathcal{U}\|_{\infty}=\sup_{m=0,\dots,n}\|U^{m}\|_{\infty}$ for $\mathcal{U}=\{U^{m}\}_{m=0}^{n}\in\mathcal{G}^{n}$.

For $n\in\{0,\dots,N-1\},$, we define a map $G^{n}:\mathcal{G}^{n}\to\mathcal{G}$ by

[TABLE]

Hence, the scheme (3.6) can be rewritten in the equivalent iterative form

[TABLE]

Definition 3.3**.**

We say that the scheme (3.10) is monotone if, for any $n=0,\dots,N-1$, $\mathcal{U},\mathcal{V}\in\mathcal{G}^{n}$, we have that

[TABLE]

where the previous inequalities are intended in the sense of the comparison of components.

Since the scheme (3.10) is explicit, for the monotonicity, we need some restriction on the approximation steps $h$ and ${\Delta t}$, as we will discuss later on.

Proposition 3.4**.**

Assume that the scheme (3.6) is monotone. Then, for $n=0,\dots,N-1$, we have

1. (i)

$G^{n}(\mathcal{U}+\lambda)=G^{n}(\mathcal{U})+\lambda$ for any $\lambda\in\mathbb{R}$, $\mathcal{U}\in\mathcal{G}^{n}$ (where we identify $\lambda$ both with the constant function on $\mathbb{T}^{2}_{h}$ and with the element of $\mathcal{G}^{n}$ such that $\lambda^{m}=\lambda$ is the constant function on $\mathbb{T}^{2}_{h}$ for any $m=0,\dots,n$); 2. (ii)

$\|G^{n}(\mathcal{U})-G^{n}(\mathcal{V})\|_{\infty}\leq\|\mathcal{U}-\mathcal{V}\|_{\infty}$ for any $\mathcal{U},\mathcal{V}\in\mathcal{G}^{n}$; 3. (iii)

For the constant $C\in\mathbb{R}$ in (G4),

[TABLE]

where $D_{h}\mathcal{U}=\{D_{h}U^{m}\}_{m=0}^{n}$; 4. (iv)

for any $\mathcal{U}\in\mathcal{G}^{n+1}$

[TABLE]

where $K=\sup_{x_{i,j}\in\mathbb{T}^{2}_{h},m=0,\dots,n+1}\|g(x_{i,j},[D_{h}U^{m}]_{i,j})\|_{\infty}$; 5. (v)

for any $\mathcal{U}\in\mathcal{G}^{n}$

[TABLE]

Proof.

1. (i)

By Lemma 3.2, we have

[TABLE] 2. (ii)

Let $\mathcal{U},\mathcal{V}\in\mathcal{G}^{n}$ and $\lambda=\|(\mathcal{U}-\mathcal{V})^{+}\|_{\infty}$. We have, in the sense of the comparison of components,

[TABLE]

for $m=0,\dots,n$. By monotonicity and commutativity,

[TABLE]

Hence, $G^{n}(\mathcal{U})-G^{n}(\mathcal{V})\leq\|(\mathcal{U}-\mathcal{V})^{+}\|_{\infty}$. Similarly, we obtain another inequality $G^{n}(\mathcal{U})-G^{n}(\mathcal{V})\geq-\|(\mathcal{U}-\mathcal{V})^{-}\|_{\infty}$. These two inequalities yield the result. 3. (iii)

Let $\tau$ be a translation operator in space, that is, $\tau_{l}U_{i,j}=U_{i+l_{1},j+l_{2}}$ for $l=(l_{1},l_{2})\in\mathbb{Z}^{2}$ for $U\in\mathcal{G}$ and defined in a similar way for $\mathcal{U}\in\mathcal{G}^{n}$. Then replacing the numerical Hamiltonian $g$ with $S$ defined in (3.7) in the assumption (G4) and from the part (ii) above, it follows that

[TABLE]

since $\|D_{+}^{1}\mathcal{U}\|_{\infty}\leq\|D_{h}\mathcal{U}\|_{\infty}$. We have similar estimates for the other components of $D_{h}G^{n}(\mathcal{U})$. Combining all these inequalities, we obtain the desired inequality. 4. (iv)

Using Lemma 3.2, we have

[TABLE] 5. (v)

By the consistency of scheme, it follows that $G^{n}(0)=-\rho_{\alpha}H(x_{i,j},0)$. Hence, by property (ii), we have

[TABLE]

∎

Proposition 3.5**.**

Assume that (3.6) is monotone and let $\{U^{n}\}_{n=0}^{N}$ be a sequence generated by the scheme with the initial condition (3.8). Then

[TABLE]

where $K=\sup_{x_{i,j}\in\mathbb{T}^{2}_{h},m=0,\dots,n}\|g(x,[D_{h}U^{m}]_{i,j})\|_{\infty}$.

Proof.

For $n=1$, (3.11) is true since $U_{i,j}^{1}=U_{i,j}^{0}-\rho_{\alpha}S(x_{i,j},h,U^{0}_{i,j},[U^{0}]_{i,j})$ with $S$ defined as in (3.7). Arguing by induction, assume now that (3.11) is true for $j\leq n$. Then using the property (iv) in Lemma 3.2, we get

[TABLE]

We observe that

[TABLE]

Moreover, by the inequality $(r+1)^{\beta}-r^{\beta}\geq\beta(r+1)^{\beta-1}$ for $r\geq 0$ and $\beta\in(0,1)$ (which follows by the concavity of the function $r^{\beta}$), we get

[TABLE]

Hence

[TABLE]

and replacing the previous inequality in (3.12), we get estimate (3.11). ∎

We discuss some classical examples of approximation scheme for Hamilton-Jacobi equations adapted to the fractional case. We consider the equation

[TABLE]

with periodic boundary condition.

Upwind scheme

Simple upwind schemes for the equation (3.13) are

[TABLE]

if $H$ is non-increasing, or

[TABLE]

if $H$ is non-decreasing. The numerical Hamiltonian is given by $g(q_{1},q_{2})=H(q_{1})$, in the first case, and by $g(q_{1},q_{2})=H(q_{2})$ in the second case. In both cases, $g$ is monotone, consistent and regular if $H$ is locally Lipschitz.

Now, we establish a condition under which the scheme (3.14) is monotone. By construction, the monotonicity with respect to $U^{n}_{j+1}$ is obvious. Also, since by Lemma 3.2 all the coefficients $c_{m}^{n+1}$ are positive, the map $G^{n}$ is increasing with respect to the variable $U_{j}^{m}$, $m=0,\dots,n-1$. Moreover, (3.14) is non-decreasing with respect to $U^{n}_{j}$ if $c^{n+1}_{n}+\frac{\rho_{\alpha}}{h}H^{\prime}\left(\frac{U^{n}_{j+1}-U^{n}_{j}}{h}\right)\geq 0$. Recalling that $c_{n}^{n+1}=2-2^{1-\alpha}$, we get the CFL condition

[TABLE]

for $p\in\mathbb{R}$. The same condition is necessary also for (3.15).

Lax-Friedrichs scheme

The Lax-Friedrichs scheme is given by

[TABLE]

where $\theta$ has to be chosen in order to satisfy the CFL condition. Therefore, the numerical Hamiltonian $g$ is

[TABLE]

for $\lambda=\rho_{\alpha}/h$ and $q_{1},q_{2}\in\mathbb{R}$. For the monotonicity of the scheme with respect to $U^{n}_{j}$, we need the condition

[TABLE]

and, for the monotonicity with respect to $U^{n}_{j\pm 1}$,

[TABLE]

Then, recalling (3.2) the monotonicity of the scheme is implied by the CFL condition

[TABLE]

for $p\in\mathbb{R}$.

Remark 3.6**.**

The CFL conditions (3.16) and (3.18) reduce to the classical ones for $\alpha=1$. In general, they become more and more restrictive for $\alpha$ decreasing to $0^{+}$. This phenomenum has been also observed in [12] in the study of approximation schemes for time-fractional conservation laws.

4 A convergence result for the finite difference scheme

In this section, we study the convergence of the scheme (3.6) following the classical stability argument in [4], where it is proved that a monotone, stable and consistent approximation scheme converges to the unique solution of the continuous Hamilton-Jacobi equations.

We recall the definition of the relaxed limit for a locally bounded sequence $\{u_{\sigma}\}_{\sigma>0}$. The upper relaxed limit is given by

[TABLE]

while the lower relaxed limit by $\liminf^{*}_{\sigma\to 0^{+}}u_{\sigma}=-\limsup^{*}_{\sigma\to 0^{+}}(-u_{\sigma}).$

We consider a sequence of approximation steps $({\Delta t},h({\Delta t}))$ such that $h({\Delta t})\to 0$ for ${\Delta t}\to 0$ and we denote with $u_{\Delta t}$ the piecewise constant extension to $\overline{Q_{T}}$ of the solution of the approximation scheme (3.6) corresponding to the previous parameters, i.e. $u_{\Delta t}(t,x)=U^{n}_{i,j}$ if $n=[t/{\Delta t}]$ (where $[\cdot]$ denotes the integer part) and $(i,j)$ such that $x\in((i-1/2)h,(i+1/2)h]\times((j-1/2)h,(j+1/2)h]$.

Theorem 4.1**.**

Assume that $g$ satisfies (G1)-(G4), $u_{0}$ is Lipschitz continuous and, for ${\Delta t}$ sufficiently small and $h=h({\Delta t})$, the scheme (3.6) is monotone. As ${\Delta t}\to 0^{+}$, the sequence $\{u_{\Delta t}\}_{{\Delta t}>0}$ converges locally uniformly to the unique viscosity solution $u$ of (1.1), (1.2).

We need a preliminary result about the convergence of the fractional derivative for a test function.

Lemma 4.2**.**

Let $\varphi\in C^{1,1}((0,T]\times\mathbb{T}^{d})\cap C([0,T]\times\mathbb{T}^{d})$ be a test function and let $(t_{\Delta t},x_{\Delta t})\to(t,x)\in(0,T)\times\mathbb{T}^{d}$ for ${\Delta t}\to 0$. Then, defined for $n_{\Delta t}+1=[t_{\Delta t}/{\Delta t}]$ the discrete fractional derivative

[TABLE]

we have

[TABLE]

Proof.

Convergence properties of the discrete fractional derivative to the continuous one are already studied in [10, 11]. Here we give a different proof of the convergence result for reader’s convenience. To simplify the notation, since in this argument only the time variable is involved, we omit the dependence of $\varphi$ on $x$. Because of the continuity of the Caputo derivative of $\varphi$ with respect to $t$ (see [15, Prop. 2.1]), it is sufficient to prove that

[TABLE]

Moreover, for a test function $\varphi$, the Caputo derivative can be defined in the standard way, see (2.1). In the rest of the proof, we omit the index ${\Delta t}$ and we write $t$, $n$ and in place of $t_{\Delta t}$, $n_{\Delta t}$. Fix $\eta>0$ such that $t>2\eta$ and let ${\bar{n}}<n$ be the greatest integer such that ${\bar{n}}{\Delta t}\leq\eta$. Then we write

[TABLE]

We estimate the two sums (multiplied by $\Gamma(1-\alpha)$) in (4.2) separately. For $0\leq j\leq\bar{n}-1$, use the integration by parts to get

[TABLE]

Hence,

[TABLE]

Observe that $\frac{1}{t-s}\leq\frac{1}{t-t_{j+1}}\leq\frac{1}{t-\eta}$ for $t_{j}\leq s\leq t_{j+1}$. For the first term of (4.3), using the inequality $t^{\alpha}-s^{\alpha}\leq(t-s)^{\alpha}$ for $0\leq s\leq t$, we have the following estimate:

[TABLE]

For the second term of (4.3), we have the following estimate:

[TABLE]

where $\omega_{\varphi}$ is a modulus of continuity of $\varphi$. Thus,

[TABLE]

Clearly, both terms converge to [math] as $\Delta t\to 0$.

Now, we estimate the second sum in (4.2). We have

[TABLE]

where $\omega_{\varphi^{\prime}}$ is the modulus of continuity of $\varphi^{\prime}$ on $[\eta,t]$. Summarizing these estimates, we conclude

[TABLE]

as $\Delta t\to 0$ . ∎

Proof of Theorem 4.1.

In order to apply the Barles-Souganidis’ convergence result (see [4]), we define for $(t,x)\in{\overline{Q_{T}}}$

[TABLE]

Note that, by definition, $\underline{u}(t,x)\leq\overline{u}(t,x)$. We claim that $\overline{u}$, $\underline{u}$ are, respectively, a viscosity subsolution and a viscosity supersolution of (1.1) such that $\overline{u}(0,x)\leq\underline{u}(0,x)$ for $x\in\mathbb{T}^{2}$. If the claim holds, then from Theorem 2.2 it follows that $\overline{u}(t,x)\leq\underline{u}(t,x)$ and therefore $u=\overline{u}\equiv\underline{u}$ is the unique viscosity solution of (1.1) in $Q_{T}$. Moreover, the definition of $\overline{u}$, $\underline{u}$ implies the uniform convergence of $\{u_{\Delta t}\}_{{\Delta t}>0}$ to $u$.

To prove the claim, we first observe that (3.11) and the continuity of $u_{0}$ implies that $\underline{u}=\overline{u}=u_{0}(x)$ for $x\in\mathbb{T}^{2}$. Clearly, by (3.8) we have $\underline{u}\leq u_{0}\leq\overline{u}$. Moreover, if $(s_{\Delta t},y_{\Delta t})\to(0,x)$ for ${\Delta t}\to 0$, then define $n_{\Delta t}=[s_{\Delta t}/{\Delta t}]$ and let $i_{\Delta t},j_{\Delta t}$ be such that $y_{\Delta t}\in((i_{\Delta t}-1/2)h,(i_{\Delta t}+1/2)h]\times((j_{\Delta t}-1/2)h,(j_{\Delta t}+1/2)h]$. We have

[TABLE]

where $L_{0}$ is the Lipschitz constant of $u_{0}$ and $K=\sup\{|g(x,q)|:x\in\mathbb{T}^{2},\,|q|\leq L_{0}\}$. Passing to the limit in the previous inequality for ${\Delta t}\to 0^{+}$, we get $\limsup_{\Delta t\to 0}u_{\Delta t}(s_{\Delta t},y_{\Delta t})\leq u_{0}(x)$ which implies, for the arbitrariness of the sequence $(s_{\Delta t},y_{\Delta t})$, $\overline{u}(0,x)\leq u_{0}(x)$. We prove similarly that $\underline{u}(0,x)\geq u_{0}(x)$.

The stability of the scheme (3.10), i.e. the boundedness of the sequence $\{u_{\Delta t}\}_{{\Delta t}>0}$ bounded uniformly in ${\Delta t}$, is clearly implied by Prop. 3.5.

To prove the consistency of the scheme, we claim that, given a test function $\varphi$ and a sequence $(t_{{\Delta t}},x_{\Delta t})=(n_{\Delta t}{\Delta t},(i_{\Delta t}h,j_{\Delta t}h))$ converging to $(t,x)\in Q_{T}$ for ${\Delta t}\to 0$, then we have

[TABLE]

where $D^{\alpha}_{\Delta t}\varphi(t_{\Delta t}+\Delta t,x_{\Delta t})$ is defined as in (4.1). The previous claim follows by Lemma 4.2 and by the equality

[TABLE]

which is consequence of the assumptions (G2)-(G3) for the numerical Hamiltonian $g$. Hence the claim (4.4) holds and we conclude that the scheme (3.10), besides monotone, it is also stable and consistent with the continuous equation (1.1).

We prove that, by the monotonicity, stability and consistency properties of the scheme, it follows that $\overline{u}$, $\underline{u}$ are, respectively, a viscosity subsolution and a viscosity supersolution of (1.1). We only show that $\overline{u}$ is a subsolution, since the argument for $\underline{u}$ is similar. First observe that the stability of the scheme implies that the sequence $\{u_{\Delta t}\}_{{\Delta t}>0}$ is bounded and therefore $\overline{u}$ is well defined. Consider a test function $\varphi$ such that $\overline{u}-\varphi$ takes a strict maximum point at $(\hat{t},\hat{x})\in\mathbb{T}^{d}\times(0,T)$. By Lemma V.I.6 in [2], there exists a sequence $\{(t_{{\Delta t}},x_{{\Delta t}})\}_{{\Delta t}>0}$ such that $u_{{\Delta t}}-\varphi$ has a maximum point at $(t_{{\Delta t}},x_{{\Delta t}})$, moreover $(t_{{\Delta t}},x_{{\Delta t}})\to(\hat{t},\hat{x})$ and $u_{\Delta t}(t_{\Delta t},x_{\Delta t})\to\overline{u}(\hat{t},\hat{x})$ for ${\Delta t}\to 0$. Set $\delta_{\Delta t}=u_{\Delta t}(t_{\Delta t}x_{\Delta t})-\varphi(t_{\Delta t},x_{\Delta t})$ and observe that

[TABLE]

Let $n_{\Delta t}$ be such that $n_{\Delta t}+1=[t_{\Delta t}/{\Delta t}]$ and $i_{\Delta t},j_{\Delta t}$ such that $x_{\Delta t}\in((i_{\Delta t}-1/2)h,(i_{\Delta t}+1/2)h]\times((j_{\Delta t}-1/2)h,(j_{\Delta t}+1/2)h]$ . Define $\Phi^{m}_{i,j}=\varphi(m{\Delta t},x_{i,j})$, for $i,j=1,\dots,1/h$, $m=0,\dots,[T/{\Delta t}]$ and let $\omega_{\varphi}$ be a modulus of continuity for $\varphi$ . By (4.5), we have

[TABLE]

for $m=0,\dots,[T/{\Delta t}]$, $i,j=0,\dots,h$. Hence, recalling that $u_{\Delta t}$ is the piecewise constant interpolation of the solution of (3.6) and the scheme can be rewritten in the equivalent form (3.10), by (4.6), Prop. 3.4.(i) and the monotonicity of the scheme, we have

[TABLE]

Hence

[TABLE]

and, using the consistency of the scheme and $((n_{\Delta t}+1){\Delta t},(i_{\Delta t}h,j_{\Delta t}h))\to(\hat{t},\hat{x})$ , we finally get

[TABLE]

∎

Remark 4.3**.**

Error estimates for the approximation of Hamilton-Jacobi equations with standard time derivative have been discussed by several authors (see for example [5], [6]). We aim to investigate this point in the future.

5 Explicit solutions and numerical tests

In this section, we implement Lax-Friedrichs scheme to test the convergence of the method. The problems we consider are in $\mathbb{R}^{d}$ and therefore the convergence is not guaranteed by the results of the previous sections, which are in the periodic case. Nevertheless, as the next examples show, we obtain a correct approximation of the exact viscosity solution of the problems.

5.1 Test 1

First, we consider the following Hamilton-Jacobi equation:

[TABLE]

It is easy to verify that, if $\alpha=1$, then the unique viscosity solution of (5.1) is given by

[TABLE]

We claim that a solution of (5.1) for $\alpha\in(0,1)$ is given by

[TABLE]

with $f(t)$ non-negative function to be determined. Replacing into the equation (5.1) for $|x|\leq\sqrt{1/f(t)}$ and taking into account the initial datum, we find that the function $f(t)$ has to satisfy the fractional differential equation

[TABLE]

We check that (5.2) gives a viscosity solution of (5.1). Note that, since $u_{\alpha}$ is defined as the minimum of two solutions of the equation, by [9, Lemma 4.1] it is a supersolution in $\mathbb{R}^{d}$. Moreover, for $|x|\neq 1/\sqrt{f(t)}$, $u_{\alpha}$ is differentiable and the equation is satisfied in point-wise sense (see [9, Prop. 2.10]), hence we only have to check the viscosity subsolution condition at $|x|=1/\sqrt{f(t)}$. Let $\varphi$ be a test function such that $u_{\alpha}-\varphi$ has a maximum point at $(t,1/\sqrt{f(t)})$ (for the other point we proceed in a similar way). Recalling (2.2) and, since

[TABLE]

we see that $\varphi$ has to satisfy

[TABLE]

Moreover, $0\leq D\varphi(t,1/\sqrt{f(t)})\leq 2f(t)/\sqrt{f(t)}$ and, therefore, substituting in (5.1), we get

[TABLE]

Nevertheless, to compute $u_{\alpha}$, we need to compute the function $f$. We look for a solution of (5.3) in the form of a power series $f(t)=\sum_{n=0}^{\infty}f_{n}t^{\alpha n}$. Replacing in the equation (5.3) and observing that $\partial_{t}^{\alpha}t^{0}=0$, we have

[TABLE]

A straightforward computation gives

[TABLE]

where $\beta_{n}=\Gamma(\alpha n+1)/\Gamma(\alpha(n-1)+1)$. Replacing the previous identity in the equation (5.4), we get

[TABLE]

Collecting the terms of the same order and recalling that, by the initial condition in (5.3), $f_{0}=1$, we find

[TABLE]

From the previous relations, we can iteratively compute the coefficients of the power series $f(t)$ and we replace in (5.2). Note that for $\alpha=1$, we get the power series of $\frac{1}{1+2t}$. As observed in [17, Section 6.2.3] where the problem (5.3) for $\alpha=1/2$ is studied via power series expansion, finding the convergence interval of (5.5) is not easy, but the series can be used for small values of $t$. Indeed we numerically find that, for each $\alpha\in(0,1]$, there is a critical time $T_{\alpha}$ for which the power series $f(t)$ converges if $t\leq T_{\alpha}$ and diverges if $t>T_{\alpha}$. The dependence of $T_{\alpha}$ on $\alpha$ is presented in Fig. 1. Hence we use the representation formula (5.2), with $f$ given by (5.3), for the exact solution of (5.1) only for $t<T_{\alpha}$. For $t\geq T_{\alpha}$, the viscosity solution of (5.1), which is defined for any $t$, cannot be anymore represented by means of (5.2) and we do not know if there is an alternative explicit formula. Alternatively, we can calculate $f$ by means of a numerical methods for (5.3), but we do not pursue this approach here.

The numerical solution at $t=0.2$ of (5.1) where $\alpha=0.8$ and $d=2$ computed by the Lax-Friedrichs scheme with $h=10^{-1}$, $\Delta t=10^{-3}$ and $\theta=1-2^{-\alpha}$ is provided in Fig. 2. We plot numerical solutions at $t=0.2$ for different values of $\alpha$ in Fig. 3 (A) for $d=1$. We observe the convergent behavior of the solutions as $\alpha\to 1$. These solutions eventually converge to the solution of the classical case.

For the convergence test, we use $l^{\infty}$ error defined by the maximum difference between the exact and numerical solutions over all nodes. From Fig. 3 (B), we determine that the convergence for the Lax-Friedrichs scheme under the CFL condition is linear.

5.2 Test 2

We consider the Hamilton-Jacobi equation of the form

[TABLE]

For $\alpha=1$, the unique viscosity solution of (5.6) is

[TABLE]

For $\alpha=1/2$, we look for a solution in the form

[TABLE]

Computing the derivatives

[TABLE]

and replacing in the equation (5.6), we get

[TABLE]

and therefore

[TABLE]

We observe that the previous formula defines a viscosity solution of the problem. Indeed, for $x\neq 0$ and $t>0$, $u_{\frac{1}{2}}$ is differentiable and, since it is a solution in point-wise sense, it is also a viscosity solution (see [9, Prop. 2.10]). For $x=0$ and $t>0$, $u_{\frac{1}{2}}$ is not differentiable with respect to $x$, hence the equation has to be verified in viscosity sense. Since $u_{\frac{1}{2}}$ behaves as $-|x|$ near $x=0$, there cannot exist a test function $\varphi$ such that $u_{\frac{1}{2}}-\varphi$ has a minimum point at $(t,0)$, hence the supersolution condition is automatically satisfied. Moreover, for any test function $\varphi$ such that $u_{\frac{1}{2}}-\varphi$ has a maximum point at $(t,0)$, recalling (2.2), we see that $\varphi$ has to satisfy

[TABLE]

and $|D\varphi(t,0)|\leq 4\sqrt{\pi}/\sqrt{t}$. Hence, also the subsolution condition at $(t,0)$ is satisfied. For $\alpha\in(0,1)$, a similar computation gives that the solution of (5.6) is given by

[TABLE]

Fig. 4 depicts the initial condition and the numerical solution at $t=0.2$ for $\alpha=0.8$ obtained using the Lax-Friedrichs scheme in 2 dimensions. The parameter $\theta$ is always chosen equal to $1-2^{-\alpha}$. The numerical solutions corresponding to different values of $\alpha$ are plotted in Fig. 5 (A) for $d=1$. We can see the same convergent behavior of the solutions as $\alpha\to 1$ as in the previous part. Moreover, from the convergence test in Fig. 5 (B), we observe that the convergence to the exact solution is linear.