The method of simplified Tikhonov regularization for a time-fractional inverse diffusion problem

Accepted Manuscript The method of simplified Tikhonov regularization for a time-fractional inverse diffusion problem Fan Yang, Chu-Li Fu, Xiao-Xiao Li

PII: DOI: Reference:

S0378-4754(17)30315-4 http://dx.doi.org/10.1016/j.matcom.2017.08.004 MATCOM 4490

To appear in:

Mathematics and Computers in Simulation

Received date : 11 December 2012 Revised date : 9 August 2017 Accepted date : 28 August 2017 Please cite this article as: F. Yang, C. Fu, X. Li, The method of simplified Tikhonov regularization for a time-fractional inverse diffusion problem, Math. Comput. Simulation (2017), http://dx.doi.org/10.1016/j.matcom.2017.08.004 This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

The method of simplified Tikhonov regularization for a time-fractional inverse diffusion problem ∗ Fan Yang1,2†, Chu-Li Fu2 , Xiao-Xiao Li1 1

School of Science, Lanzhou University of Technology, Lanzhou, Gansu, 730050, People’s Republic of China

2

School of Mathematics and Statistics, Lanzhou University, Lanzhou, Gansu, 730000, People’s Republic of China

Abstract In this paper, we consider a time-fractional inverse diffusion problem, where the data are given at x = 1 and the solution is required in the interval 0 < x < 1. This problem is typically ill-posed, i.e., the solution (if it exists) does not depend continuously on the data. The simplified Tikhonov regularization method is proposed to solve this problem. An a priori error estimate between the exact solution and its regularized approximation is obtained. Moreover, a new a posteriori parameter choice rule is proposed and the H¨older type error estimate is also obtained. Some different type examples are presented to demonstrate the feasibility and efficiency of the proposed method. keywords: Time-fractional inverse diffusion problem; Simplified Tikhonov method; A posteriori parameter choice; Error estimate

1

Introduction

Partial differential equations with fractional order arose from the studies of continuous random walk [1], Levy motion [4], and high-frequency financial data [17]. Among these studies the modeling of advection and dispersion phenomena in groundwater hydrology to simulate the transport of passive tracers carried by fluid flow in a porous medium resulted in a partial differential equation with fractional order [2]. In general, fluid flow and diffusion phenomena are governed by a fractional advection-dispersion equation. If the initial concentration distribution and boundary conditions are given, a complete recovery of the unknown solution is attainable from solving a well-posed forward problem [12]. However, in some practical problems, the boundary data can only be measured ∗

The project is supported by the National Natural Science Foundation of China (No.11561045), the Doctor Fund of Lan Zhou University of Technology. † Corresponding author, E-mail: [email protected].

1

on a portion of the boundary or some points in the solution domain. This leads to an ill-posed problem of the fractional heat diffusion equation, which means the solution does not depend continuously on the given known conditions. In this paper, we consider the following problem:   −ux (x, t) = ∂0α+ u(x, t), x > 0, t > 0,    u(x, 0) = 0, x ≥ 0, (1.1)  u(1, t) = f (t), t ≥ 0,    u(x, t)| t ≥ 0, x→∞ = 0,

where the time-fractional derivative ∂0α+ u(x, t) is the Caputo fractional derivative of order α (0 < α ≤ 1) defined by ([15]) Z t ∂u(x, s) ds 1 α , 0 < α < 1, (1.2) ∂0+ u(x, t) = Γ(1 − α) 0 ∂s (t − s)α

∂u(x, t) , α = 1. (1.3) ∂t We want to determine the temperature for 0 < x < 1. Due to the severe ill-posedness of this problem, some numerical approaches for solving such problem are needed. About the inverse problem for the fractional diffusion equation, there have been a lot of research results. In [5], the authors proved the uniqueness of an inverse problem for identifying α and diffusion coefficient p(x) in time-fractional diffusion equation. Based on space-marching mollification techniques, Muro [13, 14] gave a numerical solution to the time-fractional inverse heat conduction problem. Qian [16] gave the optimal regularization method only for α = 12 for the space-fractional inverse heat conduction problem. Xiong [19] used optimal regularization method to solve the time-fractional inverse heat conduction problem. In [23], the authors used the new mollification regularization method to solve time fractional inverse heat conduction problem. In [24], the authors used the spectral regularization method to solve the time fractional inverse heat conduction problem in a quarter plane. About identifying the unknown source for fractional diffusion equation , one can see [18, 21, 22, 27]. In [26], the authors used the new mollification regularization method to solve Cauchy problem of the time fractional diffusion equation. In [25], the authors used two regularization methods to solve a RieszFeller space-fractional backward diffusion problem. Liu and Yamamoto [11] used the quasi-reversibility regularization method to study a backward problem for timefractional diffusion equation in a bound domain. In [28], the authors used an optimal regularization method to consider space-fractional backward diffusion problem. But in [25, 11, 28], the regularization parameters which depend on the noise level and the a priori bound are selected by the a priori rule. Generally speaking, there is a defect for any a priori method, i.e., the a priori choice of the regularization parameter depends seriously on the a priori bound E of the unknown solution. However, in general, the a priori bound E cannot be known exactly in practice, and working with a wrong constant E may lead to the bad regularized solution. ∂0α+ u(x, t) =

2

In the present paper, we will consider not just the a priori choice of the regularization parameter for the simplified Tikhonov regularization method, but also the a posteriori choice of the regularization parameter will be given. To the authors’s knowledge, there are few papers for choosing the regularization parameter by the a posteriori rule for this problem. This method was first introduced by Carasso in [3], then the idea of this method has been successfully used for solving various types of ill-posed problems. Fu [9] considered the inverse heat conduction problem. Cheng et al. [6] solved a spherically symmetric inverse heat conduction problem. Feng et al. [8] solved the Cauchy problem for the Helmholtz equation. Yang et al. [20] solved the identification unknown source problem. According to the a priori and a posteriori rules for choosing regularization parameters, we obtain the H¨ older-type error estimates between the exact solution and the regularization solutions. Furthermore, the error estimates are order optimal according to [7]. The problem (1.1) is ill-posed in Hadamard’s sense, any small change on the input data can result in a dramatic change to the solution. The ill-posedness can be seen by solving the problem in the frequency domain. In order to analyze the problem in L2 (R), we define all functions to be zero for t < 0. Let Z ∞ 1 fˆ(ξ) := √ e−iξt f (t)dt (1.4) 2π −∞ be the Fourier transform of the function f (t). The solution of problem (1.1) can be easily given as the following form: α (1−x)

u ˆ(x, ξ) = e(iξ) So

where

1 u(x, t) = √ 2π

Z

∞

fˆ(ξ).

α (1−x)

(1.5) fˆ(ξ)dξ,

(1.6)

( απ |ξ|α (cos( απ 2 ) + i sin( 2 )), ξ ≥ 0, (iξ) = απ |ξ|α (cos( απ 2 ) − i sin( 2 )), ξ < 0.

(1.7)

eiξt e(iξ)

−∞

α

α

απ

|ξ| cos( 2 )(1−x) Note that (iξ)α has the positive real part |ξ|α cos( απ 2 ) and therefore the factor e is increases exponentially for 0 < x < 1 as |ξ| → ∞. So a small disturbation for the data f (t) will be amplified infinitely by this factor and lead to the integral (1.6) blow-up, therefore recovering the temperature u(x, t) from the measured data fδ (t) is severely ill-posed. Assume the exact data f (t) and the measured data fδ (t) both belong to L2 (R) and satisfy kf − fδ k ≤ δ, (1.8)

where k · k denotes L2 -norm, the constant δ > 0 represents a noise level. Moreover, assume there hold the following a priori bound: ku(0, ·)k ≤ E,

3

(1.9)

where E is a fixed positive constant. The rest of this paper is organized as follows. In Section 2, an order-optimal error estimate is obtained for the a priori parameter choice rule. The a posteriori parameter choice rule is given in Section 3, which also leads to a H¨ older-type error estimate. Numerical implement shows the effectiveness for the method in Section 4. The paper ends with a brief conclusion in Section 5.

2

The a priori simplified Tikhonov regularization method and the error estimate

Let us formulate the problem of identifying u(x, t) from (exact) data u(1, t) = f (t) as an operator equation: Au = f, (2.1) where the operator A : L2 (R) → L2 (R) acts injective. From (4.12), we know α e−(iξ) (1−x) u ˆ(x, ξ) = fˆ(ξ).

(2.2)

Obviously, (2.2) is equivalent to the following operator equation: ˆu = fˆ(ξ), Aˆ

Aˆ = F SF −1 ,

(2.3)

where F : L2 (R) → L2 (R) is the Fourier operator that maps any L2 (R) function f (t) into its fourier transform fˆ(ξ). From (2.2), we obtain ˆu = e−(iξ)α (1−x) u Aˆ ˆ(x, ξ).

(2.4)

We shall use Fourier transform to obtain a representation of approximation solution of equation (2.1). From (2.2), we know that function f (t) lies in the domain of operator A−1 provided Z ∞ α 2 ku(x, ·)k = |e(iξ) (1−x) |2 |fˆ(ξ)|2 dξ < ∞. (2.5) α |e(iξ) (1−x) |

−∞

Note that, → ∞, as |ξ| → ∞, therefore (2.5) means that fˆ(ξ) must decay rapidly at high frequencies. Such a decay is not likely to occur in the Fourier transform of measured noisy temperature history fδ (t) ∈ L2 (R). Hence we will construct an approximate solution of u(x, t) by a Tikhonov regularization method which minimizes the quantity ˆu − fˆδ k2 + βkˆ kAˆ uk2 (2.6)

over all u ˆ ∈ L2 (R). β is the regularization parameter. The following lemma will give a regularization approximate of u(x, t). Lemma 2.1. There exists a unique solution to the above minimization problem. It is given by Z ∞ α 1 e(iξ) (1−x) uδ (x, t) = √ eiξt fˆδ (ξ)dξ, (2.7) 1 + β|e(iξ)α (1−x) |2 2π −∞

4

where β is the regularization parameter. Proof: Let I denote the identity operator in L2 (R) and A∗ be the adjoint of A. Then by Theorem 2.12 in [10], the unique solution of the minimization problem (2.6) is equal to solve the following normal equation: ˆuδ (x, t) + β u Aˆ∗ Aˆ ˆδ (x, t) = Aˆ∗ fˆδ (t).

(2.8)

α Aˆ = e−(iξ) (1−x) .

(2.9)

From (2.4), we know

For the operator Aˆ is a multiplication operator, we obtain Aˆ∗ = e−(iξ)α (1−x) .

(2.10)

Due to (2.9) and (2.10), we obtain α (1−x)

(|e−(iξ)

|2 + β)ˆ uδ (x, ξ) = e−(iξ)α (1−x) fˆδ (ξ).

(2.11)

Therefore α

u ˆδ (x, ξ) =

e−(iξ)α (1−x) e(iξ) (1−x) ˆ f (ξ) = fˆδ (ξ). δ |e−(iξ)α (1−x) |2 + β 1 + β|e(iξ)α (1−x) |2

(2.12)

Finally, (2.7) holds by an inverse Fourier transform. 2 Comparing formula (2.12) with formula (1.6), we find that the Tikhonov regularization procedure consists in replacing the unknown fˆ(ξ) with an appropriately filtered Fourier transform of noisy data fˆδ (ξ). The filter in (2.12) attenuates the high frequencies in fˆδ (ξ) in a manner consistent with the goal of minimizing quantity (2.6). α By this idea, we can use a much better filter 1(1 + β|e(iξ) |2 ) to replace the filter α 1(1 + β|e(iξ) (1−x) |2 ) and give another approximation uβδ (x, t) of solution u(x, t). We define a regularization approximate solution of problem (1.1) or (2.1) for noisy data fδ (t), which is called the simplified Tikhonov regularized solution as follows: uβδ (x, t)

1 := √ 2π

Z

∞

α

eiξt

−∞

e(iξ) (1−x) ˆ fδ (ξ)dξ. 1 + β|e(iξ)α |2

(2.13)

The main conclusion of this section is Theorem 2.2. Let u(x, t) given by (1.6) be exact solution of (1.1) and uβδ (x, t) given by (2.13) be the simplified Tikhonov regularized approximation of u(x, t). Let fδ (t) be measured data at x = 1 satisfying (1.8) and the a priori condition (1.9) holds. If β is selected by δ β = ( )2 , (2.14) E then there holds the following error estimate: ku(x, ·) − uβδ (x, ·)k ≤

x 2 − x 1−x x + 1 1 − x 1−x
x 1−x . ( ) 2 +( )( ) 2 δ E 2 x 2 1+x

5

(2.15)

Proof: Due to Parseval formula, (1.6), (2.13), (1.8), (1.9) and (2.14), we have ku(x, ·) − uβδ (x, ·)k = kˆ u(x, ·) − u ˆβδ (x, ·)k = ke

(iξ)α (1−x)

≤ ke

(iξ)α (1−x)

α

e(iξ) (1−x) ˆ fˆ(ξ) − fδ (ξ)k 1 + β|e(iξ)α |2 α

e(iξ) (1−x) ˆ fˆ(ξ) − f (ξ)k 1 + β|e(iξ)α |2 α

α

e(iξ) (1−x) ˆ e(iξ) (1−x) ˆ f (ξ) − fδ (ξ)k +k α 1 + β|e(iξ) |2 1 + β|e(iξ)α |2 α

α

β|e(iξ) |2 ˆ e(iξ) (1−x) ˆ −x(iξ)α (iξ)α e k + k =k f (ξ)e (f (ξ) − fˆδ (ξ))k α 1 + β|e(iξ) |2 1 + β|e(iξ)α |2 β|e(iξ)α |2 −x(iξ)α e(iξ)α (1−x) kfˆ(ξ) − fˆδ (ξ)k ≤ sup E + sup e (iξ)α |2 (iξ)α |2 ξ∈R 1 + β|e ξ∈R 1 + β|e π

π

e(2−x)|ξ|α cos( 2 α) e(1−x)|ξ|α cos( 2 α) δ. ≤ β sup E + sup 2|ξ|α cos( π2 α) 2|ξ|α cos( π2 α) ξ∈R 1 + βe ξ∈R 1 + βe Let A(ξ) :=

e(2−x)|ξ|

α

1 + βe2|ξ|

cos( π2 α)

α

cos( π2 α)

, B(ξ) :=

e(1−x)|ξ|

α

1 + βe2|ξ|

cos( π2 α)

α

cos( π2 α)

.

(2.16)

Let s := |ξ|α cos( π2 α), so A(ξ) and B(ξ) can be rewritten as A(s) :=

e(2−x)s e(1−x)s , B(s) := . 1 + βe2s 1 + βe2s

(2.17)

Through the basic calculation, we obtain dA(s) e(2−x)s [(2 − x) − xβe2s ] = . ds (1 + βe2s )2 1

2 It is easy to find s∗ = ln( 2−x xβ ) such that

dA(s) ds

< 0, so A(s) attains its maximum at A(s) ≤ A(s∗ ) =

dA(s) ds = ∗ s . So

0. If s < s∗ ,

(2.18) dA(s) ds

> 0; if s > s∗ ,

x 2 − x 2−x 1 2−x ( ) 2 ( ) 2 . 2 x β

(2.19)

Using the same method, we obtain B(s) ≤

1 + x 1 − x 1−x 1 1−x ( ) 2 ( ) 2 . 2 1+x β

So 1 + x 1 − x 1−x 1 1−x x 2 − x 2−x 1 2−x ku(x, ·) − uβδ (x, ·)k ≤ β ( ) 2 ( ) 2 E+ ( ) 2 ( ) 2 δ 2 x β 2 1+x β

6

(2.20)

=

x 2 − x 2−x 1 + x 1 − x 1−x
x 1−x ( ) 2 + ( ) 2 δ E . 2 x 2 1+x

Remark 2.3. The coefficient function in (2.15) G(x) := is plotted in Fig. 1.

x 2−x 2−x 2 2( x )

2


1−x 1−x 2 + ( x+1 2 )( 1+x )

1.5 1.45 1.4

G(x)

1.35 1.3 1.25 1.2 1.15 1.1

0

0.2

0.4

0.6

0.8

1

x

Fig. 1. Remark 2.4. In Theorem 2.2, we don’t obtain the error estimation at x = 0. In order to obtain the error estimate at x = 0, we should use a more stronger a priori bound as follows: ku(0, ·)kH p ≤ E,

where ku(0, ·)kH p denotes the norm in Sobolev space H p (R) defined by Z ∞ 1 p ku(0, ·)kH := ( |ˆ u(0, ξ)|2 (1 + ξ 2 )p dξ) 2 .

(2.21)

(2.22)

−∞

3

The discrepancy principle and the a posteriori simplified Tikhonov regularization method

In this section, we consider the a posteriori regularization parameter choice rule. The most general a posteriori rule is the Morozov’s discrepancy principle [10]. Choose the regularization parameter β as the solution of the equation k

1 fˆδ (ξ) − fˆδ (ξ)k = τ δ, 1 + β|e(iξ)α |2

7

(3.1)

where τ > 1 is a constant and β denotes the regularization parameter. To establish existence and uniqueness of solution for equation (3.1), we need the following lemma: Lemma 3.1. Let ρ(β) := k 1+β|e1(iξ)α |2 fˆδ (ξ) − fˆδ (ξ)k, if 0 < δ < kfˆδ (ξ)k, then (a) ρ(β) is a continuous function; (b) limβ→0 ρ(β) = 0; (c) limβ→+∞ ρ(β) = kfˆδ k; (d) ρ(β) is a strictly increasing function. The proof is very easy and we omit it here. Remark 3.2. From Lemma 3.1, we know that there exists a unique solution β satisfying Eq. (3.1). Lemma 3.3. Setting ωδβ (x, ξ) := e(iξ) ing inequality holds:

α (1−x)

fˆ(ξ) −

α

e(iξ) (1−x) ˆ α f (ξ), 1+β|e(iξ) |2 δ

kωδβ (x, ·)k ≤ kωδβ (0, ·)k1−x kωδβ (1, ·)kx . Proof.

(3.2)

α

e(iξ) =e fˆδ (ξ), 1 + β|e(iξ)α |2 1 ωδβ (1, ξ) = fˆ(ξ) − fˆδ (ξ). 1 + β|e(iξ)α |2

ωδβ (0, ξ)

kωδβ (x, ·)k2

then the follow-

Z

∞

(iξ)α

fˆ(ξ) −

α

(3.3) (3.4)

e(iξ) (1−x) ˆ = |e fˆ(ξ) − fδ (ξ)|2 dξ 1 + β|e(iξ)α |2 −∞ Z ∞ 1 α fˆ (ξ)|2 dξ = |e(iξ) (1−x) |2 |fˆ(ξ) − (iξ)α |2 δ 1 + β|e Z−∞ ∞ απ 1 α = e2|ξ| (1−x) cos( 2 ) |fˆ(ξ) − fˆδ (ξ)|2 dξ 1 + β|e(iξ)α |2 −∞ Z ∞ απ 1 α = e2|ξ| (1−x) cos( 2 ) |fˆ(ξ) − fˆ (ξ)|2(1−x) (iξ)α |2 δ 1 + β|e −∞ 1 2x ˆ · |fˆ(ξ) − dξ α 2 fδ (ξ)| (iξ) 1 + β|e | Z 1−x  ∞
1 απ 1 α 2(1−x) 1−x ˆ f (ξ)) dξ ≤ e2|ξ| (1−x) cos( 2 ) (fˆ(ξ) − δ 1 + β|e(iξ)α |2 −∞ Z ∞  x

1 1 ˆδ (ξ) 2x x dξ · fˆ(ξ) − f α 1 + β|e(iξ) |2 −∞ Z ∞ 1−x απ 1 α 2 ˆ = e2|ξ| cos( 2 ) (fˆ(ξ) − f (ξ)) dξ δ 1 + β|e(iξ)α |2 −∞ Z ∞
2 x 1 ˆ · fˆ(ξ) − f (ξ) dξ δ 1 + β|e(iξ)α |2 −∞ (iξ)α (1−x)

8

= kωδβ (0, ·)k2(1−x) kωδβ (1, ·)k2x So

kωδβ (x, ·)k ≤ kωδβ (0, ·)k1−x kωδβ (1, ·)kx .

(3.5)

2 Lemma 3.4. If β is the solution of Eq. (3.1), then the following inequality also holds: δ E √ ≤ . 2(τ − 1) β

(3.6)

Proof. Due to (3.1), we obtain α

τδ = k

(iξ) |2 1 ˆδ (ξ) − fˆδ (ξ)k = k β|e f fˆδ (ξ)k α 1 + β|e(iξ) |2 1 + β|e(iξ)α |2 α

=k

α

(iξ) |2 β|e(iξ) |2 ˆδ (ξ) − fˆ(ξ)) + β|e ( f fˆ(ξ)k α 1 + β|e(iξ) |2 1 + β|e(iξ)α |2 α

α

≤k

(iξ) |2 β|e(iξ) |2 ˆδ (ξ) − fˆ(ξ))k + k β|e ( f fˆ(ξ)k α 1 + β|e(iξ) |2 1 + β|e(iξ)α |2 α

β|e(iξ) | ≤ δ + sup | √ (iξ)α |E | ξ∈R 2 β|e √ β ≤δ+ E. 2 So

δ E √ ≤ . 2(τ − 1) β

(3.7) 2

Lemma 3.5. The following inequality holds: kωδβ (1, ·)k ≤ (τ + 1)δ.

(3.8)

Proof: Due to the triangle inequality and (3.1), we obtain kωδβ (1, ·)k = kfˆ(ξ) −

1 fˆδ (ξ)k 1 + β|e(iξ)α |2

≤ kfˆ(ξ) − fˆδ (ξ)k + kfˆδ (ξ) − ≤ δ + τ δ = (τ + 1)δ.

1 fˆδ (ξ)k 1 + β|e(iξ)α |2 2

Lemma 3.6. The following inequality holds: kωδβ (0, ·)k ≤

9

4τ − 3 E. 4(τ − 1)

(3.9)

Proof: Due to the triangle inequality and (3.6), we obtain α

kωδβ (0, ·)k

(iξ)α

fˆ(ξ) −

(iξ)α

fˆ(ξ) −

= ke

e(iξ) fˆδ (ξ)k 1 + β|e(iξ)α |2 α

≤ ke

α

e(iξ) e(iξ) ˆ(ξ)k + k f (fˆ(ξ) − fˆδ (ξ))k α 1 + β|e(iξ) |2 1 + β|e(iξ)α |2 α

α

e(iξ) β|e(iξ) |2 k + sup | |δ. ≤ ke fˆ(ξ) α (iξ)α |2 1 + β|e(iξ) |2 ξ∈R 1 + β|e δ E 4τ − 3 ≤E+ √ ≤E+ = E. 4(τ − 1) 4(τ − 1) 2 β (iξ)α

2

Now we give the main result of this section. Theorem 3.7. Assume the conditions (1.8), (1.9) hold and take the solution β of Eq. (3.1) as the regularization parameter, then there holds the following error estimate: ku(x, ·) − uβδ (x, ·)k ≤ Cδ x E 1−x ,

(3.10)

4τ −3 1−x where C = (τ + 1)x ( 4(τ . −1) )

Proof: By the Parseval formula, (3.2), (3.8) and (3.9), we obtain u(x, ·) − u ˆβδ (x, ·)k = kωδβ (x, ·)k ku(x, ·) − uβδ (x, ·)k = kˆ

≤ kωδβ (0, ·)k1−x kωδβ (1, ·)kx 4τ − 3 ≤( E)1−x ((τ + 1)δ)x = Cδ x E 1−x , 4(τ − 1)

4τ −3 1−x where C = (τ + 1)x ( 4(τ . The proof of Theorem 3.7 is completed. −1) )

4

2

Several numerical examples

In this section, we present three different type examples to verify the usefulness of the proposed methods. We use the Fast Fourier and inverse Fast Fourier to complete our numerical experiments. The numerical examples are constructed in the following way: Firstly, we consider the following direct (forward) problem for the given data g(t):   −ux (x, t) = ∂0α+ u(x, t), x > 0, t > 0,    u(x, 0) = 0, x ≥ 0, (4.1)  u(0, t) = g(t), t ≥ 0,    u(x, t)| t ≥ 0. x→∞ = 0,

10

This problem is well-posed and its solution at x = 1 is given by Z ∞ 1 α f (t) = u(1, t) = √ e−iξt e−(iξ) gˆ(ξ)dξ. 2π −∞ 4

(4.2)

3

3

u(x,t) and its approximations

u(x,t) and its approximations

2 2 1 Exact solution ε=0.5 ε=0.1

0 −1 −2

1

Exact solution ε=0.5 ε=0.1

0

−1

−2 −3 −4

0

0.2

0.4

0.6

0.8

−3

1

0

0.2

0.4

t

(a)

0.6

0.8

1

t

(b) Fig.2. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.1 with Example 1: (a)x = 0.1, (b)x = 0.5.

5

3

4 2 u(x,t) and its approximations

u(x,t) and its approximations

3 2 1 Exact solution ε=0.5 ε=0.1

0 −1 −2 −3

1 Exact solution ε=0.5 ε=0.1

0

−1

−2

−4 −5

(a)

0

0.2

0.4

0.6

0.8

−3

1

0

0.2

t

0.4

0.6

0.8

1

t

(b) Fig.3. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.3 with Example 1: (a)x = 0.1, (b)x = 0.5.

Then we can obtain the perturbation data through fδ = f + ε randn(size(f )),

(4.3)

where the function “randn(·)” generates arrays of random numbers whose elements are normally distributed with mean 0, variance σ 2 = 1. The error is given by M +1 X 1 1 δ = kfδ − f kl2 = ( |fδ (ti ) − f (ti )|2 ) 2 , M +1

(4.4)

i=1

here we usually choose M = 100. Finally, we obtained the regularization solutions through the inverse problem. The bisection method is used to solve the Eq. (3.1), where we choose τ = 1.1. Example 1. Consider a smooth function g(t) = 4 sin(2πt).

11

3

3

2

2

u(x,t) and its approximations

u(x,t) and its approximations

4

1 Exact solution ε=0.5 ε=0.1

0 −1 −2

Exact solution ε=0.5 ε=0.1

0

−1

−2

−3

−3 −4

1

0

0.2

0.4

0.6

0.8

−4

1

0

0.2

0.4

t

(a)

0.6

0.8

1

t

(b) Fig.4. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.1 with Example 1: (a)x = 0.1, (b)x = 0.5.

5

4

4

3

u(x,t) and its approximations

u(x,t) and its approximations

3 2 1 0

Exact solution ε=0.5 ε=0.1

−1 −2

2 1 Exact solution ε=0.5 ε=0.1

0 −1 −2

−3 −3

−4 −5

(a)

0

0.2

0.4

0.6

0.8

−4

1

0

0.2

t

0.4

0.6

0.8

1

t

(b) Fig.5. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.3 with Example 1: (a)x = 0.1, (b)x = 0.5.

Example 2. Consider a piecewise smooth   0,    4t − 1, g(t) = 3 − 4t,    0,

function: 0 ≤ t ≤ 0.25, 0.25 < t ≤ 0.5, 0.5 < t ≤ 0.75, 0.75 < t ≤ 1.

Example 3. Consider the following discontinuous function:  1  0, 0 ≤ t ≤ 3 , g(t) = 1, 13 < t ≤ 23 ,   0, 23 < t ≤ 1.

(4.5)

(4.6)

Figs. 2-3, 6-7, 10-11 show the comparison of the numerical effects between the exact solution and its computed approximations for the a priori regularization parameter choose rule with Examples 1, 2, 3, respectively. Figs. 4-5, 8-9, 12-13 show the comparison of the numerical effects between the exact solution and its computed approximations for the a posteriori regularization parameter

12

1.2

0.8 Exact solution ε=0.5 ε=0.1

1

Exact solution ε=0.5 ε=0.1

0.7

u(x,t) and its approximations

u(x,t) and its approximations

0.6 0.8

0.6

0.4

0.2

0.5 0.4 0.3 0.2 0.1

0

−0.2

0

0

0.2

0.4

0.6

0.8

−0.1

1

0

0.2

0.4

t

(a)

0.6

0.8

1

t

(b) Fig.6. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.1 with Example 2: (a)x = 0.1, (b)x = 0.5.

2

0.8 Exact solution ε=0.5 ε=0.1

0.6 u(x,t) and its approximations

u(x,t) and its approximations

1.5

Exact solution ε=0.5 ε=0.1

0.7

1

0.5

0

0.5 0.4 0.3 0.2 0.1 0

−0.5

−0.1 −1

(a)

0

0.2

0.4

0.6

0.8

−0.2

1

t

0

0.2

0.4

0.6

0.8

1

t

(b) Fig.7. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.3 with Example 2: (a)x = 0.1, (b)x = 0.5.

choose rule with Examples 1, 2, 3, respectively. From Figs. 2-13, we firstly find that the smaller ε, the better the computed approximation is, the smaller the x is, the worse the computed approximation is, and the smaller the α is, the better the computed approximation is. Moreover, we can also easily find that the a posteriori parameter choice rule works better than the a priori parameter choice rule. Numerical test results are consistent with our theoretical results. Thirdly, from Figs. 2-13, it can be seen that the numerical solutions of Examples 2-3 are less ideal than these of Example 1. It is not difficult to see that the well-known Gibbs phenomenon and the recovered data near the non-smooth and discontinuities points are not accurate. Taking into consideration of the ill-posedness of the problem, the results presented in Figs. 2-13 are reasonable. In order to explain the effectiveness of the simplified Tikhonov regularization method for solving backward time-fractional inverse diffusion problem, we give the other regularization methods to compare with the simplified Tikhonov regularization method. We only use example 1 to computer the relative error between the regularization solution and the exact solution. The relative error is defined as follows: Error =

kregularization solution − exact solutionk . kexact solutionk

13

(4.7)

1.2

0.7 Exact solution ε=0.5 ε=0.1

0.8

0.6

0.4

0.2

0

−0.2

Exact solution ε=0.5 ε=0.1

0.6 u(x,t) and its approximations

u(x,t) and its approximations

1

0.5

0.4

0.3

0.2

0.1

0

0.2

0.4

0.6

0.8

0

1

0

0.2

0.4

t

(a)

0.6

0.8

1

t

(b) Fig.8. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.1 with Example 2: (a)x = 0.1, (b)x = 0.5.

1.5

1.2 Exact solution ε=0.5 ε=0.1

Exact solution ε=0.5 ε=0.1

1 u(x,t) and its approximations

u(x,t) and its approximations

1

0.5

0

0.8

0.6

0.4

0.2

−0.5 0

−1

0

0.2

0.4

0.6

0.8

−0.2

1

0

0.2

0.4

t

0.6

0.8

1

t

(b) Fig.9. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.3 with Example 2: (a)x = 0.1, (b)x = 0.5.

(a)

1.2

1

Exact solution ε=0.5 ε=0.1

1

Exact solution ε=0.5 ε=0.1

0.9

u(x,t) and its approximations

u(x,t) and its approximations

0.8 0.8 0.6 0.4 0.2 0

0.7 0.6 0.5 0.4 0.3 0.2

−0.2 −0.4

(a)

0.1 0

0.2

0.4

0.6

0.8

0

1

t

0

0.2

0.4

0.6

0.8

1

t

(b) Fig.10. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.1 with Example 3: (a)x = 0.1, (b)x = 0.5.

The quasi-reversibility regularization solution is obtained through solving the following problem:  α  −ux (x, t) = ∂0+ u(x, t) + µutt (x, t), x > 0, t > 0,   u(x, 0) = 0, x ≥ 0, (4.8)  u(1, t) = f (t), t ≥ 0,  δ   u(x, t)| t ≥ 0, x→∞ = 0,

14

2

1.2 Exact solution ε=0.5 ε=0.1 u(x,t) and its approximations

u(x,t) and its approximations

1.5

1

0.5

0

−0.5

−1

Exact solution ε=0.5 ε=0.1

1

0.8

0.6

0.4

0.2

0

0

0.2

0.4

0.6

0.8

−0.2

1

0

0.2

0.4

t

(a)

0.6

0.8

1

t

(b) Fig.11. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.3 with Example 3: (a)x = 0.1, (b)x = 0.5.

1.2

0.9 Exact solution ε=0.5 ε=0.1

1

Exact solution ε=0.5 ε=0.1

0.8

u(x,t) and its approximations

u(x,t) and its approximations

0.7 0.8 0.6 0.4 0.2 0

0.6 0.5 0.4 0.3 0.2 0.1

−0.2 −0.4

0 0

0.2

0.4

0.6

0.8

−0.1

1

0.2

0.4

0.6

0.8

1

t

(b) Fig.12. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.1 with Example 3: (a)x = 0.1, (b)x = 0.5.

2

1.2 Exact solution ε=0.5 ε=0.1

1

0.5

0

−0.5

−1

Exact solution ε=0.5 ε=0.1

1 u(x,t) and its approximations

u(x,t) and its approximations

1.5

(a)

0

t

(a)

0.8

0.6

0.4

0.2

0

0

0.2

0.4

0.6

0.8

−0.2

1

0

0.2

0.4

t

0.6

0.8

1

t

(b) Fig.13. The comparison of the numerical effects between the exact solution and its computed approximations for α = 0.3 with Example 3: (a)x = 0.1, (b)x = 0.5.

where µ is the regularization parameter. Using the Fourier transform about variable t, we obtain the quasi-reversibility regularization solution as follows: α

u ˆµδ (x, ξ) = So 1 uµδ (x, t) = √ 2π

Z

e(iξ) (1−x) ˆ fδ (ξ). eµ|ξ|2 (1−x) ∞

(4.9)

α

eiξt

−∞

15

e(iξ) (1−x) ˆ fδ (ξ)dξ. eµ|ξ|2 (1−x)

(4.10)

The quasi boundary value regularization solution is obtained through solving the following problem:   −ux (x, t) = ∂0α+ u(x, t), x > 0, t > 0,    u(x, 0) = 0, x ≥ 0, (4.11)  γu(0, t) + u(1, t) = fδ (t), t ≥ 0,    u(x, t)| t ≥ 0, x→∞ = 0, where γ is the regularization parameter. Using the Fourier transform about variable t, we obtain the quasi boundary value regularization solution as follows: α

u ˆγδ (x, ξ) = So uγδ (x, t)

1 =√ 2π

Z

e(iξ) (1−x) ˆ fδ (ξ). 1 + γe(iξ)α

(4.12)

α

∞

eiξt

−∞

e(iξ) (1−x) ˆ fδ (ξ)dξ. 1 + γe(iξ)α

(4.13)

Table 1 Relative errors for different ε under the a priori choice rule for three regularization methods about Example 1 (α = 0.1, x = 0.1) ε 0.1 0.05 0.01 0.005 0.001 Simplified Tikhonov method 0.0122 0.0057 0.0011 4.7743e-004 9.1554e-005 quasi-reversibility method 0.095 0.0087 0.0030 5.1033e-004 1.1106e-004 quasi boundary value method 0.094 0.0096 0.0142 4.7547e-004 1.0837e-004 Table 2 Relative errors for different ε under the a priori choice rule for three regularization methods about Example 1 (α = 0.6, x = 0.1) ε 0.1 0.05 0.01 0.005 0.001 Simplified Tikhonov method 0.8659 0.7386 0.1745 0.0774 0.0108 quasi-reversibility method 1.7425 0.7543 0.1069 0.0827 0.0152 quasi boundary value method 1.5975 0.6772 0.2052 0.0951 0.0129 From tables 1 and 2, we can see that the simplified Tikhonov regularization method is better than the other two regularization methods for solve this inverse problem.

5

Conclusions

In this paper, a simplified Tikhonov regularization method for a time-fractional inverse diffusion problem is given. The a priori and the a posteriori rules for choosing regularization parameter with strict theory analysis are presented. For a priori parameter choice rule, the regularization parameter depends on the noise level and the a priori bound. But in practice, the a priori bound is unknown exactly and the regularization

16

parameters can not be obtained exactly. Using the Morozov’s discrepancy principle, we give the a posteriori parameter choice rule which only depends on the measured data. Meanwhile, the numerical examples verified the efficiency and accuracy of the method. We only considered the temperature for 0 < x < 1 through a measurable temperature at x = 1 in this paper. This method will also be used to consider the heat flux for 0 ≤ x < 1.

References [1] E. Barkai, R. Metzler, J. Klafter, From continuous time random walks to the fractional Fokker-Planck equation, Physical Review E 61(2000)132-138. [2] A. Blumen, G. Zumofen, J. Klafter, Transport aspects in anomalous diffusion, Levy walks 40(1989)3964-3973. [3] A. Carasso, Determining surface temperature from interior observations, SIAM Journal on Applied Mathematics 42(1982)558-574. [4] A.S. Chaves, A fractional diffusion equation to descrive levy flights, Physics Letters A 239(1998)13-16. [5] J. Cheng, J. Nakagawa, M. Yamamoto, T. Yamazaki, Uniqueness in an inverse problem for one-dimensional fractional diffusion equation, Inverse Problems 16(2009)115002.(16pp). [6] W. Cheng, C.L. Fu, Z. Qian, A modified Tikhonov regularization method for a spherically symmetric three-dimensional inverse heat conduction problem, Mathematics and Computers in Simulation 75(2007)97-112. [7] H.W. Engl, M. Hanke, A. Neubauer, Regularization of Inverse Problem, Kluwer Academic, Boston, MA, 1996. [8] X.L. Feng, C.L. Fu, H. Cheng, A regularization method for solving the Cauchy problem for the Helmholtz equation, Applied Mathematical Modelling 35(2011)3301-3315. [9] C.L. Fu, Simplified Tikhonov and Fourier regularization methods on a general sideways parabolic equation, Journal of Computational and Applied Mathematics 167(2004)449-463. [10] A. Kirsch, An Introduction to the Mathematical Theory of Inverse Problems, Springer-Verlag, New York, 1996. [11] J.J. Liu, M. Yamamoto, A backward problem for the time-fractional diffusion equation, Applicable Analysis 89(11)(2010)1769-1788. [12] R. Metzler, J. Klafter, Boundary value problems for fractional diffusion equations, Physica A 278(2000)107-125.

17

[13] D.A. Murio, Stable numerical solution of a fractional-diffusion inverse heat conduction problem, Computers and Mathematics with Applications 53(10)(2007)1492-1501. [14] D.A. Murio, Time fractional IHCP with Caputo fractional derivatives, Computers and Mathematics with Applications 56(9)(2008)2371-2381. [15] I. Podlubny, Fractional Differential Equations, Acadmic Press, San Diego, 1999. [16] Z. Qian, Optimal modified method for a fractional-diffusion inverse heat conduction problem, Inverse Problems in Science and Engineering 18(4)(2010)521-533. [17] M. Raberto, E. Scalas, F. Mainardi, Waiting-times and returns in high-frequency financial data: an empirical study, Physica A 314(2002)749-755. [18] J.G. Wang, Y.B. Zhou, T. Wei, Two regularization methods to identify a spacedependent source for the time-fractional diffusion equation, Appl. Numer. Math. 68(2013)39-57. [19] X.T. Xiong, H.B. Guo, X.H. Liu, An inverse problem for a fractional diffusion equation, Journal of Computational and Applied Mathematics 236(2012)44744484. [20] F. Yang, C.L. Fu, A simplified Tikhonov regularization method for determining the heat source, Applied Mathematical Modelling 34(2010)3286-3299. [21] F. Yang, C.L. Fu, The quasi-reversibility regularization method for identifying the unknown source for time fractional diffusion equation, Applied Mathematical Modelling 39 (2015)1500-1512. [22] F. Yang, C.L. Fu, X.X. Li, A mollification regularization method for unknown source in time-fractional diffusion equation, International Journal of Computer Mathematics, 91(7)(2014)516-1534. [23] G.H. Zheng, T. Wei, A new regularization method for solving a time-fractional inverse diffusion problem, Journal of Mathematical Analysis and Applications 378(2)(2011)418-431. [24] G.H. Zheng, T. Wei, Spectral regularization method for the time fractional inverse advection-dispersion equation, Mathematics and Computers in Simulation 81(2010)37-51. [25] G.H. Zheng, T. Wei, Two regularization methods for solving a Riesz-Feller spacefractional backward diffusion problem, Inverse Problems. 26(2010)115017.(22pp). [26] G.H. Zheng, T. Wei, A new regularization method for a Cauchy problem of the time fractional diffusion equation, Advances in Computational Mathematics 36(2012)377-398.

18

[27] Z.Q. Zhang, T. Wei, Identifying an unknown source in time-fractional diffusion equation by a truncation method, Applied Mathematics and Computation 219(11)(2013)5972-5983. [28] Z.Q. Zhang, T. Wei, An optimal regularization method for space-fractional backward diffusion problem, Mathematics and Computers in Simulation. 92(2013)1427.

19