# A unified framework for the regularization of final value time-fractional diffusion equation

This paper focuses on the regularization of backward time-fractional diffusion problem on unbounded domain. This problem is well-known to be ill-posed, whence the need of a regularization method in order to recover stable approximate solution. For the problem under consideration, we present a unified framework of regularization which covers some techniques such as Fourier regularization [19], mollification [12] and approximate-inverse [7]. We investigate a regularization technique with two major advantages: the simplicity of computation of the regularized solution and the avoid of truncation of high frequency components (so as to avoid undesirable oscillation on the resulting approximate-solution). Under classical Sobolev-smoothness conditions, we derive order-optimal error estimates between the approximate solution and the exact solution in the case where both the data and the model are only approximately known. In addition, an order-optimal a-posteriori parameter choice rule based on the Morozov principle is given. Finally, via some numerical experiments in two-dimensional space, we illustrate the efficiency of our regularization approach and we numerically confirm the theoretical convergence rates established in the paper.

## Authors

• 2 publications
02/21/2022

### A variational technique of mollification applied to backward heat conduction problems

This paper addresses a backward heat conduction problem with fractional ...
04/15/2020

### Regularization of backward time-fractional parabolic equations by Sobolev equations methods

It is well-known that backward parabolic equations in which the initial ...
12/29/2020

### Estimating solution smoothness and data noise with Tikhonov regularization

A main drawback of classical Tikhonov regularization is that often the p...
09/15/2021

### Backward diffusion-wave problem: stability, regularization and approximation

We aim at the development and analysis of the numerical schemes for appr...
05/07/2021

### Identifying source term in the subdiffusion equation with L^2-TV regularization

In this paper, we consider the inverse source problem for the time-fract...
04/17/2022

### Initial state reconstruction on graphs

The presence of noise is an intrinsic problem in acquisition processes f...
10/05/2019

### Regularization of a backwards parabolic equation by fractional operators

The backwards diffusion equation is one of the classical ill-posed inver...
##### This week in AI

Get the week's most popular data science and artificial intelligence research sent straight to your inbox every Saturday.

## 1 Introduction

In this paper, we study the final value problem

 {∂γu∂tγ=Δux∈Rn,t∈(0,T)u(x,T)=g(x)x∈Rn, (1)

where we aim at recovering the initial distribution given the final distribution . In (1), and denotes the Caputo fractional derivative [kilbas2006theory] defined as

 ∂γu∂tγ=1Γ(1−γ)∫t0u′(s)(t−s)γds,

where is nothing but the Gamma function : .

Time-fractional diffusion equations usually model sub-diffusion processes such as slow and anomalous diffusion processes which failed to be described by classical diffusion models [balakrishnan1985anomalous, chechkin2005fractional, metzler2000random, podlubnv1999fractional]. Due to the high diversity of such phenomena which are not properly modeled by classical diffusion, time-fractional diffusion problems have gained much attention in last decades. Beyond applications in diffusion processes, time-fractional equation (1) has also been applied to image de-blurring [wang2013total] where the time-fractional derivative allows to capture the memory effect in image blurring.

It is well-known that the ill-posedness of equation (1) comes from the irreversibility of time of the diffusion equation, which is caused by the very smoothing property of the forward diffusion. As a result, very small perturbation of the final distribution may cause arbitrary large error in the initial distribution . Hence, a regularization method is crucial in order to recover stable approximate of the initial distribution . In this regard, many regularization methods have been applied to final value time-fractional diffusion equation. Let us mention the mollification method [van2020mollification, yang2014mollification], Fourier regularization [xiong2012inverse, yang2015fourier], the method of quasi-reversibility [liu2010backward], Tikhonov method [wang2015optimal], total variation regularization [wang2013total], boundary condition regularization [yang2013solving], non-local boundary value method [hao2019stability], truncation method [wang2012data]. Yet, the set of regularization methods applied to backward time-fractional diffusion equation still presents some sparsity, especially compared to the set regularization methods for backward classical diffusion problems.

In this paper, we describe how in the context of regularization of the final value time-fractional diffusion equation (1), the Fourier regularization [yang2015fourier] and mollification [van2020mollification] are nothing but examples of approximate-inverse [louis1990mollifier] regularization. Next, we investigate a regularization technique which yields a better trade-off between stability and accuracy compared to the Fourier regularization [yang2015fourier] and the mollification technique of Van Duc N. et al [van2020mollification]. We consider noisy setting where is approximated by a noisy data satisfying

 ||u(⋅,T)−gδ||L2(Rn)≤δ,

and we derive order-optimal convergence rates between our approximate solution and the exact solution under classical Sobolev smoothness condition

 u(⋅,0)∈Hp(Rn)with||u(⋅,0)||Hp(Rn)≤E,p>0. (2)

We also provide error estimates under the more realistic setting where both the data and the forward diffusion operator are only approximately known. The motivation here being that, in practice, the Mittag-Leffler function which plays a major role in the resolution of equation (1) can only be approximated in practice.

In Section 2, we discuss existence of solution of equation (1) and reformulate the equation into an operator equation of the form in which is a bounded linear operator on . We present key estimates necessary for the regularization analysis and illustrate the ill-posedness of recovering from . Next, we introduce the framework of regularization which includes Fourier regularization, mollification, and approximate inverse. At last, we introduce our regularization approach.

Section 3 deals with error estimates and order-optimality of our regularization technique under the smoothness condition (2). In this section, we derive error estimates between the approximate solution and the exact solution in Sobolev spaces with . We also give error estimates for the approximation of early distribution with . We end the section by presenting analogous error estimates for the case where both the data and the forward diffusion operator are only approximately known.

Section 4 is devoted to parameter selection rules which is a critical step in the application of a regularization method. Here we propose a Morozov-like a-posteriori parameter choice rule leading to order-optimal convergence rates under smoothness condition (2). We also present analogous error estimates as that obtained in Section 3 corresponding to the a-posteriori rule prescribed.

Finally, we study four numerical examples in Section 5 to illustrate the effectiveness of the regularization approach coupled with the parameter choice rule described in Section 4. Moreover, in this Section, we also carry out a numerical convergence rates analysis in order to confirm the theoretical convergence rates given in Section 3.

In the sequel, or always refers to the -norm of the function on , denotes the Sobolev norm of on and denotes operator norm of a bounded linear mapping. Throughout the paper, or (resp.

) denotes the Fourier (resp. inverse Fourier) transform of the function

defined as

 ˆf(ξ)=F(f)(ξ)=1√2πn∫Rnf(x)e−ix⋅ξdx,F−1(f)(x)=1√2πn∫Rnf(ξ)eix⋅ξdξ,ξ,x∈Rn.

## 2 Regularization

Let us start by the following result about the existence and uniqueness of solution of equation (1).

###### Proposition 1.

For all and , Problem (1) admits a unique weak solution . That is, the first equation in (1) holds in for all and for all with

 limt→T||u(⋅,t)−g||H2=0.

Proposition 1 is merely generalization of [yang2015fourier, Lemma 2.2] where only the case is considered. The idea of the proof is merely to check that the formal solution defined by (5) is the weak solution.

Now, let us define the framework that we will consider for the regularization of problem (1). Consider a data , by applying the Fourier transform in (1) with respect to variable , we get

 {∂γˆu(ξ,t)∂tγ=−|ξ|2ˆu(ξ,t)ξ∈Rn,t∈(0,T)ˆu(ξ,T)=ˆg(ξ)ξ∈Rn. (3)

By applying the Laplace transform with respect to variable in (3), one gets

 {ˆu(ξ,t)=ˆu(ξ,0)Eγ,1(−|ξ|2tγ)ξ∈Rnˆu(ξ,T)=ˆg(ξ)ξ∈Rn, (4)

where is the Mittag Leffler function [podlubnv1999fractional] defined as

 Eγ,1(z)=+∞∑k=0zkΓ(γk+1),z∈C.

From (4), we can deduce the following relation between the solution and the data

in the frequency domain:

 ˆu(ξ,0)=ˆg(ξ)Eγ,1(−|ξ|2Tγ). (5)

Moreover, we can also derive the following relation for early distribution with

 ∀t∈(0,T),ˆu(ξ,t)=Eγ,1(−|ξ|2tγ)Eγ,1(−|ξ|2Tγ)ˆg(ξ). (6)

From equations (5) and (6), we can see that the Mittag Leffler function plays an important role in time-fractional diffusion equation (1). Hence, let us recall some key estimates about the function that will be repeatedly used in the sequel.

###### Lemma 1.

Let , there exists constants and depending only on and such that

 ∀x≤0,C1Γ(1−γ)11−x≤Eγ,1(x)≤C2Γ(1−γ)11−x. (7)

For a proof of Lemma 1, see [yang2015fourier, Lemma 2.1]. From Lemma 1, we can easily derive the next Lemma.

###### Lemma 2.

Assume , then for every and ,

 1(1∨tγ)(C1Γ(1−γ)11+|ξ|2)≤Eγ,1(−|ξ|2tγ)≤1(1∧tγ)(C2Γ(1−γ)11+|ξ|2). (8)

and

 C1C2≤Eγ,1(−|ξ|2tγ)Eγ,1(−|ξ|2Tγ)≤C2C1(Tt)γ. (9)

In (8), denotes the maximum while denotes the minimum, that is, and .

From (6) and (9), we get that for every , which implies that the problem of recovering from is actually well posed. However, that of recovering is ill-posed. Indeed, from (5) and (8), we get that

 ||u(ξ,0)||≥Cγ(1+|ξ|2)|ˆg(ξ)|,withCγ=Γ(1−γ)(1∧Tγ)/C2. (10)

From (10), we see that very small perturbations in high frequencies in the data leads to arbitrary large errors in the solution . Therefore, one needs a regularization method to recover stable estimates of .

###### Remark 1.

From (5) and (8), we can also derive that

 ||u(⋅,t)||L2≤Dγ(1+|ξ|2)|ˆg(ξ)|,withDγ=Γ(1−γ)(1∨Tγ)/C1. (11)

Estimate (11) illustrates the fact that the backward time-fractional diffusion equation is less ill-posed (mildly ill-posed) on the contrary to the classical backward diffusion equation which is exponentially ill-posed. This is actually due to the asymptotic slow decay of the Mittag Leffler function compared to .

From (5), we can reformulate equation (1) into an operator equation

 Au(⋅,0)=g. (12)

where is the linear forward diffusion operator operator which maps the initial distribution to the final distribution , that is,

 A=F−1(Eγ,1(−|ξ|2Tγ))F. (13)

In the sequel, given a data , we aim at recovering . Let be a smooth real-valued function in satisfying . It is well-known that the family of functions defined by

 ∀x∈Rn,φα(x):=1αnφ(xα),

satisfies

 ∀f∈L2(Rn),φα⋆f→finL2(Rn)asα↓0, (14)

where is nothing but the convolution of the functions and defined as For , let be the mollifier operator defined by

 ∀α>0,∀f∈L2(Rn),Mαf=φα⋆f. (15)

From (14), we see that the family of operators is an approximation of unity in , that is,

 ∀f∈L2(Rn),Mαf→fin% L2(Rn)asα↓0. (16)

Let be the solution of the equation

 {∂γu∂tγ=Δux∈Rn,t∈(0,T)u(x,T)=(Mαg)(x)x∈Rn. (17)

From (5) and (6), by replacing by , and using the fact that the , one gets

 ⎧⎪ ⎪ ⎪⎨⎪ ⎪ ⎪⎩ˆuα(ξ,0)=√2πnEγ,1(−|ξ|2Tγ)ˆφα(ξ)ˆg(ξ)=√2πnˆφα(ξ)ˆu(ξ,0)ˆuα(ξ,t)=√2πnEγ,1(−|ξ|2tγ)Eγ,1(−|ξ|2Tγ)ˆφα(ξ)ˆg(ξ)=√2πnˆφα(ξ)ˆu(ξ,t),t∈(0,T). (18)

which yields that

 ∀t∈[0,T]uα(⋅,t)=Mαu(⋅,t). (19)
###### Proposition 2.

Assume that there exists a function such that the mollifier kernel verifies

 ∀α>0,ˆφ(αξ)Eγ,1(−|ξ|2Tγ)≤ϱ(α), (20)

Then the family solution of equation (17) defines a regularization method for problem (1) in .

Proof. From (18), noticing that and using the fact that the function is increasing on with , we can see that (20) implies that for every and , the mapping which maps the data to is bounded with . Moreover, from (16) and (19), we deduce that for every , converges to in as goes to .
From Proposition 2, we can see that, choosing a kernel which satisfies condition (20) allows to defines a regularization method for equation (1). Now, let us show that the family of regularization methods defined in this way actually coincides with approximate-inverse introduced by Louis and Maass [louis1990mollifier].

From the first equation in (18), we can derive that

 uα(x,0)=∫Rneixξˆφα(ξ)ˆg(ξ)Eγ,1(−|ξ|2Tγ)dξ. (21)

From (19) and (21), we can reformulate as follows

 ⎧⎪⎨⎪⎩∀x∈Rn,uα(x,0)=⟨eα(x,⋅),u(⋅,0)⟩L2,witheα(x,y)=φα(x−y)∀x∈Rn,uα(x,0)=⟨vx,α,g⟩L2,withvx,α=F−1(eixξˆφα(ξ)Eγ,1(−|ξ|2Tγ)). (22)

Moreover, one can easily check that is nothing but the solution of the adjoint equation

 A∗f=eα(x,⋅).

Hence, we deduce that the solution of equation (17) actually corresponds to the approximate-inverse [louis1990mollifier] regularized solution of equation (12), the pair being what is usually called reconstruction kernel - mollifier.

This setting of regularization we just described encompasses many regularization methods that appears distinctively in the literature of the regularization of the final value time-fractional diffusion equation. Each regularization method being a particular choice of the mollifier kernel .

For the Fourier regularization [yang2015fourier] (where ), we have

 ˆuξmax(ξ,t)=χ[−ξmax,ξmax](ξ)Eγ,1(−|ξ|2tγ)Eγ,1(−|ξ|2Tγ)ˆg(ξ), (23)

where

denotes the characteristic function of the set

equal to on and elsewhere. By comparing (23) and (18), we readily get that

 α=1/ξmax,andφ(x)=F−1(1√2πχ[−1,1])=sin(x)πx. (24)

In this case, condition (20) merely reads

 χ[−1,1](ξ/ξmax)Eγ,1(−|ξ|2Tγ)≤ϱ(1ξmax).

Using (8), we can see that this condition is satisfies with where .

For the mollification method of N. Van Duc et al. [van2020mollification], the mollifier operator is denoted where is the convolution by the so called Dirichlet kernel defined as

 Dν(x)=1πnn∏j=1sin(νxj)xj,withν>0andx∈Rn.

Hence we can deduce that, this merely corresponds to

 α=1/ν,andφ(x)=n∏j=1sin(xj)πxj. (25)

Given that the Fourier transform of the kernel is given by

 F(Dν)(ξ)=χΛ(ξ),withΛ={x∈Rn:|xj|≤ν,j=1,...,n}, (26)

we see that condition (20) merely reads

 χΛ(ξ/ν)Eγ,1(−|ξ|2Tγ)≤ϱ(1ν),

which is fulfilled with where .

By the way, from (24) and (25), we can see that the Fourier regularization and the mollification approach of N. Van Duc et al. actually coincides, the latter approach being a generalization of the former to dimensional case. From (26), we can conclude that both regularization approaches are nothing but truncation methods. That is, the regularization is done by merely throwing away high frequency components of the data, which are responsible of the ill-posedness, and conserving unchanged the remaining frequency components. In order word, these two methods can be regarded as spectral cut-off methods. It is important to notice that though high frequency components are responsible of ill-posedness, nevertheless, the still carry non-negligible information on the sought solution. Therefore, it is desirable to apply a regularization which do not suppress high frequency components but which applies much regularization to those components compared to low frequency components. Let us point out that mere truncation of high frequency components usually entails Gibbs phenomena and oscillation of the approximate solution which should be avoided as far as possible. This is actually possible by choosing a kernel whose Fourier transform is supported on the whole domain .

Now on, let us consider a mollifier kernel defined by

 ˆφ(ξ)=1√2πnexp(−τ|ξ|s),τ>0,s>0,i.e.φ=1√2πnF−1(exp(−τ|ξ|s)), (27)

where and are two free positive parameters. From (27), we can see that and satisfies .

###### Lemma 3.

Let and be two positive numbers, consider the function . Then there exists a constant depending only on such that

 supx≥0fb,d(x)≤Cb1/dasb↓0. (28)

The proof of Lemma 3 is deferred to appendix. Lemma 3 will help us to prove that the kernel given by (27) allows us to define a regularization method for equation (1).

###### Proposition 3.

Let be the mollifier operator defined by (15) with the kernel given in (27) with and being two positive numbers. Then the family of function solution of equation (17) defines a regularization method for equation (1).

Proof. In view of Proposition 2, it suffices to prove that the kernel given in (27) verifies (20). By considering (27) and estimate (8), we have

 ˆφ(αξ)Eγ,1(−|ξ|2Tγ)≤C(1+|ξ|2)exp(−ταs|ξ|s),withC=Γ(1−γ)(1∨Tγ)/√2πnC1. (29)

The right hand side in (29) is nothing but with and . Hence from (28), we deduce that there exists a constant independent on such that

 ∀ξ∈Rn,ˆφ(αξ)Eγ,1(−|ξ|2Tγ)≤Cα2asα→0,

whence (20) with .

###### Remark 2.

By defining the mollifier kernel as in (27), we can see that the regularization technique induces a more suitable treatment of frequency components. Indeed, with our choice of mollifier kernel, the amount of regularization smoothly depends on the magnitude of the frequency components: The higher the frequency, the stronger the regularization applied, and similarly, the lower the frequency, the lower the regularization applied. This is actually desirable for a regularization method given that as the frequency gets higher, the noise in the frequency components gets much more amplified, and as the frequency gets lower, the noise in the frequency components gets less and less amplified.

###### Remark 3.

From Proposition 2, we can see that the family of function defined by (27) allows to define a family of regularization methods, each regularization method being determined by the choice of the free parameter . For instance, the choice means considering a Cauchy convolution kernel while means taking a Gaussian convolution kernel.

Let us end this section by the following Lemma which gives rates of convergence of the mollifier operator corresponding to kernel defined by (27) on Sobolev spaces with .

###### Lemma 4.

Let and be the mollifier operator defined by (15) with the mollifier kernel given by (27). Then

 ∀f∈Hp(Rn),||f−Mαf||L2≤τp∧ssαp∧s||f||Hp. (30)

Proof. Let . If , using Parseval identity, we have

 ∥f−Mαf∥L2 = ∥∥[1−√2πnˆφ(αξ)]ˆf(ξ)∥∥L2 = ∥∥[1−exp(−τ(α|ξ|)s)]p/sˆf(ξ)×[1−exp(−τ(α|ξ|)s)]1−p/s∥∥L2 ≤ ∥∥[1−exp(−τ(α|ξ|)s)]p/sˆf(ξ)∥∥L2 ≤ (ταs)p/s∥∥(|ξ|s)p/sˆf(ξ)∥∥L2≤τp/sαp∥f∥Hp.

If ,

## 3 Error estimates

Henceforth, denotes the mollifier kernel defined by (27) and denotes a noisy data satisfying the noise level condition

 ||g−gδ||L2≤δ, (31)

where is the exact final distribution. Let us introduce the regularized solution corresponding to the noisy data as the solution of equation

 {∂γu∂tγ=Δux∈Rn,t∈(0,T)u(x,T)=(Mαgδ)(x)x∈Rn, (32)

Equivalently, we can define in the frequency domain by

 ˆuδα(ξ,t)=√2πnEγ,1(−|ξ|2tγ)Eγ,1(−|ξ|2Tγ)ˆφ(αξ)ˆgδ(ξ),t∈[0,T]. (33)

It is well known that without assuming a smoothness condition on the exact solution (or on the exact data ), it is impossible to exhibit a rate of convergence of regularized solution towards the exact solution [schock1985approximate]. Henceforth, we consider the following classical Sobolev smoothness condition:

 u(⋅,0)∈Hp(Rn),∥u(⋅,0)∥Hp≤E,withp>0,E>0. (34)

Before presenting the main results of this section, let us state some lemmas which will be useful in the sequel.

###### Lemma 5.

Let , and be a solution of equation on . If , then for every , and

 ||v(⋅,t)||Hp+2≤C1∧tγ||v(⋅,0)||Hp,withC=C2Γ(1−γ). (35)

Proof. The proof follows readily by applying Parseval identity and estimate (8) to equation which links and in the frequency domain.
The next lemma illustrates the fact that the Sobolev smoothness condition (34) is nothing but a Hölder source condition.

###### Lemma 6.

Let , be the solution of problem (1) and being the forward diffusion operator defined in (13). The smoothness condition is equivalent to the Hölder source condition with , satisfying

 (Γ(1−γ)(1∧Tγ)C2)p/2||u(⋅,0)||Hp≤||w||L2≤(Γ(1−γ)(1∨Tγ)C1)p/2||u(⋅,0)||Hp (36)

Proof. For , formally define in the frequency domain by

 ˆw(ξ)=Eγ,1(−|ξ|2Tγ)−p/2ˆu(ξ,0),⟺ˆu(ξ,0)=(Eγ,1(−|ξ|2Tγ)2)p/4ˆw(ξ)

From (13), we can verify that the above definition of from is merely reformulation of the equation in the frequency domain. Next, we can check that is well defined and belongs to . Finally, estimate (36) is deduced from (8).

###### Remark 4.

From Lemma 6, we can deduce that the order optimal convergence rate under smoothness condition (34) is nothing but with independent of and .

The next Lemma which generalizes [van2020mollification, Lemma 3] will be useful in the sequel for establishing Sobolev norm error estimates.

###### Lemma 7.

Let and be a solution of equation on . If then

 ∀l∈[0,p],∀t∈(0,T],∥v(⋅,t)∥Hl+2≤C(γ)(1∧tγ)∥v(⋅,0)∥2+lp+2Hp∥v(⋅,T)∥p−lp+2L2. (37)

where . Moreover,

 ∀l∈[0,p],∀t∈[0,T],∥v(⋅,t)∥Hl≤¯C(γ)∥v(⋅,0)∥2+lp+2Hp∥v(⋅,T)∥p−lp+2L2, (38)

where .

Proof. Let and be a solution of equation with . For , using Hölder inequality, we have

 ∥v(⋅,t)∥2Hl+2 = ∫Rn(1+|ξ|2)l+2∣∣Eγ,1(−|ξ|2tγ)ˆv(ξ,0)∣∣2dξ ≤ C2γ(1∧tγ)2∫Rn(1+|ξ|2)l|ˆv(ξ,0)|2dξwithCγ=C2Γ(1−γ)using(???) = C2γ(1∧tγ)2∫Rn((1+|ξ|2)p(l+2)p+2|ˆv(ξ,0)|2(l+2)p+2)((1+|ξ|2)2(l−p)p+2|ˆv(ξ,0)|2(p−l)p+2)dξ ≤ C2γ(1∧tγ)2(∫Rn(1+|ξ|2)p|ˆv(ξ,0)|2dξ)l+2p+2(∫Rn(1+|ξ|2)−2|ˆv(ξ,0)|2dξ)p−lp+2 ≤ (C2C11∨Tγ<