Enhanced shifted Jacobi operational matrices of integrals: spectral algorithm for solving some types of ordinary and fractional differential equations

We provide here a novel approach for solving IVPs in ODEs and MTFDEs numerically by means of a class of MSJPs. Using the SCM, we build OMs for RIs and RLFI for MSJPs as part of our process. These architectures guarantee accurate and eﬃcient numerical computations. We provide theoretical assurances for the eﬃcacy of an algorithm by establishing its convergence and error analysis features. We oﬀer ﬁve numerical examples to prove that our method is accurate and applicable. Through these examples, we demonstrate the greater accuracy and eﬃciency of our approach by comparing our results with previously published ﬁndings. Tables and graphs show that the method produces exact and approximate solutions that agree quite well with each other.


Introduction
A subfield of mathematics known as fractional calculus has recently attracted a lot of interest due to involved integrals and derivatives of noninteger order.Complex systems displaying long-term memory effects and anomalous diffusion phenomena, such as heat transport, can be effectively modeled and analyzed using this mathematical technique [1][2][3].Financing, biology, engineering, physics, and many branches of applied calculus are all included [4][5][6][7][8].
There has been a lot of research on numerical methods for solving IVPs and BVPs in ordinary differential equations and partial differential equations (e.g., [9][10][11][12][13][14][15][16][17][18][19][20][21][22][23][24]).To numerically solve various types of DEs, OMs constructed from orthogonal and nonorthogonal polynomials have been extensively used [25][26][27][28][29][30][31][32].As far as accuracy and computing efficiency are concerned, many algorithms have shown promise.Nonetheless, new avenues should be investigated to improve the numerical solutions in terms of accuracy and efficiency.Our new method for numerically solving ODEs takes the form y (n) (z) + n-1 q=0 η q y (q) (z) = f 1 (z), z ∈ [0, L], n = 1, 2, 3, . . ., (1.1) and, for solving MTFDEs, subject to the ICs y (j) (0) = α j , j = 0, 1, . . ., n -1, ( where η i , α i (i = 0, 1, . . ., n -1), β i , γ i (i = 0, 1, . . ., k), γ k+1 , and ν are constants such that n - The functions f 1 and f 2 are supposed to be continuous.First, we build OMs for RIs and RLFI for MSJPs using our technique.Then we apply the SCM.We are able to acquire very close approximations of the solutions since these architectures guarantee efficient and accurate numerical computations.Following these basic stages, the suggested method may solve ODEs (1.1) and MTFDEs (1.2) susceptible to ICs (1.3): (i) We transform equations (1.1) and (1.2) together with ICs (1.3) into an equivalent form with homogeneous conditions.(ii) Fully integrated forms of (1.1) and (1.2) are obtained by applying RIs and RLFI, respectively.This conversion allows for a more comprehensive representation of the problem.(iii) In the integrated forms of (1.1) and (1.2), the solution and all of its RIs and RLFI are approximated by using the constructed OMs to write them as linear combinations of MSJPs, followed by the application of the SCM.(iv) Using a suitable numerical method, the systems of algebraic equations obtained in (iii) can be solved, which provides the required numerical solutions.We examine the convergence properties and perform a comprehensive error analysis to prove that our suggested algorithm works.To prove that our method is accurate and efficient, we offer theoretical guarantees.We also include five numerical examples that show a variety of IVPs, from (1.1) to (1.3).By contrasting our findings with those of other researchers we demonstrate how our method is more precise and efficient.The provided graphs and tables show that the exact and approximate solutions correspond very well.
One notable aspect of our suggested approach is using MSJPs.An innovative method for solving the aforementioned IVPs is presented by making use of these polynomials and the OMs that are linked with them.Advantages of MSJPs over current approaches include better accuracy, faster convergence, and lower computing cost.
The paper is outlined as follows.The definitions and properties of RFI and RIs are provided in Sect. 2. Section 3 details several features of SJPs and MSJPs.The main emphasis is placed on developing new OMs for the RFI and RIs of MSJPs in Sect. 4. To construct the algorithms that are to be used to resolve the IVPs (1.1)-(1.3),this act is carried out.
Here in Sect.5, we lay out our approach to solving IVPs in ODEs and MTFDEs with the help of the built OMs and the SCM.The suggested method is the subject of the theoretical examination in Sect.6.As part of our investigation of the convergence characteristics of the method, we run an error analysis.Theoretical assurances regarding the precision and performance of the algorithm are deduced and examined.We provide five numerical examples in Sect.7 to verify that our method is accurate and applicable.To evaluate the suggested algorithm, we may look at these examples that span a variety of IVPs, from (1.1) to (1.3).We show that our method is more accurate and efficient in comparison with those published earlier.There is a very close match between the exact and approximate solutions, as seen in the tables and graphs.We highlight the merits, limits, and prospective enhancements of our algorithm in Sect.8, where we also summarize the key findings and offer conclusions based on our study.

Preliminaries and notation
This section introduces the key ideas and tools needed to construct the suggested approach.These ideas and technologies underpin our approach, helping us solve the challenge.In this context, the Riemann-Liouville definition of a fractional integration of order ν > 0 is defined as follows [8].

Definition 2.1
and , where m -1 ≤ ν < m, and m ∈ N is the smallest integer greater than ν.
For μ, ν ≥ 0 and γ > -1, the following properties of are satisfied: x γ +μ , (2.4) The RLFD of order ν > 0, denoted by R D ν , is defined as follows: On the other hand, the CFD of order ν, denoted by C D ν , is defined as follows: which can be written in the form The CFD satisfies the following properties: ) ) x γ -ν , (2.11) (2.12) The relation between the RLFI and CFD is given by [8, Eq. (2.4.6)] To accomplish the proposed algorithm, we must define the q times repeated integral of f (x) as follows.
Definition 2.2 Let f be a continuous function on the real line.Then the qth repeated integral of f , J q f , is defined as follows: which is known as the q-fold integral and has the form [33, Eq. (2.16)] (2.16) According to integral expressions (2.1) and (2.16), it is shown that (2.17) Lemma 2.1 x q-r+j , r ≤ q, (2.18) and for the case r = q, we have Accordingly, if f (j) (0) = 0, j = 0, 1, . . ., r -1, then J q f (r) (x) = J q-r f (x), r ≤ q.
Proof It is easy to prove this lemma using induction on r.
Remark 2.1 Riemann's modified form of Liouville's fractional integral operator is a direct generalization of Cauchy's formula for a q-fold integral.Moreover, in view of formula (2.16), we can see that formula (2.19) coincides with Taylor's formula with remainder.

An overview on SJPs and MSJPs
The main objective of this section is to present the fundamental characteristics of JPs and their shifted form.Furthermore, we will introduce a set of MSJPs.

An overview on SJPs
The orthogonal JPs, J (a,b) n (x), a, b > -1, satisfy the following relationship [34]: The SJPs, denoted as where w a,b L (z) = (Lz) a z b .The expansions that will serve as the foundation in this paper are the following fundamental ones [35,Sect. 11.3.4]: 1.The power form representation of J (a,b) L,n (z) is as follows: where 2. Alternatively, the expression for z k in relation to J (a,b) L,r (z) has the form where

Presenting MSJP
In this section, we define the polynomials {K (a,b) n,j (z)} j≥0 as follows: They are needed to satisfy the homogeneous form of the given ICs (1.3) for a suitable choice of p. Subsequently, these polynomials satisfy the orthogonality relation: 4 OM for RIs and RLFI for K (a,b) n,i (z) In this section, we prove Theorems 4.1 and 4.2, which give the qth integrals for all q ≥ 1 and fractional integrals of K n,i (z) in terms of the same polynomials. with where Consequently, J q K (a,b) n,N (z), q = 1, 2, . . ., n, have the form where and Proof The following formula can be obtained by combining integrating operations q times and relation (3.1): Now using formula (3.3), we obtain Expanding and collecting similar terms, after some algebra, we get Then substituting formulae (3.2) and (3.4) into (4.9),after some manipulation, yields (4.1), which can be expressed as follows: i,0 (n, q), P (a,b) i,1 (n, q), . . ., P (a,b) i,i+q (n, q), 0, . . ., 0 K (a,b) n,N+q (z), (4.10) and this expression leads to the proof of (4.3).
Theorem 4.2 I μ K (a,b) n,i (z), i ≥ 0, can be written in the form and, consequently, I μ K (a,b) n,N (z) has the form where ) is a matrix of order (N + 1) × (N + 1), which can be expressed explicitly as where and Proof Considering (3.1) and utilizing (2.4), we obtain By utilizing (3.3), (4.16) may be reformulated as (4.11), which can be represented as and this expression leads to the proof of (4.12).
Note 4.1 It is easy to see that J (0) n = I (0) n = I N+1 , where I N+1 is the identity matrix of size N +1, and hence (4.3) and (4.12) are satisfied for q = 0 and μ = 0, respectively.Remark 4.1 It is worth stating that the forms of P (a,b)  i,j (n, q) and F (μ) i,j (n) in Theorems 4.1 and 4.2 have a closed form for certain a, b.These include particular cases of Jacobi polynomials: the Chebyshev polynomials of the first and second kinds, Legendre polynomials, and ultraspherical polynomials.(4.13).This structure reduces the complexity of the algorithm, allowing it to handle larger problem sizes without excessive computational demands.

Numerical algorithm for solving ODE (1.1) and MTFDE (1.2) subject to ICs (1.3)
In this section, we propose a numerical solution for Eq.(1.1) when we apply the homogeneous form of ICs, specifically, when α j = 0 for all j = 0, 1, 2, . . ., n -1.With this respect, the basis K (a,b) n,i (z) is used to satisfy this homogeneous form of ICs.On the other hand, creating the suggested procedure requires converting (1.1) and (1.2) into equivalent forms with homogeneous conditions, also taking into account the nonhomogeneous conditions (1.3).
Note 5.1 There are a number of aspects to think about deciding which collocation points to use, such as the nature of the problem and the desired numerical solution properties.We should perform a comparison analysis with the numerical solutions computed to determine which of these choices is better.Applying both sets of collocation points to the investigated problems and evaluating them according to accuracy, convergence, and computing efficiency will shed light on the relative merits of the two options for the problem class in consideration.
A computer system with 3.60 GHz Intel(R) Core(TM) i9-10850 CPU, 10 cores, and 20 logical processors ran the calculations using Mathematica 13.3.The algorithmic steps for solving the ODE and MTFDE using MSJCOMIM are expressed in Algorithms 1 and 2, respectively:

Convergence and error analysis
In this section, we examine the convergence and error estimates of suggested method.The space S n,N is defined as follows: Additionally, we define the error between y(z) and its approximation y N (z) as In the paper, the error of the numerical scheme is analyzed by using the estimate of the L 2 norm error, and the estimate of the L ∞ norm error, 3) Theorem 6.1 [36] suppose that y(z) = z n u(z) and y N (z) is presented by (5.3) and represents the best possible approximation for y(z) out of S n,N .In that case, there exists a constant K such that and ) The following conclusion demonstrates that the obtained error converges at a fairly rapid rate.Corollary 6.1 For all N > s -1, we have the estimates and An estimate for error propagation is the focus of the next theorem, which stresses the stability of error.

Theorem 6.2 For two iterative approaches to y(z), we have
where indicates that there exists a generic constant d such that |y N+1y N | ≤ d(eL/4) N × N s-N-1 .Note 6.1 As eL/4 goes down, the error estimates in this section show that the rate of convergence changes from an inverse polynomial to an exponential one.

Numerical simulations
To demonstrate that the method given in Sect. 5 is effective and efficient, we give several examples.To measure precision, we display the MAE between the exact and approximate solutions.In particular, we demonstrate in Examples 7.1 and 7.4 that the suggested method MSJCOMIM produces the exact solution for problems with a polynomial solution of degree N .We also show the calculated errors for numerical solutions y N (z) obtained with MSJCOMIM for N = 1, . . ., 20.We can see the excellent computational accuracy in the findings summarized in Tables 1, 2, 4, 6, and 7.In addition, Tables 3, 5, and 7 compare MSJCOMIM with other techniques provided in [37][38][39][40][41].The results show that MSJCOMIM is the best method, giving more accurate predictions than the others.where g(z) is chosen such that y(z) = z 4 + 2. The application of the proposed method where y(z) = z 2 e -2zz 2 + 3.This solution agrees perfectly with the numerical solutions obtained with accuracy of 10 -16 when L = 1, 4 and N = 11, 20, respectively, as shown in Tables 1 and 2.
where g 1 (z) and g 2 (z) are chosen such that the exact solutions are y(z) = z 3 and y(z) = z 4 (z -1), respectively.The application of MSJCOMIM give the exact solution y(z) = y 1 (z) = z 3 in the form where g(z) is chosen such that y(z) = sin(αz).
The numerical solutions at N = 16 agree precisely with this solution and are presented in Tables 6 and 7, and their accuracy is 10 -16 .In Table 8 the comparison results of MAE using MSJCOMIM significantly outperform those of [37,41].Additionally, the CPU (seconds) of MSJCOMIM was found to be faster compared to the corresponding method [37, Table 1].

Conclusions
In this research, we have introduced a kind of shifted JPs that satisfy homogeneous ICs.We have also developed a new method for approximating the ODE and MTFDE solutions specified in Sect. 4 using the SCM in conjunction with the derived OMs.Using five different cases, MSJCOMIM has proven to be incredibly accurate and efficient in resolving these issues.Based on the promising results obtained in this research, we envision several potential directions for future work.Firstly, an interesting avenue would be to investigate the extension of MSJCOMIM to handle higher-dimensional problems, such as systems of ODEs and MTFDEs.This expansion would require the development of suitable multidimensional OMs and the adaptation of the spectral collocation framework.Furthermore, it could be valuable to explore the applicability of MSJCOMIM to other classes of FDEs beyond those considered in this study.Various types of FDEs exist in different scientific and engineering fields, and investigating their solutions using MSJCOMIM could provide valuable insights and contribute to advancing the field.Additionally, the theoretical findings presented in this paper open up possibilities for further research in the area of numerical methods for DEs.Exploring alternative modifications of shifted JPs or investigating the use of different OPs could lead to the development of even more accurate and efficient approximation techniques.In conclusion, the introduced MSJCOMIM has shown great potential in solving ODEs and FDEs with high accuracy.We believe that the knowledge and techniques presented in this work can serve as a foundation for addressing a broader range of DEs and inspire further advancements in the field of numerical methods for DEs.

Remark 4 . 2
The utilization of formula(3.3)  in relation(4.16)leads to the lower triangular structure of matrix

Table 4
Errors obtained for Example 7.3