Triangular Function Method is Adopted to Solve Nonlinear Stochastic It o^–Volterra Integral Equations

Abstract

This article presents the numerical solutions of nonlinear stochastic It $\widehat{o}$–Volterra integral equations by using the basis function method under the global Lipschitz condition. Integral operator matrixes of triangular functions are used to convert the nonlinear stochastic integral equations into a system of algebraic equations. Meanwhile, we gain the error of the current method, and it is demonstrated that the error accuracy of this method is higher than that of the BPFs. In the end, the feasibility, accuracy, and validity of the current method are demonstrated by numerical results.

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Stochastic Volterra integral equations are widely used in medicine, physics, biology, finance, and other fields. Babei et al. [1] examined a stochastic model for the transmission of coronavirus. Khodabin et al. [2] presented the generalized stochastic exponential population growth model. Ait-Sahalia [3] proposed the derivative securities pricing model. Since these systems rely on Gaussian white noise, various stochastic integral equations are demanded to simulate these phenomena. As everyone knows, it is difficult to obtain explicit solutions of the stochastic Volterra integral equation even if analytic solutions exist; this greatly limits the application of the equation in practical cases. Thus, it is very crucial to discuss the numerical solutions of these equations.

In the past two decades, some researchers used the usual iteration methods such as Euler and Milstain [4–7]. Others focused on basis function methods such as triangular functions, block pulse functions, Fourier series, Chebyshev cardinal wavelet, Legendre multiwavelets, Haar wavelets, Chebyshev polynomials, and Bell polynomials; see references [8–21] for details.

In the literature [13], Maleknejad et al. applied block pulse functions to explore the numerical solution of the following linear stochastic Volterra integral equation:$$\begin{array}{c}\text{(1)}& zt=z_{0}t+{{\displaystyle \int }}_0^t{a}_{1}s,tzs\mathrm{d}s+{{\displaystyle \int }}_0^t{a}_{2}s,tzs\mathrm{d}Bs,t\in 0,1,\end{array}$$where $zt$ is unknown function, $z_{0}t$ is known determinate function, and $a_{1}s,t$ and $a_{2}s,t$ are kernel functions, for $0\le s<t\le 1$. They rely on the probability space $\Omega ,\mathcal{F},\mathrm{P}$ and field ${\mathcal{F}}_t,t\ge 0$. $Bs$ is a standard Brownian motion. The second integral is in the It $\widehat{o}$ sense. Integral operator matrixes of block pulse functions are used to convert the stochastic integral equations into a system of algebraic equations. Deng et al. [14] studied the following nonlinear form of (1) applying block pulse functions:$$\begin{array}{c}\text{(2)}& zt=z_{0}t+{{\displaystyle \int }}_0^t\widehat{a}s,t\rho zs\mathrm{d}s+{{\displaystyle \int }}_0^t\tilde{a}s,t\sigma zs\mathrm{d}{B}_s^H,t\in 0,1.\end{array}$$where $\rho \cdot$ and $\sigma \cdot$ are analytic functions and satisfy some bounded and Lipschitz conditions. $B_s^H$ is the fractional Brownian motion (FBM) with Hurst parameter $H\in 1/2,1$.

Deb et al. [22] recommended triangular functions and applied them to dynamic system. Some scholars applied triangular functions to study different equations (20)–(23), for example, linear stochastic integral equation, multiorder fractional stochastic Volterra integral equation, and Volterra–Fredholm integro-differential equations.

In particular, Khodabin et al. [23] applied triangular functions to explore the numerical solution of equation (1). Compared with reference [13], the advantage of reference [23] is that it greatly reduces the computational cost and improves the accuracy of the numerical solution.

To the best of our knowledge, there are probably few studies to discuss the numerical solutions of the following nonlinear stochastic It $\widehat{o}$–Volterra integral equations:$$\begin{array}{c}\text{(3)}& zt=z_{0}t+{{\displaystyle \int }}_0^t\widehat{a}s,t\rho s,zs\mathrm{d}s+{{\displaystyle \int }}_0^t\tilde{a}s,t\sigma s,zs\mathrm{d}Bs,t\in 0,1,\end{array}$$where $zt$ is unknown function, $z_{0}t$ is known determinate function, and $\widehat{a}s,t$ and $\tilde{a}s,t$ are kernel functions, for $0\le s<t\le 1$. They rely on the probability space $\Omega ,\mathcal{F},\mathrm{P}$ and field ${\mathcal{F}}_t,t\ge 0$. $Bs$ is a standard Brownian motion. ${{\displaystyle \int }}_0^t\tilde{a}s,t\sigma s,zs\mathrm{d}Bs$ is an integral in the It $\widehat{o}$ sense. $\rho \cdot ,\cdot$ and $\sigma \cdot ,\cdot$ are analytic functions and satisfy some bounded and Lipschitz conditions. In this article, we will research nonlinear equation (3) based on triangular functions. Meanwhile, we will present Lemma 5 to deal with analytic functions $\rho \cdot ,\cdot$ and $\sigma \cdot ,\cdot$.

Most nonlinear stochastic It $\widehat{o}$–Volterra integral equations cannot be explicitly solved. This paper proposes a global approximation method to nonlinear stochastic It $\widehat{o}$–Volterra integral equations by the triangular functions. The trajectory of Brownian motion is represented as a linear combination of triangular functions. The nonlinear stochastic It $\widehat{o}$–Volterra integral equations are transformed into algebraic equations by using the operation matrix of integration, and the coefficients in this method do not need integration. In addition, compared to traditional iterative methods [4, 6, 7], the triangular function method is a global approximation method that can reduce the accumulation of errors. The method has good stability and is not easily affected by iterative parameters. Therefore, the triangular function method has higher error accuracy, lower computational complexity, and shorter computation time. This indicates that the method is effective. However, it is more difficult to establish an integration operator calculation program.

The composition of this article is as follows. We review the relevant contents of triangular functions in Section 2. We exhibit the integral operator matrixes and stochastic integral operator matrixes about triangular functions, respectively, in Sections 3 and 4. Section 5 presents the error of the current method. An effectual numerical method for (3) is proposed in Section 6. The accuracy and validity of the method are demonstrated by some numerical simulations in Section 7. Conclusions are provided in Section 8.

2. Preliminaries

2.1. Triangular Functions (TFs) and Properties [22, 23]

In the following, we define two sets of TFs on the interval [0, 1):$$\begin{array}{c}\text{(4)}& T_j^{1}t=\begin{array}{cc}1-\frac{t-j\mathcal{l}}{\mathcal{l}},& j\mathcal{l}\le t<j+1\mathcal{l};\\ 0,& \text{others},\end{array}{T}_j^{2}t=\begin{array}{cc}\frac{t-j\mathcal{l}}{\mathcal{l}},& j\mathcal{l}\le t<j+1\mathcal{l};\\ 0,& \text{others},\end{array}\end{array}$$where $j=0,\dots ,n-1$, $n=2^a$, $a$ is a positive integer, and $\mathcal{l}=1/n$. Moreover,$$\begin{array}{}\text{(5)}& T_j^{1}t+T_j^{2}t={\varphi }_{j}t,\end{array}$$where ${\varphi }_{j}t$ is the j-th block pulse function:$$\begin{array}{}\text{(6)}& {\varphi }_{j}t=\begin{array}{cc}1,& j\mathcal{l}\le t<j+1\mathcal{l};\\ 0,& \text{others}.\end{array}\end{array}$$

From the definition of TFs, it is clear that TFs are disjoint, orthogonal, and complete [22]. Therefore, we can obtain$$\begin{array}{}\text{(7)}& {{\displaystyle \int }}_0^T{T}_i^{p}x{T}_j^{p}x\mathrm{d}x={\Delta }_{ij},p=1,2,\end{array}$$where$$\begin{array}{}\text{(8)}& {\Delta }_{ij}=\begin{array}{cc}\frac{\mathcal{l}}{3},& i=j;\\ 0,& i\ne j.\end{array}\end{array}$$

We can easily define their vectors:$$\begin{array}{c}\text{(9)}& \widehat{T}t={T_0^{1}t,T_1^{1}t,\dots ,T_{n-1}^{1}t}^{\mathrm{T}},\tilde{T}t={T_0^{2}t,T_1^{2}t,\dots ,T_{n-1}^{2}t}^{\mathrm{T}},\end{array}$$and$$\begin{array}{c}\text{(10)}& \widehat{T}t{\widehat{T}}^{\mathrm{T}}t\approx \begin{array}{cccc}{T}_0^{1}t& 0& \cdots & 0\\ 0& T_1^{1}t& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & T_{n-1}^{1}t\end{array},\tilde{T}t{\tilde{T}}^{\mathrm{T}}t\approx \begin{array}{cccc}{T}_0^{2}t& 0& \cdots & 0\\ 0& T_1^{2}t& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & T_{n-1}^{2}t\end{array}.\end{array}$$

The TF vector $Tt$ is defined below:$$\begin{array}{}\text{(11)}& Tt={T_0^{1}t,T_1^{1}t,\dots ,T_{n-1}^{1}t,T_0^{2}t,T_1^{2}t,\dots ,T_{n-1}^{2}t}^T.\end{array}$$

On the basis of [25], we gain$$\begin{array}{}\text{(12)}& \tilde{T}t{\widehat{T}}^{\mathrm{T}}t=\widehat{T}t{\tilde{T}}^{\mathrm{T}}t\approx 0_{n\times n},\end{array}$$and thus$$\begin{array}{}\text{(13)}& Tt{T}^{\mathrm{T}}t\approx \text{diag}Tt.\end{array}$$

2.2. Approximation of Function

The approximation of any function $gt\in L^{2}0,1$ can be expressed as$$\begin{array}{c}\text{(14)}& gt\simeq g_{n}t={\displaystyle \sum }_{j=0}^{n-1}{\widehat{c}}_j{T}_j^{1}t+{\displaystyle \sum }_{j=0}^{n-1}{\tilde{c}}_j{T}_j^{2}t={\widehat{C}}^{\mathrm{T}}\widehat{T}t+{\tilde{C}}^{\mathrm{T}}\tilde{T}t\\ =C_n^{\mathrm{T}}Tt,\end{array}$$where$$\begin{array}{c}\text{(15)}& \widehat{C}={{\widehat{c}}_0,{\widehat{c}}_1,\dots ,{\widehat{c}}_{n-1}}_n^{\mathrm{T}},\tilde{C}={{\tilde{c}}_0,{\tilde{c}}_1,\dots ,{\tilde{c}}_{n-1}}_n^{\mathrm{T}},\\ C_n={{\widehat{c}}_0,{\widehat{c}}_1,\dots ,{\widehat{c}}_{n-1},{\tilde{c}}_0,{\tilde{c}}_1,\dots ,{\tilde{c}}_{n-1}}_{2n}^{\mathrm{T}},\\ {\widehat{c}}_j=gj\mathcal{l},\\ {\tilde{c}}_j=gj+1\mathcal{l},j=0,1,2,\dots ,n-1,\end{array}$$and $g_{n}t$ is the n-th approximation of TFs of $gt$ and $C_n$ is a coefficient vector.

The approximation of each $gt,s\in L^{2}0,1\times 0,1$ can be extended into a triangular sequence as$$\begin{array}{}\text{(16)}& gt,s=T^{\mathrm{T}}s\mathbb{Q}Tt,\end{array}$$where $Ts$ and $Tt$ are $2n_1$ dimensional and $2n_2$ dimensional triangular vectors and $\mathbb{Q}$ is a $2n_1\times 2n_2$ coefficient matrix of TFs. For convenience, let $n_1=n_2=n$. Thus, the matrix $\mathbb{Q}$ can be written as$$\begin{array}{}\text{(17)}& \mathbb{Q}={\begin{array}{cc}{\mathbb{Q}}_{11}& {\mathbb{Q}}_{12}\\ {\mathbb{Q}}_{21}& {\mathbb{Q}}_{22}\end{array}}_{2n\times 2n},\end{array}$$where$$\begin{array}{c}\text{(18)}& {\mathbb{Q}}_{11}={\begin{array}{cccc}g0,0& g0,\mathcal{l}& \cdots & g0,n-1\mathcal{l}\\ g\mathcal{l},0& g\mathcal{l},\mathcal{l}& \cdots & g\mathcal{l},n-1\mathcal{l}\\ ⋮& ⋮& \ddots & ⋮\\ gn-1\mathcal{l},0& gn-1\mathcal{l},\mathcal{l}& \cdots & gn-1\mathcal{l},n-1\mathcal{l}\end{array}}_{n\times n},{\mathbb{Q}}_{12}={\begin{array}{cccc}g0,\mathcal{l}& g0,2\mathcal{l}& \cdots & g0,n\mathcal{l}\\ g\mathcal{l},\mathcal{l}& g\mathcal{l},2\mathcal{l}& \cdots & g\mathcal{l},n\mathcal{l}\\ ⋮& ⋮& \ddots & ⋮\\ gn-1\mathcal{l},\mathcal{l}& gn-1\mathcal{l},2\mathcal{l}& \cdots & gn-1\mathcal{l},n\mathcal{l}\end{array}}_{n\times n},\\ {\mathbb{Q}}_{21}={\begin{array}{cccc}g\mathcal{l},0& g\mathcal{l},\mathcal{l}& \cdots & g\mathcal{l},n-1\mathcal{l}\\ g2\mathcal{l},0& g2\mathcal{l},\mathcal{l}& \cdots & g2\mathcal{l},n-1\mathcal{l}\\ ⋮& ⋮& \ddots & ⋮\\ gn\mathcal{l},0& gn\mathcal{l},\mathcal{l}& \cdots & gn\mathcal{l},n-1\mathcal{l}\end{array}}_{n\times n},\\ {\mathbb{Q}}_{22}={\begin{array}{cccc}g\mathcal{l},\mathcal{l}& g\mathcal{l},2\mathcal{l}& \cdots & g\mathcal{l},n\mathcal{l}\\ g2\mathcal{l},\mathcal{l}& g2\mathcal{l},2\mathcal{l}& \cdots & g2\mathcal{l},n\mathcal{l}\\ ⋮& ⋮& \ddots & ⋮\\ gn\mathcal{l},\mathcal{l}& gn\mathcal{l},2\mathcal{l}& \cdots & gn\mathcal{l},n\mathcal{l}\end{array}}_{n\times n}.\end{array}$$

Let $\mathbb{M}$ be a $2n\times 2n$ matrix. Therefore, a similar conclusion can be drawn:$$\begin{array}{}\text{(19)}& T^{\mathrm{T}}t\mathbb{M}Tt\simeq \mathrm{M}Tt,\end{array}$$where $\mathrm{M}$ is a 2n-vector whose elements equal the diagonal term of matrix $\mathbb{M}$. Thus,$$\begin{array}{}\text{(20)}& Tt{T}^{\mathrm{T}}t\mathrm{N}\simeq \widehat{\mathbb{N}}Tt,\end{array}$$where $\widehat{\mathbb{N}}=\text{diag}\mathrm{N}$ is $2n\times 2n$ diagonal matrix.

3. Integration Operational Matrix

We need to calculate the integral operator matrix in this section. First, it is easy to obtain the integral of TFs as follows:$$\begin{array}{c}\text{(21)}& {{\displaystyle \int }}_0^t{T}_{\mathrm{j}}^{1}s\mathrm{d}s=\begin{array}{cc}0,& 0\le t<j\mathcal{l};\\ t-\frac{{t-j\mathcal{l}}^2}{2\mathcal{l}}-j\mathcal{l},& j\mathcal{l}\le t<j+1\mathcal{l};\\ \frac{\mathcal{l}}{2},& j+1\mathcal{l}\le t<1,\end{array}{{\displaystyle \int }}_0^t{T}_{\mathrm{j}}^{2}s\mathrm{d}s=\begin{array}{cc}0,& 0\le t<j\mathcal{l};\\ \frac{{t-j\mathcal{l}}^2}{2\mathcal{l}},& j\mathcal{l}\le t<j+1\mathcal{l};\\ \frac{\mathcal{l}}{2},& j+1\mathcal{l}\le t<1.\end{array}\end{array}$$

And then when $t=j\mathcal{l}+j+1\mathcal{l}/2,t-{t-j\mathcal{l}}^2/2\mathcal{l}-j\mathcal{l}=3\mathcal{l}/8$. We can set $t-{t-j\mathcal{l}}^2/2\mathcal{l}-j\mathcal{l}\simeq 3\mathcal{l}/8,j\mathcal{l}\le t<j+1\mathcal{l}$, so$$\begin{array}{}\text{(22)}& {\int }_0^t{T}_j^{1}sds\approx 0,\dots ,0,\frac{3\mathcal{l}}{8},\frac{\mathcal{l}}{2},\dots ,\frac{\mathcal{l}}{2}\widehat{T}t+\tilde{T}t,\end{array}$$where the $j+1$-th element is $3\mathcal{l}/8$. Similarly,$$\begin{array}{}\text{(23)}& {\int }_0^t{T}_j^{2}sds\approx 0,\dots ,0,\frac{\mathcal{l}}{8},\frac{\mathcal{l}}{2},\dots ,\frac{\mathcal{l}}{2}\widehat{T}t+\tilde{T}t,\end{array}$$where the $j+1$-th element is $\mathcal{l}/8$.

Therefore,$$\begin{array}{c}\text{(24)}& {{\displaystyle \int }}_0^t\widehat{T}s\mathrm{d}s={\mathbb{P}}_1\widehat{T}t+{\mathbb{P}}_1\tilde{T}t,{{\displaystyle \int }}_0^t\tilde{T}s\mathrm{d}s={\mathbb{P}}_2\widehat{T}t+{\mathbb{P}}_2\tilde{T}t,\end{array}$$where the integration operational matrix is given by$$\begin{array}{c}\text{(25)}& {\mathbb{P}}_1=\frac{\mathcal{l}}{2}{\begin{array}{ccccc}\frac{3}{4}& 1& 1& \cdots & 1\\ 0& \frac{3}{4}& 1& \cdots & 1\\ 0& 0& \frac{3}{4}& \cdots & 1\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & \frac{3}{4}\end{array}}_{n\times n},{\mathbb{P}}_2=\frac{\mathcal{l}}{2}{\begin{array}{ccccc}\frac{1}{4}& 1& 1& \cdots & 1\\ 0& \frac{1}{4}& 1& \cdots & 1\\ 0& 0& \frac{1}{4}& \cdots & 1\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & \frac{1}{4}\end{array}}_{n\times n}.\end{array}$$

Therefore,$$\begin{array}{}\text{(26)}& {{\displaystyle \int }}_0^{t}Ts\mathrm{d}s\mathbb{\simeq }\mathbb{P}Tt,\end{array}$$where $\mathbb{P}$ is the integral operator matrix given below:$$\begin{array}{}\text{(27)}& \mathbb{P}={\begin{array}{cc}{\mathbb{P}}_1& {\mathbb{P}}_1\\ {\mathbb{P}}_2& {\mathbb{P}}_2\end{array}}_{2n\times 2n}.\end{array}$$

According to (14) and (26), the approximation of every integral function $gt$ is represented by$$\begin{array}{}\text{(28)}& {{\displaystyle \int }}_0^{t}gs\mathrm{d}s\simeq {{\displaystyle \int }}_0^t{C}_n^{\mathrm{T}}Tsds\simeq C_n^{\mathrm{T}}\mathbb{P}Tt.\end{array}$$

4. Stochastic Integration Operational Matrix

This section reviews the stochastic integrals of TFs (see [23] for details). Meanwhile, two methods are applied to obtain the stochastic operator matrix. First of all, we give It $\widehat{o}$ integral of TFs:$$\begin{array}{c}\text{(29)}& {\int }_0^t{T}_{\mathrm{j}}^{1}s\mathrm{d}Bs=\begin{array}{cc}0,& 0\le t<j\mathcal{l};\\ j+1Bt-Bj\mathcal{l}-\frac{1}{\mathcal{l}}{{\displaystyle \int }}_{j\mathcal{l}}^{t}s\mathrm{d}Bs,& j\mathcal{l}\le t<j+1\mathcal{l};\\ j+1Bj+1\mathcal{l}-Bj\mathcal{l}-\frac{1}{\mathcal{l}}{{\displaystyle \int }}_{j\mathcal{l}}^{j+1\mathcal{l}}s\mathrm{d}Bs,& j+1\mathcal{l}\le t<1,\end{array}{\int }_0^t{T}_{\mathrm{j}}^{2}s\mathrm{d}Bs=\begin{array}{cc}0,& 0\le t<i\mathcal{l};\\ -jBt-Bj\mathcal{l}+\frac{1}{\mathcal{l}}{{\displaystyle \int }}_{j\mathcal{l}}^{t}s\mathrm{d}Bs,& j\mathcal{l}\le t<j+1\mathcal{l};\\ -jBj+1\mathcal{l}-Bj\mathcal{l}+\frac{1}{\mathcal{l}}{{\displaystyle \int }}_{j\mathcal{l}}^{j+1\mathcal{l}}s\mathrm{d}Bs,& j+1\mathcal{l}\le t<1,\end{array}\end{array}$$where $j=0,\dots ,n-1$. Then, the stochastic integral ${{\displaystyle \int }}_{j\mathcal{l}}^{t}s\mathrm{d}Bs$ and ${{\displaystyle \int }}_{j\mathcal{l}}^{j+1\mathcal{l}}s\mathrm{d}Bs$ can be approximated by the composite trapezium rule and integration by parts as follows:$$\begin{array}{c}\text{(30)}& {{\displaystyle \int }}_{j\mathcal{l}}^{t}s\mathrm{d}Bs=tBt-j\mathcal{l}Bj\mathcal{l}-{{\displaystyle \int }}_{j\mathcal{l}}^{t}Bsds\simeq tBt-j\mathcal{l}Bj\mathcal{l}-\frac{t-j\mathcal{l}Bt+Bj\mathcal{l}}{2},\\ \text{(31)}& {{\displaystyle \int }}_{j\mathcal{l}}^{j+1\mathcal{l}}s\mathrm{d}Bs=j+1\mathcal{l}Bj+1\mathcal{l}-j\mathcal{l}Bj\mathcal{l}-{{\displaystyle \int }}_{j\mathcal{l}}^{j+1\mathcal{l}}Bsds\simeq j+1\mathcal{l}Bj+1\mathcal{l}-j\mathcal{l}Bj\mathcal{l}-\frac{\mathcal{l}Bj+1\mathcal{l}+Bj\mathcal{l}}{2}\\ =j+\frac{1}{2}\mathcal{l}Bj+1\mathcal{l}-Bj\mathcal{l}.\end{array}$$

When $t=j\mathcal{l}+j+1\mathcal{l}/2,Bt-Bj\mathcal{l}=Bj+0.5\mathcal{l}-Bj\mathcal{l}$. Considering (30), we can set$$\begin{array}{c}\text{(32)}& j+1Bt-Bj\mathcal{l}-\frac{1}{\mathcal{l}}{{\displaystyle \int }}_{j\mathcal{l}}^{t}s\mathrm{d}Bs\simeq \frac{3}{4}Bj+0.5\mathcal{l}-Bj\mathcal{l},\\ \text{(33)}& -jBt-Bj\mathcal{l}+\frac{1}{\mathcal{l}}{{\displaystyle \int }}_{j\mathcal{l}}^{t}s\mathrm{d}Bs\simeq \frac{1}{4}Bj+0.5\mathcal{l}-Bj\mathcal{l}.\end{array}$$

Then according to (30)–(33),$$\begin{array}{}\text{(34)}& {\int }_0^t\widehat{T}s\mathrm{d}Bs=\mathbb{P}{1}_s\widehat{T}t+\mathbb{P}{1}_s\tilde{T}t,\text{(35)}& {\int }_0^t\tilde{T}s\mathrm{d}Bs=\mathbb{P}{2}_s\widehat{T}t+\mathbb{P}{2}_s\tilde{T}t,\end{array}$$where$$\begin{array}{c}\text{(36)}& \mathbb{P}{1}_s={\begin{array}{ccccc}\alpha 0& \gamma 0& \gamma 0& \cdots & \gamma 0\\ 0& \alpha 1& \gamma 1& \cdots & \gamma 1\\ 0& 0& \alpha 2& \cdots & \gamma 2\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & \alpha n-1\end{array}}_{n\times n},\mathbb{P}{2}_s={\begin{array}{ccccc}\xi 0& \kappa 0& \kappa 0& \cdots & \kappa 0\\ 0& \xi 1& \kappa 1& \cdots & \kappa 1\\ 0& 0& \xi 2& \cdots & \kappa 2\\ ⋮& ⋮& ⋮& \ddots & ⋮\\ 0& 0& 0& \cdots & \xi n-1\end{array}}_{n\times n},\end{array}$$in which$$\begin{array}{c}\text{(37)}& \alpha j=\frac{3}{4}Bj+0.5\mathcal{l}-Bj\mathcal{l},\gamma j=\frac{1}{2}Bj+1\mathcal{l}-Bj\mathcal{l},\\ \xi j=\frac{1}{4}Bj+0.5\mathcal{l}-Bj\mathcal{l},\\ \kappa j=\frac{1}{2}Bj+1\mathcal{l}-Bj\mathcal{l}.\end{array}$$

From (34) and (35), we get$$\begin{array}{}\text{(38)}& {\int }_0^{t}TsdBs={\mathbb{P}}_{S}Tt,\end{array}$$in which$$\begin{array}{}\text{(39)}& {\mathbb{P}}_S={\begin{array}{cc}\mathbb{P}{1}_s& \mathbb{P}{1}_s\\ \mathbb{P}{2}_s& \mathbb{P}{2}_s\end{array}}_{2n\times 2n}.\end{array}$$

Therefore,$$\begin{array}{}\text{(40)}& {\int }_0^{t}gs\mathrm{d}Bs\simeq {\int }_0^t{C}_n^{\mathrm{T}}Ts\mathrm{d}Bs=C_n^{\mathrm{T}}{\mathbb{P}}_{S}Tt.\end{array}$$

5. Error Analysis

The error denoted by ${\epsilon }_{n}t=zt-z_{n}t$, and $z_{n}t$ is an approximation of $zt$. In this section, the error analysis of the present method is obtained. Meanwhile, it is proved that $z_{n}t$ is convergent of order $O{\mathcal{l}}^2,\mathcal{l}=1/n$.

First of all, we look back at two significant lemmas.

Lemma 1.

Assume that

(1) $gx$ is continuous on the interval [0, 1] and there is a second derivative function in (0, 1),

(2) $g_{n}x$ is an approximation with respect to $\text{TFs}$,

(3) $g^{″}x<L$ for every $x\in 0,1$, and $L$ is constant number. Then

$$\begin{array}{}\text{(41)}& {gx-g_{n}x}^2={{\epsilon }_{n}x}^2={\int }_0^1{{\epsilon }_{n}x}^2\mathrm{d}x=\mathrm{O}{\mathcal{l}}^4.\end{array}$$

Proof.

See [23] for details.

Lemma 2.

Assume $gs,t$ is a twice differentiable function on $D=0,1\times 0,1$, ${\partial }^2\mathrm{g}s,t/\partial s\partial t<M$, and ${\epsilon }_{n}s,t=gs,t-g_{n}s,t$, where $g_{n}s,t$ is an approximation with respect to $\text{TFs}$. We can gain$$\begin{array}{}\text{(42)}& {{\epsilon }_{n}s,t}^2={\int }_0^1{\int }_0^1{{\epsilon }_{n}s,t}^{2}dtds=\mathrm{O}{\mathcal{l}}^4.\end{array}$$

Proof.

See [23] for details.

Secondly, consider $\epsilon t=zt-z_{n}t$. Then we can get$$\begin{array}{c}\text{(43)}& \epsilon t=zt-z_{n}t=z_{0}t-z_{0n}t\\ +{\int }_0^t\widehat{a}s,t\rho s,zs-{\widehat{a}}_{n}s,t\rho s,z_{n}s\mathrm{d}s\\ +{\int }_0^t\tilde{a}s,t\sigma s,zs-{\tilde{a}}_{n}s,t\sigma s,z_{n}s\mathrm{d}Bs,\end{array}$$where $z_{n}t$, $z_{0n}t$, ${\widehat{a}}_{n}s,t$, and ${\tilde{a}}_{n}s,t$ are approximations of $zt,z_{0}t$, $\widehat{a}s,t$, and $\tilde{a}s,t$ by TFs.

Lastly, we give and testify the main theorem.

Theorem 3.

Suppose that bounded analytic functions $\rho t,x,\sigma t,x$ satisfy the following conditions:

(i) $\rho t,x_1-\rho t,x_2\le q_1{x}_1-x_2$, $\sigma t,x_1-\sigma t,x_2\le q_3{x}_1-x_2$,

(ii) $\rho t,x\le q_2,\sigma t,x\le q_4$,

(iii) $\widehat{a}t,s\le q_5,\tilde{a}t,s\le q_6$,

where $x_1,x_2\mathbb{\in }\mathbb{R}$, constants $q_j>0$, $j=1,\dots ,6$. Then$$\begin{array}{c}\text{(44)}& {\int }_0^T\mathrm{E}{\epsilon t}^2\mathrm{d}t={\int }_0^T\mathrm{E}{zt-z_{n}t}^2\mathrm{d}t\le \mathrm{O}{\mathcal{l}}^4,0<T<1.\end{array}$$

Proof.

$$\begin{array}{c}\text{(45)}& \mathrm{E}{\epsilon t}^2\le 3\mathrm{E}{z_{0}t-z_{0n}t}^2+\mathrm{E}{{{\displaystyle \int }}_0^t\widehat{a}s,t\rho s,zs-{\widehat{a}}_{n}s,t\rho s,z_{n}s\mathrm{d}s}^2+\mathrm{E}{{{\displaystyle \int }}_0^t\tilde{a}s,t\sigma s,zs-{\tilde{a}}_{n}s,t\sigma s,z_{n}sds}^2.\end{array}$$

According to Lipschitz’s continuity, It $\widehat{o}$ isometry, and Cauchy–Schwarz inequality, we can gain$$\begin{array}{c}\text{(46)}& \mathrm{E}{\epsilon t}^2\le 3E{z_{0}t-z_{0n}t}^2+{{\displaystyle \int }}_0^{t}E{\widehat{a}s,t\rho s,zs-\rho s,z_{n}s+\rho s,z_{n}s\widehat{a}s,t-{\widehat{a}}_{n}s,t}^2\mathrm{d}s+{{\displaystyle \int }}_0^{t}E{\tilde{a}s,t\sigma s,zs-\sigma s,z_{n}s+\sigma s,z_{n}s\tilde{a}s,t-{\tilde{a}}_{n}s,t}^2\mathrm{d}s\le 3{z_{0}t-z_{0n}t}^2+2q_1^2{q}_5^2{{\displaystyle \int }}_0^t\mathrm{E}{\epsilon s}^{2}ds+2q_2^2{{\displaystyle \int }}_0^t{\widehat{a}s,t-{\widehat{a}}_{n}s,t}^{2}ds+2q_3^2{q}_6^2{{\displaystyle \int }}_0^t\mathrm{E}{\epsilon s}^{2}ds+2q_4^2{{\displaystyle \int }}_0^t{\tilde{a}s,t-{\tilde{a}}_{n}s,t}^{2}ds.\end{array}$$

Then,$$\begin{array}{c}\text{(47)}& \mathrm{E}{\epsilon t}^2\le 3{z_{0}t-z_{0n}t}^2+2q_2^2{\int }_0^t{\widehat{a}s,t-{\widehat{a}}_{n}s,t}^{2}ds+2q_4^2{\int }_0^t{\tilde{a}s,t-{\tilde{a}}_{n}s,t}^{2}ds+6q_1^2{q}_5^2+q_3^2{q}_6^2{\int }_0^t\mathrm{E}{\epsilon s}^{2}ds=\phi t+\lambda {\int }_0^t\mathrm{E}{\epsilon s}^2\mathrm{d}s,\end{array}$$where$$\begin{array}{c}\text{(48)}& \phi t=3{z_{0}t-z_{0n}t}^2+2q_2^2{\int }_0^t{\widehat{a}s,t-{\widehat{a}}_{n}s,t}^{2}ds+2q_4^2{\int }_0^t{\tilde{a}s,t-{\tilde{a}}_{n}s,t}^{2}ds\lambda =6q_1^2{q}_5^2+q_3^2{q}_6^2.\end{array}$$

Let $\eta t=\mathrm{E}{\epsilon t}^2$, and we get$$\begin{array}{}\text{(49)}& \eta t\le \phi t+\lambda {{\displaystyle \int }}_0^t\eta v\mathrm{d}v,v\in 0,t.\end{array}$$

By Gronwall’s inequality,$$\begin{array}{}\text{(50)}& \eta t\le \phi t+\lambda {\int }_0^t{e}^{\lambda t-v}\phi v\mathrm{d}v,0\le v<t,t\in 0,1,\end{array}$$and then$$\begin{array}{c}\text{(51)}& {{\displaystyle \int }}_0^T\eta t\mathrm{d}t={{\displaystyle \int }}_0^T\mathrm{E}{\epsilon t}^2\mathrm{d}t\le {{\displaystyle \int }}_0^T\phi t+\lambda {\int }_0^t{e}^{\lambda t-v}\phi v\mathrm{d}v\mathrm{d}t\\ ={\int }_0^T\phi tdt+{\int }_0^T\lambda {\int }_0^t{e}^{\lambda t-v}\phi v\mathrm{d}v\mathrm{d}t\\ \le {\int }_0^T\phi tdt+\lambda e^{\lambda T}{\int }_0^T{\int }_0^t\phi v\mathrm{d}v\mathrm{d}t\\ =3{\int }_0^T{z_{0}t-z_{0n}t}^2\mathrm{d}t\\ +6q_2^2{\int }_0^T{\int }_0^t{\widehat{a}s,t-{\widehat{a}}_{n}s,t}^2\mathrm{d}s\mathrm{d}t\\ +6q_4^2{\int }_0^T{\int }_0^t{\tilde{a}s,t-{\tilde{a}}_{n}s,t}^2\mathrm{d}s\mathrm{d}t\\ +\lambda e^{\lambda T}3{\int }_0^T{\int }_0^t{z_{0}v-z_{0n}v}^2\mathrm{d}v\mathrm{d}t+6q_2^2{\int }_0^T{\int }_0^t{\int }_0^v{\widehat{a}s,v-{\widehat{a}}_{n}s,v}^2\mathrm{d}s\mathrm{d}v\mathrm{d}t+6q_4^2{\int }_0^T{\int }_0^t{\int }_0^v{\tilde{a}s,v-{\tilde{a}}_{n}s,v}^2\mathrm{d}s\mathrm{d}v\mathrm{d}t.\end{array}$$

By (41) and (42), we get$$\begin{array}{c}\text{(52)}& {\int }_0^T{z_{0}t-z_{0n}t}^{2}dt=k_1{\mathcal{l}}^4,{\int }_0^T{\int }_0^t{\widehat{a}s,t-{\widehat{a}}_{n}s,t}^{2}dsdt=k_2{\mathcal{l}}^4,\\ {\int }_0^T{\int }_0^t{\tilde{a}s,t-{\tilde{a}}_{n}s,t}^{2}dsdt=k_3{\mathcal{l}}^4,\\ {\int }_0^T{\int }_0^t{z_{0}v-z_{0n}v}^{2}dvdt=k_4{\mathcal{l}}^4,\\ {\int }_0^T{\int }_0^t{\int }_0^v{\widehat{a}s,v-{\widehat{a}}_{n}s,v}^{2}dsdvdt=k_5{\mathcal{l}}^4,\\ {\int }_0^T{\int }_0^t{\int }_0^v{\tilde{a}s,v-{\tilde{a}}_{n}s,v}^{2}dsdvdt=k_6{\mathcal{l}}^4.\end{array}$$

Thus,$$\begin{array}{}\text{(53)}& {\int }_0^T\mathrm{E}{\epsilon t}^2\mathrm{d}t=M{\mathcal{l}}^4\le \mathrm{O}{\mathcal{l}}^4,\end{array}$$where$$\begin{array}{}\text{(54)}& M=3k_1+6q_2^2{k}_2+6q_4^2{k}_3+\lambda e^{\lambda T}3k_4+6q_2^2{k}_5+6q_4^2{k}_6,\end{array}$$and $k_i,i=1,2,\dots ,6$ are independent nonnegative constants.

The proof is accomplished.

In [27], the authors proved that the rate of convergence for the BPF method of solutions of stochastic It $\widehat{o}$–Volterra integral equations can be estimated by $\mathrm{O}{\mathcal{l}}^2$. Compared with reference [27], the error accuracy of the approximate solution obtained by this method may be higher.

6. Numerical Method

We present the TF approximate solution to (3) in this section. For the analytic functions, the forms of $\rho \cdot ,\cdot$ and $\sigma \cdot ,\cdot$ are polynomials. We firstly introduce two important lemmas.

Lemma 4 (see [26]).

Assume that the continuous function $gt$ is Lebesgue measure, and ${\int }_{-\infty }^{\infty }{gx}^2\mathrm{d}x<\infty$; then the m-th power of $gt$ is estimated to be$$\begin{array}{c}\text{(55)}& {gt}^m={\widehat{c}}_0^m,{\widehat{c}}_1^m,\dots ,{\widehat{c}}_{n-1}^m\widehat{T}t+{\tilde{c}}_0^m,{\tilde{c}}_1^m,\dots ,{\tilde{c}}_{n-1}^m\tilde{T}t={C_n^{\mathrm{T}}}^{m}Tt,m\in N.\end{array}$$

Lemma 5.

Suppose the analytic functions $f_{1}t={{\displaystyle \sum }}_{i=0}^{r_1}{k}_{1i}{t}^i$, $f_2{z}_{n}t={{\displaystyle \sum }}_{j=0}^{r_2}{k}_{2j}{z}_n^{j}t$, $f_{3}t,z_{n}t={{\displaystyle \sum }}_{\begin{array}{c}u=0\\ a+b=u\end{array}}^{r_3}{k}_{3u}{t}^a{z}_n^{b}t$, where $r_1,r_2,r_3,a,b$ are finite nonnegative integers. $C_{n}Tt$ and $V_{n}Tt$ are the expression of $z_{n}t$ and identity function $t$ with respect to TFs. Then we have$$\begin{array}{c}\text{(56)}& f_{1}t=f_1^{\mathrm{T}}{V}_{n}Tt,f_2{z}_{n}t=f_2^{\mathrm{T}}{C}_{n}Tt,\\ f_{3}t,z_{n}t=f_3^{\mathrm{T}}{V}_n,C_{n}Tt,\end{array}$$where$$\begin{array}{c}\text{(57)}& f_1^{\mathrm{T}}{V}_n=f_1{\widehat{v}}_0,f_1{\widehat{v}}_1,\dots ,f_1{\widehat{v}}_{n-1},f_1{\tilde{v}}_0,f_1{\tilde{v}}_1,\dots ,f_1{\tilde{v}}_{n-1}{f}_2^{\mathrm{T}}{C}_n=f_2{\widehat{c}}_0,f_2{\widehat{c}}_1,\dots ,f_2{\widehat{c}}_{n-1},f_2{\tilde{c}}_0{f}_2{\tilde{c}}_1,\dots ,f_2{\tilde{c}}_{n-1},\\ f_3^{\mathrm{T}}{V}_n,C_n=f_3{\widehat{v}}_0,{\widehat{c}}_0,f_3{\widehat{v}}_1,{\widehat{c}}_1,\dots ,f_3{\widehat{v}}_{n-1},{\widehat{c}}_{n-1},f_3{\tilde{v}}_0,{\tilde{c}}_0,f_3{\tilde{v}}_1,{\tilde{c}}_1,\dots ,f_3{\tilde{v}}_{n-1},{\tilde{c}}_{n-1}.\end{array}$$

Proof.

According to Lemma 4, we can derive$$\begin{array}{c}\text{(58)}& f_{1}t={\displaystyle \sum }_{i=0}^{r_1}{k}_{1i}{t}^i={\displaystyle \sum }_{i=0}^{r_1}{k}_{1i}{\widehat{v}}_0^i,{\widehat{v}}_1^i,\dots ,{\widehat{v}}_{n-1}^i,{\tilde{v}}_0^i,{\tilde{v}}_1^i,\dots ,{\tilde{v}}_{n-1}^{i}Tt\\ =f_1^{\mathrm{T}}{V}_{n}Tt,\end{array}$$where$$\begin{array}{}\text{(59)}& f_1^{\mathrm{T}}{V}_n=f_1{\widehat{v}}_0,f_1{\widehat{v}}_1,\dots ,f_1{\widehat{v}}_{n-1},f_1{\tilde{v}}_0,f_1{\tilde{v}}_1,\dots ,f_1{\tilde{v}}_{n-1},\end{array}$$and$$\begin{array}{c}\text{(60)}& f_2{z}_{n}t={\displaystyle \sum }_{j=0}^{r_2}{k}_{2j}{z}_n^{j}t={\displaystyle \sum }_{j=0}^{r_2}{k}_{2j}{\widehat{c}}_0^j,{\widehat{c}}_1^j,\dots ,{\widehat{c}}_{n-1}^j,{\tilde{c}}_0^j,{\tilde{c}}_1^j,\dots ,{\tilde{c}}_{n-1}^{j}Tt\\ =f_2^{\mathrm{T}}{C}_{n}Tt,\end{array}$$where$$\begin{array}{}\text{(61)}& f_2^{\mathrm{T}}{C}_n=f_2{\widehat{c}}_0,f_2{\widehat{c}}_1,\dots ,f_2{\widehat{c}}_{n-1},f_2{\tilde{c}}_0,f_2{\tilde{c}}_1,\dots ,f_2{\tilde{c}}_{n-1}.\end{array}$$

Moreover, we have$$\begin{array}{c}\text{(62)}& f_{3}t,z_{n}t={\displaystyle \sum }_{\begin{array}{c}u=0\\ a+b=u\end{array}}^{r_3}{k}_{3u}{t}^a{z}_n^{b}t={\displaystyle \sum }_{\begin{array}{c}u=0\\ a+b=u\end{array}}^{r_3}{k}_{3u}{\widehat{v}}_0^a,{\widehat{v}}_1^a,\dots ,{\widehat{v}}_{n-1}^a,{\tilde{v}}_0^a,{\tilde{v}}_1^a,\dots ,{\tilde{v}}_{n-1}^{a}Tt\times {\widehat{c}}_0^b,{\widehat{c}}_1^b,\dots ,{\widehat{c}}_{n-1}^b,{\tilde{c}}_0^b,{\tilde{c}}_1^b,\dots ,{\tilde{c}}_{n-1}^{b}Tt\\ ={\displaystyle \sum }_{\begin{array}{c}u=0\\ a+b=u\end{array}}^{r_3}{k}_{3u}{T}^{\mathrm{T}}t{{\widehat{v}}_0^a,{\widehat{v}}_1^a,\dots ,{\widehat{v}}_{n-1}^a,{\tilde{v}}_0^a,{\tilde{v}}_1^a,\dots ,{\tilde{v}}_{n-1}^a}^{\mathrm{T}}\times {\widehat{c}}_0^b,{\widehat{c}}_1^b,\dots ,{\widehat{c}}_{n-1}^b,{\tilde{c}}_0^b,{\tilde{c}}_1^b,\dots ,{\tilde{c}}_{n-1}^{b}Tt\\ ={\displaystyle \sum }_{\begin{array}{c}u=0\\ a+b=u\end{array}}^{r_3}{k}_{3u}{\widehat{v}}_0^a{\widehat{c}}_0^b,{\widehat{v}}_1^a{\widehat{c}}_1^b,\dots ,{\widehat{v}}_{n-1}^a{\widehat{c}}_{n-1}^b{\tilde{v}}_0^a{\tilde{c}}_0^b,{\tilde{v}}_1^a{\tilde{c}}_1^b,\dots ,{\tilde{v}}_{n-1}^a{\tilde{c}}_{n-1}^{b}Tt\\ =f_3^{\mathrm{T}}{V}_n,C_{n}Tt,\end{array}$$where$$\begin{array}{}\text{(63)}& f_3^{\mathrm{T}}{V}_n,C_n=f_3{\widehat{v}}_0,{\widehat{c}}_0,f_3{\widehat{v}}_1,{\widehat{c}}_1,\dots ,f_3{\widehat{v}}_{n-1},{\widehat{c}}_{n-1},f_3{\tilde{v}}_0,{\tilde{c}}_0,f_3{\tilde{v}}_1,{\tilde{c}}_1,\dots ,f_3{\tilde{v}}_{n-1},{\tilde{c}}_{n-1}.\end{array}$$

In virtue of (14) and (16), we can obtain the approximate forms of $zt$, $z_{0}t$, $\widehat{a}s,t$, and $\tilde{a}s,t$ below.$$\begin{array}{}\text{(64)}& zt\simeq z_{n}t=C_n^{\mathrm{T}}Tt=T^{\mathrm{T}}t{C}_n,\text{(65)}& z_{0}t\simeq z_{0n}t=C_{0n}^{\mathrm{T}}Tt=T^{\mathrm{T}}t{C}_{0n},\text{(66)}& \widehat{a}s,t\simeq {\widehat{a}}_{n}s,t=T^{\mathrm{T}}t\widehat{\mathbb{Q}}Ts,\text{(67)}& \tilde{a}s,t\simeq {\tilde{a}}_{n}s,t=T^{\mathrm{T}}t\tilde{\mathbb{Q}}Ts,\end{array}$$where $C_n$ and $C_{0n}$ are TF coefficient vectors $\widehat{\mathbb{Q}}$ and $\tilde{\mathbb{Q}}$ are TF coefficient matrixes.

By using Lemma 5, we assume analytic functions$$\begin{array}{c}\text{(68)}& \rho t,z_{n}t=f_{1}t+f_2{z}_{n}t+f_{3}t,z_{n}t=f_1^{\mathrm{T}}{V}_{n}Tt+f_2^{\mathrm{T}}{C}_{n}Tt+f_3^{\mathrm{T}}{V}_n,C_{n}Tt,\\ \text{(69)}& \sigma t,z_{n}t=g_{1}t+g_2{z}_{n}t+g_{3}t,z_{n}t=g_1^{\mathrm{T}}{V}_{n}Tt+g_2^{\mathrm{T}}{C}_{n}Tt+g_3^{\mathrm{T}}{V}_n,C_{n}Tt,\end{array}$$where the expressions of functions $f_1\cdot$, $f_2\cdot$, $f_3\cdot$ and $g_1\cdot$, $g_2\cdot$, $g_3\cdot$ are similar.

Substituting (64)–(69) into (3), it yields$$\begin{array}{c}\text{(70)}& C_n^{\mathrm{T}}Tt=C_{0n}^{\mathrm{T}}Tt\\ +{\int }_0^t{T}^{\mathrm{T}}t\widehat{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{f}_1{V}_n+T^{\mathrm{T}}s{f}_2{C}_n+T^{\mathrm{T}}s{f}_3{V}_n,C_n\mathrm{d}s\\ +{\int }_0^t{T}^{\mathrm{T}}t\tilde{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{g}_1{V}_n+T^{\mathrm{T}}s{g}_2{C}_n+T^{\mathrm{T}}s{g}_3{V}_n,C_{n}dBs\\ =C_{0n}^{\mathrm{T}}Tt+{\int }_0^t{T}^{\mathrm{T}}t\widehat{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{f}_1{V}_n+T^{\mathrm{T}}t\widehat{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{f}_2{C}_n+T^{\mathrm{T}}t\widehat{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{f}_3{V}_n,C_{n}ds\\ +{\int }_0^t{T}^{\mathrm{T}}t\tilde{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{g}_1{V}_n+T^{\mathrm{T}}t\tilde{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{g}_2{C}_n+T^{\mathrm{T}}t\tilde{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{g}_3{V}_n,C_{n}dBs\\ =C_{0n}^{\mathrm{T}}Tt+T^{\mathrm{T}}t\widehat{\mathbb{Q}}{\int }_0^t{f}_1{\widehat{\mathbb{V}}}_n+f_2{\widehat{\mathbb{C}}}_n+f_3{\widehat{\mathbb{V}}}_n,{\widehat{\mathbb{C}}}_{n}Tsds+T^{\mathrm{T}}t\tilde{\mathbb{Q}}{\int }_0^t{g}_1{\widehat{\mathbb{V}}}_n+g_2{\widehat{\mathbb{C}}}_n+g_3{\widehat{\mathbb{V}}}_n,{\widehat{\mathbb{C}}}_{n}TsdBs\\ =C_{0n}^{\mathrm{T}}Tt+T^{\mathrm{T}}t\widehat{\mathbb{Q}}{f}_1{\widehat{\mathbb{V}}}_n+f_2{\widehat{\mathbb{C}}}_n+f_3{\widehat{\mathbb{V}}}_n,{\widehat{\mathbb{C}}}_n\mathbb{P}Tt+T^{\mathrm{T}}t\tilde{\mathbb{Q}}{g}_1{\widehat{\mathbb{V}}}_n+g_2{\widehat{\mathbb{C}}}_n+g_3{\widehat{\mathbb{V}}}_n,{\widehat{\mathbb{C}}}_n{\mathbb{P}}_{S}Tt\\ =C_{0n}^{\mathrm{T}}Tt+T^{\mathrm{T}}t{\mathbb{A}}_{1}Tt+T^{\mathrm{T}}t{\mathbb{A}}_{2}Tt\\ =C_{0n}^{\mathrm{T}}Tt+A_{1}Tt+A_{2}Tt,\end{array}$$where ${\widehat{\mathbb{C}}}_n=\text{diag}{C}_n$ and ${\widehat{\mathbb{V}}}_n=\text{diag}{V}_n$ are matrixes, ${\mathbb{A}}_1=\widehat{\mathbb{Q}}{f}_1{\widehat{\mathbb{V}}}_n+f_2{\widehat{\mathbb{C}}}_n+f_3{\widehat{\mathbb{V}}}_n,{\widehat{\mathbb{C}}}_n\mathbb{P}$ and ${\mathbb{A}}_2=\tilde{\mathbb{Q}}{g}_1{\widehat{\mathbb{V}}}_n+g_2{\widehat{\mathbb{C}}}_n+g_3{\widehat{\mathbb{V}}}_n,{\widehat{\mathbb{C}}}_n{\mathbb{P}}_S$, and $A_1$ and $A_2$ are 2n-vectors whose elements are equal to the diagonal elements of matrixes ${\mathbb{A}}_1$ and ${\mathbb{A}}_2$; then$$\begin{array}{}\text{(71)}& C_n^T=C_{0n}^T+A_1+A_2.\end{array}$$

Remark 6.

Analytic functions $\rho \cdot ,\cdot$ and $\sigma \cdot ,\cdot$ are expressed as follows:$$\begin{array}{c}\text{(72)}& \rho 1,z_{n}t=f_2{z}_{n}t=f_2^{\mathrm{T}}{C}_{n}Tt,\sigma 1,z_{n}t=g_2{z}_{n}t=g_2^{\mathrm{T}}{C}_{n}Tt.\end{array}$$

Thus, we can obtain the numerical solution of (3) as follows:$$\begin{array}{c}\text{(73)}& C_n^{\mathrm{T}}Tt=C_{0n}^{\mathrm{T}}Tt+{\int }_0^t{T}^{\mathrm{T}}t\widehat{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{f}_2{C}_{n}ds+{\int }_0^t{T}^{\mathrm{T}}t\tilde{\mathbb{Q}}Ts{T}^{\mathrm{T}}s{g}_2{C}_{n}dBs=C_{0n}^{\mathrm{T}}Tt+T^{\mathrm{T}}t\widehat{\mathbb{Q}}{\int }_0^t{f}_2{\widehat{\mathbb{C}}}_{n}Tsds+T^{\mathrm{T}}t\tilde{\mathbb{Q}}{\int }_0^t{g}_2{\widehat{\mathbb{C}}}_{n}TsdBs\\ =C_{0n}^{\mathrm{T}}Tt+T^{\mathrm{T}}t\widehat{\mathbb{Q}}{f}_2{\widehat{\mathbb{C}}}_n\mathbb{P}Tt+T^{\mathrm{T}}t\tilde{\mathbb{Q}}{g}_2{\widehat{\mathbb{C}}}_n{\mathbb{P}}_{S}Tt\\ =C_{0n}^{\mathrm{T}}Tt+T^{\mathrm{T}}t{\mathbb{Y}}_{1}Tt+T^{\mathrm{T}}t{\mathbb{Y}}_{2}Tt\\ =C_{0n}^{\mathrm{T}}Tt+Y_{1}Tt+Y_{2}Tt,\end{array}$$where ${\widehat{\mathbb{C}}}_n=\text{diag}{C}_n$ is a matrix, ${\mathbb{Y}}_1=\widehat{\mathbb{Q}}{f}_2{\widehat{\mathbb{C}}}_n\mathbb{P}$ and ${\mathbb{Y}}_2=\tilde{\mathbb{Q}}{g}_2{\widehat{\mathbb{C}}}_n{\mathbb{P}}_S$, and $Y_1$ and $Y_2$ are 2n-vectors whose elements are equal to the diagonal elements of matrixes ${\mathbb{Y}}_1$ and ${\mathbb{Y}}_2$; then$$\begin{array}{}\text{(74)}& C_n^T=C_{0n}^T+Y_1+Y_2.\end{array}$$

The present text adopts the fsolve function in MATLAB to seek the solution $C_n$ of the nonlinear equation (71) or linear equation(74). According to formula (14), we can obtain the approximate solution $z_n$ of (3).

7. Numerical Simulations

We used MATLAB software for numerical techniques. Firstly, the integration operator $\mathbb{P}$ and the integration operator ${\mathbb{P}}_s$ with Brownian motion are established by a loop function. Then, the nonlinear equation is solved by the fsolve function. Finally, the error is given by the program. We can illustrate that the above method is feasible and accurate through two numerical simulations in the present section.

Example 1.

Take into account the underlying stochastic Volterra integral equations [27]:$$\begin{array}{c}\text{(75)}& zt=z_0-{\lambda }^2{\int }_0^{t}zv1-z^{2}vdv+\lambda {\int }_0^{t}1-z^{2}vdBv,t\in 0,1,\end{array}$$and the exact solution of $zt$ is as follows:$$\begin{array}{}\text{(76)}& zt=\tanh \frac{1}{30}Bt+\mathrm{arctanh}{z}_0,\end{array}$$where $z_0=1/10$, $\lambda =1/30$. When n = 16, as shown in Figure 1, it displays the curves of approximate solutions and exact solutions calculated by TFs. Additionally, when n = 32, the curves are shown in Figure 2. The absolute error $X_e$, error standard deviations $X_s$, and error mean confidence intervals are listed in Tables 1 and 2. The average error value in Table 1 shows that the accuracy of the method of this paper is not lower than or even slightly higher than that of [27] where the numerical solution is obtained by BPF function.

[figure(s) omitted; refer to PDF]

Table 1

$X_e$, $X_s$, and 95% error mean confidence interval of Example 1, n = 16.

$t$	$X_e$	$X_s$	95% error mean confidence interval
			Lower	Upper
$1/16$	4.56535481 $\times {10}^{-6}$	2.28267741 $\times {10}^{-6}$	9.13070963 $\times {10}^{-8}$	9.03940253 $\times {10}^{-6}$
$3/16$	1.32782282 $\times {10}^{-5}$	6.63911410 $\times {10}^{-6}$	2.65564564 $\times {10}^{-7}$	2.62908918 $\times {10}^{-5}$
$5/16$	2.95558696 $\times {10}^{-5}$	1.47779348 $\times {10}^{-5}$	5.91117392 $\times {10}^{-7}$	5.85206218 $\times {10}^{-5}$
$7/16$	3.28145865 $\times {10}^{-5}$	1.64072932 $\times {10}^{-5}$	6.56291730 $\times {10}^{-7}$	6.49728812 $\times {10}^{-5}$
$9/16$	5.51699992 $\times {10}^{-5}$	2.75849996 $\times {10}^{-5}$	1.10339998 $\times {10}^{-6}$	1.09236598 $\times {10}^{-4}$

Table 2

$X_e$, $X_s$, and 95% error mean confidence interval of Example 1, n = 32.

$t$	$X_e$	$X_s$	95% error mean confidence interval
			Lower	Upper
$1/32$	2.25905255 $\times {10}^{-6}$	1.12952628 $\times {10}^{-6}$	4.51810510 $\times {10}^{-8}$	4.472924 $\times {10}^{-6}$
$5/32$	3.93200515 $\times {10}^{-5}$	1.96600258 $\times {10}^{-5}$	7.86401030 $\times {10}^{-7}$	7.785370 $\times {10}^{-5}$
$9/32$	2.79234314 $\times {10}^{-5}$	1.39617157 $\times {10}^{-5}$	5.58468627 $\times {10}^{-7}$	5.528839 $\times {10}^{-5}$
$13/32$	3.52202610 $\times {10}^{-5}$	1.76101305 $\times {10}^{-5}$	7.04405221 $\times {10}^{-7}$	6.973612 $\times {10}^{-5}$
$17/32$	5.64767604 $\times {10}^{-5}$	2.82383802 $\times {10}^{-5}$	1.12953521 $\times {10}^{-6}$	1.118240 $\times {10}^{-4}$

Example 2.

Take into account the underlying stochastic Volterra integral equations [15]:$$\begin{array}{c}\text{(77)}& zt=\frac{\pi }{32}+\frac{1}{64}{\int }_0^t\cos zs{\sin }^{3}zsds-\frac{1}{8}{{\displaystyle \int }}_0^t{\sin }^{2}zs\mathrm{d}Bs,t\in 0,1,\end{array}$$and the exact solution of $zt$ is as follows:$$\begin{array}{}\text{(78)}& zt=\mathrm{arccot}\frac{1}{8}Bt+\cot \frac{\pi }{32}.\end{array}$$

When n = 16, Figure 3 displays the curves of approximate and exact solutions calculated by this method in Example 2, and their errors are very small. Additionally, when n = 32, the curves are shown in Figure 4. From Tables 3 and 4, we can see the absolute error $X_e$, error standard deviations $X_s$, and confidence intervals for error mean. Compared with [15], the method of this paper has better simulation effect.

[figure(s) omitted; refer to PDF]

Table 3

$X_e$, $X_s$, and 95% confidence interval for error mean of Example 2, n = 16.

$t$	$X_e$	$X_s$	95% confidence interval for error mean
			Lower	Upper
$1/16$	7.92998982 $\times {10}^{-7}$	3.96499491 $\times {10}^{-7}$	1.58599796 $\times {10}^{-8}$	1.57013798 $\times {10}^{-6}$
$3/16$	1.72232704 $\times {10}^{-6}$	8.61163520 $\times {10}^{-7}$	3.44465408 $\times {10}^{-8}$	3.41020754 $\times {10}^{-6}$
$5/16$	2.97356464 $\times {10}^{-6}$	1.48678232 $\times {10}^{-6}$	5.94712927 $\times {10}^{-8}$	5.88765798 $\times {10}^{-6}$
$7/16$	5.51948609 $\times {10}^{-6}$	2.75974305 $\times {10}^{-6}$	1.10389722 $\times {10}^{-7}$	1.09285825 $\times {10}^{-5}$
$9/16$	1.17303229 $\times {10}^{-5}$	5.86516146 $\times {10}^{-6}$	2.34606458 $\times {10}^{-7}$	2.32260394 $\times {10}^{-5}$

Table 4

$X_e$, $X_s$, and 95% confidence interval for error mean of Example 2, n = 32.

$t$	$X_e$	$X_s$	95% confidence interval for error mean
			Lower	Upper
$1/32$	4.93348608 $\times {10}^{-7}$	2.46674304 $\times {10}^{-7}$	9.86697216 $\times {10}^{-9}$	9.76830243 $\times {10}^{-7}$
$5/32$	3.39374157 $\times {10}^{-6}$	1.69687079 $\times {10}^{-6}$	6.78748315 $\times {10}^{-8}$	6.71960835 $\times {10}^{-6}$
$9/32$	5.60005472 $\times {10}^{-6}$	2.80002736 $\times {10}^{-6}$	1.12001094 $\times {10}^{-7}$	1.10881082 $\times {10}^{-5}$
$13/32$	7.01784998 $\times {10}^{-6}$	3.50892499 $\times {10}^{-6}$	1.40356999 $\times {10}^{-7}$	1.38953430 $\times {10}^{-5}$
$17/32$	9.10622234 $\times {10}^{-6}$	4.55311117 $\times {10}^{-6}$	1.82124447 $\times {10}^{-7}$	1.80303202 $\times {10}^{-5}$

8. Applications in Mathematical Finance

In this section, we present a practical application example of general stock model in the financial field [13]. The market consists of a riskless cash bond, ${Bt}_{t\ge 0}$, and a single risky asset with price process ${St}_{t\ge 0}$. We can consider the following general stock model:$$\begin{array}{}\text{(79)}& \begin{array}{cc}dBt=\sin tBtdt,B_0=1,& t\in 0,0.8,\\ St=\frac{1}{10}+{\int }_0^t\ln 1+uSudu+{\int }_0^{t}uSudWu& u,t\in 0,0.8.\end{array}\end{array}$$

The solution of this model is $Bt=e^{1-\cos t}$, and $St=1/10e^{1+t\ln 1+t-t-t^3/6+{\int }_0^{t}udWu}$, where ${Wt}_{t\ge 0}$ is a Brownian motion generating the filtration ${\mathcal{F}t}_{t\ge 0}$.

${\overline{X}}_e$ is mean of error and $S_e$ is standard deviation of error. The numerical results are shown in Tables 5 and 6. It can be seen that the mean of error obtained by the triangular functions is smaller than that of the literature [13], so the triangular function is more effective and has better application in practical problems.

Table 5

${\overline{X}}_e$, $S_e$, and mean confidence interval for error in block pulse functions.

$n$	${\overline{X}}_e$	$S_e$	95% confidence interval for error mean
			Lower	Upper
8	0.02607134096	0.01053117415	0.02145585039	0.03068683153
16	0.03559123141	0.00982903575	0.03128346656	0.03989899626
32	0.05520441602	0.01451123788	0.04884458561	0.06156424644
64	0.09662823642	0.01970281731	0.08799309601	0.10526337684

Table 6

${\overline{X}}_e$, $S_e$, and mean confidence interval for error in triangular functions.

$n$	${\overline{X}}_e$	$S_e$	95% confidence interval for error mean
			Lower	Upper
8	0.01174809815	0.00539360303	0.001176636213	0.02231956009
16	0.01404862614	0.00651627133	0.001276734329	0.02682051796
32	0.00996974476	0.00443011722	0.001286715001	0.01865277453
64	0.01015361632	0.00479242698	0.000760459434	0.01954677321

9. Conclusion

Stochastic Volterra integral equations usually do not have explicit solutions, so it is very essential to discuss the numerical solutions of these equations. This article presents the numerical analysis of (3) with Brownian motion by using the basis function method under the global Lipschitz condition. Compared with the literature [9], the nonlinear stochastic It $\widehat{o}$–Volterra integral equation is transformed into an algebraic equation, and the coefficients in this method do not need integration, so the calculation is simpler, and the error accuracy of the approximate solution obtained by this method is higher. Through the above numerical simulation, it demonstrates the accuracy and validity of the current method. We think that there are further research and improvement of the triangular function method to solve stochastic Volterra integral equations. For example, we will research the nonlinear stochastic It $\widehat{o}$–Volterra integral equation driven by fractional Brownian motion and explore potential combinations of the triangular functions with other numerical methods.

Acknowledgments

This research was supported by the NSF of Hubei Province in China (no. 2023AFD013).

References

[1] A. Babaei, H. Jafari, S. Banihashemi, M. Ahmadi, "Mathematical analysis of a stochastic model for spread of Coronavirus," Chaos, Solitons & Fractals, vol. 145,DOI: 10.1016/j.chaos.2021.110788, 2021.

[2] M. Khodabin, K. Maleknejad, M. Rostami, M. Nouri, "Interpolation solution in generalized stochastic exponential population growth model," Applied Mathematical Modelling, vol. 36 no. 3, pp. 1023-1033, DOI: 10.1016/j.apm.2011.07.061, 2012.

[3] Y. Ait-Sahalia, "Testing continuous-time models of the spot interest rate," Review of Financial Studies, vol. 9 no. 2, pp. 385-426, DOI: 10.1093/rfs/9.2.385, 1996.

[4] W. Zhang, "Strong convergence of the Euler–Maruyama method for a class of stochastic Volterra integral equations," Journal of Computational Mathematics, vol. 40 no. 4, pp. 607-623, DOI: 10.4208/jcm.2101-m2020-0070, 2022.

[5] M. Li, C. Huang, J. Wen, "A two-parameter Milstein method for stochastic Volterra integral equations," Journal of Computational and Applied Mathematics, vol. 404,DOI: 10.1016/j.cam.2021.113870, 2022.

[6] G. H. Zhao, M. H. Song, Z. W. Yang, "Mean-square stability of analytic solution and Euler–Maruyama method for impulsive stochastic differential equations," Applied Mathematics and Computation, vol. 251, pp. 527-538, DOI: 10.1016/j.amc.2014.11.098, 2015.

[7] G. Lan, M. Zhao, S. Qi, "Exponential stability of Θ-EM method for nonlinear stochastic volterra integro-differential equations," Applied Numerical Mathematics, vol. 172, pp. 279-291, DOI: 10.1016/j.apnum.2021.10.012, 2022.

[8] F. Usta, M. Akyiit, F. Say, K. J. Ansari, "Bernstein operator method for approximate solution of singularly perturbed Volterra integral equations," Journal of Mathematical Analysis and Applications, vol. 507 no. 2,DOI: 10.1016/j.jmaa.2021.125828, 2022.

[9] M. S. Barikbin, A. R. Vahidi, T. Damercheli, E. Babolian, "An iterative shifted Chebyshev method for nonlinear stochastic it o-Volterra integral equations-ScienceDirect," Journal of Computational and Applied Mathematics, vol. 378, 2020.

[10] F. S. Md Nasrudin, C. Phang, "Numerical solution via operational matrix for solving prabhakar fractional differential equations," Journal of Mathematics, vol. 2022,DOI: 10.1155/2022/7220433, 2022.

[11] S. zbasi, "Fractional Bell collocation method for solving linear fractional integro-differential equations," Mathematical Sciences, vol. 18, 2022.

[12] M. Asgari, E. Hashemizadeh, M. Khodabin, K. Maleknejad, "Numerical solution of nonlinear stochastic integral equation by stochastic operational matrix based on Bernstein polynomials," Bulletin math e´matique de la Soci e´t e´des Sciences Math me´atiques de Roumanie, vol. 57105 no. 1, 2014.

[13] K. Maleknejad, M. Khodabin, M. Rostami, "Numerical solution of stochastic Volterra integral equations by a stochastic operational matrix based on block pulse functions," Mathematical and Computer Modelling, vol. 55 no. 3-4, pp. 791-800, DOI: 10.1016/j.mcm.2011.08.053, 2012.

[14] M. T. Deng, G. Jiang, T. Ke, "Numerical solution of nonlinear stochastic it o-volterra integral equations driven by fractional brownian motion using block pulse functions," Discrete Dynamics in Nature and Society, vol. 2021,DOI: 10.1155/2021/4934658, 2021.

[15] A. Pourdarvish, K. Sayevand, I. Masti, S. Kumar, "Orthonormal Bernoulli polynomials for solving a class of two dimensional stochastic volterra-fredholm integral equations," International Journal of Algorithms, Computing and Mathematics, vol. 8 no. 1, 2022.

[16] S. C. Shiralashetti, L. Lamani, "Fibonacci wavelet based numerical method for the solution of nonlinear Stratonovich Volterra integral equations," Scientific African, vol. 10,DOI: 10.1016/j.sciaf.2020.e00594, 2020.

[17] M. Asif, I. Khan, N. Haider, Q. Al-Mdallal, "Legendre multi-wavelets collocation method for numerical solution of linear and nonlinear integral equations," Alexandria Engineering Journal, vol. 59 no. 6, pp. 5099-5109, DOI: 10.1016/j.aej.2020.09.040, 2020.

[18] Y. Rostami, "A new wavelet method for solving a class of nonlinear partial integro-differential equations with weakly singular kernels," Mathematical Sciences, vol. 16 no. 3, pp. 225-235, DOI: 10.1007/s40096-021-00414-4, 2022.

[19] A. R. Yaghoobnia, M. Khodabin, R. Ezzati, "Numerical solution of stochastic it o ^ -Volterra integral equations based on Bernstein multi-scaling polynomials," Applied Mathematics-A Journal of Chinese Universities, vol. 36 no. 3, pp. 317-329, DOI: 10.1007/s11766-021-3694-9, 2021.

[20] S. N. Kiaee, M. Khodabin, R. Ezzati, A. M. Lopes, "A new approach to approximate solutions of single time-delayed stochastic integral equations via orthogonal functions," Symmetry, vol. 14 no. 10,DOI: 10.3390/sym14102085, 2022.

[21] T. A. A. Cheraghi, M. Khodabin, R. Ezzati, "Numerical solution of linear stochastic Volterra integral equations via new basis functions," Filomat, vol. 33 no. 18, pp. 5959-5966, DOI: 10.2298/fil1918959c, 2019.

[22] A. Deb, A. Dasgupta, G. Sarkar, "A new set of orthogonal functions and its application to the analysis of dynamic systems," Journal of the Franklin Institute, vol. 343 no. 1,DOI: 10.1016/j.jfranklin.2005.06.005, 2006.

[23] M. Khodabin, K. Maleknejad, F. Hossoini Shckarabi, "Application of triangular functions to numerical solution of stochastic volterra integral equations," IAENG International Journal of Applied Mathematics, vol. 43 no. 1, 2013.

[24] S. Nemati, A. H. Anas, "Numerical solution of multi-order fractional differential equations using generalized sine-cosine wavelets," Universal Journal of Mathematics and Applications, vol. 1 no. 4, pp. 215-225, DOI: 10.32323/ujma.427381, 2018.

[25] E. Babolian, Z. Masouri, S. Hatamzadeh-Varmazyar, "Numerical solution of nonlinear Volterra Fredholm integro-differential equations via direct method using triangular functions," Computers & Mathematics with Applications, vol. 58 no. 2, pp. 239-247, DOI: 10.1016/j.camwa.2009.03.087, 2009.

[26] S. K. Damarla, M. Kundu, "Numerical solution of multi-order fractional differential equations using generalized triangular function operational matrices," Applied Mathematics and Computation, vol. 263, pp. 189-203, DOI: 10.1016/j.amc.2015.04.051, 2015.

[27] X. Y. Sang, G. Jiang, J. H. Wu, "Numerical solution of nonlinear stochastic it o ^ -Volterra integral equations by block pulse functions," Mathematica Applicata, vol. 32, pp. 935-946, 2019.