Orthonormal Bernoulli Polynomials for Solving a Class of Two Dimensional Stochastic Volterra-Fredholm Integral Equations.

Analyzing the mathematical models involving Itô integral, in particular in science and engineering has received much attention, and the reason for this issue is the randomness and lack of access to the exact answer of this type of models. For this purpose, in this paper, the approximate solution of two dimensional (2-D) stochastic Volterra-Fredholm integral equations (SVFIEs) based on the operational matrix method and orthonormal Bernoulli polynomials (OBP) is investigated. Some results and convergence analysis are also presented. Finally, by presenting three examples and reviewing the results and numerical comparisons, we showed that the proposed method has an excellent performance.

Introduction and Basic Definitions

Mathematical models can greatly assist researchers around the world. Many researchers employed theoretical frameworks and numerical simulations to survey the manner of transmission of different infectious diseases. Recently, Babei et al. [1] examined a stochastic model for the transmission of Coronavirus, which is a recent major crisis in many communities. Ganji et al. [2] presented a mathematical model of a brain tumor. This model was an extension of a simple two-dimensional mathematical model of glioma growth and diffusion. In [3] the Klein-Gordon equation (KGE) which has many scientific applications such as quantum field theory, nonlinear optics and solid state physics is investigated. Jafari et al. [4], examined the population dynamics model, including the hunting-hunting problem where the logistics equation, is generalized using the fractional operator.

In many cases, real problems are dependent on some random factors which were ignored due to poor computational power. By increasing the computational power in recent years, these problems are modeled by various types of stochastic equations. Stochastic processes occur in many real issues such as control systems [5], biological population growth [6], biology and medicine [7]. In recent decades, due to the importance of stochastic differential equations (SDE) and stochastic integral equations (SIE) in modeling programs where there is considerable uncertainty, scientists have studied the stochastic process and its applications [8-12].

The 2D Volterra–Fredholm integral equations are an important class of multi dimensional integral equations which arise in various physical and biological models. These equations can be rarely solved exactly and computational complexity of mathematical operations is important obstacle for solving high dimensional stochastic integral equations. Therefore obtaining their numerical solutions have attracted the attention of researchers in numerical analysis branch and some numerical approaches have been presented to solve stochastic integral equations. Such as: operational matrix method based on hat functions [13] and block-pulse [14]. Also, some studies have been performed on different types of stochastic differential integral equations [15-17].

Until now, the number of published papers on numerical solution of multi dimensional stochastic integral equations are very few. Therefore, in this framework this study deals with the numerical solution of 2-D SVFIEs in the following formwhere . In this equation, g(u, v) and for are known functions and z(u, v) are unknown functions that we try to determine with our numerical method and is a Brownian process. Also in this equation is 2-D Itô integral.

This study applies the operational matrix method based on the OBP to obtain the approximate solution of the SVFIEs and the analysis of convergence of the proposed method is discussed.

Itô Integral

In this section, we briefly introduce the Itô integral and its properties. The reader can refer to references [8, 16, 17] for more information.

Definition 1.1

[8, 18] For , is called the Brownian motion, if it satisfies the following properties:

For , is a continuous function of t.

For , , is independent of the past.

has normal distribution i.e., and .

Remark 1.2

By using property (), we consider [16].

Definition 1.3

Assume now that and is a sequence of elementary functions which satisfies in the following formthen, the Itô integral of g [17] is defined as follows:

Property 1.4

One of the valuable properties of Itô integral is as follows [17]:

2-D Orthonormal Bernoulli Polynomials

Definition 1.5

The Bernoulli polynomial of order m is defined as follows [19]:

One of the features of polynomials in the finite-dimensional, that makes them effective in numerical methods is their orthogonality. But Bernoulli polynomials are not orthogonal, despite their very useful properties [19]. Using Gram-Schmidt orthogonalization process, these polynomials can be orthogonalized as defined below:

Definition 1.6

The orthonormal Bernoulli polynomial of order m is defined as follows:

Therefore, according to the orthogonal property, we will have:In this case, we express a set of OBP as follows:

Definition 1.7

2-D orthonormal Bernoulli polynomial of order m, n is defined as follows:which applies to the following property:

Definition 1.8

The 2-D OBP vector is defined as:orwhere and are one-dimensional OBP vector and is the Kronecker product.

Assume now that . Then, by using Eq. (1.6) we obtainwhereand is signified by Since , therefore:If and , then:whereWe can approximate the any function of two variables in the interval as follows:where andBy doing calculations, it can easily be shown that:where C is defined as follows:Similarly, an arbitrary function k(r, s, e, w) is estimated via 2-D OBP as:where K is a matrix of order .

Some Properties of Operational Matrices

Assume that . Then, we conclude thattherefore:By performing calculations from Eqs. (1.12) and (2.2), we deduce that

Theorem 2.1

Assume now that be the OBP vector. Consequently, the operational matrix of integration has the following form:where . Here, A is given in Eq. (1.15) and M is defined by:

Proof

By performing calculations, the result can be easily obtained.

Theorem 2.2

The operational matrix of integration based of integer order on the 2-D OBP vector defined in Eq. (1.11) for mixed variables is as follows:where .

From Eq. (2.7) we have

Theorem 2.3

Assume now that be OBP vector. Consequently, the stochastic operational matrix of integration based on the OBP can be defined in the following form:where , A is given in Eq. (1.15) and E is defined as:

We have:whereNow, by approximating the integral in Eq. (2.9) in sense of trapezoidal rule for , with one recoversTherefore:We approximate and , , by and , respectively. In other words By using Eq. (1.13) we can writeFurthermore, by inserting Eq. (2.12) into Eq. (2.13) it results:

To obtain and , in Eq. (2.12) we apply the Definition (1.1), where has normal distribution such that:Now, assume that , for some and write for denoting with . Remark (1.2) says that and Definition (1.1) statesso that are independent random variable of the form [20].

Theorem 2.4

The stochastic operational matrix of integration based of integer order on the 2-D OBP vector defined in Eq. (1.11) is as follows:where .

From Eq. (2.14) we have

Theorem 2.5

Assume now that is an arbitrary vector. The operational matrix of product using 2-D OBP can be given as follows:

For the first one, from Eq. (1.17) it followsNow, we approximate all functions for by 2-D OBP:whereFurthermore, for and we deduce thatTherefore:where is an matrix, which is as follows:Now, we define where:Consequently, by defining , we elicit

Method of Solution

Hereunder, we employ the operational matrix method based on 2-D OBP for solving Eq. (1.1). For this purpose we expand z(u, v), g(v, u) and for with 2-D OBP so that:with Z, G, , , and as defined in the previous section. We will now examine the integrals in Eq. (1.1), respectively. In this way we will have:For the second integral in Eq. (1.1) we have:Similarly:By placing Eq. (3.1) to Eq. (3.6) in Eq. (1.1), we have:To calculate the unknown coefficients Z from the above equation, we use suitable collocation points. Thus Eq. (3.7) becomes a linear system with equations and unknowns. After calculating the unknown coefficients with the help of Equation , we can get the approximate answer of the equation.

Some Results and Convergence Analysis

Theorem 4.1

Suppose that where , is an approximation of a continuous function z(u) on [0, 1] by the 1-D OBP. Then, the coefficients for are bounded as follows:whereand is an arbitrary constant such that .

Using the OBP, z(u) can be approximated aswhere can be determined bySince the function z(u) is continuous on [0, 1], thenThus, we will have:

Theorem 4.2

Suppose that is an approximation of a continuous function z(u) on [0, 1] by the 1-D OBP. Then, the error bound is as follows:which is presented in Eq. (4.2).

Suppose that and . Therefore:Since , then:

Theorem 4.3

Suppose is the best approximation for the function f(u, v) by the 2-D OBP. In this case, we havewhich means that if , then . Here, is a big- notation.

See [21].

Note 4.4

By an argument similar to Theorem (4.3) , we have:

Theorem 4.5

Suppose . In this case, we havewhere , P is a probability function and is cumulative standard normal distribution function.

See [22].

Theorem 4.6

Let z(u, v) and are the exact and 2-D OBP approximate solution of Eq. (1.1). Also let and are the 2-D OBP of g(u, v) and for , respectively. Also assume that: Then, we have:

We haveThe approximate equation of Eq. (4.12) is:Therefore, from Eqs. (4.12) and (4.13) we have:We now estimate Eq. (4.14) by 2D-integrals mean value theorem for every and . Therefore:On the other hand, for By assumptions (i) and (ii) for and Note (4.4) we obtain:Substituting Eq. (4.17) in Eq. (4.15) and using Note (4.4) implies that:Therefore, we deduce thatSo:

Note 4.7

By an argument similar to Theorem (4.6), we have:

Numerical Feedbacks

In this section, we show the efficiency and accuracy of the suggested method by giving three examples. All the numerical calculations are performed on an Intel Core i7 laptop by running some programming codes written in MATLAB R2016b.

Example 5.1

Assuming , , and in Eq. (1.1), consider the SVFIE that was studied by [13, 14]:In this case, assuming , the exact solution of the equation is .

In Table 1, the values of absolute error at different points for and are reported which indicates the accuracy of the suggested scheme. Table 2, provides a comparison of maximum absolute error for different values of m in points (0.1, 0.2) and (0.3, 0.7) between the presented scheme and previously published algorithms by Fallahpour et al. [14].

Table 1

Absolute error of Example 5.1 for , and

(u,v)	m=2	m=3
(0,0)	4.3509×10^-12	1.3028×10^-13
(0.1,0.1)	1.0652×10^-12	2.4236×10^-13
(0.2,0.2)	3.7256×10^-11	1.9625×10^-11
(0.3,0.3)	2.5477×10^-11	1.2239×10^-12
(0.4,0.4)	2.2593×10^-10	1.1144×10^-11
(0.5,0.5)	3.3341×10^-10	1.5870×10^-12
(0.6,0.6)	1.4257×10^-11	1.6791×10^-13
(0.7,0.7)	2.3569×10^-11	2.9614×10^-12
(0.8,0.8)	2.2910×10^-11	4.0209×10^-12
(0.9,0.9)	1.8456×10^-10	3.4962×10^-12
(1,1)	1.2958×10^-9	5.7024×10^-10

Table 2

Comparing maximum absolute error of Example 5.1

	Maximum error in (0.1, 0.2)
m	[14]	[13]	Our  method
2	0.0478382	0.117×10^-11	0.2796×10^-12
3	0.0072947	-	0.1285×10^-13
4	0.0033901	0.222×10^-14	0.5035×10^-16

Also, absolute error for and are plotted in Figs. 1 and 2 , respectively. By examining the absolute error presented in Figs. 1 and 2 for different values of m, we will find that by increasing m, the error decreases substantially.

Fig. 1

The absolute errors with for Example 5.1

Fig. 2

The absolute errors with for Example 5.1

Example 5.2

Consider the SVFIE that was studied by [14]:whereThe exact solution of (5.2) is .

In Table 3, the values of absolute error at different points for and are reported which indicates the accuracy of the suggested scheme. Table 4, provides a comparison of maximum absolute error for different values of m in points (0.1, 0.2) and (0, 0.6) between the presented scheme and previously published algorithms by Fallahpour et al. [14].

Table 3

Absolute error of Example 5.2 for , and

(u,v)	m=2	m=3
(0,0)	3.2398×10^-5	0.7164×10^-9
(0.1,0.1)	1.6148×10^-5	3.6453×10^-8
(0.2,0.2)	2.4195×10^-4	1.9021×10^-8
(0.3,0.3)	2.936×10^-5	1.2293×10^-7
(0.4,0.4)	1.2864×10^-4	2.5051×10^-8
(0.5,0.5)	2.3285×10^-4	2.7438×10^-7
(0.6,0.6)	1.7705×10^-4	2.6781×10^-7
(0.7,0.7)	3.2308×10^-4	2.9154×10^-8
(0.8,0.8)	6.1029×10^-5	4.9020×10^-8
(0.9,0.9)	4.6548×10^-5	4.4286×10^-8
(1,1)	4.5927×10^-4	5.3481×10^-7

Table 4

Comparing maximum absolute error of Example 5.2

m	Maximum error in (0.1, 0.2)
	[14]	Our method
2	0.1737880	1.5267×10^-4
3	0.0105896	2.9784×10^-7
4	0.0665286	3.7765×10^-9

Also, absolute error for and are plotted in Figs. 3 and 4, respectively. Three dimensional graphs of the exact and approximation solution is shown in Fig. 5. It can easily be seen from the figures that as the number of m increases the approximate solution converges towards the exact solution. Figure 6 shows the approximate solution and the exact solution at the top and the absolute error at the bottom, for and .

Fig. 3

The absolute errors with for Example 5.2

Fig. 4

The absolute errors with for Example 5.2

Fig. 5

The absolute errors with for Example 5.2

Fig. 6

The approximate solution and the exact solution and absolute errors of Example 5.2 for and

Example 5.3

Consider the SVFIE that was studied by [14]:whereThe exact solution of (5.4) is .

In Table 5, the values of absolute error at different points for and are reported which indicates the accuracy of the suggested scheme. Table 6, provides a comparison of maximum absolute error for different values of m in points (0.6, 0.8) and (0.4, 0.9) between the presented scheme and previously published algorithms in [13, 14].

Table 5

Absolute error of Example 5.3 for , and

(u,v)	m=2	m=3
(0,0)	1.294×10^-5	0.617×10^-7
(0.1,0.1)	1.469×10^-4	4.806×10^-6
(0.2,0.2)	6.635×10^-4	1.436×10^-6
(0.3,0.3)	9.891×10^-4	1.292×10^-6
(0.4,0.4)	1.918×10^-3	2.970×10^-6
(0.5,0.5)	2.056×10^-3	9.892×10^-6
(0.6,0.6)	2.870×10^-3	1.523×10^-5
(0.7,0.7)	2.106×10^-3	1.119×10^-5
(0.8,0.8)	1.875×10^-3	5.973×10^-6
(0.9,0.9)	3.042×10^-5	2.807×10^-6
(1,1)	1.597×10^-3	4.742×10^-6

Table 6

Comparing maximum absolute error of Example 5.3

m	Maximum error in (0.6, 0.8)
	[14]	Our method
2	0.7683170	6.537×10^-5
4	0.156307	5.668×10^-8

Also, absolute error for , and are plotted in Figs. 7, 8 and 9, respectively. Figure 10 shows the approximate solution and the exact solution for and at the top and the 3-D graphs of the approximation solution at the bottom.

Fig. 7

The absolute errors with for Example 5.3

Fig. 8

The absolute errors with for Example 5.3

Fig. 9

The absolute errors with for Example 5.3

Fig. 10

The approximate solution and the exact solution for and and 3-D approximate solution of Example 5.3

Conclusion

Obtaining approximate answers for Itô integral with their applications to modelling in science and engineering is very important offers. Also it is a comprehensive examination to the most important issues of stochastic differential equations. Due to computational complexity, in this study, we develop a new matrix approach based on 2D-OBP to solve two dimensional stochastic Volterra–Fredholm integral equations with the given boundary conditions. Based on this, 2-D operational matrix of integration and 2-D stochastic operational matrix of integration is obtained, which leads to the converting the problem into a linear system. Some results and convergence analysis are also presented. Finally, by presenting three examples and reviewing the results and numerical comparisons, we showed that the proposed method has an excellent performance. In other words, numerical examples show that the suggested method is very accurate, simple and powerful tool to solve various problems. Another considerable advantage of this method in case of linear problems is that the exact solution can be obtained from mth-order approximation when the solution is a polynomial of degree equal to or less than m.