arXiv:1810.04836v2 [math.AP] 5 Feb 2019 WELL-POSEDNESS OF TIME-FRACTIONAL ADVECTION-DIFFUSION-REACTION EQUATIONS IN “FCAA” JOURNAL William McLean1 , Kassem Mustapha2 , Raed Ali2 , Omar Knio3 Abstract We establish the well-posedness of an initial-boundary value problem for a general class of linear time-fractional, advection-diffusion-reaction equations, allowing space- and time-dependent coefficients as well as initial data that may have low regularity. Our analysis relies on novel energy methods in combination with a fractional Gronwall inequality and properties of fractional integrals. MSC 2010 : Primary 26A33; Secondary 35A01, 35A02, 35B45, 35D30, 35K57, 35Q84, 35R11. Key Words and Phrases: Fractional PDE, weak solution, Volterra integral equation, fractional Gronwall inequality, Galerkin method. 1. Introduction The main scope of this paper is to investigate the existence and uniquesness of the weak solution of a linear, time-fractional problem of the form
~ + a∂ 1−α u + bu = g ~ ∂ 1−α u − Gu (1.1) ∂t u − ∇ · κ∇∂ 1−α u − F t t t for x ∈ Ω and 0 < t ≤ T . The parameter α in the fractional derivative lies in the range 0 < α < 1, and the spatial domain Ω ⊆ Rd (d ≥ 1) ~ a and b, as well as the is bounded and Lipschitz. The coefficients F~ , G, source term g, are assumed to be known functions of x and t, whereas the generalized diffusivity κ = κ(x) may depend only on x but is permitted to be a real, symmetric positive-definite matrix. We impose homogeneous Dirichlet boundary conditions, u(x, t) = 0 for x ∈ ∂Ω and 0 ≤ t ≤ T , c Year Diogenes Co., Sofia pp. xxx–xxx, DOI: ...................... (1.2) 2 W. McLean, K. Mustapha, R. Ali, O. Knio and the initial condition u(x, 0) = u0 (x) for x ∈ Ω. (1.3) The Riemann–Liouville fractional derivative of order 1 − α is defined via the fractional integral of order α: with ωα (t) = tα−1 /Γ(α) we have Z t ωα (t − s)v(x, s) ds. ∂t1−α v(x, t) = ∂t I α v(x, t) where I α v(x, t) = 0 Wpk (Ω) We denote by the usual Sobolev space of functions whose partial derivatives of order k or less belong to Lp (Ω). The following regularity assumptions on the coefficients will be used:
~ ∈ C 2 [0, T ]; W 1 (Ω)d , κ ∈ L∞ (Ω)d×d , F~ , G ∞ (1.4)
a, b ∈ C 1 [0, T ]; L∞ (Ω) . In addition, to ensure that the spatial operator v 7→ −∇·(κ∇v) is uniformly elliptic on Ω, we assume that the minimal eigenvalue of κ(x) is bounded away from zero, uniformly for x ∈ Ω. Different classes of time-fractional models (typically considered only for scalar κ) arise as special cases of (1.1), including ~ = 0, • fractional Fokker–Planck equations [3, 9, 14, 27], when G a = b = 0 and g = 0; ~ =G ~ = 0; • fractional reaction-diffusion equations [10, 11], when F ~ ~ • fractional cable equations [17], when F = G = 0; • fractional advection-dispersion (or fractional convection-diffusion) ~ = F~ (x), G ~ = 0 and a = b = 0. equations [23], when F Consider the simplest non-trivial case, when κ is the identity matrix with ~ = 0, a = b = 0 and g = 0, so that (1.1) reduces to the fracF~ = G tional subdiffusion equation: ∂t u − ∇2 ∂t1−α u = 0. Let ϕ denote a Dirichlet eigenfunction of the Laplacian on Ω, with corresponding eigenvalue λ > 0, that is, −∇2 ϕ = λϕ in Ω with ϕ|∂Ω = 0. For the special choice of initial data u0 = ϕ(x), the solution of the initial-boundary value problem (1.1)–(1.3) the separable form u(x, t) = Eα (−λtα )ϕ(x), where P∞ has n Eα (z) = n=0 z /Γ(1 + nα) is the Mittag–Leffler function. Notice that ∂tm u = O(tα−m ) as t → 0. Moreover, we can extend the classical method of separation of variables for the heat equation to construct a series solution for arbitrary initial data u0 ∈ L2 (Ω), and the regularity properties of the solution u follow from this representation [25]. Such an explicit construction is no longer possible for the solution of the general equation (1.1). Instead, we proceed by formally integrating (1.1) in time, multiplying both sides by a test function v, and applying the first WELL-POSEDNESS . . . 3 Green identity over Ω to arrive at the weak formulation Z t ~ (s)∂s1−α u(s) − G(s)u(s), ~ κ∇∂s1−α u(s) − F ∇v ds hu(t), vi + 0 Z t Z t hg(s), vi ds (1.5) + a(s)∂s1−α u(s) + b(s)u(s), v ds = hu0 , vi + 0 0 H01 (Ω), where we have suppressed the dependence of the funcfor all v ∈ tions on x, and where h·, ·i denotes the inner product in L2 (Ω) or L2 (Ω)d . Numerical methods for particular cases of (1.1) were extensively studied over the last two decades [7, 12, 13, 16, 18, 20, 21, 28, 31, 33]. However, due to various types of mathematical difficulties, proof of the well-posedness of the continuous problem is almost missing despite its importance, apart from ~ = 0 and a = b = 0. In this paper, we address the case [25] when F~ = G these fundamental questions. A related preprint [19] treats the fractional ~ = 0 and a = b = 0) via a Fokker–Planck equation (that is, the case G different, and somewhat simpler, chain of estimates that, for instance, does not use the quadratic operator Qµ1 defined below in Section 2. ~ = 0 and b = 0, If the coefficients F~ and a are independent of t, and if G 1−α then by applying the fractional integration operator I to both sides of (1.1) we obtain C α ∂t u − ∇ · (κ∇u − F~ u) + au = g̃, (1.6) where C ∂tα u = I 1−α ∂t u denotes the Caputo fractional derivative and where g̃ = I 1−α g. Existence and uniqueness results for (1.6) were studied by several authors, including Zacher [34], Alikhanov [1], Sakamoto and Yamamoto [30] and Kubica and Yamamoto [15]. Some of these papers include results for time-dependent coefficients, but in that case (1.6) is no longer equivalent to (1.1). To recast the weak formulation (1.5) as a Volterra integral equation, we introduce two bounded linear operators, firstly K1 (t) : H01 (Ω) → H −1 (Ω) defined by hK1 (t)v, wi = hκ∇v, ∇wi − hF~ (t)v, ∇wi + ha(t)v, wi for v, w ∈ H01 (Ω), and secondly K2 (t) : L2 (Ω) → H −1 (Ω) by ~ hK2 (t)v, wi = hb(t)v, wi − hG(t)v, ∇wi for v ∈ L2 (Ω) and w ∈ H01 (Ω). The variational problem (1.5), subject to the initial condition (1.3), can then be written more succinctly as Z t Z t   1−α g(s) ds. (1.7) u(t) + K1 (s)∂s u(s) + K2 (s)u(s) ds = f (t) ≡ u0 + 0 0 4 W. McLean, K. Mustapha, R. Ali, O. Knio Integrating by parts, and using a dash to indicate a derivative in time, leads to Z t Z t 1−α α K1′ (s)I α u(s) ds K1 (s)∂s u(s) ds = K1 (t)I u(t) − 0 0  Z t Z t ′ = ωα (z − s)K1 (z) dz u(s) ds, ωα (t − s)K1 (t) − with K1′ (t) : 0 s 1 −1 H0 (Ω) → H (Ω) given by hK1′ (t)v, wi = −hF~ ′ (t)v, ∇wi Thus, u satisfies u(t) + Z + ha′ (t)v, wi. t K(t, s)u(s) ds = f (t) for 0 ≤ t ≤ T , (1.8) 0 where K(t, s) : H01 (Ω) → H −1 (Ω) is the weakly-singular, operator-valued kernel Z t ωα (z − s)K1′ (z) dz. (1.9) K(t, s) = ωα (t − s)K1 (t) + K2 (s) − s Following some technical preliminaries in Section 2, we apply the Galerkin method in Section 3 to project the problem (1.8) to a finite dimensional subspace X ⊆ H01 (Ω), thereby obtaining an approximate solution uX : [0, T ] → X. Using delicate energy arguments and a fractional Gronwall inequality, we prove a priori estimates for uX that are uniform with respect to the dimension of X, allowing us in Section 4 (Theorems 4.1 and 4.2) to establish the existence and uniqueness of a weak solution u to the original problem (1.1)–(1.3), provided (1.4) holds. The regularity of the weak solution u will be studied in a companion paper [26]. 2. Preliminaries and notations Our subsequent analysis makes frequent use of two quadratic operators defined, for µ ≥ 0 and 0 ≤ t ≤ T , by Z t Z t µ µ µ kI µ φk2 ds. hφ, I φi ds and Q2 (φ, t) = Q1 (φ, t) = 0 0 I 0φ These operators coincide when µ = 0 because = φ, and so we write Q0 = Q01 = Q02 . If we put φ(t) = 0 for t > T , then the Laplace transRT µ φ(z) = z −µ φ̂(z), so form φ̂(z) = 0 e−zt φ(t) dt is an entire function and Id it follows by the Plancherel Theorem that Z cos(πµ/2) ∞ −µ µ y kφ̂(iy)k2 dy ≥ 0, (2.10) Q1 (φ, T ) = π 0 WELL-POSEDNESS . . . 5 assuming that φ is real-valued; see also [29, Theorem 2]. Note that because ωµ ∈ L1 (0, T ), the fractional integral defines a bounded linear operator

I µ : Lp (0, T ), L2 (Ω) → Lp (0, T ), L2 (Ω) for 1 ≤ p ≤ ∞. (2.11) Also, I µ+ν = I µ I ν because ωµ ∗ ων = ωµ+ν for µ > 0 and ν > 0; here, ∗ denotes the Laplace convolution. The next four lemmas establish key inequalities satisfied by Qµ1 and Qµ2 . Lemma 2.1. If 0 < α < 1 and ǫ > 0, then Z t Qα1 (φ, t) + ǫ Qα1 (ψ, t), hφ, I α ψi ds ≤ 4ǫ(1 − α)2 0 2tα Qα (φ, t), Qα2 (φ, t) ≤ 1−α 1 Qα1 (φ, t) ≤ 2tα Q0 (φ, t), Z t tα Q0 (φ, t) + ǫ Qα1 (ψ, t). hφ, I α ψi ds ≤ 2ǫ(1 − α)2 0 (2.12) (2.13) (2.14) (2.15) P r o o f. The first three inequalities are proved by Le, McLean and Mustapha [18, Lemma 3.2]. The fourth inequality follows from (2.12) and (2.14). ✷
For the next result, note that if φ ∈ W11 (0, T ); X for a normed space X, then φ : [0, T ] → X is absolutely continuous and (∂t I αφ − I α∂t φ)(t) = φ(0)ωα (t) for 0 < t ≤ T . (2.16)
Lemma 2.2. If 0 < α ≤ 1, then for φ ∈ L2 (0, t), L2 (Ω) , Z t ωα (t − s)Qα1 (φ, s) ds. Qα2 (φ, t) ≤ 2 0
P r o o f. Assume first that φ ∈ W11 (0, T ), L2 (Ω) and let ψ = I α φ. Since ψ(0) = 0, the Caputo fractional derivative of ψ is C α ∂t ψ = I 1−α (ψ ′ ) = (I 1−α ψ)′ − ψ(0)ω1−α = (I 1 φ)′ = φ. Recalling an identity of Alikhanov [2, Corollary1],
2 ψ(t), C ∂tα ψ(t) = C ∂tα kψk2 (t) Z s ′  Z t ψ (q) dq 2 1 α ds, + 2Γ(1 − α) 0 (t − s)1−α 0 (t − q)α 6 W. McLean, K. Mustapha, R. Ali, O. Knio we see that
2hφ, I α φi = 2hC ∂tα ψ, ψi ≥ C ∂tα kψk2 = I 1−α (kI α φk2 )′ , and thus (2.17)

′
′ I 1 kI α φk2 = I 2 kI α φk2 = I 1+α I 1−α kI α φk2

≤ 2I 1+α hφ, I α φi = 2I α I 1 hφ, I α φi , which is equivalent to the desired

inequality. Now let φ ∈ L2 (0, T ), L2 (Ω) , and choose φn ∈ W11 (0, T ), L2 (Ω) such RT that 0 kφn (t) − φ(t)k2 dt → 0 as n → ∞. Using (2.11) with µ = α and p = 2, it follows that Qα1 (φn , t) → Qα1 (φ, t) and Qα2 (φn , t) → Qα2 (φ, t), uniformly for t ∈ [0, T ], which implies the result in the general case. ✷ The next lemma will eventually enable us to establish pointwise (in time) estimates for u(t).
Lemma 2.3. Let 0 < α ≤ 1. If the function φ ∈ W11 (0, t); L2 (Ω) satisfies φ(0) = I α φ′ (0) = 0, then kφ(t)k2 ≤ 2ω2−α (t) Qα1 (φ′ , t). P r o o f. For α = 1, equality holds: Z t 1 ′ 2ω1 (t)Q1 (φ , t) = 2 hφ′ , φi ds = kφ(t)k2 . 0 For 0 < α < 1, applying the operator I 1 to both sides of (2.17) with φ′ in place of φ, and using I α φ′ (0) = 0, we observe that,
I 1−α kI α φ′ k2 (t) ≤ 2Qα1 (φ′ , t). (2.18) Put ψ(t) = I α φ′ . Since φ = I 1 φ′ = I 1−α ψ, 2 Z t 2 ω1−α (t − s)kψ(s)k ds kφ(tk ≤ 0 Z t Z t ω1−α (t − s)kψ(s)k2 ds ω1−α (t − s) ds ≤ 0 0   = ω2−α (t) I 1−α kI α φ′ k2 (t), and hence the desired result follows immediately after using (2.18). ✷ Lemma 2.4. If 0 ≤ µ ≤ ν ≤ 1, then Qν2 (φ, t) ≤ 2t2(ν−µ) Qµ2 (φ, t). P r o o f. See Le, McLean and Mustapha [18, Lemma 3.1]. ✷ WELL-POSEDNESS . . . 7 We will make essential use of the following fractional Gronwall inequality. Lemma 2.5. Let β > 0 and T > 0. Assume that a and b are nonnegative, non-decreasing functions on the interval [0, T ]. If q : [0, T ] → R is an integrable function satisfying Z t ωβ (t − s)q(s) ds for 0 ≤ t ≤ T , 0 ≤ q(t) ≤ a(t) + b(t) 0 then q(t) ≤ a(t)Eβ b(t)tβ
for 0 ≤ t ≤ T . P r o o f. See Dixon and McKee [8, Theorem 3.1]. ✷ Let M denote the operator of pointwise multiplication by t, that is, (Mφ)(t) = tφ(t), and note the commutator property MI µ − I µM = µI µ+1 , (2.19) for any real µ ≥ 0. We will need the following estimates involving the linear operator Bψµ defined (for suitable ψ and φ) by Z t µ µ ψ ′ (s) I µ φ(s) ds. (2.20) (Bψ φ)(t) = ψ(t) I φ(t) − 0

1 (0, T ); L (Ω)d and φ ∈ W 1 (0, T ); L (Ω) , Lemma 2.6. If ψ ∈ W∞ 2 ∞ 1 then there is a constant C (depending only on ψ, µ and T ) such that for 0 ≤ t ≤ T , Q0 (Bψµ φ, t) ≤ CQµ2 (φ, t), (2.21) Q0 (MBψµ φ, t) + Q0 (I 1 Bψµ φ, t) ≤ Ct2 Qµ2 (φ, t),

Q0 (MBψµ φ)′ , t ≤ CQµ2 (Mφ)′ , t + CQµ2 (Mφ, t) + CQµ2 (φ, t). (2.22) (2.23) P r o o f. The assumption on ψ implies that Z t µ 2 µ 2 k(I µ φ)(s)k2 ds, k(Bψ φ)(t)k ≤ Ck(I φ)(t)k + C 0 and (2.21) follows after integrating in time. By the Cauchy–Schwarz inequality, Z t µ µ 2 1 µ 2 2 2 k(Bψµ φ)(s)k2 ds, k(MBψ φ)(t)k + k(I Bψ φ)(t)k ≤ t k(Bψ φ)(t)k + t 0 8 W. McLean, K. Mustapha, R. Ali, O. Knio and (2.22) follows after integrating in time. The third identity in (2.19) implies that
MBψµ φ = ψ I µ Mφ + µI µ+1 φ − MI 1 (ψ ′ I µ φ) and therefore, differentiating with respect to t,

(MBψµ φ)′ = ψ ′ I µ Mφ+µI µ+1 φ +ψ (I µ Mφ)′ +µI µφ −(I 1 +M)(ψ ′ I µ φ). Thus, noting that (I µ Mφ)′ = I µ(Mφ)′ by (2.16), with kI µ+1 φ(t)k2 = kI 1 (I µ φ)(t)k2 ≤ tQµ2 (φ, t) and kI 1 (ψ ′ I µ φ)(t)k2 ≤ CtQµ2 (φ, t), we have k(MBψµ φ)′ (t)k2 ≤ CkI µ (Mφ)(t)k2 + CkI µ (Mφ)′ (t)k2 + Ck(I µ φ)(t)k2 + CtQµ2 (φ, t), so (2.23) follows after integrating in time. ✷ 3. The projected equation Suppose that X is a finite-dimensional subspace of H01 (Ω), equipped with the induced norm: kvkX = kvkH01 (Ω) . We define a bounded linear operator KX (t, s) : X → X in terms of K(t, s) in (1.9) by hKX (t, s)v, wi = hK(t, s)v, wi for v, w ∈ X and 0 ≤ s ≤ t ≤ T , and let fX (t) denote the L2 -projection onto X of f (t) from (1.7), that is, hfX (t), wi = hf (t), wi for w ∈ X and 0 ≤ t ≤ T . In this way, we arrive at a finite dimensional reduction of the Volterra equation (1.8), Z t KX (t, s)uX (s) ds = fX (t) for 0 ≤ t ≤ T . (3.24) uX (t) + 0 In the next theorem, we outline a self-contained proof of existence and uniqueness under relaxed assumptions on the coefficients in the fractional PDE (1.1). Similar results for scalar-valued kernels are shown by Linz [22, §3.4], Becker [4], and Brunner [5]. Henceforth, C will denote a generic constant that may depend on the coefficients in (1.1), the spatial domain Ω, the time interval [0, T ], the fractional exponent α, the parameter η, and the integer m in (1.4). However, any dependence on the subspace X is indicated explicitly by writing CX . We let Y = C([0, T ]; X) with the norm kvkY = max0≤t≤T kv(t)kX . WELL-POSEDNESS . . . 9 Theorem 3.1. Assume that the coefficients in (1.1) satisfy

1 ~ ∈ L∞ (0, T ); L∞ (Ω)d , κ ∈ L∞ (Ω)d×d , F~ ∈ W∞ (0, T ); L∞ (Ω)d , G

1 a ∈ W∞ (0, T ); L∞ (Ω) , b ∈ L∞ (0, T ); L∞ (Ω) . Assume, in addition, that the source term g : (0, T ] → L2 (Ω) is a measurable function satisfying kg(t)k ≤ M tη−1 for 0 < t ≤ T , (3.25) where M and η are positive constants, and that the initial data u0 ∈ L2 (Ω). Then, the weakly-singular Volterra integral equation (3.24) has a unique solution uX ∈ Y , and moreover kuX kY ≤ CX kfX kY ≤ CX (ku0 k + M ). P r o o f. Our assumptions on u0 and g ensure that fX ∈ Y . The kernel (1.9) has the form K(t, s) = ωα (t − s)G(t, s) + H(t, s), where G(t, s) = K1 (t) − Γ(α)(t − s) Z 1 0
ωα (y)K1′ s + (t − s)y dy and H(t, s) = K2 (s) for 0 ≤ s ≤ t ≤ T . Our assumptions on the coefficients of the fractional PDE (1.1) ensure that G and H are continuous mappings from the closed triangle △ = { (t, s) : 0 ≤ s ≤ t ≤ T } into the space of bounded linear operators H01 (Ω) → H −1 (Ω). Likewise, KX (t, s) = ωα (t − s)GX (t, s) + HX (t, s), where GX (t, s) : X → X and HX (t, s) : X → X are defined by hGX (t, s)v, wi = hG(t, s)v, wi and hHX (t, s)v, wi = hH(t, s)v, wi for (t, s) ∈ △ and v, w ∈ X. Since X is finite dimensional, GX and HX are continuous functions from △ into the space of bounded linear operators X → X. Hence, there is a positive constant γX such that kKX (t, s)vkX ≤ γX ωα (t − s)kvkX for (t, s) ∈ △ and v ∈ X, so we can define the Volterra operator KX : Y → Y by Z t KX (t, s)v(s) ds for 0 ≤ t ≤ T and v ∈ Y . KX v(t) = 0 We see that kKX vkY ≤ γX ω1+α (T )kvkY . In fact, using the semigroup property, Z t ωα (t − s)ωβ (s) ds = ωα+β (t), 0 10 W. McLean, K. Mustapha, R. Ali, O. Knio we obtain the following estimate for the operator norm of the nth power of KX , Z t n n n ωnα (t − s) ds = γX ω1+nα (T ) for n ≥ 1. kKX kY →Y ≤ γX max 0≤t≤T 0 P∞ n n It follows that the sum RX = n=1 (−1) KX defines a bounded linear operator with ∞ X n ω1+nα (T )γX = Eα (γX T α ) − 1. kRX kY →Y ≤ n=1 This operator is the resolvent for KX , that is, uX + KX uX = fX if and only if uX = fX − RX fX , implying the existence and uniqueness of uX ∈ Y , as well as the a priori estimate claimed in the theorem. ✷ For a scalar, weakly-singular, second-kind Volterra equation, it is known that if fX admits an expansion in powers of t and tα , then so does the solution uX ; see Lubich [24, Corollary 3], and also Brunner, Pedas and Vainikko [6, Theorem 2.1] (with ν = 1−α). To outline a proof that a similar
result holds for systems of Volterra equations, let Cαm = Cαm [0, T ]; X denote the space of continuous functions v : [0, T ] → X that are C m on the half-open interval (0, T ] and for which the seminorm |v|j,α = sup tj−α kv (j) (t)kX is finite for 1 ≤ j ≤ m. 0<t≤T We make Cαm into a Banach space by defining the obvious norm: kvkm,α = kvkY + m X |v|j,α . j=1 Theorem 3.2. Let m ≥ 1, and strengthen the assumptions (1.4) by requiring

1 ~,G ~ ∈ C m+1 [0, T ]; W∞ F (Ω)d and a, b ∈ C m [0, T ]; L∞ (Ω) . If u0 ∈ L2 (Ω) and g : (0, T ] → X is C m with kg(i−1) (t)k ≤ M tα−i for 1 ≤ i ≤ m, then uX ∈ Cαm and kuX km,α ≤ CX kfX km,α ≤ CX (ku0 k + M ). P r o o f. Our assumptions on u0 and g imply that fX ∈ Cαm . Using the substitution z = s + (t − s)y in (1.9), we find that if j + k ≤ m and 0 ≤ s < t ≤ T , then ∂tk (∂t + ∂s )j K(t, s)v H −1 (Ω) ≤ CX (t − s)α−1−k kvkH01 (Ω) for v ∈ H01 (Ω), WELL-POSEDNESS . . . 11 and, since X is finite dimensional, ∂tk (∂t + ∂s )j KX (t, s)v X ≤ CX (t − s)α−1−k kvkX for v ∈ X. Hence, the Volterra operator KX : Cαm → Cαm is compact [32, Theorem 6.1]. Theorem 3.1 implies that the homogeneous equation, uX + KX uX = 0, has only the trivial solution uX = 0, and therefore the inhomogeneous equation uX + KX uX = fX is well-posed not only in Y but also in Cαm . ✷ Our goal in the remainder of this section is to obtain bounds for kuX (t)k and k∇uX (t)k with constants that are independent of X. Our proof relies on a sequence of technical lemmas. To simplify our estimates, we rescale the time variable, if necessary, so that the minimal eigenvalue of κ is bounded below by unity:
(3.26) λmin κ(x) ≥ 1 for x ∈ Ω. In this way, hκ∇v, ∇vi ≥ k∇vk2 for v ∈
H01 (Ω), and we see from (2.10) that for (real-valued) φ ∈ C [0, T ]; H01 (Ω) , Z Z t cos(πµ/2) ∞ −µ µ y hκ∇φ̂(iy), ∇φ̂(iy)i dy hκI ∇φ, ∇φi ds = π 0 Z0 ∞ cos(πµ/2) ≥ y −µ k∇φ̂(iy)k2 dy, π 0 so Z Z t t hκI µ ∇φ, ∇φi ds ≥ 0 0 hI µ ∇φ, ∇φi ds = Qµ1 (∇φ, t). (3.27) Since (1.7) is equivalent to (1.8), if v ∈ X then  Z t KX (t, s)uX (s) ds, v 0 = Z t 0 K1 (s)∂s1−α uX , v ds + Z t K2 (s)uX (s), v ds 0 = κ(I α ∇uX )(t), ∇v − (B1 uX )(t), ∇v + (B2 uX )(t), v , where Z t  ~ F~ (s)∂s1−α φ(s) + G(s)φ(s) ds, 0 (3.28) Z t  1−α a(s)∂s φ(s) + b(s)φ(s) ds. B2 φ(t) = 0
Assuming φ ∈ Cα1 ([0, T ]; X , we may integrate by parts and use the notation (2.20) to write ~ 1 φ(t) = B ~ 1 = Bα + B1 B ~ ~ F G and B2 = Baα + Bb1 . (3.29) 12 W. McLean, K. Mustapha, R. Ali, O. Knio Thus, the solution of (3.24) satisfies ~ 1 uX )(t), ∇v + (B2 uX )(t), v huX (t), vi + hκ∇I α uX (t), ∇vi − (B = hfX (t), vi for v ∈ X, (3.30) which yields the following estimates (with C independent of X). Lemma 3.1. For 0 ≤ t ≤ T , the solution uX of the Volterra equation (3.24) satisfies the a priori estimates Qα1 (uX , t) + Qα2 (∇uX , t) ≤ CtαQ0 (fX , t) and Q0 (uX , t) + Qα1 (∇uX , t) ≤ CQ0 (fX , t). P r o o f. From (3.30), ~ 1 uX (t)k2 + 1 kB2 uX (t)k2 uX (t), v + hκ∇I α uX (t), ∇vi ≤ 21 k∇vk2 + 12 kB 2 + 12 kvk2 + fX (t), v . Choosing v = I αuX (t) we have hκ∇I α uX (t), ∇vi = hκ∇v, ∇vi ≥ k∇vk2 because of (3.26). Thus, after canceling the term 12 k∇vk2 and integrating in time, we see that ~ 1 uX , t) + 1 Q0 (B2 uX , t) + 1 Qα2 (uX , t) Qα1 (uX , t) + 12 Qα2 (∇uX , t) ≤ 12 Q0 (B 2 2 Z t + fX (s), I α uX (s) ds. (3.31) 0 Using the representation (3.29) and the achieved estimate (2.21), ~ 1 uX , t) ≤ 2Q0 (B α uX , t) + 2Q0 (B 1 uX , t) Q0 (B ~ ~ F G ≤ CQα2 (uX , t) + CQ12 (uX , t) ≤ CQα2 (uX , t), where, in the final step, we used Lemma 2.4. In the same way, Q0 (B2 uX , t) ≤ CQα2 (uX , t). Using (2.15) with φ = fX , ψ = uX and ǫ = 1/2, we deduce that Qα1 (uX , t) + 21 Qα2 (∇uX , t) ≤ CQα2 (uX , t) + Ctα Q0 (fX , t) + 21 Qα1 (uX , t). Hence, applying Lemma 2.2 with φ = uX , we can show that the function q(t) = Qα1 (uX , t) + Qα2 (∇uX , t) satisfies Z t ωα (t − s)Qα1 (uX , s) ds. q(t) ≤ Ctα Q0 (fX , t) + C 0 Since Qα1 (uX , s) ≤ q(s), Lemma 2.5 implies the first estimate. WELL-POSEDNESS . . . 13 ~ 1 uX )(t), ∇v = ∇ · B ~ 1 uX (t), v To show the second estimate, use − (B in (3.30) to obtain ~ 1 uX )(t)k2 huX (t), vi + κ∇I α uX (t), ∇v ≤ 21 kvk2 + 32 k∇ · (B + 32 k(B2 uX )(t)k2 + 23 kfX (t)k2 . Choosing v = uX (t), integrating in time, and using (3.27), we have α 1 0 2 Q (uX , t)+Q1 (∇uX , t) ~ 1 uX , t)+CQ0 (B2 uX , t)+CQ0 (fX , t). ≤ CQ0 (∇· B Since
∇ · (BFα~ uX )(t) = ∇ · F~ (t) I α uX (t) + F~ (t) · I α ∇uX (t) Z t 
∇ · F~ ′ (s) I αuX (s) + F~ ′ (s) · I α∇uX (s) ds (3.32) − 0 it follows that k∇ · (BFα~ uX )(t)k2 ≤ CkI α uX (t)k2 + CkI α∇uX (t)k2 Z t  kI α uX (s)k2 + kI α ∇uX (s)k2 ds, +C 0 ≤ CQα2 (uX , t) + CQα2 (∇uX , t). In the implying that · same way, Q0 (∇ · B 1~ uX , t) ≤ CQ12 (uX , t) + CQ12 (∇uX , t) and therefore, G by Lemma 2.4, Q0 (∇ B α~ uX , t) F ~ 1 uX , t) ≤ CQα2 (uX , t) + CQα2 (∇uX , t). Q0 (∇ · B Recall Q0 (B2 uX , t) ≤ CQα2 (uX , t) and let q(t) = Q0 (uX , t) + Qα1 (∇uX , t). It follows using Lemma 2.2 and (2.14) that q(t) ≤ CQα2 (uX , t) + CQα2 (∇uX , t) + CQ0 (fX , t) Z t   0 ωα (t − s) Qα1 (uX , s) + Qα1 (∇uX , s) ds ≤ CQ (fX , t) + C 0 Z t 0 α ≤ CQ (fX , t) + Ct ωα (t − s)q(s) ds. 0 We may now apply Lemma 2.5 to complete the proof. ✷ The function MuX (t) = tuX (t) satisfies a similar estimate to the first one in Lemma 3.1, but with an additional factor t2 on the right-hand side. Lemma 3.2. The solution uX of (3.24) satisfies Qα1 (MuX , t) + Qα2 (M∇uX , t) ≤ Ct2+α Q0 (fX , t) for 0 ≤ t ≤ T . 14 W. McLean, K. Mustapha, R. Ali, O. Knio P r o o f. Multiplying both sides of (3.30) by t, and applying the third identity in (2.19), we find that (since κ is independent of t) hMuX , vi + κ(I α M + αI α+1 )∇uX , ∇v ~ 1 uX , ∇vi + hM(fX − B2 uX ), vi, (3.33) = hMB whereas integrating (3.30) in time gives ~ 1 uX , ∇vi + I 1 (fX − uX − B2 uX ), v , hκI α+1 ∇uX , ∇vi = hI 1 B so, after eliminating hκI α+1 ∇uX , ∇vi, ~ 1 uX , ∇vi hMuX , vi + hκI α M∇uX , ∇vi = h(M − αI 1 )B + h(M − αI 1 )(fX − B2 uX ) + αI 1 uX , vi ~ 3 uX k2 + 1 kB4 uX k2 + 1 kvk2 + (M−αI 1 )fX +αI 1 uX , v , ≤ 12 k∇vk2 + 12 kB 2 2 ~ 3 φ = (M − αI 1 )B ~ 1 φ and B4 φ = (M − αI 1 )B2 . By choosing where B α v = I MuX , we have hκI α M∇uX , ∇vi = hκ∇v, ∇vi ≥ k∇vk2 so, after canceling the term 21 k∇vk2 and integrating in time, Qα1 (MuX , t) + 12 Qα2 (M∇uX , t) ≤ 12 Q0 (B3 uX , t) + 12 Q0 (B4 uX , t) + 21 Qα2 (MuX , t) Z t Z t 1 α I 1 uX , I α MuX ds. + (M − αI )fX , I MuX ds + α 0 0 Using (2.15), we find that Z t (M − αI 1 )fX , I α MuX ds ≤ Ctα Q0 (M − αI 1 )fX , t) + 41 Qα1 (MuX , t) 0 and Z t 0 I 1 uX , I α MuX ds ≤ Ctα Q0 (I 1 uX , t) + 41 Qα1 (MuX , t), so Qα1 (MuX , t) + Qα2 (M∇uX , t) ≤ Q0 (B3 uX , t) + Q0 (B4 uX , t) + 2Qα2 (MuX , t) + Ctα Q0 (M − αI 1 )fX , t) + Ctα Q0 (I 1 uX , t). Since ~ 3 = (M − αI 1 )B α + (M − αI 1 )B 1 B ~ ~ F G and B4 = (M − αI 1 )Baα + (M − αI 1 )Bb1 , WELL-POSEDNESS . . . 15 the estimate (2.22) gives ~ 3 uX , t) + Q0 (B4 uX , t) ≤ Ct2 Qα2 (uX , t) + Ct2 Q12 (uX , t) Q0 (B ≤ Ct2 Qα2 (uX , t), where, in the last step, we used Lemma 2.4 with µ = α and ν = 1. We easily verify that Q0 (M − αI 1 )fX , t) ≤ Ct2 Q0 (fX , t), and by Lemma 2.4 with µ = 0 and ν = 1, Q0 (I 1 uX , t) = Q12 (uX , t) ≤ t2 Q0 (uX , t). Thus, the function q(t) = Qα1 (MuX , t) + Qα2 (M∇uX , t) satisfies q(t) ≤ Ct2 Qα2 (uX , t) + 2Qα2 (MuX , t) + Ct2+α Q0 (fX , t) + Ct2+α Q0 (uX , t). By (2.13) and Lemma 3.1, t2 Qα2 (uX , t) + t2+α Q0 (uX , t) ≤ Ct2+α Q0 (uX , t) ≤ Ct2+α Q(fX , t), and therefore, using Lemma 2.2 with φ = MuX , Z t ωα (t − s)q(s) ds, q(t) ≤ Ct2+α Q0 (fX , t) + C 0 The result now follows by applying Lemma 2.5. ✷ Lemma 3.3. The solution uX of (3.24) satisfies, for 0 ≤ t ≤ T ,


Qα1 (MuX )′ , t + Qα2 (M∇uX )′ , t ≤ Ctα Q0 (fX , t) + Ctα Q0 (MfX )′ , t . P r o o f. By differentiating (3.33) with respect to t, we have ~ 5 uX − ακI α ∇uX , ∇v (MuX )′ , v + κ∇(I α MuX )′ , ∇v = B + (MfX )′ − B6 uX , v , (3.34) ~ 5 φ = (MB ~ 1 φ)′ and B6 φ = (MB2 φ)′ . Hence, where B ~ 5 uX k2 + 1 kB6 uX k2 (MuX )′ , v + κ∇(I α MuX )′ , ∇v ≤ 12 k∇vk2 + kB 2 + 21 kvk2 + CkI α∇uX k2 + (MfX )′ , v . Putting v = I α (MuX )′ , we can cancel 12 k∇vk2 because v = (I α MuX )′ by (2.16). Thus, by integrating in time and using (2.15) to show Z t

(MfX )′ , I α (MuX )′ ds ≤ CtαQ0 (MfX )′ , t + 21 Qα1 (MuX )′ , t , 0 16 W. McLean, K. Mustapha, R. Ali, O. Knio and using (3.27), we arrive at the estimate

~ 5 uX , t) + Q0 (B6 uX , t) Qα1 (MuX )′ , t + Qα2 (M∇uX )′ , t ≤ 2Q0 (B

+ Qα2 (MuX )′ , t + CQα2 (∇uX , t) + Ctα Q0 (MfX )′ , t . Since ~ 5 uX = (MB α uX )′ + (MB 1 uX )′ B ~ ~ G F and B6 uX = (MBaα uX )′ + (MBb1 uX )′ , it follows from (2.23) that
~ 5 uX , t) + Q0 (B6 uX , t) ≤ CQα (MuX )′ , t + CQα (MuX , t) Q0 (B 2 2 + CQα2 (uX , t). By Lemmas 2.4, 3.1 and 3.2, Qα2 (MuX , t) + Qα2 (uX , t) ≤ Ctα Qα1 (MuX , t) + Ctα Qα1 (uX , t) ≤ C(t2+2α + t2α )Q0 (fX , t) and Qα2 (∇uX , t) ≤ CtαQ0 (fX , t). Hence, the function

q(t) = Qα1 (MuX )′ , t + Qα2 (M∇uX )′ , t satisfies

q(t) ≤ Ctα Q0 (fX , t) + Ctα Q0 (MfX )′ , t + CQα2 (MuX )′ , t . Finally, by Lemma 2.2, Z t Z t

α α ′ ′ Q2 (MuX ) , t ≤ C ωα (t−s)Q1 (MuX ) , s ds ≤ C ωα (t−s)q(s) ds, 0 0 and the desired estimate follows by Lemma 2.5. ✷ Lemma 3.4. The solution uX of (3.24) satisfies, for 0 ≤ t ≤ T ,


Q0 (MuX )′ , t + Qα1 (M∇uX )′ , t ≤ CQ0 (fX , t) + CQ0 (MfX )′ , t ~ 5 uX , ∇v = ∇ · B ~ 5 uX (t), v in (3.34), we obtain P r o o f. Using − B ~ 5 uX k2 + 2kB6 uX k2 (MuX )′ , v + κI α (M∇uX )′ , ∇v ≤ 12 kvk2 + 2k∇ · B + k(MfX )′ k2 − αhκI α ∇uX , ∇vi. WELL-POSEDNESS . . . 17 Choosing v = (MuX )′ , integrating in time, and using (3.27) yields


′ α ′ 0 0 1 0 ~ 2 Q (MuX ) , t + Q1 (M∇uX ) , t ≤ 2Q ∇ · B5 uX , t + 2Q (B6 uX , t) Z t
+ Q0 (MfX )′ , t − α (M∇uX )′ (s), κI α ∇uX (s) ds. Recall from (3.32) that ∇ · the notation BFα~ · ∇φ B α~ φ F = 0 α B ~φ ∇·F ~ (t) · I α ∇φ − =F Z t + B α~ ∇φ, where we have used F· F~ ′ (s) · I α∇φ(s) ds. 0 Thus, ~ 1 uX ~ 5 uX = ∇ · MB ∇·B
′ ~ 1 uX = M∇ · B
′
′
′ 1 = M∇ · BFα~ uX + M∇ · BG ~ uX
′
′ α α + MB ∇u = MB∇· u X X ~ ~ F· F
′
′ 1 1 + MB∇· u MB ∇u + , X X ~ ~ G G· and so, by (2.23),

~ 5 uX , t + Q0 (B6 uX , t) ≤ CQα (MuX )′ , t + CQα (MuX , t) Q0 ∇ · B 2 2
+ CQα2 (uX , t) + CQα2 (M∇uX )′ , t + CQα2 (M∇uX , t) + CQα2 (∇uX , t). By (2.12), Z t
(M∇uX )′ (s), κI α ∇uX (s) ds ≤ 12 Qα1 (M∇uX )′ , t + CQα1 (∇uX , t), 0

and thus the function q(t) = Q0 (MuX )′ , t + Qα1 (M∇uX )′ , t satisfies
q(t) ≤ CQα2 (MuX )′ , t + CQα2 (MuX , t) + CQα2 (uX , t)
+ CQα2 (M∇uX )′ , t + CQα2 (M∇uX , t) + CQα2 (∇uX , t)
+ CQ0 (MfX )′ , t + CQα1 (∇uX , t)
≤ CQα2 (MuX )′ , t + CtαQα1 (MuX , t) + Ctα Qα1 (uX , t)
+ CQα2 (M∇uX )′ , t + Ct2+α Q0 (fX , t) + CtαQ0 (fX , t)
+ CQ0 (MfX )′ , t + CQ0 (fX , t), where, in the second step, we used Lemmas 2.2, 3.1 and 3.2. A further application of Lemmas 3.1 and 3.2 yields

q(t) ≤ CQ0 (MfX )′ , t + CQ0 (fX , t) + CQα2 (MuX )′ , t
+ CQα2 (M∇uX )′ , t . 18 W. McLean, K. Mustapha, R. Ali, O. Knio

Lemma 2.2 implies that Qα2 (MuX )′ , t + Qα2 (M∇uX )′ , t is bounded by Z t 

 ωα (t − s) Qα1 (MuX )′ , s + Qα1 (M∇uX )′ , s ds C 0 Z t ωα (t − s)q(s) ds, ≤C Qα1 )′ , s
Ctα Q0 0
′ ) ,s , (MuX ≤ (MuX which follows by where we used Lemma 2.4. Finally, Lemma 2.5 implies the desired estimate. ✷ The preceding lemmas yield the main result for this section. Theorem 3.3. Assume that the coefficients satisfy (1.4), that the initial data u0 ∈ L2 (Ω) and that the source term satisfies (3.25). Then, the solution uX of the projected Volterra equation (3.24) satisfies (with C independent of X)
for 0 ≤ t ≤ T . kuX (t)k2 + tα k∇uX (t)k2 ≤ C ku0 k2 + M 2 t2η P r o o f. Applying Lemma 2.3 with φ = MuX , we see that Lemma 3.3 gives
t2 kuX (t)k2 = kMuX (t)k2 ≤ Ct1−α Qα1 (MuX )′ , t
≤ CtQ0 (fX , t) + CtQ0 (MfX )′ , t . Define gX : [0, T ] → X by hgX (t), vi = hg(t), vi for v ∈ X, and observe that ′ = f + Mg . We find using fX = u0 + I 1 gX and (MfX )′ = fX + MfX X X (3.25) that  Z t
0 0 ′ 2 1 2 2 Q (fX , t) + Q (MfX ) , t ≤ C ku0 k + kI gk + kMgk ds 0
≤ Ct ku0 k2 + M 2 t2η , (3.35) so the estimate for the first term kuX (t)k2 follows at once. Similarly, applying Lemma 2.3 with φ = (M∇uX )′ followed by Lemma 3.4, we have
t2+α k∇uX (t)k = tα kM∇uX (t)k2 ≤ CtQα1 (M∇uX )′ , t
≤ CtQ0 (fX , t) + CtQ0 (MfX )′ , t , implying the estimate for the second term tα k∇uX (t)k2 . ✷ 4. The weak solution We will now establish that the weak formulation (1.5) of the initialboundary value problem (1.1)–(1.3) is well-posed. The proof relies on our WELL-POSEDNESS . . . 19 estimates from Section 3 and also the following local Hölder continuity properties of uX . Lemma 4.1. If 0 < δ ≤ t1 < t2 ≤ T , then
kuX (t2 ) − uX (t1 )k2 ≤ Cδ−2 t2 ku0 k2 + M 2 t2η 2 (t2 − t1 ) and
  kI α ∇uX (t2 )−I α∇uX (t1 )k ≤ C ku0 k+M tη2 δα−2 (t2 −t1 )+δ−α/2 (t2 −t1 )α . P r o o f. The Cauchy–Schwarz inequality implies that Z t2 Z t2 2 2 ′ kuX (t2 ) − uX (t1 )k = uX (s) ds ≤ (t2 − t1 ) ku′X (s)k2 ds, t1 t1 and by the second inequality of Lemma 3.1, together with Lemma 3.4, Z t2 Z t2 s−2 k(MuX )′ (s) − uX (s)k2 ds ku′X (s)k2 ds = t1 t1 Z t2
k(MuX )′ k2 + kuX k2 ds
 0 = 2δ−2 Q0 (MuX )′ , t2 + Q0 (uX , t2 )
  ≤ Cδ−2 Q0 MfX )′ , t2 + Q0 (fX , t2 ) . ≤ 2δ −2 The first result now follows from (3.35). To prove the second, we write Z t1 −δ/2   ωα (t2 − s) − ωα (t1 − s) ∇uX (s) ds I α ∇uX (t2 ) − I α∇uX (t1 ) = 0 Z t2 Z t1   ωα (t2 − s)∇uX (s) ds, ωα (t2 − s) − ωα (t1 − s) ∇uX (s) ds + + t1 t1 −δ/2 and deduce from Theorem 3.3 that where
kI α ∇uX (t2 ) − I α ∇uX (t1 )k ≤ C ku0 k + M tη2 I1 + I2 + I3 ), I1 = I2 = I3 = Z Z Z 0 t1 −δ/2  t1 t1 −δ/2 t2   ωα (t1 − s) − ωα (t2 − s) s−α/2 ds,  ωα (t1 − s) − ωα (t2 − s) s−α/2 ds, ωα (t2 − s)s−α/2 ds. t1 20 W. McLean, K. Mustapha, R. Ali, O. Knio By the mean value theorem, ωα (t1 − s) − ωα (t2 − s) = (t2 − t1 )|ωα−1 (ξ)| with t1 − s < ξ < t2 − s, and if 0 < s < t1 − δ/2 then t1 − s > δ/2 so Z t1 −δ/2 ds I1 ≤ (t2 − t1 )|ωα−1 (δ/2)| sα/2 0  2−α 1 − α (t1 − δ/2)1−α/2 2 (t2 − t1 ). ≤ δ 1 − α/2 Γ(α) Moreover, −α/2 I2 ≤ (δ/2) α/2 = (2/δ)  Z t1 t1 −δ/2   ωα (t1 − s) − ωα (t2 − s) ds  ωα+1 (t2 − t1 ) + ωα+1 (δ/2) − ωα+1 (t2 − t1 + δ/2) ≤ (2/δ)α/2 ωα+1 (t2 − t1 ) Rt and I3 ≤ δ−α/2 t12 ωα (t2 − s) ds = δ−α/2 ωα+1 (t2 − t1 ). ✷ Our existence theorem is stated as follows. Note the weak continuity at t = 0 asserted in part 5; we show in the companion paper [26] that the solution u is continuous on the closed interval [0, T ] provided u0 ∈ Ḣ µ (Ω) for some µ > 0. Theorem 4.1. Assume that the coefficients satisfy (1.4), that the source term satisfies (3.25), and that the initial data u0 ∈ L2 (Ω). Then, the initial-boundary value problem (1.1)–(1.3) has a weak solution u. More precisely, there exists a function u : [0, T ] → L2 (Ω) with the following properties. (1) The restriction u : (0, T ] → L2 (Ω) is continuous. (2) If 0 < t ≤ T , then u(t) ∈ H01 (Ω) with
ku(t)k + tα/2 k∇u(t)k ≤ C ku0 k + M tη . (3) The functions I α u and B2 u are continuous from the closed interval [0, T ] ~ 1 u are continuous from [0, T ] to L2 (Ω)d . to L2 (Ω). Likewise, I α ∇u and B α ~ 1 u = 0 and u(0) = u0 . (4) At t = 0 we have I u = B2 u = 0, I α ∇u = B (5) If t → 0, then hu(t), vi → hu(0), vi for each v ∈ L2 (Ω). (6) For 0 ≤ t ≤ T and v ∈ H01 (Ω), ~ 1 u)(t), ∇v + h(B2 u)(t), vi = hf (t), vi. hu(t), vi + κ(I α ∇u)(t), ∇v − (B (4.36) WELL-POSEDNESS . . . 21 P r o o f. Let ψ1 , ψ2 , ψ3 , . . . be a sequence of functions spanning a dense subspace of H01 (Ω). For each integer n ≥ 1, let Xn = span{ψ1 , ψ2 , . . . , ψn } and for brevity denote the solution of (3.30) with X = Xn by un = uX , and likewise write fn = fX , so that hun (t), vi+hκ(I α ∇un )(t), ∇vi−h(B1 un )(t), ∇vi+h(B2 un )(t), vi = hfn (t), vi (4.37) for v ∈ Xn and 0 < t ≤ T . We see from Theorem 3.3 and Lemma 4.1 that, whenever 0 < δ < T , the sequence of functions un is bounded and
equicontinuous in C [δ, T ]; L2 (Ω) . By choosing a sequence of values of δ tending to zero we can select a subsequence, again denoted by un , such that un (t) converges in L2 (Ω) for 0 < t ≤ T . We may therefore define u(t) = lim un (t) for 0 < t ≤ T , n→∞ and this function satisfies property 1 because, given any fixed δ ∈ (0, T ), the limit is uniform for t ∈ [δ, T ]. Similarly,
the functions I α∇un are bounded and equicontinuous in C [δ, T ]; L2 (Ω)d so I α ∇u : (0, T ] → L2 (Ω)d is continuous. In fact, it will follow from (4.39) below that kI α ∇u(t)k → 0 as t → 0, so I α ∇u : [0, T ] → L2 (Ω)d is continuous. By Theorem 3.3,
kun (t)k ≤ C ku0 k + M tη for 0 < t ≤ T ,
so by sending n → ∞ we conclude that ku(t)k ≤ C ku0 k + M tη . Also, for 0 < t ≤ T,
|hun (t), vi| ≤ Ckun (t)kH01 (Ω) kvkH −1 (Ω) ≤ Ct−α/2 ku0 k + M tη kvkH −1 (Ω) and sending n → ∞ it follows that
|hu(t), vi| ≤ Ct−α/2 ku0 k + M tη kvkH −1 (Ω) for all v ∈ L2 (Ω),
≤ Ct−α/2 (ku0 k + M tη , establishing so u(t) ∈ H01 (Ω) with ku(t)kH01 (Ω) property 2. Since ku(t)k is bounded, I α u is continuous on [0, T ] with Z t ωα (t − s)ku(s)k ds kI α u(t)k ≤ 0 (4.38) Z t

α−1 η η α (t − s) ku0 k + M s ds ≤ C ku0 k + M t t , ≤C 0 and similarly α kI ∇u(t)k ≤ C Z t 0

(t − s)α−1 s−α/2 ku0 k + M sη ds ≤ C(ku0 k + M tη tα/2 . (4.39) 22 W. McLean, K. Mustapha, R. Ali, O. Knio Likewise, for n ≥ 1,
kI α un (t)k ≤ C ku0 k + M tη tα
kI α ∇un (t)k ≤ C ku0 k + M tη tα/2 . (4.40) ~ 1 u and B2 u follow from (2.20) and (3.29), completing the Continuity of B proof of property 3, with Z t
~ 1 u)(t)k + k(B2 u)(t)k ≤ Ck(I α u)(t)k + C k(I α u)(s)k + ku(s)k ds k(B 0
≤ C ku0 k + M tα . (4.41) Property 4 follows from the estimates (4.38), (4.39) and (4.41). If 0 ≤ δ < t ≤ T , then Z t α α ωα (t − s)kun (s) − u(s)k ds k(I un )(t) − (I u)(t)k ≤ ≤C Z ≤ Cδ 0 α and 0 δ α−1 η Z
t ku0 k + M s ds + (t − s)α−1 kun (s) − u(s)k ds δ
ku0 k + M δη + α−1 (t − δ)α max kun (s) − u(s)k, (t − s) δ≤s≤t showing that I αun (t) → I αu(t) in L2 (Ω), uniformly for t ∈ [δ, T ]. In fact, the convergence is uniform for t ∈ [0, T ], owing to the estimates (4.38) and (4.40). Therefore, we see using (2.20) and (3.29) that, for v ∈ H01 (Ω), ~ 1 un )(t), ∇v → (B ~ 1 u)(t), ∇v (B and (B2 un )(t), v → (B2 u)(t), v . Since hfn , ψj i = hf, ψj i for j ≤ n, we have lim hfn (t), ψj i = hf (t), ψj i n→∞ for all j ≥ 1 and 0 ≤ t ≤ T , and therefore hfn (t), vi → hf (t), vi for all v ∈ L2 (Ω). Thus, by sending n → ∞ in (4.37), it follows that (4.36) holds for v ∈ H01 (Ω) and 0 < t ≤ T . In light of (4.41) and (4.39), the variational equation (4.36) is satisfied when t = 0 if and only if hu(0), vi = hu0 , vi for all v ∈ H01 (Ω), which is the case if and only if we define u(0) = u0 . Moreover, if t → 0 then hu(t), vi → hf (0), vi = hu0 , vi, for each v ∈ H01 (Ω), and hence by density for each v ∈ L2 (Ω), establishing properties 5 and 6. ✷ Remark 4.1. Since our estimates rely on Lemma 2.1, the constant C in part 2 of Theorem 4.1 becomes unbounded as α → 1. However, this behavior appears to be an artifact of our method of proof. In the limiting case when α = 1 and (1.1) reduces to a parabolic PDE, a simple energy argument combined with the classical Gronwall inequality yields the a priori WELL-POSEDNESS . . . 23 estimate   Z t kg(s)k ds ku(t)k ≤ C ku0 k + for 0 ≤ t ≤ T ; 0 see also the alternative analysis [19] of the fractional Fokker–Planck equation. Theorem 4.2. The weak solution of the initial-boundary value problem (1.1)–(1.3) is unique. More precisely, under the same assumptions as Theorem 4.1, there is at most one function u that
satisfies (4.36) and is such that u and I αu belong to L2 (0, T ); L2 (Ω) , and I α ∇u belongs
to L2 (0, T ); L2 (Ω)d . P r o o f. The problem is linear, so it suffices to show that if u0 = 0 and g(t) ≡ 0 then u(t) ≡ 0. Thus, suppose that ~ 1 u)(t), ∇v + h(B2 u)(t), vi = 0 hu(t), vi + κ(I α ∇u)(t), ∇v − (B for 0 < t ≤ T and v ∈ H01 (Ω). Proceeding as in the proof of (3.31), we have ~ 1 u, t) + 1 Q0 (B2 u, t) + 1 Qα (u, t) Qα1 (u, t) + 12 Qα2 (∇u, t) ≤ 12 Q0 (B 2 2 2 ≤ CQα2 (u, t), where the final step used (2.20), (2.21) and Lemma 2.4. Thus, applying Lemma 2.2, the function q(t) = Qα1 (u, t) + Qα2 (∇u, t) satisfies Z t ωα (t − s)q(s) ds, q(t) ≤ CQα2 (∇u, t) ≤ C 0 and therefore q(t) = 0 for 0 ≤ t ≤ T by Lemma 2.5. In particular, Qα1 (u, T ) = 0, so if we put u(t) = 0 for t > T then the Laplace transform of u satisfies û(iy) = 0 for −∞ < y < ∞ by (2.10), implying that u(t) = 0 for 0 ≤ t ≤ T . ✷ Acknowledgements The authors thank the University of New South Wales (Faculty Research Grant “Efficient numerical simulation of anomalous transport phenomena”), the King Fahd University of Petroleum and Minerals (project No. KAUST005) and the King Abdullah University of Science and Technology. 24 W. McLean, K. Mustapha, R. Ali, O. Knio References [1] A. A. Alikhanov. A priori estimates for solutions of boundary value problems for fractional-order equations. Differential Equations, 46:660–666, 2010. [2] Anatoly A. Alikhanov. Boundary value problems for the diffusion equation of the variable order in differential and difference settings. Applied Mathematics and Computation, 219:3938–3946, 2012. [3] C. N. Angstmann, B. I. Henry, B. A. Jacobs, and A. V. McGann. A time-fractional generalised advection equation from a stochastic process. Chaos, Solitons and Fractals, 102:175–183, 2017. Future Directions in Fractional Calculus Research and Applications. [4] Leigh C. Becker. Resolvents and solutions of weakly singular linear Volterra integral equations. Nonlinear Analysis, 74:1892–1912, 2011. [5] Hermann Brunner. Volterra Integral Equations: an Introduction to Theory and Applications. Cambridge University Press, 2017. [6] Hermann Brunner, Arvet Pedas, and Gennadi Vainikko. The piecewise polynomial collocation method for nonlinear weakly singular Volterra equations. Math. Comp., 68:1079–1095, 1999. [7] Eduardo Cuesta, Christian Lubich, and Cesar Palencia. Convolution quadrature time discretization of fractional diffusive-wave equations. Math. Comp., 75:673–696, 2006. [8] J. Dixon and S. McKee. Weakly singular Gronwall inequalities. ZAMM Z. Angew. Math. Mech., 66:535–544, 1986. [9] B. I. Henry, T. A. M. Langlands, and P. Straka. Fractional Fokker– Planck equations for subdiffusion with space- and time-dependent forces. Phys. Rev. Lett., 105:170602, 2010. [10] B. I. Henry, T. A. M. Langlands, and S. L. Wearne. Anomalous diffusion with linear reaction dynamics: From continuous time random walks to fractional reaction-diffusion equations. Phys. Rev. E, 74:031116, Sep 2006. [11] B. I. Henry and S. L. Wearne. Fractional reaction-diffusion. Physica A: Statistical Mechanics and its Applications, 276:448–455, 2000. [12] Bangti Jin, Buyang Li, and Zhi Zhou. Discrete maximal regularity of time-stepping schemes for fractional evolution equations. Numer. Math., 138:101–131, 2018. [13] Samir Karaa and Amiya K. Pani. Error analysis of a FVEM for fractional order evolution equations with nonsmooth initial data. ESAIM: M2AN, 52(2):773–801, 2018. [14] J. Klafter and I. M. Sokolov. First Steps in Random Walks. Oxford University Press, 2011. WELL-POSEDNESS . . . 25 [15] Adam Kubica and Masahiro Yamamoto. Initial-boundary value problems for fractional diffusion equations with time-dependent coefficients. Fract. Calc. Appl. Anal., 21:276–311, 2018. [16] T. A. M. Langlands and B. I. Henry. The accuracy and stability of an implicit solution method for the fractional diffusion equation. Journal of Computational Physics, 205(2):719–736, 2005. [17] T. A. M. Langlands, B. I. Henry, and S. L. Wearne. Fractional cable equation models for anomalous electrodiffusion in nerve cells. SIAM J. Appl. Math., 71:1168–1203, 2011. [18] Kim Ngan Le, William McLean, and Kassem Mustapha. A semidiscrete finite element approximation of a time-fractional Fokker–Planck equation with non-smooth initial data. SIAM J. Sci. Computing. To appear. [19] Kim Ngan Le, William McLean, and Martin Stynes. Existence, uniqueness and regularity of the solution of the time-fractional Fokker–Planck equation with general forcing. Comm. Pure Appl. Anal., 2018. to appear. [20] Hong-Lin Liao, Dongfang Li, and Jiwei Zhang. Sharp error estimate of the nonuniform l1 formula for linear reaction-subdiffusion equations. SIAM J. Numer. Anal., 66(2):1112–1133, 2018. [21] Yumin Lin and Chuanju Xu. Finite difference/spectral approximations for the time-fractional diffusion equation. Journal of Computational Physics, 225(2):1533–1552, 007. [22] Peter Linz. Analytical and Numerical Methods for Volterra Equations. Studies in Applied and Numerical Mathematics. SIAM, Philadelphia, 1985. [23] F. Liu, V. V. Anh, I. Turner, and P. Zhuang. Time fractional advectiondispersion equation. J. Appl. Math. Computing, 13:233–245, 2003. [24] Ch. Lubich. Runge–Kutta theory for Volterra and Abel integral equations of the second kind. Math. Comp., 41:87–102, 1983. [25] William McLean. Regularity of solutions to a time-fractional diffusion equation. ANZIAM J., 52:123–138, 2010. [26] William McLean, Kassem Mustapha, Raed Ali, and Omar Knio. Regularity theory for time-fractional advection-diffusion-reaction equations. In preparation. [27] R. Metzler, E. Barkai, and J. Klafter. Deriving fractional Fokker– Planck equations from a generalised master equation. Europhys. Lett., 46:431–436, 1999. [28] Kassem Mustapha. Time-stepping discontinuous Galerkin methods for fractional diffusion problems. Numer. Math., 130(3):497–516, 2015. 26 W. McLean, K. Mustapha, R. Ali, O. Knio [29] J. A. Nohel and D. F. Shea. Frequency domain methods for volterra equations. Advances in Mathematics, 22:278–304, 1976. [30] Kenichi Sakamoto and Masahiro Yamamoto. Initial value/boundary value problems for fractional diffusion-wave equations and applications to some inverse problems. Journal of Mathematical Analysis and Applications, 382(1):426–447, 2011. [31] Martin Stynes, Eugene O’Riordan, and José Luis Gracia. Error analysis of a finite difference method on graded meshes for a time-fractional diffusion equation. SIAM J. Numer. Anal., 55(2):1057–1079, 2017. [32] Gennadi Vainikko. Weakly singular integral equations. Lecture Notes, University of Tartu, Helsinki University of Technology, 2006–2007. [33] S. B. Yuste and L. Acedo. An explicit finite difference method and a new von Neumann stability analysis for fractional diffusion equations. SIAM. J. Numer. Anal., 42(5):1862–1874, 2005. [34] Rico Zacher. Weak solutions of abstract evolutionary integrodifferential equations in Hilbert spaces. Funkcialaj Ekvacioj, 52:1–18, 2009. 1 School of Mathematics and Statistics The University of New South Wales Sydney 2052, Australia e-mail: w.mclean@unsw.edu.eu 2 Department of Mathematics and Statistics KFUPM, Dhahran 31261, Saudi Arabia e-mail: kassem@kfupm.edu.su, g201305090@kfupm.edu.sa 3 Computer, Electrical, Mathematical Sciences and Engineering Division KAUST, Thuwal 23955, Saudi Arabia e-mail: Omar.Knio@kaust.edu.sa