Reconstruction of a Space-Time Dependent Source in Subdiffusion Models via a Perturbation Approach
RReconstruction of a Space-Time Dependent Source in SubdiffusionModels via a Perturbation Approach ∗ Bangti Jin † Yavar Kian ‡ Zhi Zhou § February 8, 2021
Abstract
In this article we study inverse problems of recovering a space-time dependent source componentfrom the lateral boundary observation in a subidffusion model. The mathematical model involves aDjrbashian-Caputo fractional derivative of order α ∈ (0 ,
1) in time, and a second-order elliptic operatorwith time-dependent coefficients. We establish a well-posedness and a conditional stability result for theinverse problems using a novel perturbation argument and refined regularity estimates of the associateddirect problem. Further, we present an algorithm for efficiently and accurately reconstructing the sourcecomponent, and provide several two-dimensional numerical results showing the feasibility of the recovery.
Key words : inverse source problem, subdiffusion, time-dependent coefficient, conditional stability, re-construction
This work is concerned with inverse source problems (ISPs) of identifying a space-time dependent componentof the source in the subdiffusion model in a cylindrical domain from the lateral Cauchy data on a part ofthe boundary. Let d ≥
2, Ω = ω × ( − (cid:96), (cid:96) ), ω ⊂ R d − be an open bounded domain with a C boundary, andfix T > x ∈ Ω, we write x = ( x (cid:48) , x d ), with x (cid:48) ∈ ω and x d ∈ ( − (cid:96), (cid:96) ). For m = 0 , u : ∂ αt u + A ( t ) u = F, in Ω × (0 , T ] ,u ( x,
0) = 0 , in Ω ,∂ mx d u ( x (cid:48) , (cid:96), t ) = 0 , on ω × (0 , T ) ,∂ x d u ( x (cid:48) , − (cid:96), t ) = 0 , on ω × (0 , T ) ,u ( x, t ) = 0 , on ∂ω × ( − (cid:96), (cid:96) ) × (0 , T ) . (1.1)In the model (1.1), α ∈ (0 ,
1) is fixed, and the notation ∂ αt u denotes the so-called Djrbashian-Caputofractional derivative of order α in time, which, for α ∈ (0 , ∂ αt u ( t ) = 1Γ(1 − α ) (cid:90) t ( t − s ) − α u (cid:48) ( s )d s, where Γ( z ) = (cid:82) ∞ s z − e − s d s for (cid:60) ( z ) > α approaches 1 − , thefractional derivative ∂ αt u recovers the usual first-order derivative u (cid:48) ( t ), and accordingly, the model coincides ∗ The work of B.J. is partially supported by UK EPSRC grant EP/T000864/1, that of Y.K. by the French National ResearchAgency ANR (project MultiOnde) grant ANR-17-CE40-0029, and that of Z.Z. by Hong Kong RGC grant (No. 15304420). † Department of Computer Science, University College London, Gower Street, London WC1E 6BT, UK ( [email protected] ) ‡ Centre de Physique Th´eorique (CPT), UMR-7332, Aix Marseille Universit´e, Campus de Luminy, Case 907, 13288 Marseillecedex 9, France ([email protected]) § Department of Applied Mathematics, The Hong Kong Polytechnic University, Hung Hom, Kowloon, Hong Kong([email protected]) a r X i v : . [ m a t h . NA ] F e b ith the standard diffusion equation. A ( t ) is a time-dependent second-order strongly elliptic operator,defined by A ( t ) u ( x ) = − d (cid:88) i,j =1 ∂ x i ( a ij ( x, t ) ∂ x j u ) + q ( t ) u, x ∈ Ω , where a = [ a ij ] di,j =1 ∈ C (Ω × [0 , T ]; R d × d ) is a symmetric matrix-valued function and satisfies suitableregularity conditions given in Assumption 2.1 below, and q ∈ C ([0 , T ]; L ∞ (Ω)) is nonnegative.The model (1.1) has received much attention in recent years, known under the name of “subdiffusion”or “time-fractional diffusion”, due to its extraordinary modeling capability for describing anomalously slowdiffusion processes arising in a wide range of applications in physics, engineering and biology. At a microscopiclevel, it can be derived from continuous time random walk with a heavy-tailed waiting time distribution (witha divergent mean) in the sense that the probability density function of the walker appearing at time t > x ∈ R d satisfies a differential equation of the form (1.1) (in the whole space R d ).The model (1.1) has been successfully employed in describing many practical applications, e.g., diffusion ofcharge carriers in amorphous photoconductors, diffusion in fractal domains [34], ion transport in columnexperiments [11], and subsurface flow [2]. We refer interested readers to the comprehensive reviews [32, 31]for physical motivations of the mathematical model and long lists of successful applications.The ISPs of interest are to determine some information of the source F from the measurement on a sub-boundary ω × { (cid:96) } ⊂ ∂ Ω of the domain Ω. Note that the boundary measurement is insufficient to uniquelydetermine a general source F (see, e.g., [22, Section 1.3.1]), and additional assumptions have to be imposedin order to restore unique recovery. Often it is formulated as the spatial or temporal component of the source F ( x, t ). In this work, the source F is assumed to be of the form F ( x, t ) = f ( x (cid:48) , t ) R ( x, t ) . (1.2)The condition (1.2) can be interpreted as that an unknown source f ( x (cid:48) , t ) depends only on the depth variable x (cid:48) and t in the case of d = 2, which corresponds to a layer structure, and on the planar location ( x , x ) and t but not on the depth in the case of d = 3, which can be a good approximation if the domain Ω is very thinin the direction of x . Note that it arises also naturally in linearizing the inverse potential problem, wherethe potential coefficient q depends only x (cid:48) and t [14]. We investigate the following two inverse problems: (i) ISPn is to recover f ( x (cid:48) , t ) from the boundary observation u | ω ×{ (cid:96) }× (0 ,T ) for m = 1 in (1.1) and (ii) ISPd isto recover f ( x (cid:48) , t ) from the flux measurement ∂ x d u | ω ×{ (cid:96) }× (0 ,T ) for m = 0 in (1.1) (i.e., n and d refer to theNeumann and Dirichlet boundary condition, respectively).This work is devoted to the theoretical analysis of ISPn and ISPd. In Theorem 3.1, we prove a well-posedness result for ISPn. This is achieved by combining the technique developed in [23], improved regularityestimates and with a novel perturbation argument from [18]. Further, in Theorem 4.1, we establish aconditional stability result under additional regularity condition on f ( x (cid:48) , t ) for ISPd. To the best of ourknowledge, this is the first work rigorously analyzing ISPs of recovering a space-time dependent sourcecomponent in subdiffusion with time-dependent coefficients. The main technical challenges in the studyinclude the nonlocality of the time-fractional derivative ∂ αt u and the time-dependence of the operator A ( t ).The nonlocality essentially limits the solution regularity pickup, and thus sharp regularity estimates forincompatible problem data are needed, which is especially delicate due to limited smoothness of the domainΩ. This is achieved in Proposition 4.1 by using a refined regularity pickup from [9, Lemma 2.4], exploitingthe cylindrical structure of the domain Ω. The time dependence of the elliptic operator A ( t ) precludesthe application of the standard separation of variable technique that has been predominant in existingstudies. This challenge is overcome by a perturbation argument and maximal L p regularity for time-fractionalproblems, which plays an important role in the analysis of ISPd. In Section 5, we derive the adjoint problemfor computing the gradient of a quadratic misfit functional and analyze the regularity of the adjoint variable.Further, we describe the conjugate gradient algorithm for recovering f ( x (cid:48) , t ), and provide extensive numericalexperiments to illustrate the feasibility of the recovery. The well-posedness, stability and reconstructionalgorithm represent the main contributions of this work.Last we situate this work in the existing literature. ISPs of recovering a part of information of the source F in a subdiffusion model from the lateral or terminal data represent an important class of applied inverseproblems, and have been extensively studied in the past decade. Most of the works devoted to this problem2ave been stated for sources F ( x, t ) = p ( t ) q ( x ) and can be divided into two groups: (i) inverse t -sourceproblem of recovering p ( t ) [36, 39, 8, 28] and (ii) inverse x -source problem of recovering q ( x ) [37, 41, 16, 35].For example, by the decay property of the Mittag-Leffler function E α,β ( z ), a two-sided stability result ofrecovering p ( t ) was shown in [36], if the observation u ( x , t ) satisfies x ∈ supp( q ). The unique recovery ofthe spatial component q ( x ) by interior observation was proved in [16] using Duhamel’s principle and uniquecontinuation principle, which also gave an iterative reconstruction algorithm. All these works are concernedwith recovering only either p ( t ) or q ( x ). The work [22] showed the simultaneous recovery of p and q undersuitable assumptions. We refer interested readers to the reviews [21, 30] for further pointers to theoreticaland numerical results. See also the work of [26] for the unique recovery of a general source F from thefull knowledge of the solution of problem (1.1), with A independent of t , on Ω × ( T , T ), with T ∈ (0 , T ).Recently the authors [23] proved a first uniqueness and stability results for the ISPs of recovering f ( x (cid:48) , t )of the subdiffusion model in a cylindrical domain. (See also Isakov [14] for relevant results for the standardparabolic problem in the half space.) The analysis in [23] relies on some representation of solutions by meanof E α,β ( z ) which are unavailable for elliptic operators with time-dependent coefficients. This work extendsthe results in [23] to the case of the time-dependent diffusion coefficients, and further, by exploiting themaximal L p regularity, we substantially relax the regularity requirement on f ( x (cid:48) , t ) for conditional stability.Inverse problems for subdiffusion with time-dependent coefficients have been scarcely studied so far, dueto a lack of mathematical tools, when compared with the time-independent counterpart. The only work weare aware of on an ISP with a time-dependent elliptic operator is [38], which showed the unique recovery of aspatial component from terminal measurement using an energy argument, which seems nontrivial to extendto the case of f ( x (cid:48) , t ). See also the works [42] for recovering a time-dependent factor in the diffusion coefficient a ( t ), where the special structure does allow applying the establish separation of variable technique. Thus,stable recovery for ISPs in the case of time-dependent coefficients remains challenging. This work presentsone promising approach to overcome the challenge (i.e., perturbation argument), and in particular it allowsestablishing the stable recovery.The rest of the paper is organized as follows. In Section 2, we state the assumptions and preliminaryestimates. Then in Sections 3 and 4, we prove the well-posedness of ISPn and conditional stability of ISPd,respectively. In Section 5, we describe a numerical algorithm for recovering f ( x (cid:48) , t ) for both ISPs, and provideseveral numerical experiments to showcase the feasibility of the recovery. Throughout, the notation c denotesa generic constant which may change from line to line. For a bivariate function g ( x, t ) or g ( x (cid:48) , t ), we oftenabbreviate it to g ( t ) as a vector-valued function by suppressing the dependence on the spatial variable. Now we collect several preliminary results. For m = 0 ,
1, we define two realizations A ( t ) and ˜ A ( t ) in L (Ω)of the elliptic operator A , with their domains respectively given by D ( A ( t )) = { v ∈ H (Ω) : A ( t ) v ∈ L (Ω) } , (2.1) D ( ˜ A ( t )) = { v ∈ H (Ω) : v | ∂ω × ( − (cid:96),(cid:96) ) = 0 , A ( t ) v ∈ L (Ω) , ∂ mx d v | x d = (cid:96) = 0 , ∂ x d v | x d = − (cid:96) = 0 } , (2.2)and let A ∗ = A ( t ∗ ) and ˜ A ∗ = ˜ A ( t ∗ ) for any t ∗ ∈ [0 , T ]. Note that we abuse the notation ˜ A ( t ) for both m = 0 and m = 1, which will be clear from the context. For any s ≥ A s ∗ and ˜ A s ∗ denote the fractionalpower of A ∗ and ˜ A ∗ via spectral decomposition, and the associated graph norms by (cid:107) · (cid:107) D ( A s ∗ ) and (cid:107) · (cid:107) D ( ˜ A s ∗ ) ,respectively. Let E ∗ ( t ) and ˜ E ∗ ( t ) be the solution operators corresponding to the source F , associated withthe elliptic operators A ∗ and ˜ A ∗ , respectively, defined by [20, Section 3.1] E ∗ ( t ) := 12 π i (cid:90) Γ θ,δ e zt ( z α + A ∗ ) − d z and ˜ E ∗ ( t ) := 12 π i (cid:90) Γ θ,δ e zt ( z α + ˜ A ∗ ) − d z, (2.3)with the contour Γ θ,δ ⊂ C (oriented with an increasing imaginary part) given byΓ θ,δ = { z ∈ C : | z | = δ, | arg z | ≤ θ } ∪ { z ∈ C : z = ρe ± i θ , ρ ≥ δ } . θ ∈ ( π , π ) so that z α ∈ Σ αθ for z ∈ Σ θ := { z ∈ C \{ } : | arg( z ) | ≤ θ } . Further, we employthe operator ˜ S ∗ ( t ) (corresponding to the initial data) defined by˜ S ∗ ( t ) := 12 π i (cid:90) Γ θ,δ e zt z α − ( z α + ˜ A ∗ ) − d z. Then it is known that [20, (3.8)] dd t ˜ S ∗ ( t ) = − ˜ A ∗ ˜ E ∗ ( t ) . (2.4)The next lemma summarizes the smoothing properties of E ∗ ( t ), ˜ E ∗ ( t ) and ˜ S ∗ ( t ). The notation (cid:107) · (cid:107) denotes the operator norm on L (Ω). Lemma 2.1 ([20, Lemma 1]) . For any β ∈ [0 , , there hold for any t ∈ (0 , T ] t α ( β − (cid:107) A β ∗ E ∗ ( t ) (cid:107) ≤ c and t α ( β − (cid:107) ˜ A β ∗ ˜ E ∗ ( t ) (cid:107) + t α (cid:107) ˜ A ∗ ˜ E ∗ ( t ) (cid:107) + t βα (cid:107) ˜ A β ∗ ˜ S ∗ ( t ) (cid:107) ≤ c. Throughout, we make the following assumption on the diffusion coefficient matrix a . The regularity a ∈ C ([0 , T ]; C (Ω; R d × d )) ∩ C ([0 , T ]; C (Ω; R d × d )) is sufficient for Lemma 2.2 and the required regularityproperties for our analysis. (ii) is a structural condition to enable unique recovery. The notation · and | · | denote standard Euclidean inner product and norm, respectively, on R d . Assumption 2.1.
The coefficient q ∈ C ([0 , T ]; L ∞ (Ω)) ∩ L ∞ (0 , T ; W , ∞ (Ω)) , and the diffusion coefficientmatrix a ∈ C ([0 , T ]; C (Ω; R d × d )) ∩ C ([0 , T ]; C (Ω; R d × d )) satisfies the following conditions. (i) There exists λ ∈ (0 , such that for any ( x, t ) ∈ Ω × [0 , T ] , λ | ξ | ≤ a ( x, t ) ξ · ξ ≤ λ − | ξ | , ∀ ξ ∈ R d . (ii) a jd ( x (cid:48) , ± (cid:96), t ) = 0 , x (cid:48) ∈ ω and j = 1 , . . . , d − , and ∂ x d a ij ( t ) = 0 , for i, j = 1 , . . . , d − . Note that the cylindrical domain Ω = ω × ( − (cid:96), (cid:96) ) is only Lipschitz continuous. Thus, some extra as-sumptions on the domain and the coefficient matrix a are needed in order to guarantee high-order Sobolevregularity of the elliptic operator A ( t ) with suitable boundary conditions. In the analysis, we need the fol-lowing elliptic regularity pickup: (i) and (ii) are sufficient for the analysis in Sections 3 and 4, respectively.(i) holds under the assumption that the domain ω is convex and a id = 0 , ∂ x j a dd = 0 , ∂ x d a ij = 0 , i, j ∈ { , . . . , d − } . (2.5)Indeed, if ω is convex, then Ω is convex and the desired assertion follows from [10, Theorems 3.2.1.2 and3.2.1.3]. This can be verified using a separation of variable argument [9, Lemma 2.4]. Besides the condition(2.5), if the domain ω is of class C , the separation of variable argument similar to [9, Lemma 2.4] impliesAssumption (cid:101) H in Definition 2.1(ii).
Definition 2.1. (i)
A tuple (Ω , A ( t )) is said to satisfy Assumption H mn , m, n = 0 , : if for any t ∈ [0 , T ] and any f ∈ L (Ω) , the following boundary value problem A ( t ) v = f, in Ω ,v = 0 , on ∂ω × ( − (cid:96), (cid:96) ) ,∂ mx d v ( x (cid:48) , (cid:96) ) = 0 , on ω,∂ nx d v ( x (cid:48) , − (cid:96) ) = 0 , on ω, admits a unique solution v ∈ H (Ω) and there exists c = c ( A , m, n, Ω) such that (cid:107) v (cid:107) H (Ω) ≤ c (cid:107) f (cid:107) L (Ω) . (ii) A tuple (Ω , A ( t )) is said to satisfy Assumption (cid:101) H , if for all v ∈ H max(1 ,s ) (Ω) satisfying A ( t ) v ∈ H s (Ω) , s ∈ [0 , , there holds v ∈ H s (Ω) and there exists a constant c = c ( A , s, Ω) such that (cid:107) v (cid:107) H s (Ω) ≤ c ( (cid:107)A v (cid:107) H s (Ω) + (cid:107) v (cid:107) H s (Ω) ) . Lemma 2.2.
Under Assumptions 2.1(i) and
H00 / H01 / H11 , for any t, s ∈ [0 , T ] and β ∈ [0 , , there hold (cid:107) A β ∗ ( I − A ( t ) − A ( s )) v (cid:107) L (Ω) ≤ c | t − s |(cid:107) A β ∗ v (cid:107) L (Ω) , ∀ v ∈ D ( A β ∗ ) , (2.6) (cid:107) ˜ A β ∗ ( I − ˜ A ( t ) − ˜ A ( s )) v (cid:107) L (Ω) ≤ c | t − s |(cid:107) ˜ A β ∗ v (cid:107) L (Ω) , ∀ v ∈ D ( ˜ A β ∗ ) . (2.7) Proof.
For the operator A ( t ), the case β = 0 is contained in [18, Corollary 3.1]. To show the estimate for β = 1, fix t, s ∈ [0 , T ], v ∈ D ( A ∗ ). From Assumption H00, we deduce D ( A ∗ ) = H (Ω) ∩ H (Ω) = D ( A ( t )) = D ( A ( s )), i.e., v ∈ D ( A ( t )) and v ∈ D ( A ( s )). Moreover, applying again Assumption H00, we get (cid:107) A ∗ ( I − A ( t ) − A ( s )) v (cid:107) L (Ω) ≤ c (cid:107) ( I − A ( t ) − A ( s )) v (cid:107) H (Ω) ≤ c (cid:107) A ( t )( I − A ( t ) − A ( s )) v (cid:107) L (Ω) = c (cid:107) A ( t ) v − A ( s ) v (cid:107) L (Ω) , with c > t and s . Combining this estimate with [18, eq. (2.6)] and AssumptionH00, we obtain (2.6) for β = 1 and q ≡
0. We can extend this result to q (cid:54)≡
0, since for q ∈ C ([0 , T ]; L (Ω)),the mean value theorem implies (cid:107) q ( t ) v − q ( s ) v (cid:107) L (Ω) ≤ (cid:107) ∂ t q (cid:107) L ∞ (0 ,T ; L ∞ (Ω)) | t − s |(cid:107) v (cid:107) L (Ω) ≤ c | t − s |(cid:107) A β ∗ v (cid:107) L (Ω) . The case β ∈ (0 ,
1) follows by interpolation. The proof of (2.7) is identical under Assumption H01 / H11.Below we need Bochner-Sobolev spaces W s,p (0 , T ; X ), for a UMD space X (see [13] for the definition ofUMD spaces, which include Sobolev spaces W s,p (Ω) with 1 < p < ∞ ). For any s ≥ ≤ p < ∞ ,we denote by W s,p (0 , T ; X ) the space of functions v : (0 , T ) → X , with the norm defined by complexinterpolation. Equivalently, the space is equipped with the quotient norm (cid:107) v (cid:107) W s,p (0 ,T ; X ) := inf (cid:101) v (cid:107) (cid:101) v (cid:107) W s,p ( R ; X ) := inf (cid:101) v (cid:107)F − [(1 + | ξ | ) s F ( (cid:101) v )( ξ )] (cid:107) L p ( R ; X ) where the infimum is taken over all possible (cid:101) v that extend v from (0 , T ) to R , and F denotes the Fouriertransform. The following norm equivalence result will be used extensively. Lemma 2.3.
Let α ∈ (0 , and p ∈ [1 , ∞ ) with αp > . If v (0) = 0 and ∂ αt v ∈ L p (0 , T ; X ) , then v ∈ W α,p (0 , T ; X ) and (cid:107) v (cid:107) W α,p (0 ,T ; X ) ≤ c (cid:107) ∂ αt v (cid:107) L p (0 ,T ; X ) . Meanwhile, if v (0) = 0 , v ∈ W α,p (0 , T ; X ) , then ∂ αt v ∈ L p (0 , T ; X ) and (cid:107) ∂ αt v (cid:107) L p (0 ,T ; X ) ≤ c (cid:107) v (cid:107) W α,p (0 ,T ; X ) . Proof.
Let g = ∂ αt v ∈ L p (0 , T ; X ). Then v ( t ) = α ) (cid:82) t ( t − s ) α − g ( s ) d s. This and Young’s convolutioninequality imply (cid:107) v (cid:107) L p (0 ,T ; X ) ≤ c (cid:107) g (cid:107) L p (0 ,T ; X ) . Let ˜ g be the extension of g from L p (0 , T ; X ) to L p ( R ; X ) byzero, i.e., ˜ g ( t ) = 0 for t ∈ ( −∞ , ∪ ( T, ∞ ) and ˜ g ( t ) = g ( t ) for t ∈ (0 , T ). Then let (cid:101) v ( t ) = α ) (cid:82) t −∞ ( t − s ) α − (cid:101) g ( s ) d s, which satisfies (cid:101) g ( t ) = α ) dd t (cid:82) t −∞ ( t − s ) − α (cid:101) v ( s ) d s. Then there holds [25, p. 90](i ξ ) α F ( (cid:101) v )( ξ ) = F ( (cid:101) g )( ξ ) , and (cid:101) v is an extension of v . Consequently, we have (cid:107) (cid:101) v (cid:107) W α,p ( R ; X ) = (cid:107)F − [(1 + | ξ | ) α F ( (cid:101) v )( ξ )] (cid:107) L p ( R ; X ) = (cid:107)F − [ K ( ξ )(1 + (i ξ ) α ) F ( (cid:101) v )( ξ )] (cid:107) L p ( R ; X ) with K ( ξ ) = (1 + | ξ | ) α (1 + (i ξ ) α ) − . Note that lim | ξ |→ | K ( ξ ) | = 1 and lim | ξ |→∞ | K ( ξ ) | = 1, so K ( ξ ) isuniformly bounded. Similarly, ξ dd ξ K ( ξ ) = α | ξ | | ξ | (1 + | ξ | ) α (1 + (i ξ ) α ) − + α (1 + | ξ | ) α (1 + (i ξ ) α ) − (i ξ ) α
5s also bounded. Therefore, vector-valued Mikhlin multiplier theorem (see, e.g. [6] or [43, Proposition 3])indicates that K ( ξ ) is a Fourier multiplier, and hence (cid:107) v (cid:107) W α,p (0 ,T ; X ) ≤ (cid:107) (cid:101) v (cid:107) W α,p ( R ; X ) ≤ c (cid:107)F − [(1 + (i ξ ) α ) F ( (cid:101) v )( ξ )] (cid:107) L p ( R ; X ) ≤ c (cid:107) (cid:101) v (cid:107) L p ( R ; X ) + c (cid:107) g (cid:107) L p ( R ; X ) ≤ c (cid:107) g (cid:107) L p ( R ; X ) = c (cid:107) g (cid:107) L p (0 ,T ; X ) = c (cid:107) ∂ αt u (cid:107) L p (0 ,T ; X ) . To prove the second assertion, let v ∈ C ∞ ([0 , T ]; X ) with v (0) = 0, and we extend v from (0 , T ) to a function (cid:101) v ∈ W α,p ( R ; X ) satisfying (cid:101) v ( t ) = 0 for all t ≤ (cid:107) (cid:101) v (cid:107) W α,p ( R ; X ) ≤ c (cid:107) v (cid:107) W α,p (0 ,T ; X ) . (2.8)Then it is direct that −∞ ∂ αt v ( t ) := α ) dd t (cid:82) t −∞ ( t − s ) − α (cid:101) v ( s ) d s = ∂ αt v ( t ) for t ∈ (0 , T ) , and (cid:107) ∂ αt v (cid:107) L p (0 ,T ; X ) = (cid:107) −∞ ∂ αt (cid:101) v (cid:107) L p (0 ,T ; X ) ≤ (cid:107) −∞ ∂ αt (cid:101) v (cid:107) L p ( R ; X ) = (cid:107)F − (i ξ ) α F ( (cid:101) v )( ξ ) (cid:107) L p ( R ; X ) = (cid:107)F − K ( ξ )(1 + | ξ | ) α F ( (cid:101) v )( ξ ) (cid:107) L p ( R ; X ) with K ( ξ ) = | ξ | α (1 + | ξ | ) − α . Note that both | K ( ξ ) | and | ξ dd ξ K ( ξ ) | are uniformly bounded, hence it is aFourier multiplier. Then we have (cid:107) ∂ αt v (cid:107) L p (0 ,T ; X ) ≤ (cid:107)F − (1 + | ξ | ) α F ( (cid:101) v )( ξ ) (cid:107) L p ( R ; X ) ≤ c (cid:107) (cid:101) v (cid:107) W α,p ( R ; X ) . This together with (2.8) and the density of C ∞ ([0 , T ]; X ) in W α,p (0 , T ; X ) leads to the second assertion.We need the following version of Gronwall’s inequality (see, e.g., [40] or [12, Exercise 3, p. 190]). Lemma 2.4.
Let c, r > and y, a ∈ L (0 , T ) be nonnegative functions satisfying y ( t ) ≤ a ( t ) + c (cid:90) t ( t − s ) r − y ( s )d s, t ∈ (0 , T ) . Then there exists c = c ( r, T ) > such that y ( t ) ≤ a ( t ) + c (cid:90) t ( t − s ) r − a ( s )d s, t ∈ (0 , T ) . This section is devoted ISPn, i.e. recovering f ( x (cid:48) , t ) in problem (1.1) with m = 1 from u | ω ×{ (cid:96) }× (0 ,T ) . Thedirect problem is given by ∂ αt u + A ( t ) u = F, in Ω × (0 , T ] ,u ( x,
0) = 0 , in Ω ,∂ x d u ( x (cid:48) , ± (cid:96), t ) = 0 , on ω × (0 , T ) ,u ( x, t ) = 0 , on ∂ω × ( − (cid:96), (cid:96) ) × (0 , T ) . (3.1)Subdiffusion with time-dependent coefficients has recently been studied in [27, 18, 20], where well-posednessand several regularity estimates have been established. Our description largely follows the approach devel-oped in [18, 20]. Throughout, for the prefactor R ( x, t ) in the source F , we make the following assumption. Assumption 3.1.
The function R ∈ L ∞ (Ω × (0 , T )) satisfies ∂ x d R ∈ L ∞ (Ω × (0 , T )) and that there exists c R > such that | R ( x (cid:48) , (cid:96), t ) | ≥ c R for any ( x (cid:48) , t ) ∈ ω × (0 , T ) . Now we give several regularity estimates for the direct problem (3.1). First we derive a representation ofthe solution u . The key step is to reformulate problem (3.1) into ∂ αt u ( t ) + ˜ A ∗ u ( t ) = F ( t ) + ( ˜ A ∗ − ˜ A ( t )) u ( t ) , ∀ t ∈ (0 , T ] . u which satisfies u ( t ) = (cid:90) t ˜ E ∗ ( t − s ) F ( s )d s + (cid:90) t ˜ E ∗ ( t − s )( ˜ A ∗ − ˜ A ( s )) u ( s )d s. Then by setting t to t ∗ , we can use Lemma 2.2 to estimate the integral (cid:82) t ∗ ˜ E ∗ ( t ∗ − s )( ˜ A ∗ − ˜ A ( s )) u ( s )d s .The next result collects a priori estimates on the solution u to problem (3.1). Lemma 3.1.
Let Assumption 2.1(i) be fulfilled. Then the solution u to problem (3.1) satisfies (cid:107) u ( t ) (cid:107) H (Ω) ≤ c (cid:90) t ( t − s ) α − (cid:107) F ( s ) (cid:107) L (Ω) d s, ∀ t ∈ (0 , T ] , (3.2) and also the following maximal L p regularity (cid:107) ∂ αt u (cid:107) L p (0 ,T ; L (Ω)) + (cid:107) u (cid:107) L p (0 ,T ; D ( ˜ A ∗ )) ≤ c (cid:107) F (cid:107) L p (0 ,T ; L (Ω)) , ∀ < p < ∞ . (3.3) Proof.
The estimate (3.2) can be found in [20, Theorem 2] (with k = 0), and (3.3) in [18, Theorem 2.1].Further, we denote by u f the solution of problem (3.1) to explicitly indicate its dependence on f . Firstwe show that the inverse problem is indeed ill-posed on the space L (0 , T ; L ( ω )). Corollary 3.1.
Under Assumptions 2.1(i) and 3.1, the map f (cid:55)→ u f | L (0 ,T ; L ( ω )) is linear and compact on L (0 , T ; L ( ω )) .Proof. The linearity is obvious. The compactness is direct from Lemma 3.1. In fact, by the maximal L p regularity in Lemmas 3.1 and 2.3 and Assumption 3.1, we have (cid:107) u f (cid:107) W α, (0 ,T ; L (Ω)) + (cid:107) u f (cid:107) L (0 ,T ; H (Ω)) ≤ c (cid:107) f (cid:107) L (0 ,T ; L ( ω )) . Thus, u f ∈ W α, (0 , T ; L (Ω)) ∩ L (0 , T ; H (Ω)). Meanwhile, by interpolation, the space W α, (0 , T ; L (Ω)) ∩ L (0 , T ; H (Ω)) embeds compactly into L (0 , T ; H (Ω)) [4, Theorem 5.2], which, by the trace theorem,embeds continuously into L (0 , T ; L ( ω )). Thus the map f (cid:55)→ u f | ω ×{ (cid:96) }× (0 ,T ) is compact on L (0 , T ; L ( ω )).Let w = ∂ x d u f . Then w satisfies ∂ αt w + A ( t ) w = − ∂ x d A ( t ) u f ( t ) + ∂ x d F ( t ) , in Ω × (0 , T ] ,w (0) = 0 , in Ω ,w = 0 , on ∂ Ω × (0 , T ) . (3.4)By applying the perturbation argument and using the operator A ( t ), the solution w to problem (3.4) can berepresented by w ( t ) = (cid:90) t E ∗ ( t − s ) (cid:0) − ∂ x d A ( s ) u f ( s ) + ∂ x d F ( s ) (cid:1) d s + (cid:90) t E ∗ ( t − s )( A ∗ − A ( s )) w ( s )d s. (3.5)Noting the definition w = ∂ x d u f and the condition ∂ x d a ij ( t ) = 0 for j, j = 1 , . . . , d − − ( ∂ x d A ( t )) u = ∂ x d ( ∂ x d a dd ( t ) ∂ x d u f ) + d − (cid:88) j =1 (cid:2) ∂ x j ( ∂ x d a jd ( t ) ∂ x d u f ) + ∂ x d ( ∂ x d a jd ( t ) ∂ x j u f ) (cid:3) + ∂ x d q ( t ) u f := B ( t ) w + B ( t ) u f , where the (time-dependent) operators B ( t ) and B ( t ) are respectively given by B ( t ) w := ∂ x d a dd ( t ) ∂ x d w + 2 d − (cid:88) j =1 ∂ x d a jd ( t ) ∂ x j w + d − (cid:88) j =1 ∂ x j ∂ x d a jd ( t ) w, ( t ) u := d (cid:88) j =1 ( ∂ x d a jd ( t )) ∂ x j u − ∂ x d q ( t ) u. The next result gives useful bounds on w := ∂ x d u f . Lemma 3.2.
Let Assumptions 2.1,
H00 , H11 and 3.1 be fulfilled. Then there exists a unique solution w ∈ L (0 , T ; H (Ω)) with A ( t ) w ( t ) , ∂ αt w ∈ L (0 , T ; L (Ω)) to problem (3.4) , and for any β ∈ [1 , , (cid:107) w ( t ) (cid:107) H β (Ω) ≤ c (cid:90) t ( t − s ) (1 − β ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s, t ∈ (0 , T ) , where the constant c depends only on R , A , β and T .Proof. By the a priori estimate (3.3) with p = 2, we have − ∂ x d A ( t ) u f , ∂ x d F ∈ L (0 , T ; L (Ω)). ThenLemma 3.1 shows that problem (3.4) has a unique solution w ∈ L (0 , T ; H (Ω)) with A ( t ) w ( t ) , ∂ αt w ∈ L (0 , T ; L (Ω)). Next, we prove the H β (Ω) bound on w ( t ). We define the operators K : L (0 , T ; H (Ω)) → L (0 , T ; H (Ω)) and K : L (0 , T ; L ( ω )) → L (0 , T ; H (Ω)), respectively, by K v ( t ) = (cid:90) t E ∗ ( t − s ) B v ( s )d s,K f ( t ) = (cid:90) t E ∗ ( t − s ) B u f ( s )d s + (cid:90) t E ∗ ( t − s ) ∂ x d R ( s ) f ( s )d s. By Lemma 2.1, we have (cid:107) K v ( t ∗ ) (cid:107) H β (Ω) ≤ c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) B v ( s ) (cid:107) L (Ω) d s ≤ c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) v ( s ) (cid:107) H (Ω) d s. (3.6)Similarly, by Lemma 2.1, under Assumption 3.1, we have (cid:107) K f ( t ∗ ) (cid:107) H β (Ω) ≤ c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s + c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) u f ( s ) (cid:107) H (Ω) d s. (3.7)Meanwhile, under Assumption 3.1 and the estimate (3.2), we deduce (cid:107) u f ( t ) (cid:107) H (Ω) ≤ (cid:90) t ( t − s ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s. Consequently, (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) u f ( s ) (cid:107) H (Ω) d s ≤ c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:90) s ( s − ξ ) α − (cid:107) f ( ξ ) (cid:107) L ( ω ) d ξ d s = c (cid:90) t ∗ (cid:107) f ( ξ ) (cid:107) L ( ω ) (cid:90) t ∗ ξ ( t ∗ − s ) (1 − β ) α − ( s − ξ ) α − d s d ξ ≤ cT α (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s. This and the estimate (3.7) imply (cid:107) K f ( t ∗ ) (cid:107) H β (Ω) ≤ c T (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s. (3.8)Next, by Lemmas 2.1 and 2.2, (2.6) and Assumption H00, we have (cid:13)(cid:13)(cid:13) (cid:90) t ∗ E ∗ ( t ∗ − s )( A ∗ − A ( s )) w ( s )d s (cid:13)(cid:13)(cid:13) H β (Ω) ≤ c (cid:90) t ∗ (cid:107) A ∗ E ∗ ( t ∗ − s ) (cid:107)(cid:107) A β ∗ ( I − A − ∗ A ( s )) w ( s ) (cid:107) L (Ω) d s c (cid:90) t ∗ ( t ∗ − s ) − ( t ∗ − s ) (cid:107) w ( s ) (cid:107) H β (Ω) d s = c (cid:90) t ∗ (cid:107) w ( s ) (cid:107) H β (Ω) d s. This estimate, (3.6), (3.8) and the solution representation (3.5) lead to (cid:107) w ( t ∗ ) (cid:107) H β (Ω) ≤ c (cid:90) t ∗ (cid:107) w ( s ) (cid:107) H β (Ω) d s + c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s. This and Gronwall’s inequality in Lemma 2.4 imply the desired H β (Ω) bound. This completes the proof.Now we can state a well-posedness result for ISPn. Below we use the notation L (0 , T ; L ( ω )) and L ( ω × (0 , T )) interchangeably since they are isomorphic by Fubini–Tonelli theorem. Theorem 3.1.
Let Assumptions 2.1,
H00 , H11 and 3.1 be fulfilled. Then for any f ∈ L ( ω × (0 , T )) , thesolution u of problem (3.1) satisfies u ∈ H ( − (cid:96), (cid:96) ; L ( ω × (0 , T ))) , ∂ αt u , A ( t ) u ∈ H ( − (cid:96), (cid:96) ; L ( ω × (0 , T ))) .Thus, the map h : ( x (cid:48) , t ) (cid:55)→ [ ∂ αt u + ( A ( t ) + a dd ( t ) ∂ x d ) u ]( x (cid:48) , (cid:96), t ) R ( x (cid:48) , (cid:96), t ) ∈ L ( ω × (0 , T )) (3.9) is well defined, and further, there exists a bounded linear operator H : L (0 , T ; L ( ω )) → L (0 , T ; L ( ω )) such that f solves f = h + H f, (3.10) which is well-posed. Finally, for every pair ( h, f ) ∈ L ( ω × (0 , T )) × L ( ω × (0 , T )) satisfying (3.10) , thesolution u of problem (3.1) satisfies (3.9) .Proof. By Lemmas 3.1 and 3.2, problem (1.1) has a solution u f ∈ L (0 , T ; H (Ω)), with A ( t ) u f , ∂ αt u f ∈ L (0 , T ; L (Ω)) and w = ∂ x d u ∈ L (0 , T ; H (Ω)), ∂ αt ∂ x d u f , A ( t ) ∂ x d u f ∈ L (0 , T ; L (Ω)). Hence, x d (cid:55)→ u f ( · , x d , · ) ∈ H ( − (cid:96), (cid:96) ; L ( ω × (0 , T ))) ∩ H ( − (cid:96), (cid:96) ; L (0 , T ; H ( ω ))) ,x d (cid:55)→ ∂ αt u f ( · , x d , · ) ∈ H ( − (cid:96), (cid:96) ; L ( ω × (0 , T ))) . By the trace theorem, we can restrict ∂ x d w = ∂ x d u f , ∂ αt u f and A ( t ) u f , to the boundary x d = (cid:96) . Thus thegoverning equation in problem (3.1) implies that for ( x (cid:48) , t ) ∈ ω × (0 , T ), a dd ( t ) ∂ x d w ( x (cid:48) , (cid:96), t ) = a dd ( t ) ∂ x d u f ( x (cid:48) , (cid:96), t )= [ ∂ αt u f + ( A ( t ) + a dd ( t ) ∂ x d ) u f ]( x (cid:48) , (cid:96), t ) − R ( x (cid:48) , (cid:96), t ) f ( x (cid:48) , t )= R ( x (cid:48) , (cid:96), t )[ h ( x (cid:48) , t ) − f ( x (cid:48) , t )] , (3.11)with the function h ( x (cid:48) , t ) given by (3.9). Let the operator H : L ( ω × (0 , T )) → L ( ω × (0 , T )) be defined by[ H φ ]( x (cid:48) , t ) = a dd ( t ) ∂ x d u φ ( x (cid:48) , (cid:96), t ) R ( x (cid:48) , (cid:96), t ) , where u φ ( x (cid:48) , x d , t ) denotes the solution to problem (3.1) with F = φR . Then it follows from (3.11) that f isthe solution to f = h + H f. Moreover, by Lemma 3.2, trace inequality and the definition w = ∂ x d u f , we deduce that for any β ∈ ( , (cid:107) f ( t ) (cid:107) L ( ω ) ≤ (cid:107) h ( t ) (cid:107) L ( ω ) + (cid:107)H f (cid:107) L ( ω ) ≤ (cid:107) h ( t ) (cid:107) L ( ω ) + c (cid:107) ∂ x d u f ( · , (cid:96), t ) (cid:107) L ( ω ) ≤ (cid:107) h ( t ) (cid:107) L ( ω ) + c (cid:107) w ( t ) (cid:107) H β (Ω) ≤ (cid:107) h ( t ) (cid:107) L ( ω ) + c (cid:90) t ( t − s ) (1 − β ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s. This and the standard Gronwall’s inequality in Lemma 2.4 yield (cid:107) f ( t ) (cid:107) L ( ω ) ≤ (cid:107) h ( t ) (cid:107) L ( ω ) + c (cid:90) t ( t − s ) (1 − β ) α − (cid:107) h ( s ) (cid:107) L ( ω ) d s, (cid:107) f (cid:107) L (0 ,T ; L ( ω )) ≤ c (cid:107) h (cid:107) L (0 ,T ; L ( ω )) . This shows the well-posedness of equation (3.10) and the recovery of f from the data h . Last, fix ( h , f ) ∈ L (0 , T ; L ( ω )) × L (0 , T ; L ( ω )) satisfying (3.10) with h = h and consider u ∈ L (0 , T ; H (Ω)) solvingproblem (1.1) with F = hR . The preceding argument shows that one can define h ∈ L (0 , T ; L ( ω )) givenby (3.9) and f solves (3.10). This implies h = f − H f = h . Therefore, we have h = h , and this completesthe proof of the theorem. Remark 3.1.
Theorem 3.1 actually gives a reconstruction algorithm for recovering f for ISPn, if the givendata g † ( x (cid:48) , t ) = u f † ( x (cid:48) , (cid:96), t ) is sufficiently accurate so that the derivatives ∂ αt g † and A ( t ) g † in (3.9) can beevaluated accurately. For noisy data g δ , one can proceed in two steps: first suitably mollify the data g δ sothat the mollified data is smooth, and then apply the fixed point iteration. In this section, we establish a conditional stability result for ISPd, i.e., recovering f ( x (cid:48) , t ) in problem (1.1)with m = 0 from the lateral flux observation ∂ x d u | ω ×{ (cid:96) }× (0 ,T ) . The direct problem is given by ∂ αt u + A ( t ) u = f R, in Ω × (0 , T ] ,u = 0 , on ∂ω × ( − (cid:96), (cid:96) ) × (0 , T ) ,u ( · , (cid:96), · ) = 0 , on ω × (0 , T ] ,∂ x d u ( · , − (cid:96), · ) = 0 , on ω × (0 , T ) ,u (0) = 0 , in ω × ( − (cid:96), (cid:96) ) . (4.1)Note that the estimates in Lemma 3.1 remain valid for problem (4.1). The next result gives an improvedregularity result, under extra regularity and compatibility assumptions on the source F . This result plays acentral role in the stability analysis. Proposition 4.1.
Let Assumptions 2.1(i),
H01 and (cid:101) H be fulfilled, γ ∈ ( , , F ∈ W , α (1 − γ ) (0 , T ; L (Ω)) ∩ L α (1 − γ ) (0 , T ; H γ (Ω)) and F (0) = 0 . Then for any β ∈ ( , γ ) , problem (4.1) has a unique weak solution u ∈ L α (1 − γ ) (0 , T ; H β (Ω)) ∩ W α, α (1 − γ ) (0 , T ; H β (Ω)) with (cid:107) u (cid:107) L α (1 − γ ) (0 ,T ; H β (Ω)) + (cid:107) u (cid:107) W α, α (1 − γ ) (0 ,T ; H β (Ω)) ≤ c (cid:16) (cid:107) F (cid:107) W , α (1 − γ ) (0 ,T ; L (Ω)) + (cid:107) F (cid:107) L α (1 − γ ) (0 ,T ; H γ (Ω)) (cid:17) . Proof.
By Sobolev embedding, F ∈ L ∞ (0 , T ; L (Ω)), and the existence and uniqueness of a weak solution u ∈ L q (0 , T ; D ( ˜ A ∗ )) for all q ∈ (1 , ∞ ) follows directly from Lemma 3.1 with (cid:107) u (cid:107) L q (0 ,T ; D ( ˜ A ∗ )) ≤ c (cid:107) f (cid:107) W , α (1 − γ ) (0 ,T ; L ( ω )) . (4.2)It suffices to show the claimed regularity. Using the operator ˜ A ( t ) in (2.2) and the perturbation argument,since u (0) = 0, the solution u can be represented by u ( t ) = (cid:90) t ˜ E ∗ ( s ) F ( t − s )d s + (cid:90) t ˜ E ∗ ( s )( ˜ A ∗ − ˜ A ( t − s )) u ( t − s )d s. (4.3)Then applying ˜ A ∗ to both sides of the identity, and using the governing equation give ∂ αt u ( t ) = − ˜ A ∗ u ( t ) + F ( t ) + ( ˜ A ∗ − ˜ A ( t )) u ( t )= − (cid:90) t ˜ A ∗ ˜ E ∗ ( s ) F ( t − s )d s − (cid:90) t ˜ A ∗ ˜ E ∗ ( s )( ˜ A ∗ − ˜ A ( t − s )) u ( t − s )d s + F ( t ) + ( ˜ A ∗ − ˜ A ( t )) u ( t ) . t at t ∗ in the identity, applying the identity (2.4) and integration by parts formula to the firstintegral and noting the condition F ( · ,
0) = 0 and the fact ˜ S ∗ (0) = I , we obtain ∂ αt u ( t ∗ ) = − (cid:90) t ∗ ˜ A ∗ ˜ E ∗ ( s ) F ( t ∗ − s )d s − (cid:90) t ∗ ˜ A ∗ ˜ E ∗ ( s )( ˜ A ∗ − ˜ A ( t ∗ − s )) u ( t ∗ − s )d s + F ( t ∗ )= − (cid:90) t ∗ ˜ S ∗ ( s ) F (cid:48) ( t ∗ − s )d s − (cid:90) t ∗ ˜ A ∗ ˜ E ∗ ( s )( ˜ A ∗ − ˜ A ( t ∗ − s )) u ( t ∗ − s )d s. (4.4)Then it follows from Lemmas 2.1 and 2.2 that (cid:107) ˜ A β ∗ ∂ αt u ( t ∗ ) (cid:107) L (Ω) ≤ (cid:90) t ∗ (cid:107) ˜ A β ∗ ˜ S ∗ ( t ∗ − s ) (cid:107)(cid:107) F (cid:48) ( s ) (cid:107) L (Ω) d s + (cid:90) t ∗ (cid:107) ˜ A β ∗ E ∗ ( t ∗ − s ) (cid:107)(cid:107) ˜ A ∗ ( I − ˜ A − ∗ ˜ A ( s )) u ( s ) (cid:107) L (Ω) d s ≤ c (cid:90) t ∗ ( t ∗ − s ) − βα (cid:107) F (cid:48) ( s ) (cid:107) L (Ω) d s + c (cid:90) t ∗ ( t ∗ − s ) − βα − ( t ∗ − s ) (cid:107) ˜ A ∗ u ( s ) (cid:107) L (Ω) d s ≤ c (cid:90) t ∗ ( t ∗ − s ) − βα (cid:107) F (cid:48) ( s ) (cid:107) L (Ω) d s + (cid:90) t ∗ ( t ∗ − s ) − βα (cid:107) ˜ A ∗ u ( s ) (cid:107) L (Ω) d s. Note that t − βα ∈ L p (0 , T ) for any p ∈ (1 , γα ] ⊂ (1 , βα ), by the choice β < γ . Now choosing r = α (1 − γ ) , p = γα and q = − γ ) α in the following Young’s convolution inequality (cid:107) f ∗ g (cid:107) L r (0 ,T ) ≤ (cid:107) f (cid:107) L p (0 ,T ) (cid:107) g (cid:107) L q (0 ,T ) , ∀ p, q, r ≥ r − + 1 = p − + q − , (4.5)we deduce (cid:107) ˜ A β ∗ ∂ αt u ( t ∗ ) (cid:107) L α (1 − γ ) (0 ,T ; L (Ω)) ≤ c (cid:0) (cid:107) F (cid:48) ( s ) (cid:107) L q (0 ,T ; L (Ω)) + (cid:107) ˜ A ∗ u ( s ) (cid:107) L q (0 ,T ; L (Ω)) (cid:1) . Further, it follows from the representation (4.3) and Lemmas 2.1 and 2.2 that (cid:107) ˜ A β ∗ u ( t ∗ ) (cid:107) L (Ω) ≤ (cid:90) t ∗ (cid:107) ˜ A β ∗ ˜ E ∗ ( s ) (cid:107)(cid:107) F ( t ∗ − s ) (cid:107) L (Ω) d s + (cid:90) t ∗ (cid:107) ˜ A ∗ ˜ E ∗ ( s ) (cid:107)(cid:107) ˜ A β ∗ ( I − ˜ A − ∗ ˜ A ( t ∗ − s )) u ( t ∗ − s ) (cid:107) L (Ω) d s ≤ (cid:90) t ∗ s (1 − β ) α − (cid:107) F ( t ∗ − s ) (cid:107) L (Ω) d s + c (cid:90) t ∗ s − s (cid:107) u ( t ∗ − s ) (cid:107) L (Ω) d s ≤ c (cid:107) F (cid:107) L ∞ (0 ,T ; L (Ω)) t (1 − β ) α ∗ + (cid:90) t ∗ (cid:107) ˜ A β ∗ u ( s ) (cid:107) L (Ω) d s. This and Gronwall’s inequality directly imply lim t → + (cid:107) ˜ A β ∗ u ( t ) (cid:107) L (Ω) = 0. Hence, from Lemma 2.3 andAssumption H01, we deduce u ∈ W α, α (1 − γ ) (0 , T ; D ( ˜ A β ∗ )) ⊂ W α, α (1 − γ ) (0 , T ; H β (Ω)). Thus, we concludethat for any fixed t ∈ (0 , T ], the solution u satisfies A ( t ) u ( t ) = F ( t ) − ∂ αt u ( t ) , in Ω u ( x (cid:48) , (cid:96), t ) = 0 , on ω,∂ x d u ( x (cid:48) , − (cid:96), t ) = 0 , on ω,u ( x, t ) = 0 , on ∂ω × ( − (cid:96), (cid:96) ) . Note that for any t ∈ (0 , T ], there holds A ( · ) u ( · ) = F ( · ) − ∂ αt u ( · ) ∈ L α (1 − γ ) (0 , T ; H β (Ω)). Then byAssumption (cid:101) H, we obtain u ∈ L α (1 − γ ) (0 , T ; H β ) (Ω)). This completes the proof.11he conditional stability analysis employs the regularity estimates in Proposition 4.1. Let u be thesolution to problem (4.1), and let v = ∂ x d u . Then v satisfies ∂ αt v + A ( t ) v = H + f ∂ x d R, in Ω × (0 , T ] ,v = 0 , on ∂ω × ( − (cid:96), (cid:96) ) × (0 , T ) ,v ( · , (cid:96), · ) = ∂ x d u ( · , (cid:96), · ) , on ω × (0 , T ] ,v ( · , − (cid:96), · ) = 0 , on ω × (0 , T ) ,v (0) = 0 , in ω × ( − (cid:96), (cid:96) ) , (4.6)with the function H given by H ( x (cid:48) , x d , t ) = − ∂ x d A ( t ) u = d (cid:88) i,j =1 ∂ x i ( ∂ x d a ij ( t ) ∂ x j u ) − ∂ x d q ( t ) u. Unlike problem (3.4) in Section 3, problem (4.6) involves a nonzero Dirichlet boundary condition, andthus requires a different analysis. We employ an extension approach to derive the requisite bound on v . For r ≥ X α,r denotes the the space W α,r (0 , T ; L ( ω )) ∩ L r (0 , T ; H ( ω )) with the norm (cid:107) v (cid:107) X α,r = (cid:107) v (cid:107) W α,r (0 ,T ; L ( ω )) + (cid:107) v (cid:107) L r (0 ,T ; H ( ω )) . Assumption 4.1.
Let ω ∈ C , γ ∈ ( , , f ∈ W , − γ ) α (0 , T ; L ( ω )) ∩ L α (1 − γ ) (0 , T ; H γ ( ω )) with f (0) = 0 , ∂ t R ∈ L ∞ (Ω × (0 , T )) and R ∈ C ([0 , T ]; W , ∞ (Ω)) . Lemma 4.1.
Let Assumptions 2.1, 3.1,
H00 , H01 and 4.1 be fulfilled, and let u be the solution of problem (4.1) . Then for any β ∈ ( , γ ) , the solution v to problem (4.6) satisfies (cid:107) ∂ x d v ( · , (cid:96), t ) (cid:107) L ( ω ) ≤ c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α, α (1 − γ ) + c (cid:90) t ( t − s ) α (1 − β ) − (cid:107) f ( · , s ) (cid:107) L ( ω ) d s, ∀ t ∈ (0 , T ) . where the constant c depend on R , Ω , T , α , β , γ and A .Proof. Let r = α (1 − γ ) . Since f ∈ W , α (1 − γ ) (0 , T ; L ( ω )) and ∂ t R ∈ L ∞ (Ω × (0 , T )), we deduce F = f R ∈ W , α (1 − γ ) (0 , T ; L (Ω)). The assumptions f ∈ L r (0 , T ; H γ ( ω )) and R ∈ C ([0 , T ]; W , ∞ (Ω)) im-ply F = L α (1 − γ ) (0 , T ; H γ (Ω)). Further, f ( · ,
0) = 0 in ω indicates F ( · ,
0) = 0 in Ω. Thus, F = f R satisfies the conditions in Proposition 4.1, and since ω ∈ C , Assumption (cid:101) H holds. By Proposition 4.1, u ∈ W α,r (0 , T ; H β (Ω)) ∩ L r (0 , T ; H β ) (Ω)), for any β ∈ ( , γ ), and by the trace theorem, there holds( x (cid:48) , t ) (cid:55)→ ∂ x d u ( x (cid:48) , (cid:96), t ) ∈ W α,r ([0 , T ]; H β − ( ω )) ∩ L r ([0 , T ]; H +2 β ( ω )) , (4.7)and ∂ x d u ( · , − (cid:96),
0) = 0. Next we split the solution v to problem (4.6) into v = v + v , with the functions v and v , respectively, solving ∂ αt v + A ( t ) v = 0 , in Ω × (0 , T ) ,v = 0 , on ∂ω × ( − (cid:96), (cid:96) ) × (0 , T ) ,v ( · , (cid:96), · ) = ∂ x d u ( · , (cid:96), · ) , on ∈ ω × (0 , T ) ,v ( · , − (cid:96), t ) = 0 , on ω × (0 , T ) ,v (0) = 0 , in Ω , and ∂ αt v + A ( t ) v = H + f ∂ x d R, in Ω × (0 , T ) ,v = 0 , on ∂ Ω × (0 , T ) ,v (0) = 0 , in Ω . v and v . To bound v , we first extend ∂ x d u ( · , (cid:96), · ) from ω × (0 , T ) to Ω × (0 , T ). Indeed,by the regularity estimate (4.7) and using the classical lifting theorem for Sobolev spaces [29, Chapter 1,Theorem 9.4] and Assumption H00, we deduce that there exists a function G ∈ L r (0 , T ; H (Ω)) satisfying − ∆ x G ( x, t ) = 0 , ( x, t ) ∈ Q (4.8) G ( x (cid:48) , x d , t ) = ∂ x d u ( x (cid:48) , (cid:96), t ) , x d = (cid:96), ( x (cid:48) , t ) ∈ ω × (0 , T ) , , x d = − (cid:96), ( x (cid:48) , t ) ∈ ω × (0 , T ) , , t = 0 , ( x (cid:48) , x d ) ∈ ω × ( − (cid:96), (cid:96) ) . (4.9)Clearly, G satisfies the following estimate (cid:107) G (cid:107) L r (0 ,T ; H (Ω)) ≤ c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) L r (0 ,T ; H ( ω )) . Moreover, by interpreting G as the solution to (4.8)-(4.9) in the transposition sense (see e.g. [29, Chapter 2,Theorem 6.3] or [5]), G ∈ W α,r (0 , T ; L (Ω)) satisfies (cid:107) G (cid:107) W α,r (0 ,T ; L (Ω)) ≤ c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) W α,r (0 ,T ; L ( ω )) . Consequently, we have (cid:107) G (cid:107) W α,r (0 ,T ; L (Ω)) + (cid:107) G (cid:107) L r (0 ,T ; H (Ω)) ≤ (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r . (4.10)Then we can decompose v into v = G + w , with the function w solving ∂ αt w + A ( t ) w = F , in Ω × (0 , T ) ,w = 0 , on ∂ Ω × (0 , T ) ,w (0) = 0 , in ω × ( − (cid:96), (cid:96) ) , with F = − ∂ αt G − A ( t ) G. Since ∂ x d u ( x (cid:48) , (cid:96),
0) = 0 for x (cid:48) ∈ ω , the uniqueness of the solution of problem(4.8)–(4.9) implies G ( · ,
0) = 0. Then, direct computation with Lemma 2.3 gives (cid:107) F (cid:107) L r (0 ,T ; L (Ω)) ≤ (cid:107) ∂ αt G (cid:107) L r (0 ,T ; L (Ω)) + (cid:107)A ( t ) G (cid:107) L r (0 ,T ; L (Ω)) ≤ c ( (cid:107) G (cid:107) W α,r (0 ,T ; L (Ω)) + (cid:107) G (cid:107) L r (0 ,T ; H (Ω)) ) . (4.11)Thus using the operator A ( t ) defined in (2.1) and the perturbation argument, we have w ( t ∗ ) = (cid:90) t ∗ E ∗ ( t ∗ − s ) F ( s )d s + (cid:90) t ∗ E ∗ ( t ∗ − s )( A ∗ − A ( s )) w ( s )d s. By Lemmas 2.1 and 2.2, (cid:107) A β ∗ w ( t ∗ ) (cid:107) L (Ω) ≤ (cid:90) t ∗ (cid:107) A β ∗ E ∗ ( t ∗ − s ) (cid:107)(cid:107) F ( s ) (cid:107) L (Ω) d s + (cid:90) t ∗ (cid:107) A ∗ E ∗ ( t ∗ − s ) (cid:107)(cid:107) A β ∗ ( I − A − ∗ A ( s )) w ( s ) (cid:107) L (Ω) d s ≤ (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) F ( s ) (cid:107) L (Ω) d s + c (cid:90) t ∗ ( t ∗ − s ) − ( t ∗ − s ) (cid:107) A β ∗ w ( s ) (cid:107) L (Ω) d s ≤ c (cid:107) F (cid:107) L r (0 ,T ; L (Ω)) + (cid:90) t ∗ (cid:107) A β ∗ w ( s ) (cid:107) L (Ω) d s. It follows from this estimate and Gronwall’s inequality that w ∈ L ∞ (0 , T ; D ( A β ∗ )) with (cid:107) w (cid:107) L ∞ (0 ,T ; D ( A β ∗ ))) ≤ c T (cid:107) F (cid:107) L r (0 ,T ; L (Ω)) . (cid:107) v (cid:107) L ∞ (0 ,T ; D ( A β ∗ )) ≤ (cid:107) w (cid:107) L ∞ (0 ,T ; D ( A β ∗ )) + (cid:107) G (cid:107) L ∞ (0 ,T ; D ( A β ∗ )) ≤ c (cid:0) (cid:107) F (cid:107) L r (0 ,T ; L (Ω)) + (cid:107) G (cid:107) L ∞ (0 ,T ; H β (Ω)) (cid:1) . Meanwhile, the condition β ∈ ( , γ ) and [4, Theorem 5.2] imply the following embedding inequality (cid:107) w (cid:107) L ∞ (0 ,T ; H β (Ω)) ≤ c ( (cid:107) w (cid:107) W α,r (0 ,T ; L (Ω)) + (cid:107) w (cid:107) L r (0 ,T ; H (Ω)) ) . The last two estimates together give (cid:107) v (cid:107) L ∞ (0 ,T ; D ( A β ∗ )) ≤ c (cid:0) (cid:107) G (cid:107) W α,r (0 ,T ; L (Ω)) + (cid:107) G (cid:107) L r (0 ,T ; H (Ω)) (cid:1) . This and the estimate (4.10) imply (cid:107) v (cid:107) L ∞ (0 ,T ; D ( A β ∗ )) ≤ c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r . (4.12)Moreover, by Assumption H00 and the trace inequality (cid:107) ∂ x d v ( · , (cid:96), · ) (cid:107) L ∞ (0 ,T ; L ( ω )) ≤ c (cid:107) v (cid:107) L ∞ (0 ,T ; D ( A β ∗ )) , weget (cid:107) ∂ x d v ( · , (cid:96), · ) (cid:107) L ∞ (0 ,T ; L ( ω )) ≤ c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r . (4.13)Next we bound v . Note that the solution v ( t ) can be represented by v ( t ) = (cid:90) t E ∗ ( t − s )[ H ( s ) + ∂ x d F ( s )]d s + (cid:90) t E ∗ ( t − s )( A ∗ − A ( s )) v ( s )d s. Thus, by Lemmas 2.1 and 2.2 and Assumption 3.1, we get (cid:107) A β ∗ v ( t ∗ ) (cid:107) L (Ω) ≤ c (cid:90) t ∗ ( t ∗ − s ) α (1 − β ) − (cid:2) (cid:107) f ( s ) (cid:107) L ( ω ) + (cid:107) H ( s ) (cid:107) L (Ω) (cid:3) d s + c (cid:90) t ∗ (cid:107) A β ∗ v ( s ) (cid:107) L (Ω) d s. (4.14)In light of Assumption 2.1(ii) and the definition v = ∂ x d u , we have H ( t ) = ∂ x d a dd ( t ) ∂ x d u + 2 d − (cid:88) j =1 ∂ x d a jd ( t ) ∂ x j ∂ x d u + d (cid:88) j =1 ∂ x j ∂ x d a jd ( t ) ∂ x d u + d − (cid:88) j =1 ∂ x d a jd ( t ) ∂ x j u − ∂ x d q ( t ) u = ∂ x d a dd ( t ) ∂ x d v + 2 d − (cid:88) j =1 ∂ x d a jd ( t ) ∂ x j v + (cid:16) d (cid:88) j =1 ∂ x j ∂ x d a jd ( t ) (cid:17) v + d − (cid:88) j =1 ∂ x d a jd ( t ) ∂ x j u − ∂ x d q ( t ) u, from which it directly follows that (cid:107) H ( t ) (cid:107) L (Ω) ≤ c (cid:0) (cid:107) v ( t ) (cid:107) H (Ω) + (cid:107) u ( t ) (cid:107) H (Ω) (cid:1) , t ∈ (0 , T ] . By Lemma 3.1 with β = (which holds also for problem (4.1)) and Assumption 3.1, we have (cid:107) u ( t ) (cid:107) H (Ω) ≤ c (cid:90) t ( t − s ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s. This and (4.12) lead to (cid:107) H ( t ) (cid:107) L (Ω) ≤ c (cid:90) t ( t − s ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s + c ( (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r + (cid:107) v ( t ) (cid:107) H (Ω) ) , which together with (4.14) yields (cid:107) A β ∗ v ( t ∗ ) (cid:107) L (Ω) ≤ c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) f ( s ) (cid:107) L ( ω ) d s + c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r c (cid:90) t ∗ ( t ∗ − s ) (1 − β ) α − (cid:107) A β ∗ v ( s ) (cid:107) L (Ω) d s. This estimate and Gronwall’s inequality in Lemma 2.4 then imply (cid:107) A β ∗ v ( t ) (cid:107) L (Ω) ≤ c (cid:90) t ( t − s ) (1 − β ) α − (cid:107) f ( · , s ) (cid:107) L ( ω ) d s + c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r . It follows from this estimate, Assumption H00 and the trace inequality that (cid:107) ∂ x d v ( · , (cid:96), t ) (cid:107) L ( ω ) ≤ c (cid:90) t ( t − s ) α (1 − β ) − (cid:107) f ( s ) (cid:107) L ( ω ) d s + c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r . Finally, combining this bound with the estimate (4.13) yields the desired assertion.Now we can state a conditional stability result for ISPd.
Theorem 4.1.
Let Assumptions 2.1, 3.1,
H00 , H01 and 4.1 be fulfilled, and u the solution of problem (4.1) .Then there exists a constant c depending on R , Ω , T , α , γ , p and A such that (cid:107) f (cid:107) L ∞ (0 ,T ; L ( ω )) ≤ c (cid:16) (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) L α (1 − γ ) (0 ,T ; H ( ω )) + (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) W α, α (1 − γ ) (0 ,T ; L ( ω )) (cid:17) . Proof.
First projecting the governing equation in (4.1) onto the lateral boundary ω × { (cid:96) } × (0 , T ) and thenusing the fact that, for all ( x (cid:48) , t ) ∈ ω × (0 , T ), we have u ( x (cid:48) , (cid:96), t ) = 0, and thus for any ( x (cid:48) , t ) ∈ ω × (0 , T ), f ( x (cid:48) , t ) R ( x (cid:48) , (cid:96), t ) = (cid:104) a dd ( t ) ∂ x d u + 2 d − (cid:88) j =1 a jd ( t ) ∂ x j ∂ x d u + d (cid:88) j =1 ∂ x j a jd ( t ) ∂ x d u (cid:105) ( x (cid:48) , (cid:96), t ) . This, Assumption 3.1, and the definition v = ∂ x d u imply that for all t ∈ (0 , T ), there holds (cid:107) f ( t ) (cid:107) L ( ω ) ≤ c − R c (cid:0) (cid:107) ∂ x d v ( · , (cid:96), t ) (cid:107) L ( ω ) + (cid:107) ∂ x d u ( · , (cid:96), t ) (cid:107) H ( ω ) (cid:1) . Under the condition γ ∈ ( , r = α (1 − γ ) , [4, Theorem 5.2] implies (cid:107) ∂ x d u ( · , (cid:96), t ) (cid:107) L ∞ (0 ,T ; H ( ω )) ≤ c ( (cid:107) ∂ x d u ( · , (cid:96), t ) (cid:107) W α,r (0 ,T ; L ( ω )) + (cid:107) ∂ x d u ( · , (cid:96), t ) (cid:107) L r (0 ,T ; H ( ω )) ) . The last two estimates and Lemma 4.1 imply (cid:107) f ( t ) (cid:107) L ( ω ) ≤ c (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) X α,r + c (cid:90) t ( t − s ) α (1 − β ) − (cid:107) f ( s ) (cid:107) L ( ω ) d s. Then Gronwall’s inequality in Lemma 2.4 implies the desired assertion, and this completes the proof.
Remark 4.1.
Theorem 4.1 shows the influence of the order α on the stability: the larger is the order α , thestronger temporal regularity (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) W α, α (1 − γ ) (0 ,T ; L ( ω )) on the data u | ω ×{ (cid:96) }× (0 ,T ) the stability needs.This agrees with the smoothing property of the solution operator, and shows also the beneficial influence ofanomalous diffusion. Theorem 4.1 improves the corresponding result in [23, Theorem 1.4] ( with δ > ) : (cid:107) f (cid:107) L ∞ (0 ,T ; L ( ω )) ≤ c ( (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) L ∞ (0 ,T ; H ( ω )) + (cid:107) ∂ x d u ( · , (cid:96), · ) (cid:107) W , ∞ (0 ,T ; H δ ( ω )) ) , This improvement is achieved by the maximal L p regularity and the interpolation inequality. Remark 4.2.
In the spirit of [23, Corollary 1.5], Theorem 4.1 allows proving the stable recovery of a class ofthe zeroth order coefficients q from the flux data ∂ x d u | ω ×{ (cid:96) }× (0 ,T ) . This requires the existence of a solution toproblem (4.1) in W , ∞ (0 , T ; W , ∞ (Ω)) ∩ L ∞ (0 , T ; W , ∞ (Ω)) . The latter can be achieved using the argumentof Proposition 4.1, and we leave the details to future investigation. Numerical experiments and discussions
In this section, we present several numerical experiments to illustrate the feasibility of recovering the space-time dependent f from lateral boundary observation. First we describe a numerical algorithm for recovering f for ISPn (and the algorithm for ISPd is similar).We employ an iterative regularization technique, which approximately minimizes J ( f ) := (cid:107) u f − g δ (cid:107) L (0 ,T ; L ( ω )) , (5.1)where u f denotes the solution to the direct problem in (3.1) with F = f R . By Corollary 3.1, the map u f : L (0 , T ; L ( ω )) → L (0 , T ; L ( ω )) is linear and compact, and thus standard regularization theory [7, 15]can be applied to justify the reconstruction technique. In particular, when equipped with an appropriatestopping criterion, the approximate minimizer obtained by gradient type methods, e.g., gradient descentand conjugate gradient method, will converge to the exact source component f † as the noise level tends tozero, and further it will converge at a certain rate dependent of the “regularity” of f † , when equipped witha suitable stopping criterion.To (approximately) minimize the functional J ( f ), we employ the conjugate gradient (CG) method [3].The main computational effort is to compute the gradient, which can be done efficiently using the adjointtechnique. Specifically, let v be the solution to the following adjoint problem tR ∂ αT v + A ( t ) v = 0 , ( x (cid:48) , x d , t ) ∈ Ω × (0 , T ] , t I − αT v ( x, T ) = 0 , in Ω ,∂ x d v ( x (cid:48) , (cid:96), t ) = u f − g δ , on ω × (0 , T ) ,∂ x d v ( x (cid:48) , − (cid:96), t ) = 0 , on ω × (0 , T ) ,v ( x, t ) = 0 , on ∂ω × ( − (cid:96), (cid:96) ) × (0 , T ) , (5.2)where the notation t I − αT v ( t ) and tR ∂ αT v denotes the right-sided Riemann-Liouville fractional integral andderivative of v , defined respectively by [25] t I − αT v ( t ) = 1Γ(1 − α ) (cid:90) Tt ( s − t ) − α v ( s )d s and tR ∂ αT v ( t ) = − − α ) dd t (cid:90) Tt ( s − t ) − α u ( s )d s. Then we have the following representation of the gradient J (cid:48) ( f ) of J ( f ). Proposition 5.1.
The gradient J (cid:48) ( f ) of the functional J ( f ) is given by J (cid:48) ( f ) = (cid:90) (cid:96) − (cid:96) Rv d x d , (5.3) where v is the solution to the adjoint problem (5.2) .Proof. The derivation follows a standard procedure. The directional derivative J (cid:48) ( f )[ h ] of the functional J with respect to f in the direction h ∈ L (0 , T ; L ( ω )) is given by J (cid:48) ( f )[ h ] = ( u h , u f − g δ ) L (0 ,T ; L ( ω )) , where u h is the solution to problem (3.1) with h in place of f (or the source F = hR ). Multiplying theequation for u h with a test function φ ( x, t ) and then integration by parts yield (cid:90) T (cid:90) Ω ( φ ∂ αt u h + a ∇ u h · ∇ φ )d x d t = (cid:90) T (cid:90) Ω Rhφ d x d t. (5.4)Meanwhile, the weak formulation for the adjoint solution v is given by (cid:90) T (cid:90) Ω ( φ tR ∂ αT v + a ∇ v · ∇ φ )d x d t = (cid:90) T (cid:90) ω ( u f − g δ ) φ d x (cid:48) d t. (5.5)16hen taking φ = v in (5.4) and φ = u h in (5.5) and the identity (cid:82) T (cid:82) Ω v∂ αt u h d x d t = (cid:82) T (cid:82) Ω u htR ∂ αT v d x d t [25,p. 76, Lemma 2.7] (in view of the zero initial / terminal conditions) and subtracting the two identities give (cid:90) T (cid:90) Ω Rhv d x d t = (cid:90) T (cid:90) ω ( u f − g δ ) u h d x (cid:48) d t. This and the definition of the derivative J (cid:48) ( f ) shows the assertion.The next result gives the regularity of the adjoint variable v using the perturbation argument. Theorem 5.1.
Let g δ ∈ L (0 , T ; L ( ω )) , and Assumption H11 be fulfilled. Then there exists a uniquesolution v ∈ L (0 , T ; H − (cid:15) (Ω)) ∩ W α − (cid:15), (0 , T ; L (Ω)) , for any (cid:15) > , to the adjoint problem (5.2) .Proof. For any fixed t ∈ [0 , T ], let N ( t ) be the Neumann map defined by φ = N ( t ) ψ , with φ solving A ( t ) φ = 0 , in Ω ,∂ x d φ ( x (cid:48) , − (cid:96) ) = 0 , on ω,∂ x d φ ( x (cid:48) , (cid:96) ) = ψ, on ω,φ ( x ) = 0 , on ∂ω × ( − (cid:96), (cid:96) ) . (5.6)It is known that (cid:107) N ( t ) ψ (cid:107) H (Ω) ≤ c (cid:107) ψ (cid:107) L ( ω ) with a range R ( N ( t )) = D ( ˜ A − (cid:15) ), for any small (cid:15) > ψ = u f − g δ ∈ L (0 , T ; L ( ω )). Case i: A ( t ) ≡ A ∗ and N ( t ) ≡ N ∗ . Note that the solution v to problem (5.2) can be represented by usingthe operators ˜ E ( t ) and N ( t ) [24, Section 2.2]: v ( t ) = (cid:90) Tt ˜ A ∗ ˜ E ∗ ( s − t ) N ∗ ψ ( s ) d s. (5.7)Thus, by Lemma 2.1 and Young’s inequality, for any θ ∈ [0 , ) and (cid:15) ∈ (0 , − θ − ), there holds (cid:107) v (cid:107) L (0 ,T ; D ( ˜ A θ ∗ )) ≤ (cid:90) T (cid:107) ˜ A θ + + (cid:15) ∗ E ∗ ( t ) (cid:107) d t (cid:16) (cid:90) T (cid:107) ˜ A − (cid:15) ∗ N ∗ ψ ( t ) (cid:107) L (Ω) d t (cid:17) ≤ c (cid:15) (cid:107) ψ (cid:107) L (0 ,T ; L ( ω )) . It follows from equation (5.2) that tR ∂ αT v + ˜ A ∗ ( v − N ∗ ψ ) = 0 , which implies tR ∂ αT v ∈ L (0 , T ; D ( ˜ A θ − ∗ )). Then by a similar argument as in the proof of Lemma 2.3,we deduce v ∈ W α, (0 , T ; D ( ˜ A θ − ∗ )) (see also [19, Theorem 2.1]). Then by interpolation, we derive v ∈ W α − (cid:15), (0 , T ; L (Ω)) for any (cid:15) > t I − αT v ( t ) = (cid:90) Tt ˜ A ∗ ˜ S ∗ ( s − t ) N ∗ ψ ( s ) d s, and Young’s inequality, we deduce for any θ ∈ (1 − α , ) (cid:107) t I − αT v ( t ) (cid:107) L (Ω) ≤ (cid:90) Tt (cid:107) ˜ A − θ ∗ ˜ S ∗ ( s − t ) (cid:107) (cid:107) ˜ A θ ∗ N ∗ ψ ( s ) (cid:107) L (Ω) d s ≤ c (cid:90) Tt ( s − t ) − (1 − θ ) α (cid:107) ψ ( s ) (cid:107) L ( ω ) d s ≤ c (cid:16) (cid:90) Tt ( s − t ) − − θ ) α d s (cid:17) (cid:107) ψ (cid:107) L (0 ,T ; L ( ω )) ≤ c ( T − t ) − (1 − θ ) α (cid:107) ψ (cid:107) L (0 ,T ; L ( ω )) . Therefore the terminal condition t I − αT v ( T ) = 0 holds. Case ii: time-dependent elliptic operator A ( t ). We rewrite the adjoint problem (5.2) as17 R ∂ αT v ( t ) + ˜ A ∗ ( v − N ψ )( t ) = (cid:0) ˜ A ∗ − ˜ A ( t ) (cid:1) ( v − N ψ )( t ) , (5.8)with t I − αT v ( T ) = 0. Then the solution v ( t ) can be represented as v ( t ) = (cid:90) Tt ˜ A ∗ ˜ E ∗ ( s − t ) N ( t ) ψ ( s ) d s + (cid:90) Tt ˜ A ∗ ˜ E ∗ ( s − t ) (cid:0) ˜ A ∗ − ˜ A ( s ) (cid:1) ( v ( s ) − N ( s ) ψ ( s )) d s. (5.9)By Lemmas 2.1 and 2.2, we deduce that for any θ ∈ [0 , ) and (cid:15) ∈ (0 , − θ − ), (cid:107) v ( t ∗ ) (cid:107) D ( ˜ A θ ) ≤ c (cid:90) Tt ∗ ( s − t ∗ ) − ( θ + + (cid:15) ) α (cid:107) ψ ( s ) (cid:107) L ( ω ) d s + c (cid:90) Tt ∗ (cid:107) v ( s ) (cid:107) D ( ˜ A θ ) d s + c (cid:90) Tt ∗ (cid:107) ψ ( s ) (cid:107) L ( ω ) d s. (5.10)Then squaring both sides of (5.10) and integrating over [ t , T ] lead to (cid:107) v (cid:107) L ( t ,T ; D ( ˜ A θ )) ≤ c (cid:15) (cid:107) ψ (cid:107) L (0 ,T ; L ( ω )) + c (cid:90) Tt (cid:107) v (cid:107) L ( t,T ; D ( ˜ A θ )) d t. This together with Gronwall’s inequality implies that for any θ ∈ [0 , ) (cid:107) v (cid:107) L (0 ,T ; D ( ˜ A θ )) ≤ c (cid:107) ψ (cid:107) L (0 ,T ; L ( ω )) . Then it follows from (5.8) that tR ∂ αT v ∈ L (0 , T ; D ( ˜ A θ − ∗ )), and by Lemma 2.3, v ∈ W α, (0 , T ; D ( ˜ A θ − ∗ )).Then by interpolation, we derive u ∈ W α − (cid:15), (0 , T ; L (Ω)) for any (cid:15) > Remark 5.1.
It follows from Theorem 5.1 that for data g δ ∈ L (0 , T ; L ( ω )) , the gradient J (cid:48) ( f ) belongsto L (0 , T ; H − (cid:15) ( ω )) ∩ W α − (cid:15), (0 , T ; L ( ω )) for any small (cid:15) > , if the factor R is smooth, and further, t I − αT J (cid:48) ( f )( x (cid:48) , T ) = 0 for x (cid:48) ∈ ω and J (cid:48) ( f )( x (cid:48) , t ) = 0 for ( x (cid:48) , t ) ∈ ∂ω × (0 , T ) . These conditions will impactthe convergence behavior of the conjugate gradient method, dependent of the regularity of f † . Now we can describe the conjugate gradient method [3] for minimizing J . The complete procedureis listed in Algorithm 1. In the algorithm, Steps 6-7 compute the conjugate descent direction, and Step 8computes the optimal step size using the sensitivity problem. In general, the algorithm converges within tensof iterations; see the numerical experiments below. At each iteration, the algorithm involves three forwardsolves (direct problem, adjoint problem and sensitivity problem), which represent the main computationaleffort. For the stopping criterion at Step 11, we employ the discrepancy principle [33, 7, 15], i.e., k ∗ = arg min { k ∈ N : (cid:107) u f k − g δ (cid:107) L (0 ,T ; L ( ω )) ≤ cδ } , (5.11)where c > δ = (cid:107) g † − g δ (cid:107) L (0 ,T ; L ( ω )) is the noise level of the data g δ .Algorithm 1 can also be applied to ISPd, by viewing the zero Dirichlet data on ω × { (cid:96) } × (0 , T ) as themeasurement, and then the measurement ∂ x d u | ω ×{ (cid:96) }× (0 ,T ) as the Neumann data on ω × { (cid:96) } × (0 , T ) forproblem (3.1). However, the discrepancy principle (5.11) cannot be applied directly, due to a lack of thenoise level for the Dirichlet boundary data. Now we present several examples to illustrate the feasibility of recovering f . The domain Ω is taken to bethe unit square Ω = ( − , ) , with ω = ( − , ), q ≡
0, and the final time T = 1. The direct and adjointproblems are all discretized by the standard continuous piecewise linear Galerkin method in space, andbackward Euler convolution quadrature in time; see [18] for the error analysis for relevant direct problemsand the review [17] for various numerical schemes. The domain Ω is first divided into M small squares eachof width 1 /M , and then further divided into triangles by connecting the upper right with the lower left vertexof each square to obtain a uniform triangulation. For the inversion step, we take M = 100 and N = 1000.The same spatial and temporal mesh is used for approximating f . The factor R ( x, t ) is fixed at R ≡
1. Theexact data g † on the lateral boundary ω × { (cid:96) } × (0 , T ) is obtained by solving the direct problem (1.1) withthe exact f † on a finer mesh. The noisy boundary data g δ is generated from the exact data g † by g δ ( x (cid:48) , t ) = g † ( x (cid:48) , t ) + ε (cid:107) g † (cid:107) L ∞ ( ω × (0 ,T )) ξ ( x (cid:48) , t ) , ∀ ( x (cid:48) , t ) ∈ ω × (0 , T ) , lgorithm 1 Conjugate gradient method for minimizing the functional J in (5.1). Initialize f , and set k = 0. for k = 0 , . . . , K do Solve for u k from problem (3.1) with F = f k R , and compute the residual r k = u k | ω × (cid:96) × (0 ,T ) − g δ . Solve for v k from problem (5.2) with r k . Compute the gradient J (cid:48) ( f k ) by (5.3). Compute the conjugate coefficient γ k by γ k = , k = 0 , (cid:107) J (cid:48) ( f k ) (cid:107) L (0 ,T ; L ( ω )) (cid:107) J (cid:48) ( f k − ) (cid:107) L (0 ,T ; L ( ω )) , k ≥ . Compute the conjugate direction d k by d k = − J (cid:48) ( f k ) + γ k d k − . Solve for u d k from problem (3.1) with F = d k R . Compute the step size s k by s k = − ( u d k , r k ) L (0 ,T ; L ( ω )) (cid:107) u d k (cid:107) L (0 ,T ; L ( ω )) . Update the source component f k +1 = f k + s k d k . Check the stopping criterion. end for where ξ ( x (cid:48) , t ) follows the standard normal distribution, and ε denotes the relative noise level. In Algorithm1, the maximum number of CG iterations is fixed at 50, and the constant c in (5.11) is taken to be c = 1 . f by the L error E ( ˆ f ) defined by E ( ˆ f ) = (cid:107) ˆ f − f † (cid:107) L (0 ,T ; L ( ω )) . First we illustrate the case of time-independent coefficients.
Example 5.1. a ( x , x ) = 1 + sin( πx ) x (1 − x ) and f † ( x , t ) = (0 . − x ) t ( T − t ) e t . The numerical results for Example 5.1 are presented in Table 1, where the numbers in the bracket denotethe stopping index determined by the discrepancy principle (5.11). For noisy data g δ , the method reachesconvergence within ten iterations, and thus it is fairly efficient. It is observed that as the relative noiselevel (cid:15) increases from zero to 5e-2, the error e also increases, whereas the required number of CG iterationsdecreases. For a fixed noise level ε , the reconstruction error E tends to decrease with the order α , and allthe reconstructions are fairly accurate; see Fig. 1 for typical reconstructions and the associated pointwiseerrors e = ˆ f − f † . These results clearly show the feasibility of recovering f from the lateral boundary data,corroborating the theoretical results in [23].Table 1: The reconstruction errors for Example 5.1. α \ ε χ S denotes the characteristicfunction of a set S . 19a) exact (b) ε =1e-2 (c) ε =5e-2Figure 1: Reconstructions and the pointwise errors for Example 5.1 with ε =1e-2 and ε =5e-2. Example 5.2.
The diffusion coefficient a is given by a ( x , x , t ) = (1 + sin( πx ) x (1 − x ))(1 + sin t ) , andconsider two different source components. (i) f † ( x , t ) = sin( x + ) πt ( T − t ) e t . (ii) f † ( x , t ) = sin( x + ) πt ( T − t ) e t χ [0 , . ( t ) . In case (i), f is smooth in time, but it is discontinuous for case (ii). The results for Example 5.2 are shownin Table 2. The results for case (i) are largely comparable with that for Example 5.1, and all the observationsremain valid; see also Fig. 2. The behavior of ISPn is largely independent of the fractional order α , due tothe good regularity and compatibility of f † . In sharp contrast, the results for case (ii) exhibit a differenttrend: for a fixed noise level (cid:15) , the reconstruction error E increases with the order α , and also it takes moreCG iterations to reach the convergence (see also Fig. 4). This is attributed to the discontinuity in time of f † and the regularity of the adjoint v in problem (5.2): the temporal regularity of the adjoint v increasessteadily with α , cf. Theorem 5.1, which makes it increasingly harder to approximate a discontinuous f † .This is clearly visible from the error plots in Fig. 3, where the errors around the discontinuity dominate.This is especially pronounced for α = 0 .
50 and α = 0 . α \ ε ε =1e-2 (c) ε =5e-2Figure 2: Reconstructions and the pointwise errors for Example 5.2(i) with ε =1e-2 and ε =5e-2.(a) exact (b) α = 0 .
25 (c) α = 0 .
50 (d) α = 0 . ε =1e-2. Now we present two examples for ISPd, with the setting similar to Example 5.2.
Example 5.3.
The diffusion coefficient a is given by a ( x , x , t ) = (1 + sin( π ( x + 0 . . − x ))(1 + sin t ) ,and consider two different source components. (i) f † ( x , t ) = sin( x + ) πt ( T − t ) e t . (ii) f † ( x , t ) = sin( x + ) πt ( T − t ) e t χ [0 , . ( t ) . Note that case (ii) does not satisfy the condition of Theorem 4.1. The numerical results for Example 5.3are shown in Table 3, where the stopping index is taken so that the reconstruction error E is smallest (since21a) α = 0 .
25 (b) α = 0 . α \ ε (cid:15) , the error E increases with α , and also it takes more CG iterations to reach convergence,due to the mismatch between the temporal regularity of f and J (cid:48) ( f ). This is also clear from the error plotsin Fig. 6, where the errors around the discontinuity become increasingly dominating as α increases.These numerical results indicate that indeed it is feasible to recover a space-time dependent source fromthe lateral boundary observation in a cylindrical domain for both time-independent and dependent diffusioncoefficients, and standard regularization techniques, e.g., conjugate gradient method (when equipped withthe discrepancy principle (5.11)), can deliver accurate reconstructions for both exact and noisy data. Thisprovides numerical evidences to the theoretical results in Theorems 3.1 and 4.1. References [1] P. Acquistapace, F. Flandoli, and B. Terreni. Initial-boundary value problems and optimal control fornonautonomous parabolic systems.
SIAM J. Control Optim. , 29(1):89–118, 1991.[2] E. E. Adams and L. W. Gelhar. Field study of dispersion in a heterogeneous aquifer: 2. spatial momentsanalysis.
Water Res. Research , 28(12):3293–3307, 1992.[3] O. M. Alifanov, E. A. Artyukhin, and S. V. Rumyantsev.
Extreme Methods for Solving Ill-PosedProblems with Applications to Inverse Heat Transfer Problems . Begell House, New York, 1995.22a) exact (b) ε =1e-2 (c) ε =5e-2Figure 5: Reconstructions and the pointwise errors for Example 5.3(i) with ε =1e-2 and ε =5e-2.(a) exact (b) α = 0 .
25 (c) α = 0 .
50 (d) α = 0 . ε =1e-2.[4] H. Amann. Compact embeddings of vector-valued Sobolev and Besov spaces. Glas. Mat. Ser. III ,35(55)(1):161–177, 2000.[5] M. Berggren. Approximations of very weak solutions to boundary-value problems.
SIAM J. Numer.Anal. , 42(2):860–877, 2004.[6] J. Bourgain. Vector-valued singular integrals and the H -BMO duality. In Probability Theory andHarmonic Analysis (Cleveland, Ohio, 1983) , pages 1–19. Dekker, New York, 1986.[7] H. W. Engl, M. Hanke, and A. Neubauer.
Regularization of Inverse Problems . Kluwer Academic,Dordrecht, 1996. 238] K. Fujishiro and Y. Kian. Determination of time dependent factors of coefficients in fractional diffusionequations.
Math. Control Relat. Fields , 6(2):251–269, 2016.[9] P. Gaitan and Y. Kian. A stability result for a time-dependent potential in a cylindrical domain.
InverseProblems , 29(6):065006, 18, 2013.[10] P. Grisvard.
Elliptic Problems in Nonsmooth Domains . Pitman, Boston, MA, 1985.[11] Y. Hatano and N. Hatano. Dispersive transport of ions in column experiments: An explanation oflong-tailed profiles.
Water Res. Research , 34(5):1027–1033, 1998.[12] D. Henry.
Geometric Theory of Semilinear Parabolic Equations . Springer-Verlag, Berlin-New York,1981.[13] T. Hyt¨onen, J. van Neerven, M. Veraar, and L. Weis.
Analysis in Banach Spaces. Vol. I. Martingalesand Littlewood-Paley Theory . Springer, Cham, 2016.[14] V. Isakov.
Inverse Source Problems . AMS, Providence, RI, 1990.[15] K. Ito and B. Jin.
Inverse Problems: Tikhonov Theory and Algorithms . World Scientific Publishing Co.Pte. Ltd., Hackensack, NJ, 2015.[16] D. Jiang, Z. Li, Y. Liu, and M. Yamamoto. Weak unique continuation property and a related inversesource problem for time-fractional diffusion-advection equations.
Inverse Problems , 33(5):055013, 22,2017.[17] B. Jin, R. Lazarov, and Z. Zhou. Numerical methods for time-fractional evolution equations withnonsmooth data: a concise overview.
Comput. Methods Appl. Mech. Engrg. , 346:332–358, 2019.[18] B. Jin, B. Li, and Z. Zhou. Subdiffusion with a time-dependent coefficient: analysis and numericalsolution.
Math. Comp. , 88(319):2157–2186, 2019.[19] B. Jin, B. Li, and Z. Zhou. Pointwise-in-time error estimates for an optimal control problem withsubdiffusion constraint.
IMA J. Numer. Anal. , 40(1):377–404, 2020.[20] B. Jin, B. Li, and Z. Zhou. Subdiffusion with time-dependent coefficients: improved regularity andsecond-order time stepping.
Numer. Math. , 145(4):883–913, 2020.[21] B. Jin and W. Rundell. A tutorial on inverse problems for anomalous diffusion processes.
InverseProblems , 31(3):035003, 40, 2015.[22] Y. Kian, E. Soccorsi, Q. Xue, and M. Yamamoto. Identification of time-varying source term in time-fractional diffusion equations. Preprint, arXiv:1911.09951, 2019.[23] Y. Kian and M. Yamamoto. Reconstruction and stable recovery of source terms and coefficients ap-pearing in diffusion equations.
Inverse Problems , 35(11):115006, 24, 2019.[24] Y. Kian and M. Yamamoto. Well-posedness for weak and strong solutions of non-homogeneous initialboundary value problems for fractional diffusion equations. Preprint, arXiv:2004.14305, 2020.[25] A. A. Kilbas, H. M. Srivastava, and J. J. Trujillo.
Theory and Applications of Fractional DifferentialEquations . Elsevier Science B.V., Amsterdam, 2006.[26] N. Kinash and J. Janno. An inverse problem for a generalized fractional derivative with an applica-tion in reconstruction of time- and space-dependent sources in fractional diffusion and wave equations.
Mathematics , 7:1138, 2019.[27] A. Kubica and M. Yamamoto. Initial-boundary value problems for fractional diffusion equations withtime-dependent coefficients.
Fract. Calc. Appl. Anal. , 21(2):276–311, 2018.[28] Z. Li and Z. Zhang. Unique determination of fractional order and source term in a fractional diffusionequation from sparse boundary data.
Inverse Problems , 36:115013, 2020.2429] J.-L. Lions and E. Magenes.
Non-homogeneous Boundary Value Problems and Applications. Vol. I .Springer-Verlag, New York-Heidelberg, 1972.[30] Y. Liu, Z. Li, and M. Yamamoto. Inverse problems of determining sources of the fractional partialdifferential equations. In
Handbook of Fractional Calculus with Applications. Vol. 2 , pages 411–429. DeGruyter, Berlin, 2019.[31] R. Metzler, J. H. Jeon, A. G. Cherstvy, and E. Barkai. Anomalous diffusion models and their properties:non-stationarity, non-ergodicity, and ageing at the centenary of single particle tracking.
Phys. Chem.Chem. Phys. , 16(44):24128–24164, 2014.[32] R. Metzler and J. Klafter. The random walk’s guide to anomalous diffusion: a fractional dynamicsapproach.
Phys. Rep. , 339(1):77, 2000.[33] V. A. Morozov. On the solution of functional equations by the method of regularization.
Soviet Math.Dokl. , 7:414–417, 1966.[34] R. R. Nigmatullin. The realization of the generalized transfer equation in a medium with fractalgeometry.
Phys. Stat. Sol. B , 133:425–430, 1986.[35] W. Rundell and Z. Zhang. Recovering an unknown source in a fractional diffusion problem.
J. Comput.Phys. , 368:299–314, 2018.[36] K. Sakamoto and M. Yamamoto. Initial value/boundary value problems for fractional diffusion-waveequations and applications to some inverse problems.
J. Math. Anal. Appl. , 382(1):426–447, 2011.[37] K. Sakamoto and M. Yamamoto. Inverse source problem with a final overdetermination for a fractionaldiffusion equation.
Math. Control Relat. Fields , 1(4):509–518, 2011.[38] M. Slodiˇcka. Uniqueness for an inverse source problem of determining a space-dependent source in anon-autonomous time-fractional diffusion equation.
Frac. Calc. Appl. Anal. , 23(6):1702–1711, 2020.[39] T. Wei, X. L. Li, and Y. S. Li. An inverse time-dependent source problem for a time-fractional diffusionequation.
Inverse Problems , 32(8):085003, 24, 2016.[40] H. Ye, J. Gao, and Y. Ding. A generalized Gronwall inequality and its application to a fractionaldifferential equation.
J. Math. Anal. Appl. , 328(2):1075–1081, 2007.[41] Y. Zhang and X. Xu. Inverse source problem for a fractional diffusion equation.
Inverse Problems ,27(3):035010, 12, 2011.[42] Z. Zhang. An undetermined coefficient problem for a fractional diffusion equation.
Inverse Problems ,32(1):015011, 21, 2016.[43] F. Zimmermann. On vector-valued Fourier multiplier theorems.