AN UNDETERMINED TIME-DEPENDENT COEFFICIENT IN A FRACTIONAL DIFFUSION EQUATION

. In this work, we consider a FDE (fractional diﬀusion equation) with a time-dependent diﬀusion coeﬃcient a ( t ). This is an extension of [13], which deals with this FDE in one-dimensional space. For the direct problem, given an a ( t ) , we establish the existence, uniqueness and some regularity properties with a more general domain Ω and right-hand side F ( x,t ). For the inverse problem–recovering a ( t ) , we introduce an operator K one of whose ﬁxed points is a ( t ) and show its monotonicity, uniqueness and existence of its ﬁxed points. With these properties, a reconstruction algorithm for a ( t ) is created and some numerical results are provided to illustrate the theories.


1.
Introduction. This paper considers the fractional diffusion equation (FDE) with a continuous and positive coefficient function a(t) : (1) C D α t u(x, t) − a(t)Lu(x, t) = F (x, t), x ∈ Ω, t ∈ (0, T ]; u(x, t) = 0, (x, t) ∈ ∂Ω × (0, T ]; u(x, 0) = u 0 (x), where Ω is a bounded and smooth subset of R n , n = 1, 2, 3, −L is a symmetric uniformly elliptic operator defined as (a ij (x)u xi ) xj + c(x)u with conditions (2) a ij , c ∈ C 2 (Ω) (i, j = 1, . . . , n), ∂Ω is C 3 , and C D α t is the left-sided Djrbashian-Caputo α-th order derivative with respect to time t. The definition for C D α t is with Gamma function Γ(·) and the nearest integer n with α ≤ n. In this paper, we are assuming a subdiffusion process, i.e. α ∈ (0, 1). This simplifies the definition of This work is an extension of [13] from a simple space domain Ω to R n , considers the more general analysis for the direct problem and contains an existence argument for the inverse problem of recovering a(t). This paper consists of two parts; the direct problem and the inverse problem. For the direct problem, we build the spectral representation of the weak solution u(x, t; a). The notation u(x, t; a) is used for displaying the dependence of the solution u on the diffusivity a(t). Then the existence, uniqueness and regularity results are proved with several assumptions on the coefficient function a(t). Unlike [13], the right hand side function F (x, t) is not of the form f (x)g(t), so that the proof of regularity is more delicate. For the inverse problem, we use the single point flux data a(t) ∂u ∂ − → n (x 0 , t; a) = g(t), x 0 ∈ ∂Ω to recover the coefficient a(t) (We choose the data a(t) ∂u ∂ − → n (x 0 , t; a) = g(t) instead of the classical flux ∂u ∂ − → n (x 0 , t; a) because in practice, a(t) ∂u ∂ − → n (x 0 , t; a) is usually measured as the flux). For the reconstruction, we only consider to recover a continuous and positive a(t) to match the assumptions set in the direct problem. Acting a flux data, we introduce an operator K one of whose fixed points is the coefficient a(t). Using the weak maximum principle [7], we establish the monotonicity and uniqueness of the fixed points of operator K, and the proof of uniqueness leads to a numerical reconstruction algorithm. Since we consider a multidimensional domain Ω here, the Sobolev Embedding Theorem yields that we need to add the condition (2) on the operator −L to ensure the C 1 -regularity of the series representation of u. Then the operator K is well-defined, where the proofs can be seen in section 4. This is a significant difference from [13]. Furthermore, an existence argument of the fixed points of K is included by this paper, which [13] does not contain.
The rest of this paper follows the following structure. In section 2, we collect some preliminary results about fractional calculus and the eigensystem of −L. The direct problem is discussed in section 3, i.e. we establish the existence, uniqueness and some regularity results of the weak solution for FDE (1). Then section 4 deals with the inverse problem of recovering a(t). Specifically, an operator K is introduced at the beginning of this section, then its monotonicity and uniqueness of its fixed points give an algorithm to recover the coefficient a(t). In particular, the existence argument of the fixed points of K is included by this section. In section 5, some numerical results are presented to illustrate the theoretical basis.
The following lemma about the composition between C D α t and the fractional integral I α t is presented in [11].

3.2.
Existence and uniqueness. In order to show the existence and uniqueness of the weak solution (4), we state the following lemma [5,Theorem 3.25].
then there exists a unique solution y(t) for the Cauchy-type problem, which satisfies C D α t y ∈ C[0, T ]. The theorem of existence and uniqueness for u(x, t; a) follows from Lemma 3.4.
Theorem 3.5 (Existence and Uniqueness). Suppose Assumption 3.1 holds. Under Definition 3.2, there exists a unique weak solution u(x, t; a) of FDE (1) with the spectral representation (4) and for each n ∈ N + , u n (t; a) ∈ C[0, T ] is the unique solution of the fractional ODE (5) with C D α t u n (t; a) ∈ C[0, T ]. Proof. From the spectral representation (4), it suffices to show the existence and uniqueness of u n (t; a), n ∈ N + . Fix n ∈ N + , Assumption 3.1 (a) and (b) yield that the fractional ODE (5) satisfies the conditions of Lemma 3.4. Hence the existence and uniqueness for u n (t; a) hold.

3.3.
Sign of u n (t; a). In this part, we state two properties of u n (t; a) which play important roles in building the regularity of u(x, t; a).
The following corollary, which concerns the sign of u n (t; a), follows from Lemma 3.6 directly. Proof. Assumption 3.1 gives that λ n a(t) ∈ C + [0, T ]. Then the proof is completed by applying Lemma 3.6 to the fractional ODE (5).
Denote the weak solutions of FDEs (7) and (8) by u r (x, t; a) and u i (x, t; a), respectively ("r" and "i" denote the initials of "right-hand side" and "initial condition"). The following lemma about u r (x, t; a) and u i (x, t; a) follows from Lemma 3.3 and Theorem 3.5.
Lemma 3.7. Suppose Assumption 3.1 holds. Then u r (x, t; a) and u i (x, t; a) are the unique solutions for FDEs (7) and (8), respectively, with the spectral representations as where u r n (t; a), u i n (t; a) satisfy the following fractional ODEs (10) C D α t u r n (t; a) + λ n a(t)u r n (t; a) = F n (t), u r n (0; a) = 0, n ∈ N + ; (11) C D α t u i n (t; a) + λ n a(t)u i n (t; a) = 0, u i n (0; a) = b n , n ∈ N + . Moreover, Theorem 3.5 ensures the weak solution u(x, t; a) of FDE (1) can be written as u(x, t; a) = u r (x, t; a) + u i (x, t; a), i.e. u n (t; a) = u r n (t; a) + u i n (t; a), n ∈ N + .
From Lemma 3.9, we obtain are the unique solutions of fractional ODEs (13) and (14) respectively with a(t) ≡ q a on [0, T ]. The next two lemmas concern the regularity of u r,+ (x, t; a) and C D α t u r,+ (x, t; a), respectively.
where the last inequality is obtained from (18). By the Monotone Convergence Theorem, we have For each n ∈ N + , [10] gives the explicit representation of u r,+ n (t; q a ) which together with Young's inequality leads to

Lemmas 2.2, 2.3 and 2.4 give the bound of
which together with (19) and the completeness of {φ n (x) : n ∈ N + } in L 2 (Ω) gives where the constant C only depends on a(t). This completes the proof.
The following corollary follows immediately from the proofs of Lemmas 3.10 and 3.11.
. From Lemmas 3.10, 3.11, Corollary 2 and (15), we are able to deduce the regularity for u r (x, t; a) and C D α t u r (x, t; a). Lemma 3.12 (Regularity of u r ).
If we impose a higher regularity on F, we can obtain the regularity estimate of u r C([0,T ];H 2 (Ω)) .
, where C depends on Ω, −L and a(t).
Proof. For each t ∈ [0, T ], we have Hence, Inverse Problems and Imaging Volume 11, No. 5 (2017), 875-900 By [10, Lemma 3.4], we have , and the constant C depends on Ω, −L and a(t). Similarly, we can show From the above proof for u r, . The estimates of u r,+ , u r,− , C D α t u r,+ and C D α t u r,− yield the desired result and complete this proof.

3.4.2.
Regularity of u i . In this part we consider the regularity of u i . Just as in the regularity results for u r , we first state two lemmas which concern the positivity and monotonicity of u i , respectively. Lemma 3.13. With the representation (9) and the fractional ODE (11), for each Proof. This is a directly result of Corollary 1.
Proof. Fix n ∈ N + , from the fractional ODE (11), the functions u i n (t; a 1 ) and u i n (t; a 2 ) satisfy the following system If b n ≥ 0, Corollary 1 shows that u i n (t; a 1 ), u i n (t; a 2 ) ≥ 0. Also, Lemma 3.13 and a 1 ≤ a 2 ensures the right side of (21) is nonnegative, which together with Proof. Given t ∈ [0, T ], the direct calculation and Lemma 3.14 yield that Recall that [10] established the representation as u i n (t; q a ) = b n E α,1 (−λ n q a t α ), n ∈ N + . Hence, by Lemma 2.1, which leads to u i 2 L 2 (0,T ;H 2 (Ω)) ≤ C For the estimate of C D α t u i (x, t; a), (9) and (11) yield which together with (17) gives , where the last inequality follows from (22). This result implies that , which together with (23) completes the proof.
Moreover, with a stronger condition on u 0 , such as assuming u 0 ∈ H 2 (Ω)∩H 1 0 (Ω), we can deduce the C-regularity estimate of u i .
meanwhile, the following estimates have been shown in the proof of Theorem 3.15  Moreover, if the conditions u 0 ∈ H 2 (Ω) ∩ H 1 0 (Ω) and F ∈ C θ ([0, T ]; L 2 (Ω)), 0 < θ < 1 are added, we have: 4. Inverse problem-reconstruction of the diffusion coefficient a(t). In this section, we discuss how to recover the coefficient a(t) through the output flux data All cross the inverse problem work, the operator −L is assumed to satisfy the condition (2), then the expression ∂φn ∂ − → n (x 0 ) makes sense. We only consider this reconstruction in the space C + [0, T ], which can be regarded as the admissible set for a(t). To this end, we introduce an operator K, which will be shown to have a fixed point consisting of the desired coefficient a(t).

4.1.
Operator K. The operator K is defined as To analyze K, we make the following assumptions.
Assumption 4.1. u 0 , F and g should satisfy the following restrictions: The next remark shows that the equality in the definition of K is valid.

Remark 2.
For the inverse problem, the right-hand side function F (x, t) and the initial condition u 0 (x) are input data, which, at least in some circumstance, can be assumed to be controlled. Even though Assumption 4.1 (a), (b) and (c) appear restrictive, it is not hard to construct functions that satisfy them. For example, in (a) if u 0 = cφ k for some c > 0, then Assumption 4.1 (a) will be satisfied. This will also be true if is also a linear combination of {φ n : n ∈ N + } with positive coefficients. For (c), by the completeness of {φ n : n ∈ N + } in L 2 (Ω), there should exist N ∈ N + s.t.
This together with a ∈ C + [0, T ] gives that g > 0. The continuity of g follows from the ones of a and u n (t; a), n ∈ N + , which are derived from the admissible set C + [0, T ] and Theorem 3.5, respectively. Therefore, Assumption 4.1 (d) is reasonable and can be attained.
Recall that the numerator g > 0, so that the lower bound g(t) ∂u0 > 0, which gives that D(K) is a subspace of C + [0, T ]. Also, The strict positivity of the dominator avoids D(K) degenerating to an empty set. In order to show the well-definedness of K, Assumption 4.1 (a), (b) and (c) will be used. Furthermore, Assumption 4.1 (a) and (b) are crucial to build the monotonicity of operator K; meanwhile, Assumption 4.1 (c) is stated for the uniqueness of fixed points of K.
For the operator K, we have the following lemmas. Proof. For each ψ ∈ D(K), Theorem 3.5 ensures that there exists a unique u n (t; ψ) for n ∈ N + , which implies the existence and uniqueness of Kψ.
Then it is suffice to show the dominator which completes the proof. Proof. Given ψ ∈ D(K). The continuity of Kψ follows from the continuity of u n (t; ψ) for each n ∈ N + and the continuity of g, which are established by Theorem 3.5 and Assumption 4.1 (d) respectively. For each n ∈ N + , (5) ensures u n (t; ψ) satisfies C D α t u n (t; ψ) + λ n ψ(t)u n (t; ψ) = F n (t), u n (0; ψ) = b n . Taking I α t on both sides of the above ODE and using Lemma 2.6 yield that u n (t; ψ) + λ n I α t [ψ(t)u n (t; ψ)] = I α t F n + b n . From the proof of Lemma 4.2, we have u n (t; ψ) ≥ 0 on [0, T ], which together with λ n > 0, the positivity of ψ and the definition of I α t yields that λ n I α t [ψ(t)u n (t; ψ)] ≥ 0. Since u n (t; ψ) ≥ 0 and λ n I α t [ψ(t)u n (t; ψ)] ≥ 0, we deduce that 0 ≤ u n (t; ψ) ≤ which together with g > 0 yields that where the last inequality follows from Remark 3. The above result and the continuity of Kψ lead to Kψ ∈ D(K), which is the expected result.

4.3.
Uniqueness. In order to show the uniqueness, we state two lemmas.
Lemma 4.5. If a 1 , a 2 ∈ D(K) are both fixed points of K with a 1 ≤ a 2 , then a 1 ≡ a 2 .
Proof. Pick a fixed point a(t), then which gives by taking I α t on both sides. Similarly, taking I α t on the both sides of (5) and applying Lemma 2.6 yield that ; a), n ∈ N + , which together with (25) generates In (26), the convergence of the two series in C[0, T ] is supported by Assumption 4.1, Remark 1 and the fact that 0 < λ 1 ≤ λ 2 ≤ · · · . Given two fixed points a 1 , a 2 with a 1 ≤ a 2 , then a 1 and a 2 should satisfy (26) simultaneously, which gives In the proof of Theorem 4.4, we have shown that u n (t; a 1 ) ≥ u n (t; a 2 ) ≥ 0. Also recall that λ −1 a 1 (t)) ≡ 0 on [0, T ]; while the proof of Lemma 4.2 yields that u N (t; a 2 ) > 0. Hence, we have a 1 = a 2 on [0, T ], which completes the proof.
Before showing uniqueness, we introduce a successive iteration procedure which will generate a sequence converging to a fixed point if it exists. Set , a n+1 = Ka n , n ∈ N.
Then this iteration reproduces a sequence {a n : n ∈ N} which is contained by D(K) due to Lemma 4.3.
Lemma 4.6. If there exists a fixed point a(t) ∈ D(K) of operator K, then the sequence {a n : n ∈ N} will converge to a(t).
Proof. a 0 is the lower bound of D(K) and {a n : n ∈ N} ⊂ D(K) yield that a 0 ≤ a 1 . Using Theorem 4.4, we have a 1 = Ka 0 ≤ Ka 1 = a 2 , i.e. a 1 ≤ a 2 . The same argument gives a 2 = Ka 1 ≤ Ka 2 = a 3 . Continue this process, we can deduce a 0 ≤ a 1 ≤ a 2 ≤ . . . , which means {a n : n ∈ N} is increasing. Since the results that a 0 is the lower bound of D(K) and a(t) ∈ D(K), it holds a 0 ≤ a. Applying Theorem 4.4 to this inequality, we obtain a 1 = Ka 0 ≤ Ka = a, i.e. a 1 ≤ a. This argument generates a n ≤ a, n ∈ N, which means a(t) is an upper bound of {a n : n ∈ N}.
We have proved {a n : n ∈ N} is an increasing sequence in D(K) with an upper bound a(t), which leads to {a n : n ∈ N} is convergent in D(K) and the limit is smaller than a(t). Denote the limit of {a n : n ∈ N} by a. We have a ∈ D(K), a ≤ a and a is a fixed point of K in D(K). Hence, Lemma 4.5 yields a = a, which is the desired result. Now, we are able to prove the uniqueness of fixed points of K. Proof. Let a 1 , a 2 ∈ D(K) be both fixed points of K. Lemma 4.6 implies that a n → a 1 and a n → a 2 , which leads to a 1 = a 2 and completes this proof.

4.4.
Existence. Assumption 4.1 is not sufficient to deduce the existence of the fixed points of K since D(K) has no upper bound so that an increasing sequence in D(K) may not be convergent. In this part, we discuss the existence of fixed points, by providing some extra conditions. Assumption 4.8. Additional assumptions on u 0 , F and g: Define the subspace D(K) of D(K) as  Proof. Given ψ ∈ D(K) , we have proved Kψ ∈ C + [0, T ] and For each n ∈ N + , let w n (t; ψ) = u n (t; ψ) − b n , (5) yields the following ODE by direct calculation and this proof is complete.
The existence conclusion is derived from Lemmas 4.6 and 4.9. Proof. Lemma 4.6 yields the sequence {a n : n ∈ N} is increasing, while Lemma 4.9 gives {a n : n ∈ N} ⊂ D(K) . Then {a n : n ∈ N} is an increasing sequence with an upper bound g(t) ∂u0 ∂ − → n (x 0 ) −1 , which implies the convergence of {a n : n ∈ N}.
Denote the limit by a, clearly a is a fixed point of K. Also, the closedness of D(K) yields that a ∈ D(K) . Therefore, a is a fixed point of K in D(K) , which confirms the existence.

4.5.
Main theorem for the inverse problem and reconstruction algorithm.  The following reconstruction algorithm for a(t) is based on Theorem 4.11.

Table 1. Numerical Algorithm
Iteration algorithm to recover the coefficient a(t) 1: Set up the right-hand side function F (x, t) and the initial condition u 0 (x), then measure the output flux data g(t). F , u 0 and g should satisfy Assumption 4.1; 2: Set the initial guess as a 0 (t) = g(t) ∂u0 ; 3: for k = 1,...,N do 4: Using the L1 time-stepping [3] to compute u(x, t; a k−1 ), which is the weak solution of FDE (1) with coefficient function a k−1 ; 5: Update the coefficient a k−1 by a k = Ka k−1 ; 6: Check stopping criterion a k − a k−1 L 2 [0,T ] ≤ 0 for some 0 > 0; 7: end for 8: output the approximate coefficient function a N .

5.
Numerical results for inverse problem.
The fourth step of Algorithm 1 includes solving the direct problem of FDE (1) numerically. To this end, we choose L1 time stepping [3,6] to discretize the term C D α t u(x, t) :

5.2.
Numerical results for noise free data. In this part, we set Ω = (0, 1), x 0 = 0, T = 1, Lu = u xx , pick u 0 (x) = − sin πx, F (x, t) = −(t + 1) sin πx and consider the following two coefficients: (a1) smooth coefficient: a(t) = sin 5πt + 1.3; (a2) nonsmooth coefficient ("smile" function): In experiment (a1), the exact coefficient we pick is a smooth function. Figure 1 shows the initial guess and the first three iterations, while Figure 2 presents the exact and approximate coefficients. From these two figures, we observe that {a n : n ∈ N} converges to a(t) monotonically, which illustrates Theorems 4.4 and 4.11. Moreover, the L 2 error of the approximation in Figure 2 is a−a N L 2 [0,T ] = 1.04×10 −6 , which implies us the L 2 error of this approximation may be bounded by the stopping criterion number 0 . This guess is confirmed by Figure 4 and can be expressed as Several attempts of experiment (a1) for different α ∈ (0, 1) are taken to find the dependence of the convergence rate of Algorithm 1 on the fractional order α, which is shown in Figure 3. This figure shows the amounts of iterations required, i.e. N, corresponding to different α, which imply that restricted α ∈ (0, 1), the larger α is, the faster the convergence rate of Algorithm 1 is. This phenomenon is explained in [4] by a property of the Mittag-Leffler function; for α ∈ (0, 1), the larger α is, the faster the decay rate of E α,1 (−z) is as z → ∞.
The definition of D(K) restricts the coefficient a(t) in the space C + [0, T ], however, the results of experiment (a2) indicate that Algorithm 1 still works for nonsmooth a(t), which means the numerical restriction on a(t) can possibly be extended from a(t) ∈ C + [0, T ] to a(t) ∈ L ∞ [0, T ]. For discontinuous a(t),

5.3.
Numerical results for noisy data. In this subsection, we will consider data polluted by noise. Set g be the exact data and denote the noisy data by g δ with relative noise level δ, i.e. (g − g δ )/g L ∞ [0,T ] ≤ δ. Then the perturbed operator K δ is  Also, the sequence {a δ,n : n ∈ N} can be obtained from the iteration , a δ,n+1 = K δ a δ,n , n ∈ N.  Since δ is a small positive number and g is a strictly positive function, we can assume g δ is still positive, which means Theorem 4.11 still holds for K δ . Hence, if there exists a fixed point a δ ∈ D(K δ ), the sequence {a δ,n : n ∈ N} will converge to a δ monotonically and we denote the limit by a δ . Algorithm 1 is still able to be used to recover a δ after a slightly modification−replacing g and K by g δ and K δ , respectively. We take the experiments (a1) and (a2) with noise level δ > 0. Figures 7 and 8 present the exact and approximate coefficients under δ = 3% for experiments (a1) and (a2) respectively. From figures 7 and 8, we observe that the smaller |a(t)| is, the better the approximation is. This can be explained by δ means the relatively noise level, i.e. we pick g δ = (1 + ζδ)g in the codes, where ζ follows a uniform distribution on [−1, 1]. Figure 9 illustrates that