A Biorthogonal Hermite Cubic Spline Galerkin Method for Solving Fractional Riccati Equation

: This paper is devoted to the wavelet Galerkin method to solve the Fractional Riccati equation. To this end, biorthogonal Hermite cubic Spline scaling bases and their properties are introduced, and the fractional integral is represented based on these bases as an operational matrix. Firstly, we obtain the Volterra integral equation with a weakly singular kernel corresponding to the desired equation. Then, using the operational matrix of fractional integration and the Galerkin method, the corresponding integral equation is reduced to a system of algebraic equations. Solving this system via Newton’s iterative method gives the unknown solution. The convergence analysis is investigated and shows that the convergence rate is O ( 2 − s ) . To demonstrate the efﬁciency and accuracy of the method, some numerical simulations are provided.


Introduction
One of the most important classes of nonlinear ordinary differential equations (ODEs) that plays a remarkable role in engineering, mathematics, and science is the Riccati equation. Count Riccati has studied the particular version of the Riccati equation for the first time in 1724. Since there is a close relationship between the homogeneous differential equation of the second-order and the Riccati equation, we can imagine many applications for this equation. This equation is closely related to the one-dimensional static Schrödinger equation and the solitary wave solution of nonlinear PDEs [1,2]. Furthermore, this equation also plays a vital role in modeling classical and modern dynamical systems [3,4].
Because of the importance of this type of differential equation, several analytical and numerical methods have been used to solve it. In [5], the authors used new fractional bases based on the classical Legendre wavelet. In this work, the desired equation is solved using the operational matrix for Caputo fractional derivative and applying the Tau method. Rabiei et al. [6] introduced Boubaker wavelets of the fractional-order and used Definition 1. Given β ∈ R + , let Γ(β) is the Gamma function. The Riemann-Liouville fractional integral operator I β a of order β is determined by where [a, b] is a finite interval on R.
It can be easy to directly verify that the fractional integration from the power functions is a yield power function of the same form, via It follows from [25] that the fractional integral operator I β a is bounded. To this end, we have the following Lemma. [25]). The operator I β a is bounded in L p ([a, b]), i.e., Definition 2. The Riemann-Liouville operator of the fractional derivative is defined by where α ∈ R + , [α] + 1 := n ∈ N and D n := d n dx n .

Biorthogonal Hermite Cubic Spline Scaling Bases
The biorthogonal Hermite cubic Spline scaling bases (BHCSSb) ψ 1 and ψ 2 are defined via and It follows from [15] that ψ 1 , ψ 2 ∈ C 1 (R) and fulfill Hermite interpolation where δ i,j denotes the Kronecker delta. Assume that the subspace V s ⊂ L 2 (R) is spanned by where s ∈ Z + ∪ {0}, B := {1, . . . , 2 s − 1} and ψ k s,b := ψ k (2 s . − b). Motivated by the multiresolution properties [27], we know that these spaces are nested V s ⊂ V s+1 . Thus, considering ψ = (ψ 1 , ψ 2 ) as a vector function of the scaling function, it is easy to show that this vector satisfies the matrix refinement equation via, in which and ∀b ∈ {−1, 0, 1}, H b = O (O is the zero matrix). The vector function Ψ satisfies the following symmetry properties where Due to this relation, one can say that ψ 1 is symmetric and ψ 2 is antisymmetric. Using (13) and (14), we can write Because the Hermite cubic spline multiwavelet system is biorthogonal, there exists a dual multi-generatorψ = (ψ 1 ,ψ 2 ) that satisfies the biorthogonality condition, i.e., where I 2 is the identity matrix of size two and ., . denotes the L 2 -inner product. This dual multi-generator generates another multiresolution spaceṼ s ⊂ L 2 (R), which is biorthogonal to V s . In order to construct the dual scaling functionsψ 1 ,ψ 2 , we utilize the refinement relation for primal and dual scaling functions and insert them into the biorthogonal relation (16). This gives rise to the discrete duality relation [15] ∑ l∈Z H lH T l+2b = 2δ 0,b I 2 , b ∈ Z. (17) In which the refinement maskH is chosen to bẽ , By reindexing the scaling functions via the set ϕ = {ϕ 1 , ϕ 2 , . . . , ϕ 2 s+1 }, whose elements are equal to 1] . Now, we introduce the operator P s that is based on multi-scaling functions, which allows us to approximate any function u ∈ L 2 (R) as follows where the coefficients u l = u, ϕ l for l = 1, . . . , 2 s+1 are computed by using the Hermit type interpolation property of BHCSSb, Now, for additional simplification, assume that Ψ s is a vector function of dimension 2 s+1 whose ith element is ϕ i (x). Similarly, the vector U is chosen to be a vector of the same dimension of Ψ s for which the ith element is u i . According to this introduction, (18) can be rewritten via It follows from Theorem 2 in [14] that the error of approximation (18) can be bounded via the following theorem.

Representation of Fractional Integral Operator in BHCSSb
The fractional integration of the vector function Ψ s (x) can be expressed by where I β is the Riemann-Liouville fractional integral operational matrix of dimension To find the elements of matrix I β , we continue the following process. Given β ∈ R + , the Riemann-Liouville fractional integral operator I β 0 , acting on ψ k (2 s x − b) for k = 1, 2, can be represented by To evaluate this integral, we check out the four cases due to the support of

4.
If x ≥ b+1 2 s then for k = 1, 2, we get The above integrals can be evaluated explicitly in terms of β, s, b for all values of b ∈ B for given s ∈ R + . We use a library function "int" available in Maple to evaluate the above integrals analytically. Thus, using the above-obtained integrals, the Riemann-Liouville It follows from (22) that the fractional integration of vector function Ψ s (x) takes the form where Γ(x) is a vector function whose elements are obtained via and Now, we can find the entries of matrix I β through expanding each of the components of the vector function Γ(x) by Biorthogonal Hermite cubic spline multi-scaling functions [14] as where M = 6 Γ(β + 4) 2 −sβ (β + 1), and the 2 × 2 block matrices Finally, we introduce the matrices ∆ i , i = 1, . . . , 2 s − 1 as follows.

Wavelet Galerkin Method
In the present section, we utilize the wavelet Galerkin method based on BHCSSb to solve the Riccati Equation (1). To derive the approximate solution, we suppose that the unknown solution can be approximated by where U is a vector of dimension N that should be determined. Assume that β ∈ R + , n = −[−α], and f , g, and h are continuous functions. Then, it is easy to show that the function u(x) is a solution of the Riccati Equation (1), if, and only if, it satisfies the integral equation Using (31), we can approximate all terms in the right side of (1) as follows Inserting Equation (33) into (32) and using the operational matrix of fractional integration I β , we have the residual as follows . We would like to reduce the residual to zero. There are several methods to do this. However, in this work, we use the wavelet Galerkin method. The biorthogonality of BHCSSb ( Ψ s ,Ψ s ) yields the linear or nonlinear system where F is a vector function of U. This function may be linear or nonlinear, and it depends on the function h. To find the unknown vector U , we utilize Newton's method in the nonlinear type and the generalized minimal residual method (GMRES method) [28] in the linear type. Proof. If f (x) is a continuous function, we can directly find the following error via Theorem 1

Convergence Analysis
Since the functions g and h are continuous, then there exist a constant C 2 such that max{g, h} ≤ C 2 . It follows from Lemma 1 that I β 0 (gu) − I β 0 P s (gu) ≤ and where we have Taking the norm from both sides of (40) and using the triangle inequality, it follows from Theorem 1 that

Numerical Experiments
Example 1. Consider the fractional Riccati equation subject to the initial condition u(0) = 0. The exact solution is reported in [6] and is u(x) = x β+1 Γ(β+2) . Table 1 shows a comparison between our proposed method and the Bernoulli wavelets method [29]. We observe that the wavelet Galerkin method based on BHCSSb gives better results than the Bernoulli wavelets method. To illustrate the effect of refinement level s on L 2 -errors, Table 2 is reported. It is worth emphasizing that these results verify our convergence analysis, and by increasing this parameter, the L 2 -errors decrease. To show the accuracy of the method, Figure 1 is plotted. In this figure, we can see a compare between the exact and approximate solutions. Figure 2 demonstrates the approximate solutions for different values of β on the left side and corresponding absolute errors on the right.

Example 2.
The second example is dedicated to the fractional Riccati equation subject to the initial condition u(0) = 0. There is no exact solution to the problem here. However, in the case of β = 1, the exact solution would form u(x) = e 2x −1 e 2x +1 [5,6]. Figure 3 displays the approximate solution for different values of β. As we expect, when β = 1, the corresponding solutions tend to the solution at it. Table 3 shows a comparison of the proposed method and the fractional-order Legendre wavelet method [5]. Table 3. Comparison of the absolute value of residual between the proposed method and fractionalorder Legendre wavelet method [5] for Example 2.
subject to the initial condition u(0) = 0. The exact solution is reported in [6] and is u(x) = x 2 . Figure 4 illustrates a comparison between the exact and approximate solution. The absolute errors are reported in Table 4.

Conclusions
In this paper, we applied the wavelet Galerkin method to solve the fractional Riccati equation. To this end, we utilized the Biorthogonal cubic Hermite spline multiwavelets and the operational matrix for fractional integration to reduce the desired equation to a set of nonlinear algebraic systems. The convergence analysis is investigated and shows that the convergence rate is O(2 −s ). Some numerical simulations and results demonstrate the ability and efficiency of the method.

Data Availability Statement:
The data presented in this study is available on request from the corresponding author.

Conflicts of Interest:
The writers state that they have no known personal relationships or competing financial interest that could have appeared to affect the work reported in this work.

Abbreviations
The following abbreviations are used in this manuscript:

R
The real numbers R + The positive real number N The natural numbers Z + The positive integers C The space of continuous functions C n The space of functions which are n times continuously differentiable L p The spaces of p-integrable functions ODE Ordinary differential equations BHCSSb Biorthogonal Hermite cubic Spline scaling bases