Galerkin Approximation for Stochastic Volterra Integral Equations with Doubly Singular Kernels

This paper is concerned with the more general nonlinear stochastic Volterra integral equations with doubly singular kernels, whose singular points include both s = t and s = 0. We propose a Galerkin approximate scheme to solve the equation numerically, and we obtain the strong convergence rate for the Galerkin method in the mean square sense. The rate is min{2− 2(α1 + β1), 1− 2(α2 + β2)} (where α1, α2, β1, β2 are positive numbers satisfying 0 < α1 + β1 < 1, 0 < α2 + β2 < 2 ), which improves the results of some numerical schemes for the stochastic Volterra integral equations with regular or weakly singular kernels. Moreover, numerical examples are given to support the theoretical result and explain the priority of the Galerkin method.


Introduction
Stochastic Volterra integral equations (SVIEs) were firstly introduced by Berger and Mizel [1,2] and can be regarded as generalizations of stochastic differential equations or deterministic Volterra integral equations. Due to the advantages of the SVIEs in describing memory, heritability and ubiquitous noise perturbations, these kinds of equations have extensive applications in many fields (e.g., finance [3], control science [4] and mathematical biology [5]) and have spurred great research enthusiasm. An area of particular interest in the study of SVIEs has been numerical analysis because the analytic solution to SVIEs is rarely known and the numerical approximations provide a powerful tool for understanding the behavior of the solution. Up to now, most numerical methods are developed to deal with SVIEs with regular kernels. For example, Tudor [6] first proposed a one-step numerical approximation for Itô-Volterra equations and obtained a basic convergence theorem. Wen and Zhang [7,8] studied the rectangular method. Methods based on the operational matrix were introduced in [9][10][11]. Xiao et al. [12] introduced a split-step collocation method for SVIEs. The Euler-Maruyama (EM) methods were discussed in [13][14][15]. Recently, Conte et al. [16] introduced improved stochastic θ-methods, having better stability properties, for the numerical integration of stochastic Volterra integral equations.
However, stochastic Volterra integral equations with singular kernels arise while dealing with some problems in stochastic partial differential equations [17][18][19], such as certain heat conduction problems with mixed boundary conditions, and in the analysis of the fractional Brownian motion [20]. It is therefore of great value and significance to consider SVIEs involving the kernels with singularities. Wang [21] established the existence and uniqueness theorem for stochastic Volterra equations with singular kernels. Zhang [22] investigated the convergence of the EM method for this kind of equation. In particular, Li et al. [23] recently proposed a θ-Euler-Maruyama scheme and a Milstein scheme to solve SVIEs with weakly singular kernels numerically. It is observed that the singular SVIEs have been relatively less studied. It is probably because the singularity of the integrand kernel brings us more difficulties. In this case, the key Itô formula, which is a powerful and necessary tool in the study of stochastic differential equations (SDEs), is not available and we have to seek other tools and techniques.
As we can see, SVIE (1) involves both diagonal and boundary singular kernels, that is, the kernel functions are singular when s = t or s = 0, which are more general than the case considered in [23][24][25], whose singular point only includes s = t. The existence and uniqueness theorem of the true solution as well as the strong convergence of the Euler-Maruyama method for Equation (1) have been developed in [26]. Li et al. [27] obtained the asymptotic separation for Equation (1).
Our aim here is to establish the Galerkin approximations for Equation (1). This method is high-order accurate, easy-to-handle complicated geometries and side conditions, highly parallelizable, nonlinear stable and has the ability to capture discontinuities without spurious oscillations [28]. Hence, it has been extensively used for stochastic partial differential equations [29][30][31] and stochastic fractional differential equations [32,33], while it has not been analyzed to singular SVIEs.
The rest of this paper is organized as follows.
In Section 2, we present some preliminaries. Section 3 provides the main result and Section 4 covers numerical experiments to illustrate the efficiency of the numerical method. We conclude this work in Section 5.
(1) If u satisfies (2) If u satisfies Throughout this paper, we impose the following conditions on the coefficients f and g.

Assumption 1.
Assume that there exists a constant L > 0 such that for any x, y ∈ R, we have

A Numerical Approximation of White Noise for Stochastic Integral
Next, we introduce the approximation of the noise. Following the method discussed in [29], let where We substitute dB(t) with dB(t) in Equation (1) to obtain the following equation Letx(t) be the solution of the approximate SVIEs

The Discontinuous Galerkin Method
Because (6) is a system of ODEs, we apply the standard discontinuous Galerkin method for deterministic ODEs [28]. Multiplying (6) by v, integrating over I n and using integration by parts, we obtain the discontinuous Galerkin formulation We define the piecewise polynomial space V p h = v : v| I n ∈ P p (I n ), n = 0, 1, . . . , N − 1 as the space of polynomials of degree at most p in each interval I n , where P p (I n ) is the set of all polynomials of degree less or equal than p on I n . Because polynomials in V p h are allowed to have discontinuities across element boundaries, we use v(t ± n ) = lim s→0 ± v(t n + s) to denote the left limit and the right limit of v at t n . Next, we approximate eachx i (t) by a piecewise polynomial where we used the classical upwind numerical flux. We will refer to this discontinuous Galerkin scheme as the stochastic discontinuous Galerkin scheme.

Main Result
Theorem 1 ([26]). Under Assumption 1, there exists a unique strong solution x(t) to Equation (1). Moreover, there exists a constant C > 0 such that The main result of this paper is the following theorem.
Theorem 2. Under Assumption 1, there exists a constant C > 0 such that where, for simplicity of notation, C is a positive constant which depends on α 1 , α 2 , β 1 , β 2 , T, L and may change from line to line during the rest of the paper, while its specific form is of unimportance.
In order to prove Theorem 2, we need the following Lemma.
Lemma 3. Let x(t) be the solution of (1), then for any 0 ≤ t 1 < t 2 ≤ T, we have Proof. Observe that It follows that On the one hand, Cauchy-Schwarz inequality and basic calculus imply that On the other hand, using Itô isometry, we obtain The proof is completed.
Proof of Theorem 2. Let e(t) = x(t) −x(t), then we have For the term J 1 , it follows from Assumption 1 that Next, let us estimate the term J 2 . Observe that In order to estimate J 3 , we assume there exists an integer N 0 ≤ N such that t N 0 ≤ t < t N 0 +1 , then According to (4), we obtain Based on the fact On the one hand, it is easy to see that On the other hand, Consequently, Similarly, we can obtain Combining (9) with (10), we obtain Next, we estimate the term J 4 . For simplicity, let t = t N 0 +1 , then Combining this with the above estimates, we can derive that A simple application of the generalized Gronwall's inequality in Lemma 2 yields Finally, Lemma 1 helps us derive the required assertion. We complete the proof.

Numerical Simulation
In this section, we will verify the numerical solution of the stochastic Volterra integral equations with doubly singular kernels and the strong convergence rate of the Galerkin method. Example 1. Let us consider the following stochastic Volterra integral equations with doubly singular kernels When the terminal time t N = T = 1, the positive arguments α i and β i (i = 1, 2) take the following four cases: Obviously, functions f and g satisfy the Assumption 1. By Theorem 2, the convergence rate of the Galerkin scheme is min{2 − 2(α 1 + β 1 ), 1 − 2(α 2 + β 2 )}, we obtain the convergence rate of Table 1. In order to verify the numerical solution of Equation (11) and the strong convergence rate of the Galerkin method, we adopt five different step sizes h = 2 −4 , 2 −5 , 2 −6 , 2 −7 , 2 −8 on the same Brownian motion to compute the four cases, respectively. The numerical solution by the Galerkin method with small step size h * = 2 −14 is used to represent the unknown true solution. Figure 1 draws the calculated error results, which are consistent with Theorems 1 and 2. Compared with the classic Euler-Maruyama method or Milstein method, the Galerkin method is more accurate in algorithm estimation by space discretization, for example, when β 1 = β 2 = 0, Equation (1) became SVIE with weakly singular kernels, revealing the convergence rate of the Euler-Maruyama method is min{1 − α 1 , 1 2 − α 2 } in [23].
Due to the appearance of the singularity in the above stochastic integral, it is difficult for us to illustrate the convergence rate of the Galerkin scheme numerically. We set α 1 = 0.3, α 1 = 0.5, β 1 = 0.2, β 2 = 0.3. We regard the numerical solution yielded by small step size h * = 2 −13 as the 'exact' solution. Moreover, the corresponding numerical solutions are generated by four different step sizes h = 2 −4 , 2 −5 , 2 −6 and 2 −7 , respectively. The mean square errors of the Galerkin scheme are calculated at the terminal time t N = T = 1 by where the expectation is approximated by averaging over M = 1000 Brownian sample paths. The mean square errors are plotted in Figure 2 in a loglog scale. In these plots, the reference lines and error lines are parallel to each other, revealing the convergence rate of the Galerkin scheme is min{2 − 2(α 1 + β 1 ), 1 − 2(α 2 + β 2 )}. Therefore, the convergence rate is 0.4.

Conclusions
In this work, we investigated the Galerkin approximation for a class of nonlinear stochastic Volterra integral equations with doubly singular kernels. These kinds of equations are more general, for example, when β 1 = β 2 = 0, Equation (1) becomes the SVIE with weakly singular kernels studied in [23]. Because the Itô formula is not avail-able, the classical proof techniques can no longer be used. With the help of new tools such as classical fractional calculus, the Mittag-Leffler-type function and generalized Gronwall inequalities with singular kernels, we obtained the strong convergence rate min{2 − 2(α 1 + β 1 ), 1 − 2(α 2 + β 2 )} for the Galerkin method, which improves the corresponding result of the θ-Euler-Maruyama scheme in [23]. In forthcoming works, we will consider other numerical methods for SVIEs with doubly singular kernels and verify whether the order of convergence is optimal. Moreover, we will try to analyze the stability of the numerical methods.