New exact and numerical solutions with their stability for Ito integro-differential equation via Riccati–Bernoulli sub-ODE method

This paper points out several exact travelling wave solutions for -dimensional Ito integro-differential equation via the Riccati–Bernoulli sub-ODE approach. We also aim to develop a numerical solution of the respective equation using central finite difference formulas. The stability of the presented exact and numerical solutions is also deduced using the Hamiltonian system and Von Neumann's concept, respectively. Moreover, the numerical schemes are studied in terms of their accuracy. The relative error arising from executing the numerical method is exhibited. We compare our results with others published in some articles. The accomplished numerical solutions are successfully compared with the analytical ones. The used processes can be extended to solve more integrable problems as well as non-integrable ones.


Introduction
Most natural events are modelled and analysed using partial differential equations (PDEs). Nonlinear evolution equations usually contribute in giving efficient, acceptable and practical results when they are used to describe some phenomena. Therefore, such equations are commonly utilized to deal with some sophisticated problems arising in various fields such as biology, chemistry, fluid dynamics, fibres, condensed matter physics, heat transfer, plasma physics, optics, etc. Many experts have established a considerable number of techniques to discover and construct the exact solutions of the developed model. Among the existing ones, we reveal the following: the Adomian decomposition process [1,2], the sine-cosine principal [3,4], the rank analysis approach [5], the tanh-sech principal [6,7], the extended tanh process [8,9], the complex hyperbolic function concept [10,11], the exp(−f (ζ ))expansion approach [12,13]. For more methods, one can see [14][15][16][17][18][19][20][21][22][23].
Ito established a (1 + 1)-dimensional Ito equation in 1980. This equation is a generalization of the bilinear KdV problem [24]. The (1 + 1)-dimensional Ito equation reads as follows: According to [25], Equation (1) is primarily related to the Fokker-Planck model. Ito equation has received much attention from several researchers. For instance, the authors in [26] have extracted the exact solutions of Equation (1) using three different techniques, namely, the generalized direct algebraic approach, the long simple equation process, and the modified F-expansion technique. In [25], the authors have applied the extended homoclinic test approach to construct some exact soliton solutions, including two soliton solutions, breather type of soliton, and periodic type of soliton. Moreover, the Kudryashov method is employed in [27] to obtain the exact solution of the (1 + 1)-dimensional form of the generalized Ito integrodifferential equation.
Since our understanding of building exact travelling wave solutions for Ito equation is extensively based on minimal techniques, we employ the Riccati-Bernoulli sub-ODE method on the governing equation to discover some exact travelling wave solutions in the form of trigonometric and hyperbolic functions. The proposed technique, which allows us to achieve complicated algebraic calculations, is used to construct peaked wave solutions, solitary wave solutions, and exact wave solutions for nonlinear evolution equations [28]. Once the Riccati-Bernoulli equation approach is used, the proposed equation is changed to be a system of algebraic equations that can be easily solved. The central finite difference methods are also invoked to seek numerical results. Verily, the achieved results effectively contribute to explain some practically physical models and some complex problems.

Summary of the Riccati-Bernoulli sub-ODE process
This section is mainly devoted to briefly summarize the Riccati-Bernoulli Sub-ODE approach, as presented in [28]. Consider the following NPDE: where x and t are independent variables. S(x, t) plays the role of the solution of Equation (2).
• Here, we present some forms for the travelling wave solutions of Equation (4) depending on the value of n and the obtained parameters.

Travelling wave solution
This section concentrates on constructing some exact travelling wave solutions for (1 + 1)-dimensional form of the generalized Ito integro-differential equation given by Here, α is a user-specified parameter. Equation (16) is firstly simplified by assuming that U = x . Thus, Equation (16) becomes where x and t are the spatial and time variables, respectively. Substituting the transform into Equation (17) and integrating twice w.r.t. η and taking the integration constants by zero yield Here, k is a real number and w indicates the wave speed. Plugging Equation (5) and Equation (7) into Equation (19) and taking n = 0 give Equating the coefficients of ψ m to zero and solving the obtained algebraic system lead to In both cases, we observe that a 2 1 − 4a 0 a 2 = w k 3 . It is worth noting that U i,j and i,j are the exact solutions of Equation (16) and Equation (17), respectively. The index i indicates the solution number while j shows the case number. Hence, the exact solutions of Equation (17) and Equation (16) can be formed as follows.
(1) If n = 1, a 1 = 0 and a 2 = 0. Then, the solution of Equation (17) under the values of the second case is given by Hence, the exact solution of Equation (16), which we obtain by taking the first derivative for 1,2 with respect to x, is (2) If n = 1, a 0 = 0 and w k < 0. Then, the solutions are given byIn case 1: In case 2: Thus, the exact solutions of Equation (16), in both cases j = 1 and 2, are expressed as and (3) If n = 1, a 0 = 0 and w k > 0. Then, the solutions are shown asIn case 1: In case 2: where μ 0 is an arbitrary constant. As a result, the exact solutions of Equation (16), in both cases j = 1 and 2, are and

Hamiltonian system
Hamiltonian system [29,30] which is given by is employed in this article to figure out the interval in which the achieved travelling wave solutions are stable. Here, w denotes the wave speed, θ plays the role of the momentum, a, b are the endpoints of the specified domain, and (ζ ) is the relevant solution. The required condition for the stability is given by

Stability of the travelling wave solution
In this subsection, the Hamiltonian system is used to seek the stability of the obtained exact solutions of Equation (17). Applying Equation (36) on Equation (17) over the domain [0, 60] leads to Therefore, Differentiating Equation (39) with respect to w which indicates wave speed yields In the subsequent section, we consider the numerical solution of Equation (17) by employing the finite differences. The spatial derivatives in Equation (17) are semi-discretized using central finite difference formulas with second-order accurate. Temporal derivatives are kept continuous. Once we have the result, we will obtain the numerical solution of Equation (16) by using the second-order finite difference formulas of the spatial direction in , where U = x .

Numerical solution
The central finite differences are employed in this section to establish the numerical results of Equation (17) over the computational domain [x L , x R ]. Here, x L and x R indicate the endpoints of the domain. The numerical solution of Equation (16) is also manifested. The adaptive moving mesh, line and normal techniques [31][32][33] can be used to build the numerical results of the governing equation. In order to develop the numerical solution of Equation (17), we convert Equation (17) into a firstorder differential system in time. We assume that t = V x . Then, Equation (17) becomes Equation (42) can be converted into the following system: Here, x and t play the role of the spatial and time variables, respectively. To gain a much faster and more stable scheme, we modify the second equation of system (43). The considered boundary conditions, which we can deduce from Figure 3(a), are given by where x L and x R are the boundaries of the computational domain. The boundary conditions (44) lead us to have fictitious points that we use to obtain the spatial derivatives at the endpoints of the domain. We need to generate initial conditions for system (43). Hence, we consider the function V which is given by Achieving the above integral leads to Therefore, the initial conditions are given by We begin with splitting the interval [x L , x R ] into subintervals so that where h denotes the space increment. As a result, the numerical solutions at the mesh points are defined by n and V n . The spatial derivatives are semi-discretized while the temporal derivatives are kept continuous. Hence, the resulting scheme is shown as follows: where the discretization of V xxx | n and ( V x ) x | n is given by subject to the boundary conditions The initial conditions are given by evaluating It is worth mentioning that we invoke the method of lines which depends on a standard ODE solver in FORTRAN software, DASPK solver [34]. This solver applies a standard backward differentiation formulas to approximate time derivatives. DASPK uses an iterative numerical method and also uses LU factorization to approximate the Jacobian matrix of the linearized system. This leads to successful and powerful performance. To have a less bandwidth for the matrix which results from the performance of the method, we use a unique system numbering for the unknowns 1 , V 1 , 2 , V 2 , . . . , N−1 , V N−1 , N and V N .

Accuracy of the scheme
This subsection deals with the accuracy of the numerical scheme using Taylor expansion. The finite differences are used to develop the accuracy of the schemes. Plugging Taylor expansion into the numerical schemes and simplifying the results lead to the accuracy of the schemes which are O( t 2 , h 2 ). Thus, the accuracy of the schemes are from the second-order in time and space.

Stability of the scheme
This subsection is assigned to test the stability of the obtained schemes using Von Neumann's concept. We first convert system (43) to a linear system by taking c 0 = x | m n . Von Neumann's definition is given by where γ 0,1 and ξ are constants. Plugging Equations (51) into Equations (49) gives (52) Hence, system (52) can be easily expressed as follows: Then, where A = r 1 0 The values of r 1 and r 2 are where c 0 = x | m n , τ = t/h, and t is the time increment. Note that Von Neumann's condition of the stability is that the maximum eigenvalue of A must be less than or equal one. Thus, the eigenvalues of the matrix A can be obtained using Maple or Mathematica software to have It can be simply observed that the maximum value of the eigenvalues is one. Consequently, the respective numerical scheme is unconditionally stable.

Comparison and discussion
Throughout this article, we have constructed novel and more general travelling wave solutions than the exact solutions obtained in various published papers such as those written in [27]. The Kudryashov method was used in [27] to present one rational exact solution based on the exponential function. However, we have utilized the Riccati-Bernoulli Sub-ODE process to construct several solitary wave solutions. The achieved solutions in this article are expressed in the form of trigonometric and hyperbolic functions. The Riccati-Bernoulli Sub-ODE technique is more effective than the Kudryashov       method. Furthermore, the used numerical technique is reliable and effective. The agreement between the analytical and numerical results is highly the same. As can be observed in Figure 5, both exact and numerical results coincide and are identical. Figure 6 exhibits the exact travelling wave solutions for U(x, t) in Equation (16) and (x, t) in Equation (17), which are almost the same as the numerical solutions presented in Figure 7 for the same equations. Figures 1, 4, 5, 6 and 7 are plotted under the value of the parameters k = 2, μ 0 = −18, α = 0.5, w = 3. Figure 2 is sketched under the parameter values k = −10, μ 0 = 15, α = 1, w = 1, to satisfy the condition w k < 0, while Figure 3 is depicted under the parameter values k = 2, μ 0 = 15, α = 1, w = 1, to satisfy the condition w k > 0. Regarding the stability, the Hamiltonian system shows that the travelling wave solution is stable over the selected domain [0, 60] with the parameter values w = 3.0, k = 2.0, α = 0.50, μ 0 = −18.0, and t = 0 → 30. We utilize Von Neumann's concept to testify the relevant numerical scheme's stability and we find that the scheme is unconditionally stable. Taylor expansion is used to examine the accuracy of the numerical scheme. We obtain that the accuracy is from second order in time and space. That is O( t 2 , h 2 ). Further, Table 1 illustrates an essential tool for recognizing the validity of the applied numerical procedure. As can be easily observed, L 2 error reaches 1.0700e−06 in acceptable CPU time which is 1.48 × 10 −1 s. Nevertheless, the error becomes better and smaller 5.3700e−07, but the CPU time which is 3.32 × 10 −1 s increases slightly. It can be concluded that the error declines when a tiny space size is used, as shown in Figure 8.  Table 1.

Conclusion
To sum up, the Riccati-Bernoulli sub-ODE method has been successfully employed on the generalized Ito equation to establish its exact travelling wave and numerical solutions. Some exact travelling wave solutions have been expressed in forms of trigonometric and hyperbolic functions. The accomplished numerical solutions have been precisely compared with the obtained exact solutions, and we have recognized that both solutions almost coincide with each other. This can be easily observed in the presented figures. It can be pointed out that the suggested numerical technique provides reliable, successful and practical results. We have proved that the obtained exact solutions are stable on a specific interval with certain parameters. The Von Neumann stability has emphasized that the scheme is always stable. The accuracy of the schemes is the second order in time and space. The L 2 error dramatically diminishes for a tiny space size, as depicted in Figure 8. Ultimately, it can be deduced that the used techniques can be applied to more sophisticated nonlinear integrodifferential equations.