Next Article in Journal
Analytical Solution and Numerical Simulation of Heat Transfer in Cylindrical- and Spherical-Shaped Bodies
Next Article in Special Issue
Spillover Effects of Green Finance on Attaining Sustainable Development: Spatial Durbin Model
Previous Article in Journal
Computation of the Exact Forms of Waves for a Set of Differential Equations Associated with the SEIR Model of Epidemics
Previous Article in Special Issue
Informal Sector, ICT Dynamics, and the Sovereign Cost of Debt: A Machine Learning Approach
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Applications of Modified Bessel Polynomials to Solve a Nonlinear Chaotic Fractional-Order System in the Financial Market: Domain-Splitting Collocation Techniques

by
Mohammad Izadi
1 and
Hari Mohan Srivastava
2,3,4,5,*
1
Department of Applied Mathematics, Faculty of Mathematics and Computer, Shahid Bahonar University of Kerman, Kerman 76169-14111, Iran
2
Department of Mathematics and Statistics, University of Victoria, Victoria, BC V8W 3R4, Canada
3
Department of Medical Research, China Medical University Hospital, China Medical University, Taichung 40402, Taiwan
4
Department of Mathematics and Informatics, Azerbaijan University, 71 Jeyhun Hajibeyli Street, AZ1007 Baku, Azerbaijan
5
Center for Converging Humanities, Kyung Hee University, 26 Kyungheedae-ro, Dongdaemun-gu, Seoul 02447, Republic of Korea
*
Author to whom correspondence should be addressed.
Computation 2023, 11(7), 130; https://doi.org/10.3390/computation11070130
Submission received: 29 May 2023 / Revised: 27 June 2023 / Accepted: 30 June 2023 / Published: 3 July 2023
(This article belongs to the Special Issue Quantitative Finance and Risk Management Research)

Abstract

:
We propose two accurate and efficient spectral collocation techniques based on a (novel) domain-splitting strategy to handle a nonlinear fractional system consisting of three ODEs arising in financial modeling and with chaotic behavior. One of the major numerical difficulties in designing traditional spectral methods is in the handling of model problems on a long computational domain, which usually yields to loss of accuracy. One remedy is to split the underlying domain and apply the spectral method locally in each subdomain rather than on the global domain of interest. To treat the chaotic financial system numerically, we use the generalized version of modified Bessel polynomials (GMBPs) in the collocation matrix approaches along with the domain-splitting strategy. Whereas the first matrix collocation scheme is directly applied to the financial model problem, the second one is a combination of the quasilinearization method and the direct first numerical matrix method. In the former approach, we arrive at nonlinear algebraic matrix equations while the resulting systems are linear in the latter method and can be solved more efficiently. A convergence theorem related to GMBPs is proved and an upper bound for the error is derived. Several simulation outcomes are provided to show the utility and applicability of the presented matrix collocation procedures.

1. Introduction

There has been a great deal of interest in working with nonlinear ordinary differential equations (ODEs) and dynamical systems over the past decades. Diverse physical as well as nonphysical phenomena including economical, epidemiological, and biological models belong to these category. Among others, nonlinear systems with chaotic behavior are of particularly of interest. Chaos theory of nonlinear dynamics is a branch of mathematical science that deals with highly complex systems, where small changes now in their inputs cause large changes in their outputs later [1]. The notion of chaos was first introduced to the modern world by Edward Lorenz in 1972 with conceptualization of “Butterfly Effect”. By understanding this theory we are able to predict the behavioral representation of a complex system more accurately. Chaotic systems are unstable since they tend not to resist any outside disturbances but instead react in significant manner. As a meteorologist, Lorenz developed a mathematical model for the air movement in the atmosphere. It caused various differences in the outcome of the given model. In this manner he discovered the principle of sensitive dependence on initial conditions, which is now considered as a key element in any chaotic system [2]. Indeed, systems with chaotic behavior are common in nature and diverse real-life phenomena can also be classified in the framework of a chaotic system. Beside of their appearance in meteorology, they can be found also in solar system, chemical reactions, climate dynamic, population dynamic, economics, and so on.
In the economic setting, Ma and Chen [3,4] were one of the pioneering who studied and analyzed the chaotic economic systems. The study of memory effects and fractional-order derivative on a chaotic financial system was given by Chen [5]. The concept of fractional calculus provides more insight in predicting the behavior of complex dynamical systems especially those with chaotic characteristics, see [6]. Various fractional order systems with chaos have been considered in literature such as Chua [7], Chen [8], Rössler [9], Lorentz [10] to name a few. Finding the exact or analytical solutions of the most fractional-order systems cannot be possible. Alternately, the researchers have been tried to develop effective semi-analytical as well as numerical procedures to deal with such systems of fractional order. Among diverse existing proposed methods, we refer to the recent published works including the spectral collocation approaches [11,12,13], the fractional Adams-Bashforth approach [14], the residual power series method and Laplace technique [15], the variational iteration scheme [16], and the wavelet technique [17] to name a few.
The main objective of the current research is to investigate a real-life model arising in the financial market. The model equation consists of three differential equations with quadratic nonlinearities given by [5]
L C D t ψ u ( t ) = λ 1 u ( t ) + u ( t ) v ( t ) + w ( t ) , L C D t ψ v ( t ) = u 2 ( t ) λ 2 v ( t ) + 1 , L C D t ψ w ( t ) = u ( t ) λ 3 w ( t ) ,
where t > 0 and λ i for i = 1 , 2 , 3 are constants. Here, λ 1 denotes the saving amount rate, λ 2 represents the costs per investment, and λ 3 shows the elasticity of demand of commercial markets. In the model (1) by L C D t ψ we denote the fractional derivative operator of order ψ in the sense of Liouville–Caputo (LC). Let u 0 , v 0 , and w 0 be given real numbers. The supplemented initial conditions are as follow
u ( 0 ) = u 0 , v ( 0 ) = v 0 , w ( 0 ) = w 0 .
It should be emphasized that the dynamic analysis and chaotic intrinsic nature of the financial system (1) have already been investigated by researchers, see cf. [5].
Based on our knowledge, only a few works have been proposed to solve the system (1) numerically. The Adams-Bashforth-Moulton scheme used in [18] to deal with the variable-order fractional counterpart of (1). The integer-order of the above model solved by a robust spectral Chebyshev method [19]. The fractional-order of model (1) in the sense of both Caputo and Atangana-Baleanu-Caputo (ABC) derivatives considered in [20] and solved by the modified homotopy analysis transform method. They also proved the existence and uniqueness of the solutions based on Picard-Lindelof’s approach. The variable-order of the underlying model in the sense of ABC derivative invsetiagted in [21]. A computational stochastic approach based on the artificial neural networks introduced recenlty in [22]. Furthermore, the authors in [23] proposed a new robust nonlinear controller that stabilizes the chaotic finance system (1). Some other financial models related to (1) with some proposed approximate methods were given in [24,25,26,27].
In this work, we solve the nonlinear system (1) with two spectral collocation approaches based on (novel) modified Bessel functions. To preserve the high accuracy of the proposed methods on a long computational domain, we further employ the domain-splitting methodology. In the first and direct collocation procedure, we convert the financial model into an algebraic system of nonlinear equations. On the other hand, in the next and more efficient algorithm, we utilize first the methodology of quasilinearization to (1) and then employ the former direct technique to the resulting sequence of the linearized system of equations. In the latter technique we arrive at a family of linear matrix equations to be solved easily. Using the domain-splitting methodology enables us to keep the accuracy high on a long time domain. For this purpose, it is sufficient to divide the domain into several subdomains and then employ the traditional collocation schemes on each subintervals more accurately.
The balance of the current work is given next. A review of fundamental facts on fractional integral and derivatives is provided in Section 2. Section 3 contains of introducing the fundamental definitions and main properties of the modified Bessel functions. The generalized form of these polynomials is introduced next. The error estimate of the given basis functions is also established. In Section 4, we first develop the direct GMBPs collocation based on domain-splitting strategy. Afterward, the hybrid QLM-GMBPs procedure relied on the quasilinearization process is implemented in detail. We also testify the accuracy of two proposed collocation technique via defining the technique of residual error functions in this section. Experiments and computational results are carried out in Section 5. A concise conclusion of the current manuscript is discussed in Section 6.

2. Basic Preliminaries on Fractional Calculus

Here and in this section, our aim is to provide some facts associated to fractional calculus. The definition of Liouville–Caputo (LC) is first given and some of its properties will be mentioned. At the end, we provide a theorem for the generalized Taylor’s formula with LC operator.
Definition 1
([6]). Let suppose that e ( t ) is a continuously differentiable function r-times. The LC fractional operator L C D t ψ of order ψ > 0 of e ( t ) is defined as follows
L C D t ψ ( e ( t ) ) = 0 I t r ψ ( e ( r ) ( t ) ) , i f r 1 < ψ < r , e ( r ) ( t ) , i f ψ = r , r N ,
where
0 I t ψ ( e ( t ) ) = 1 Γ ( ψ ) 0 t e ( z ) ( t z ) 1 ψ d z , t > 0 .
Let c 1 and c 2 be two constants. The LC operator has linear property in the sense that
L C D t ψ ( c 1 e 1 ( t ) + c 2 e 2 ( t ) ) = c 1 L C D t ψ ( e 1 ( t ) ) + c 2 L C D t ψ ( e 2 ( t ) ) ,
where two functions e i ( t ) for i = 1 , 2 are r-times continuously differentiable. In addition to the above linear property, we utilize the next properties of LC operator will be used frequently throughout the paper as follow
L C D t ψ ( c 1 ) = 0 ,
L C D t ψ t ω = Γ ( 1 + ω ) Γ ( 1 + ω ψ ) t ω ψ , i f ω N 0 a n d ω ψ , o r ω N 0 a n d ω > ψ , 0 , i f ω N 0 a n d ω < ψ ,
where N 0 : = N { 0 } , ψ and ψ denote the floor and ceiling functions of ψ respectively. The next Theorem states the generalized form of Taylor series expansion when we have the LC fractional derivative operator. To see a proof, we refer the readers to [28]. Towards this end, let us have the convention that
L C D t r ψ = L C D t ψ L C D t ψ L C D t ψ r times .
Theorem 1.
Let 0 < ψ 1 and the functions L C D t r ψ e ( t ) are continuous on ( 0 , T ] for r = 0 , 1 , , p . Then, the series representation of e ( t ) is written as follows
e ( t ) = r = 0 p t r ψ Γ ( 1 + r ψ ) L C D t r ψ e ( 0 + ) + t ψ ( p + 1 ) Γ ( 1 + ψ ( p + 1 ) ) L C D t ( p + 1 ) ψ e ( ξ ) , t [ 0 , T ] ,
for some ξ ( 0 , T ) .
The immediate consequence of the former Theorem is given next.
Corollary 1
(Fractional Taylor inequality). Let assume that the hypotheses of Theorem 1 hold and we have C max : = sup t [ 0 , T ] | L C D t ψ ( p + 1 ) | . Then, we get the upper bound for remainder of Taylor’s formula as
e ( t ) r = 0 p t r ψ Γ ( 1 + r ψ ) L C D t r ψ e ( 0 + ) t ψ ( p + 1 ) Γ ( 1 + ψ ( p + 1 ) ) C max , t [ 0 , T ] .

3. The Modified Bessel Polynomials: An Overview

Let us begin first by introducing the modified Bessel polynomials (MBPs) were studied in [29] to define a set of totally positive polynomials. The original Bessel polynomials were first given by Krall and Frank [30] systematically. Due to vast applications, this set of polynomials has been recently used to solve various model problems in science and applied engineering ranging from integer to fractional models. see [31,32,33,34,35,36]. See also [37] for a recent introductory overview of these polynomials.

3.1. The Basic Facts on MBPs

Let assume that B p ( s ) denotes the Bessel polynomial of order p. The modified Bessel polynomial MB p ( s ) of this order is defined by
MB p ( s ) : = p ! ( 2 p ) ! B p ( 2 s ) , p 0 .
In the explicit format, the MBPs can be written as
MB p ( s ) = p ! ( 2 p ) ! j = 0 p ( p + j ) ! j ! ( p j ) ! s j , s 0 .
By using p = 0 , 1 , one can easily check that MB 0 ( s ) = 1 and MB 1 ( s ) = 1 2 + s . A simple calculation is also indicated that MB 2 ( s ) = 1 12 + 1 2 s + s 2 . The next recursive relation will be used to get the remaining MBPs for p = 2 , 3 , as
MB p + 1 ( s ) = s MB p ( s ) + 1 4 ( 4 p 2 1 ) MB p 1 ( s ) .
Therefore, the next a few MBPs are obtained via the former formula given as
MB 3 ( s ) = 1 120 + 1 10 s + 1 2 s 2 + s 3 , MB 4 ( s ) = 1 1680 + 1 84 s + 3 28 s 2 + 1 2 s 3 + s 4 , MB 5 ( s ) = 1 30240 + 1 1008 s + 1 72 s 2 + 1 9 s 3 + 1 2 s 4 + s 5 .
It is evident the leading coefficients are all equal to one and MB p ( 0 ) = p ! / ( 2 p ) ! .
The additional property of these polynomials is that they are polynomial solutions of the next second-order differential equations
s 2 MB p ( s ) + ( 1 + 2 s ) MB p ( s ) = p ( p + 1 ) MB p ( s ) , p N 0 .
The another interesting characteristic of these MBPs is the orthogonality relation. This condition can be deduced on the circle | s | = 1 / 2 . The orthogonality condition was shown with regard to weight function w ( s ) e 1 / s . Thus, we have
1 2 π i | s | = 1 2 MB p ( s ) MB q ( s ) w ( s ) d s = 0 , p q , ( 1 ) p 2 p + 1 p ! ( 2 p ) ! 2 , p = q .
We lastly emphasize that the locations of roots of the MBPs MB p ( s ) are inside the circle | s | = 1 / 2 . The exception is MB 1 ( s ) , which is zero at s = 1 / 2 .
In practical investigations, we shall utilize the MBPs on the interval [ 0 , T ] . Consequently, the fractional form or the generalized MBPs are taken into consideration and will be defined next:
Definition 2.
The generalized MBPs (GMBPs) on [ 0 , T ] signify by MB p δ ( s ) is defined via relation
MB p δ ( s ) = MB p s δ T , s [ 0 , T ] ,
where 0 < δ 1 is a real number.
Note, with δ = 1 , we recover the integer-order form of these polynomials on [ 0 , T ] . From this new transformation, we can represent the series form (7)
MB p δ ( s ) = j = 0 p m p , j s δ j , m p , j : = p ! ( 2 p ) ! ( p + j ) ! T j j ! ( p j ) ! ,
for p = 0 , 1 , With the above given change of variable, the orthogonality condition (10) will be modified accordingly. The orthogonality of the set { MB p δ ( s ) } p = 0 is derived with respect to the weight function w δ ( s ) = s δ 1 e T / s δ . We have indeed the next formula
1 2 π i | s | = 1 2 MB p δ ( s ) MB q δ ( s ) w δ ( s ) d s = T ( 1 ) p ( 2 p + 1 ) δ p ! ( 2 p ) ! 2 δ p q ,
where δ q p stands for the Kronecker delta function.

3.2. Function Approximation: Convergence and Error Analysis of GMBPs

It is known that a function e ( t ) belongs to the weighted L 2 space on [ 0 , T ] can be represented as a summation of GMBPs. It follows that
e ( t ) = p = 0 ξ p MB p δ ( t ) , t [ 0 , T ] .
To make known the unknown coefficients ξ p , p 0 , one can utilize the orthogonality property (13). However, to deal with our model problem (1) in practical situations, we need to cut the above series form into a finite series, say with ( P + 1 ) basis functions, where P N . Therefore, we can approximate e ( t ) in the form
e ( t ) e P , δ ( t ) : = p = 0 P ξ p MB p δ ( t ) = Λ δ ( t ) Z P , t [ 0 , T ] ,
where we have used the following vectors to rewrite it in a concise representation format
Λ δ ( t ) : = MB 0 δ ( t ) MB 1 δ ( t ) MB P δ ( t ) , Z P : = ξ 0 ξ 1 ξ P T .
So, we restrict ourselves to a finite-dimensional subspace W P L w δ 2 [ 0 , T ] defined by
W P : = Span MB 0 δ ( t ) , MB 1 δ ( t ) , , MB P δ ( t ) .
We now investigate the norm of the error term between e ( t ) and its approximation e P , δ ( t ) defined by the formula
E P ( t ) : = e ( t ) e P , δ ( t ) .
The following results indicates that the norm E P 2 , w δ converges to zero if one tends P to infinity.
Theorem 2.
Let 0 < δ 1 and assume that L C D t r δ ( t ) , ( r = 0 , 1 , , P ) belong to the space of continuous funtions over ( 0 , T ] . If e P , δ ( t ) = Λ δ ( t ) Z P be the closet possible approximation to e ( t ) in the space W P , then, e P , δ ( t ) converges to e ( t ) as P , i.e.,
E P 2 , w δ 0 .
Proof. 
We first set
A P ( t ) = e ( 0 + ) + t δ Γ ( δ + 1 ) L C D t δ e ( 0 + ) + + t P δ Γ ( P δ + 1 ) L C D t P δ e ( 0 + ) .
Since 0 < δ 1 , we can apply Theorem 1 to e ( t ) and its generalized Taylor expansion series A P ( t ) of e ( t ) . The immediate consequence of Corollary 1 is that
| A P ( t ) e ( t ) | t ( P + 1 ) δ Γ ( ( P + 1 ) δ + 1 ) C max , 0 < t < T .
We have e P , δ ( t ) W P and we know that e P , δ ( t ) is the finest possible approximation to e ( t ) . The conclusion is
e ( t ) e P , δ ( t ) 2 , w δ e ( t ) h ( t ) 2 , w δ , h W P .
If utilize h ( t ) = A P ( t ) W P in the preceding relation, we find that
e ( t ) e P , δ ( t ) 2 , w δ 2 e ( t ) A P ( t ) 2 , w δ 2 = 0 T | e ( s ) A P ( s ) | 2 w δ ( s ) d s ,
where w δ ( s ) = s δ 1 e T / s δ . The definition of the error term followed by using (16) we get
E P 2 , w 2 C max Γ ( ( P + 1 ) δ + 1 ) 2 0 T s δ ( 2 P + 3 ) 1 e T / s δ d s .
We then use the inequality e T / s δ < 1 / e holds for all s > 0 . By a simple integration, we find that
E P 2 , w 2 C max Γ ( ( P + 1 ) δ + 1 ) 2 T δ ( 2 P + 3 ) δ ( 2 P + 3 ) e .
Our proof is established by performing the square roots and tending P to infinity. □

4. Domain-Splitting QLM-GMBPs/GMBPs Collocation Techniques

As previously mentioned in the Section 1, the spectral collocation procedure may be failed to converge if the computational domain of interest is large. To overcome this inaccuracy, the idea is to split the given (large) domain into a sequence of sub-domains and then apply the spectral approach in each subdomain. To this end, we consider a (uniform) partitioning of the interval [ 0 , T ] into M 1 subintervals
I m : = [ t m 1 , t m ] , m = 1 , 2 , , M ,
where t 0 : = 0 and t M : = T . The length of each I m will denote by | I m | and is equal to = T / M . When M = 1 , we recover the traditional spectral collocation strategy, see cf. [31,32,33,34,35,36]. In each subdomain I m we now consider the approximate solution to be in the form (15) as
e P , δ m ( t ) : = p = 0 P ζ p m MB p δ ( t ) = Λ δ ( t ) Z P m , t I m ,
where Z P m : = ζ 0 m ζ 1 m ζ P m T is the vector of unknown coefficients and Λ δ ( t ) is defined in (15). After obtaining the approximate solutions for m = 0 , 1 , , M 1 , the desired approximation on the whole domain [ 0 , T ] can be expressed as
e P , δ ( t ) = m = 0 M χ m ( t ) e P , δ m ( t ) ,
where the characteristic function χ m ( t ) is given by
χ m ( t ) : = 1 , t I m , 0 , otherwise .
Furthermore, at each subinterval I m , we utilize the followinng collocation nodes
t m , k : = t m 1 + t m t m 1 P k , k = 0 , 1 , , P .
Finally, we emphasize that the initial conditions on the first subinterval I 0 are taken as the original initial conditions for our model problem (1) whilst on the remaining subintervals I m , m 1 we will use the approximate solutions obtained on the previous subintervals I m 1 . The fundamental idea of the spectral collocation for the model (1) will be illustrated next. Towards this end, we only consider the proposed approach on an arbitrary domain I m for m = 0 , 1 , , M 1 .

4.1. The GMBPs Collocation Approach

Let us first represent the financial system (1) into a matrix format on the subinterval I m . It can be expressed as
L C D t ψ y ( t ) + π y ( t ) + π 1 n 1 ( t ) + π 2 n 2 ( t ) = g ( t ) , t I m .
Here, we have utilized the following notations
y ( t ) = u ( t ) v ( t ) w ( t ) , π = λ 1 0 0 0 λ 2 0 1 0 λ 3 , π 1 = 1 0 0 , π 2 = 0 1 0 , g = π 2 .
Furthermore, the nonlinear terms n i ( t ) for i = 1 , 2 are
n 1 ( t ) = [ u ( t ) · v ( t ) ] , n 2 ( t ) = [ u 2 ( t ) ] .
Now, the task is that to state the unknown solutions u ( t ) , v ( t ) , and w ( t ) as a finite combination of cutted series form (17) using ( P + 1 ) GMBPs bases. It follows that
u ( t ) u P , δ m ( t ) = p = 0 P ζ p , 1 m MB p δ ( t ) = Λ δ ( t ) Z P , 1 m , v ( t ) v P , δ m ( t ) = p = 0 P ζ p , 2 m MB p δ ( t ) = Λ δ ( t ) Z P , 2 m , w ( t ) w P , δ m ( t ) = p = 0 P ζ p , 3 m MB p δ ( t ) = Λ δ ( t ) Z P , 3 m , t I m .
To find the coefficients { ζ p , q m } p = 0 P for q = 1 , 2 , 3 , we apply the spectral matrix collocation procedure by using the GMBPs. A decomposition of the vector Λ δ ( t ) is provided next:
Lemma 1.
We can write the vector Λ δ ( t ) as the following
Λ δ ( t ) = X P δ ( t ) E P ,
where
X P δ ( t ) = 1 t t 2 δ t P δ ,
represents the vector of standard bases and the second matrix E P is a non-singular matrix whose components are given in (12). It has a triangular structure and det ( E P ) = 1 given by
E P = 1 m 1 , 0 m 2 , 0 m P 1 , 0 m P , 0 0 1 m 2 , 1 m P 1 , 1 m P , 1 0 0 1 m P 1 , 2 m P , 2 0 0 0 1 m P , P 1 0 0 0 0 1 .
Proof. 
By expanding the relation (12) in terms of powers of t followed by using the definitions of X P δ ( t ) and E P we conclude the result easily. □
By inserting the expression Λ δ ( t ) in (21) into (20) we render
u P , δ m ( t ) = Λ δ ( t ) Z P , 1 m = X P δ ( t ) E P Z P , 1 m , v P , δ m ( t ) = Λ δ ( t ) Z P , 2 m = X P δ ( t ) E P Z P , 2 m , w P , δ m ( t ) = Λ δ ( t ) Z P , 3 m = X P δ ( t ) E P Z P , 3 m . t I m .
The next job is to compute the ψ -derivative of the aforementioned approximate solutions. It is evident that we need to calculate only the ψ -derivative of X P δ ( t ) . Towards this end, we exploit the property (4) as well property (5) to calculate L C D t ψ X P δ ( t ) , which is denoted by X P ψ , δ ( t ) . Note that we can calculate it algorithmically with linear complexity O ( P ) , see Algorithm 3.1 in [38]. Therefore, we get
L C D t ψ u P , δ m ( t ) = X P ψ , δ ( t ) E P Z P , 1 m , L C D t ψ v P , δ m ( t ) = X P ψ , δ ( t ) E P Z P , 2 m , L C D t ψ w P , δ m ( t ) = X P ψ , δ ( t ) E P Z P , 3 m . t I m .
In the vector format, we can approximate the true solutions of system (19) in the form
y ( t ) y P ( t ) : = u P , δ m ( t ) v P , δ m ( t ) w P , δ m ( t ) , L C D t ψ y ( t ) L C D t ψ y P ( t ) : = L C D t ψ u P , δ m ( t ) L C D t ψ v P , δ m ( t ) L C D t ψ w P , δ m ( t ) .
By introducing the following matrix notations
X ˜ ( t ) = X P δ ( t ) 0 0 0 X P δ ( t ) 0 0 0 X P δ ( t ) , E ˜ = E P 0 0 0 E P 0 0 0 E P , W ˜ ( t ) = X P ψ , δ ( t ) 0 0 0 X P ψ , δ ( t ) 0 0 0 X P ψ , δ ( t ) ,
we are able to rewrite y P ( t ) and L C D t ψ y P ( t ) in the forms of matrices. It follows that
y P ( t ) = X ˜ ( t ) E ˜ Z m , L C D t ψ y P ( t ) = W ˜ ( t ) E ˜ Z m ,
where we have used the notation Z m = Z P , 1 m Z P , 2 m Z P , 3 m T .
We now insert the collocation nodes (18) of the subinterval I m into the matrix represenation form (19). The next set of equations are obtained
L C D t ψ y ( t m , k ) + π y ( t m , k ) + π 1 n 1 ( t m , k ) + π 2 n 2 ( t m , k ) = g ( t m , k ) , k = 0 , 1 , , P .
If we introduce the next notations
Y m ( ψ ) = L C D t ψ y P ( t m , 0 ) L C D t ψ y P ( t m , 1 ) L C D t ψ y P ( t m , P ) , Y m = y P ( t m , 0 ) y P ( t m , 1 ) y P ( t m , P ) , Π = π 0 0 0 π 0 0 0 π , g m = g ( t m , 0 ) g ( t m , 1 ) g ( t m , P ) , n i , m = n i ( t m , 0 ) n i ( t m , 1 ) n i ( t m , P ) , Π i = π i 0 0 0 π i 0 0 0 π i , i = 1 , 2 ,
we are able to represent the system of equations in a concise form as
Y m ( ψ ) + Π Y m + Π 1 N 1 , m + Π 2 N 2 , m = G m .
Our goal is to write the involved vectors in (27) in a matrix format.
Lemma 2.
In matrix representation form, two vectors y P ( t ) and L C D t ψ y P ( t ) in (25) at the collocation nodes (18) can be written as
Y m = X ˜ ¯ E ˜ Z m , Y m ( ψ ) = W ˜ ¯ E ˜ Z m .
Here, two matrices X ˜ ¯ and W ˜ ¯ are defined as
X ˜ ¯ = [ X ˜ ( t m , 0 ) X ˜ ( t m , 1 ) X ˜ ( t m , P ) ] T , W ˜ ¯ = [ W ˜ ( t m , 0 ) W ˜ ( t m , 1 ) W ˜ ( t m , P ) ] T ,
where three matrices X ˜ , E ˜ , W ˜ as well as the vector Z m are introduced in (25).
The proof of the foregoing Lemma can be done easily by just considering relations (25) followed by inserting the nodes (18) into them. Similarly, two nonlinear terms N i , m for i = 1 , 2 still need to be written in the matrix expressions. We have
Lemma 3.
We can express N i , m in (27) for i = 1 , 2 in the forms
N 1 , m = U ^ 1 V ^ 1 , N 2 , m = U ^ 1 U ^ 2 .
Here, we have the representations U ^ 1 = X ^ E ^ Z ^ 1 m , and V ^ 1 = X ¯ ^ E ¯ 1 Z m so that
X ^ = X P δ ( t m , 0 ) 0 0 0 X P δ ( t m , 1 ) 0 0 0 X P δ ( t m , P ) , E ^ = E P 0 0 0 E P 0 0 0 E P , Z ^ 1 m = Z P , 1 m 0 0 0 Z P , 1 m 0 0 0 Z P , 1 m , X ¯ ^ = X P δ ( t m , 0 ) X P δ ( t m , 1 ) X P δ ( t m , P ) , E ¯ 1 = R E P R , R = 0 0 0 0 0 0 0 0 0 ( P + 1 ) × ( P + 1 ) , Z m = Z P , 1 m Z P , 2 m Z P , 3 m .
Also, we have U ^ 2 = X ¯ ^ E ¯ 2 Z m , where E ¯ 2 = E P R R .
Proof. 
To prove the first assertion, let us note the following relation
N 1 , m = u P , δ m ( t m , 0 ) · v P , δ m ( t m , 0 ) u P , δ m ( t m , 1 ) · v P , δ m ( t m , 1 ) u P , δ m ( t m , P ) · v P , δ m ( t m , P ) = u P , δ m ( t m , 0 ) 0 0 0 u P , δ m ( t m , 1 ) 0 0 0 u P , δ m ( t m , P ) v P , δ m ( t m , 0 ) v P , δ m ( t m , 1 ) v P , δ m ( t m , P ) .
In accordance to relations (22), it is no difficult to show that the first diagonal matrix is written as U ^ 1 and the second column matrix can be expressed as V ^ 1 . Similar arguments hold for N 2 , m , where the second column vector is terms of u P , δ m instead of v P , δ m . □
We now turn our attention to (27) and insert all obtained matrix forms into it. To be more precise, we replace Y m ( ψ ) , Y m , N 1 , m , and N 2 , m by the corresponding representation form in relations (28) and (29). We reach at the fundamental matrix equation, which has the representation
A m Z m = G m , o r A m ; G m ,
where
A m : = W ˜ ¯ E ˜ + Π X ˜ E ˜ + Π 1 X ^ E ^ Z ^ 1 m X ¯ ^ E ¯ 1 + Π 2 X ¯ ^ E ¯ 2 .
It should be noted that the resultant matrix Equation (30) is a nonlinear system with 3 P + 3 unknowns ζ p , q m for p = 0 , 1 , , P and q = 1 , 2 , 3 to be determined by our proposed algorithm.
Next, we incorporate the intial conditions (2) into (30). Towards this end, we provide the matrix form of (2). By tending t 0 in the first relation in (25) we reach at the next form
A m , 0 Z m = G m , 0 , A m , 0 : = X ˜ ( 0 ) E ˜ , G m , 0 = u 0 v 0 w 0 .
Let us emphasize that the three constants u 0 , v 0 , and w 0 are at hand by (2). The rows one to third of the augmented matrix [ A m ; G m ] will now be substituted by the above three rows [ A m , 0 ; G m , 0 ] . Thus, we reach at the modified and final version of fundamental matrix equation as
A m ˜ Z m = G m ˜ , o r A m ˜ ; G m ˜ .
We solve (31) to acquire the approximate solution of the financial model (1) on the subintervals I m for m = 1 , 2 , , M . As stated before, for the first subinterval I 0 , we use the original initial conditions (2). By iterating on m, the approximate solutions on I m are obtained by using the previously approximate solutions on I m 1 as the initial conditions evaluated at the starting point of I m 1 .

4.2. Computing Errors via Residual Error Function Technique

The true exact form of the solutions of the financial system (1) are not yet known. To check the quality of approximations, we compute the residual error functions (REFs) rather than the difference between the exact and obtained numerical solutions. Towards this end, we first calculate the three approximations u P , δ m ( t ) , v P , δ m ( t ) , and w P , δ m ( t ) via the proposed GMBPs matrix collocation algorithm as described above. Afterwards, we insert these approximations into the financial system (1) followed by defining the REFs as the difference between the left and right hand sides of the involved equations. In the other words, we set the REFs for (1) on the subintervals I m as
R 1 , P m ( t ) : = L C D t ψ u P , δ m ( t ) + λ 1 u P , δ m ( t ) u P , δ m ( t ) v P , δ m ( t ) w P , δ m ( t ) 0 , R 2 , P m ( t ) : = L C D t ψ v P , δ m ( t ) + u P , δ m ( t ) 2 + λ 2 v P , δ m ( t ) 1 0 , R 3 , P m ( t ) : = L C D t ψ w P , δ m ( t ) + u P , δ m ( t ) + λ 3 w P , δ m ( t ) 0 .

4.3. The QLM-GMBPs Collocation Approach

The direct matrix collocation procedure based on the GMBPs is illustrated in the last part for the nonlinear chaotic system (1). We expect that a combination of the proposed approach with the domain-splitting leads to a better accuracy of the obtained approximate solutions. However, using a large value for the number of basis function (P) may cause of some inefficiencies during the execution of the algorithm, which comes from the intrinsic nonlinearity of the model under consideration. The idea of linearization has usually been worked to accelerate the running time of proposed algorithms. Applying the technique of quasi-lineariztion method (QLM) has been successfully proved by many published works, see cf. [39,40,41]. The QLM aimed to not only linearize the model but also keep the accuracy of the proposed algorithm in the same level as the directly employed algorithm.
Let us we rewrite the original model (1) in the form
L C D t ψ y ( t ) = F ( t , y ( t ) ) ,
where
y ( t ) = u ( t ) v ( t ) w ( t ) , F ( t , y ( t ) ) = f 1 ( t ) f 2 ( t ) f 3 ( t ) = λ 1 u ( t ) + u ( t ) v ( t ) + w ( t ) u 2 ( t ) λ 2 v ( t ) + 1 u ( t ) λ 3 w ( t ) .
By starting the initial rough guess y 0 ( t ) for the solution y ( t ) of system (33), we propose the following iterative QLM as follows
L C D t ψ y e ( t ) F ( t , y e 1 ( t ) ) + F y ( t , y e 1 ( t ) ) y e ( t ) y e 1 ( t ) , e = 1 , 2 , ,
where the matrix F y represents the Jacobian matrix of size 3 × 3 . Applying the linearization, we get
L C D t ψ y e ( t ) + J e 1 ( t ) y e ( t ) = g e 1 ( t ) , e = 1 , 2 , .
Here, we have utilized the notations
y e ( t ) = u e ( t ) v e ( t ) w e ( t ) , J e 1 ( t ) = v e 1 ( t ) + λ 1 u e 1 ( t ) 1 2 u e 1 ( t ) λ 2 0 1 0 λ 3 , g e 1 ( t ) = u e 1 ( t ) v e 1 ( t ) 1 + u e 1 2 ( t ) 0 .
The above system has also the initial conditions y e ( 0 ) = u 0 v 0 w 0 T in accordance to (2).
Now, the direct GMBPs matrix approach can be applied to the set of linearized equations (34) without any difficulties stem from the nonlinear terms in (19). Supposedly, we have already the approximate solution y P ( e ) ( t ) available for the solution y e ( t ) of (34) in the iteration number e 1 . Then, we take the next approximation on I m as
u P , δ , m ( e ) ( t ) = p = 0 P ζ p , m , 1 ( e ) MB p δ ( t ) = Λ δ ( t ) Z P , m , 1 ( e ) = X P δ ( t ) E P Z P , m , 1 ( e ) , v P , δ , m ( e ) ( t ) = p = 0 P ζ p , m , 2 ( e ) MB p δ ( t ) = Λ δ ( t ) Z P , m , 2 ( e ) = X P δ ( t ) E P Z P , m , 2 ( e ) , w P , δ , m ( e ) ( t ) = p = 0 P ζ p , m , 3 ( e ) MB p δ ( t ) = Λ δ ( t ) Z P , m , 3 ( e ) = X P δ ( t ) E P Z P , m , 3 ( e ) ,
where Z P , m , q ( e ) = ζ 0 , m , q ( e ) ζ 1 , m , q ( e ) ζ P , m , q ( e ) T for q = 1 , 2 , 3 are the vectors of unknowns in the iteration e = 1 , 2 , . Moreover, the vector X P δ ( t ) and matrix E P are defined in (21). We continue by inserting the collocation nods (18) into (34) followed by introducing the following matrix and vector notations
J ^ e 1 = J e 1 ( t m , 0 ) 0 0 0 J e 1 ( t m , 1 ) 0 0 0 J e 1 ( t m , P ) , Z m ( e ) = Z P , m , 1 ( e ) Z P , m , 2 ( e ) Z P , m , 3 ( e ) , g e 1 = g e 1 ( t m , 0 ) g e 1 ( t m , 1 ) g e 1 ( t m , P ) .
By employing the relations (28), we reach at the next fundamental matrix equation
A m ( e ) Z m ( e ) = G e 1 , or [ A m ( e ) ; G e 1 ] ,
where A m ( e ) : = W ˜ ¯ + J ^ e 1 X ˜ ¯ E ˜ for e 1 on I m . Here, the matrices X ˜ ¯ , W ˜ ¯ , and E ˜ are given in (28). In analogue to the direct approach, we enter the initial y e ( 0 ) = u 0 v 0 w 0 T into the matrix equation (36). By solving the alternate matrix equation we arrive at the unknown coefficients Z m ( e ) efficiently and our job is finished.
Finally, let us mention that the REFs related to the QLM-GMBPs approach are obtained as (32). We denote them by R P , m , q ( e ) ( t ) for q = 1 , 2 , 3 defined on I m and on the iteration e 1 . For a fixed value of e, we further calculate the numerical order of convergence related to the achieved REFs in the infinity norm as
L , q L , q ( P ) : = max 1 m M max t I m R P , m , q ( e ) ( t ) , q = 1 , 2 , 3 .
Hence, we set
Ord P , q : = log 2 L , q ( P ) L , q ( 2 P ) , q = 1 , 2 , 3 .

5. Experimental and Graphical Results

The financial system (1) will now be solved by using the GMBPs/QLM-GMBPs collocation techniques as described in the last parts. In order to validate our results numerically, different model parameters will be considered to show the performance of the proposed spectral collocation approaches. Moreover, various values of fractional order ψ are used in the simulations to show the applicability of the presented algorithms. All the experimental simulations are performed by utilizing Matlab 2021a executed on a personnel laptop. Also, note that employing δ = 5 gives us sufficiently accurate results in the QLM-GMBPs collocation procedure. The starting vector function y 0 ( t ) is taken as the zero function.
Example 1.
To begin, we consider the parameters related to the chaotic financial system (1) taking from [22] as follow
λ 1 = 0.3 , λ 2 = 0.2 , λ 3 = 0.1 .
The initial conditions are given by
u 0 = v 0 = w 0 = 0.2 .
We set P = 5 and take ψ , δ = 1 in the computations. We first focus on the performance of the novel QLM-GMBPs collocation technique based on domain-splitting strategy compared to traditional version of this algorithm. For this purpose, we consider T = 3 and M = 3 and M = 1 in these two algorithms respectively. The result of approximations are depicted in Figure 1 and Figure 2 for the three solutions u ( t ) , v ( t ) , and w ( t ) . For validation, we also plot the results obtained via ode45 routine in Matlab. It is evident that the outcomes of domain-splitting are accurate than those obtained by using M = 1 especially for the interest rate u ( t ) . The gap between the results of both approaches will increase if one takes a longer domain of computation.
We next show that the results of direct and indirect collocation schemes i.e., the GMBPs and QLM-GMBPs matrix algorithms are coincided. Let us take T = 1 and set the parameters as above. We only report the results of the direct GMBPs collocation matrix approach for t [ 0 , 1 ] as follow
u 5 , 1 1 ( t ) = 0.00038622572 t 5 0.0056827876 t 4 + 0.018763555 t 3 0.026138427 t 2 + 0.17985295 t + 0.2 , v 5 , 1 1 ( t ) = 0.00013208249 t 5 0.0010136806 t 4 + 0.0021680058 t 3 0.12829803 t 2 + 0.92005533 t + 0.2 , w 5 , 1 1 ( t ) = 0.00093962517 t 5 0.0045590298 t 4 + 0.011052039 t 3 0.078820297 t 2 0.22001844 t + 0.2 .
It should be noted that the outcomes of QLM-GMBPs technique are very close to the above approximations and we omit them to save space. However, the achieved REFs related to both approaches are visualized in Figure 3 graphically. On the left plot, we show the REFs of direct GMBPs collocation method using P = 5 , 10 , 15 while on the right picture we depict the REFs utilizing P = 5 , 10 , 15 , 20 in the QLM-GMBPs matrix technique. It is clearly seen that in the direct approach we have problem with the convergence of the method if we increase P. The reason may be the fact that the used nonlinear solver fails to give us a reasonable approximation. On the other hand, we can utilize larger values of P in the QLM-GMBPs technique to reach at more accurate results. Moreover, the latter iterative method has more efficacy than the direct method. For instance, for P = 10 , the required time to solve the fundamental matrix Equation (31) is about 23.523 s while the elapsed time for the corresponding matrix equation via QLM-GMBPs procedure is approximately 2.001 s. Consequently, in the following experiments, we usually run the latter numerical approach due to its efficacy.
We further show the utility of QLM-GMBPs approach by using P = 5 and a relatively large interval T = 50 . The results of approximations on the first and last subibtervals I 1 and I 50 are obtained as follow
u 5 , 1 , 1 5 ( t ) = 0.00038622573 t 5 0.0056827876 t 4 + 0.018763555 t 3 0.026138427 t 2 + 0.17985295 t + 0.2 , v 5 , 1 , 1 5 ( t ) = 0.00013208242 t 5 0.0010136805 t 4 + 0.0021680057 t 3 0.12829803 t 2 + 0.92005533 t + 0.2 , w 5 , 1 , 1 5 ( t ) = 0.00093962517 t 5 0.0045590298 t 4 + 0.011052039 t 3 0.078820297 t 2 0.22001844 t + 0.2 ,
for t 0 τ 1 and
u 5 , 1 , 50 5 ( t ) = 0.10225624 t 5 25.288941 t 4 + 2501.5109 t 3 123712.94 t 2 + 3058923.9 t 30251900.0 , v 5 , 1 , 50 5 ( t ) = 0.056917986 t 5 + 14.097621 t 4 1397.1477 t 3 + 69254.485 t 2 1716963.7 t + 17032175.0 , w 5 , 1 , 50 5 ( t ) = 0.014757577 t 5 + 3.6916847 t 4 369.36648 t 3 + 18477.276 t 2 462146.33 t + 4623649.6 ,
for t 49 τ 50 . In Figure 4 and Figure 5 we depict all approximate solutions on the domains I m for m = 1 , 2 , , 50 . The results of the ode45 are also plotted to show the validity of our approximations. By looking at these figures we find that the results of both numerical strategies are overlapped. We get a higher order accuracy if we increase the number of bases or use a lager value of subdivisions M. We next take P = 10 . The 3D pictorial representation of uvw-plot as well as the 2D depiction of uv-plot are displayed in Figure 6, which show the chaotic dynamical behavior of the financial system (1). Furthermore, the similar results related to 2D representations of uw/vw-plot are shown in the next Figure 7. For validation, the numerical approximations using the ode45 are also displayed in this figure for comparison purposes. It is clear that the outcomes of both approaches are coincided.
We now pay attention to the non-integer cases and consider ψ = 0.75 . Using the domain-splitting approach along with the QLM-GMBPs strategy we plot the approximate solutions in Figure 8. We utilize the values of parameter δ to be both 1 and equal to ψ = 0.75 . As shown in this figure, the results obtained to δ = 1 and δ = ψ are very close together. However, the REFs achieved by using δ = ψ are smaller in magnitude compared to those obtained by δ = 1 as presented in Figure 9. Note that one gets more accurate results by increasing P. Towards this end, we employ different values of P in the next experiments and report the maximum values of REFs and their associated order of convergence, i.e., Ord P , q for q = 1 , 2 , 3 . These simulation results are tabulated in Table 1. We use M = 1 , ψ = 0.5 , and P = 4 , 8 , 16 , 32 . The parameter δ is taken as δ = 1 , ψ . Clearly, the exponential order of convergence is visible for the results presented in Table 1. Of course, a slightly higher order of accuracy is seen for the case δ = 1 compared to δ = 0.5 .
Finally, we employ diverse values of fractional orders ψ = 0.25 , 0.5 , 0.75 and visualize the numerical simulations for the three approximated solutions of system (1). We utilize P = 10 , M = 3 , and the parameter δ is equal to one. These experimental results are depicted in Figure 10 and Figure 11. We also plotted the results related to ψ = 1 in these figures.

6. Conclusions

Two spectral matrix collocation methods based on (novel) domain-splitting strategy developed to handle more accurately the nonlinear fractional chaotic system arsing in financial market. The considered fractional derivative is Liouville–Caputo operator and a modified and generalized version of Bessel polynomials (GMBPs) is used in the implementation of the spectral approaches. In the first and direct collocation method, we converted the underlying model into a nonlinear algebraic system whereas in the second and efficient version of the first method, we arrived at a family of algebraic system of equations to be solve iteratively. To solve the model on a long domain, we employed the technique of domain-splitting together with aforementioned spectral schemes, which lead to a higher-order accuracy compared to traditional spectral approaches. A convergence theorem related to GMBPs established and an upper bound for the error was given. To support the theoretical findings, diverse computational experiments were conducted to show the utility of the proposed collocation approaches. The results presented through figures and tables indicate the improvement of the spectral methods based on domain-splitting strategy.

Author Contributions

Conceptualization, M.I. and H.M.S.; methodology, M.I. and H.M.S.; software, M.I.; validation, M.I. and H.M.S.; formal analysis, M.I. and H.M.S.; funding acquisition, H.M.S.; investigation, M.I. and H.M.S.; writing-original draft preparation, M.I.; writing-review and editing, M.I. and H.M.S. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Data Availability Statement

Not Applicable.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Devaney, R.L.; Keen, L. Chaos and Fractals: The Mathematics behind the Computer Graphics; American Mathematical Society: Providence, RI, USA, 1989. [Google Scholar]
  2. Biswas, H.R.; Hasan, M.M.; Bala, S.K. Chaos theory and its applications in our real life. Barishal Univ. J. Part 1 2018, 5, 123–140. [Google Scholar]
  3. Ma, Y.C.J. Study for the bifurcation topological structure and the global complicated character of a kind of non-linear finance system (I). Appl. Math. Mech. 2001, 22, 1240–1251. [Google Scholar] [CrossRef]
  4. Ma, Y.C.J. Study for the bifurcation topological structure and the global complicated character of a kind of non-linear finance system (II). Appl. Math. Mech. 2001, 22, 1375–1382. [Google Scholar] [CrossRef]
  5. Chen, W.C. Nonlinear dynamics and chaos in a fractional-order financial system. Chaos Solitons Fractals 2008, 36, 1305–1314. [Google Scholar] [CrossRef]
  6. Kilbas, A.A.; Srivastava, H.M.; Trujillo, J.J. Theory and Application of Fractional Differential Equations. In North-Holland Mathematics Studies; Elsevier: Amsterdam, The Netherlands, 2006; Volume 204. [Google Scholar]
  7. Hartle, T.T.; Lorenzo, C.F.; Qammer, H.K. Chaos in a fractional order Chua’s system. IEEE Trans CAS-I 1995, 42, 485–490. [Google Scholar] [CrossRef]
  8. Li, C.; Peng, G. Chaos in Chen’s system with a fractional order. Chaos Solitons Fractals 2004, 22, 443–450. [Google Scholar] [CrossRef]
  9. Li, C.; Chen, G. Chaos and hyperchaos in fractional order Rossler equation. Physica A 2004, 341, 55–61. [Google Scholar] [CrossRef]
  10. Grigorenko, I.; Grigirenko, E. Chaotic dynamics of the fractional Lorenz system. Phys. Rev. Lett. 2003, 91, 034101. [Google Scholar] [CrossRef]
  11. Izadi, M.; Srivastava, H.M. Fractional clique collocation technique for numerical simulations of fractional-order Brusselator chemical model. Axioms 2022, 11, 654. [Google Scholar] [CrossRef]
  12. Sabermahani, S.; Ordokhani, Y. Numerical solution of a fractional epidemic model via general Lagrange scaling functions with bibliometric analysis. In Mathematical Analysis of Infectious Diseases; Academic Press: Cambridge, MA, USA, 2022; pp. 305–320. [Google Scholar]
  13. Bidarian, M.; Saeedi, H.; Baloochshahryari, M.R. A Legendre Tau method for numerical solution of multi-order fractional mathematical model for COVID-19 disease. Comput. Methods Differ. Equ. 2023. [Google Scholar] [CrossRef]
  14. Owolabi, K.M.; Atangana, A.; Gómez-Aguilar, J.F. Fractional Adams-Bashforth scheme with the Liouville-Caputo derivative and application to chaotic systems. Discret. Contin. Dyn. Syst. S 2021, 14, 2455. [Google Scholar] [CrossRef]
  15. Alquran, M.; Alsukhour, M.; Ali, M.; Jaradat, I. Combination of Laplace transform and residual power series techniques to solve autonomous n-dimensional fractional nonlinear systems. Nonlinear Eng. 2021, 10, 282–292. [Google Scholar] [CrossRef]
  16. Adel, M.; Khader, M.M.; Ahmad, H.; Assiri, T.A. Approximate analytical solutions for the blood ethanol concentration system and predator-prey equations by using variational iteration method. AIMS Math. 2023, 8, 19083–19096. [Google Scholar] [CrossRef]
  17. Sunthrayuth, P.; Dutt, H.M.; Ghani, F.; Arefin, M.A. The analysis of fractional-order system delay differential equations using a numerical Method. Complexity 2022, 2022, 3570667. [Google Scholar] [CrossRef]
  18. Ma, S.; Xu, Y.; Yue, W. Numerical solutions of a variable-order fractional financial system. J. Appl. Math. 2012, 2012, 417942. [Google Scholar] [CrossRef] [Green Version]
  19. Moutsinga, C.R.B.; Pindza, E.; Mare, E. A robust spectral integral method for solving chaotic finance systems. Alex. Eng. J. 2020, 59, 601–611. [Google Scholar] [CrossRef]
  20. Farman, M.; Akgül, A.; Saleem, M.U.; Imtiaz, S.; Ahmad, A. Dynamical behaviour of fractional-order finance system. Pramana 2020, 164, 94. [Google Scholar] [CrossRef]
  21. Shabani, A.; Sheikhani, A.R.; Aminikhah, H. Robust Control for variable order time fractional financial system. Dyn. Syst. Appl. 2020, 29, 111–122. [Google Scholar] [CrossRef]
  22. Junswang, P.; Sabir, Z.; Raja, M.A.Z.; Adel, W.; Botmart, T.; Weera, W. Intelligent networks for chaotic fractional-order nonlinear financial model. Comput. Mater. Contin. 2022, 72, 5015–5030. [Google Scholar] [CrossRef]
  23. Ahmad, I.; Ouannas, A.; Shafiq, M.; Pham, V.T.; Baleanu, D. Finite-time stabilization of a perturbed chaotic finance model. J. Adv. Res. 2022, 32, 1–14. [Google Scholar] [CrossRef]
  24. Gao, W.; Veeresha, P.; Baskonus, H.M. Dynamical analysis fractional-order financial system using efficient numerical methods. Appl. Math. Sci. Eng. 2023, 31, 2155152. [Google Scholar] [CrossRef]
  25. Xu, C.J.; Liu, Z.X.; Liao, M.X.; Yao, L.Y. Theoretical analysis and computer simulations of a fractional order bank data model incorporating two unequal time delays. Expert Syst. Appl. 2022, 199, 116859. [Google Scholar] [CrossRef]
  26. Srivastava, H.M.; Raghavan, D.; Nagarajan, S. A comparative study of the stability of some fractional-order cobweb economic models. Rev. Real Acad. Cienc. Exactas Fís. Natur. Ser. A Mat. (RACSAM) 2022, 116, 98. [Google Scholar] [CrossRef]
  27. Edeki, S.O.; Okoli, D.C.; Ahmad, H.; Wong, W.K. Approximate series solutions of a one-factor term structure model for bond pricing. Ann. Financ. Econ. 2021, 16, 1–22. [Google Scholar] [CrossRef]
  28. Odibat, Z.M.; Shawagfeh, N.T. Generalized Taylor’s formula. Appl. Math. Comput. 2007, 186, 286–293. [Google Scholar] [CrossRef]
  29. Van Rossum, H. Totally positive polynomials. Proc. Kon. Nederl. Akud. Wetensch. Ser. A 1965, 68; Zndag. Math. 1965, 27, 305–315. [Google Scholar] [CrossRef]
  30. Krall, H.L.; Frink, O. A new class of orthogonal polynomials: The Bessel polynomials. Trans. Am. Math. Soc. 1949, 65, 100–115. [Google Scholar] [CrossRef]
  31. Izadi, M.; Cattani, C. Generalized Bessel polynomial for multi-order fractional differential equations. Symmetry 2020, 12, 1260. [Google Scholar] [CrossRef]
  32. Ahmed, S.S.; MohammedFaeq, S.J. Bessel collocation method for solving Fredholm-Volterra integro-fractional differential equations of multi-high order in the Caputo sense. Symmetry 2021, 13, 2354. [Google Scholar] [CrossRef]
  33. Izadi, M. Numerical approximation of Hunter-Saxton equation by an efficient accurate approach on long time domains. U.P.B. Sci. Bull. Ser. A 2021, 83, 291–300. [Google Scholar]
  34. Izadi, M.; Yüzbası, S.; Noeiaghdam, S. Approximating solutions of non-linear Troesch’s problem via an efficient quasi-linearization Bessel approach. Mathematics 2021, 9, 1841. [Google Scholar] [CrossRef]
  35. Yüzbası, S.; Izadi, M. Bessel-quasilinearization technique to solve the fractional-order HIV-1 infection of CD4+ T-cells considering the impact of antiviral drug treatment. Appl. Math. Comput. 2022, 431, 127319. [Google Scholar] [CrossRef]
  36. Izadi, M.; Yüzbası, S.; Cattani, C. Approximating solutions to fractional-order Bagley-Torvik equation via generalized Bessel polynomial on large domains. Ric. Mat. 2023, 72, 235–261. [Google Scholar] [CrossRef]
  37. Srivastava, H.M. An introductory overview of Bessel polynomials, the generalized Bessel polynomials and the q-Bessel polynomials. Symmetry 2023, 15, 822. [Google Scholar] [CrossRef]
  38. Izadi, M.; Yüzbası, S.; Adel, W. A new Chelyshkov matrix method to solve linear and nonlinear fractional delay differential equations with error analysis. Math. Sci. 2022. [Google Scholar] [CrossRef]
  39. Aznam, S.M.; Ghani, N.A.; Chowdhury, M.S. A numerical solution for nonlinear heat transfer of fin problems using the Haar wavelet quasilinearization method. Results Phys. 2019, 14, 102393. [Google Scholar] [CrossRef]
  40. Izadi, M. A combined approximation method for nonlinear foam drainage equation. Sci. Iran. 2022, 29, 70–78. [Google Scholar] [CrossRef]
  41. Ahmed, S.; Jahan, S.; Nisar, K.S. Hybrid Fibonacci wavelet method to solve fractional-order logistic growth model. Math. Meth. Appl. Sci. 2023. [Google Scholar] [CrossRef]
Figure 1. Approximate solutions of u ( t ) (left) and v ( t ) (right) obtained via QLM-GMBPs based on traditional ( M = 1 ) and domain-splitting ( M = 3 ) in Example 1 with P = 5 , ψ , δ = 1 .
Figure 1. Approximate solutions of u ( t ) (left) and v ( t ) (right) obtained via QLM-GMBPs based on traditional ( M = 1 ) and domain-splitting ( M = 3 ) in Example 1 with P = 5 , ψ , δ = 1 .
Computation 11 00130 g001
Figure 2. Approximate solutions of w ( t ) obtained via QLM-GMBPs based on traditional ( M = 1 ) and domain-splitting ( M = 3 ) in Example 1 with P = 5 , ψ , δ = 1 .
Figure 2. Approximate solutions of w ( t ) obtained via QLM-GMBPs based on traditional ( M = 1 ) and domain-splitting ( M = 3 ) in Example 1 with P = 5 , ψ , δ = 1 .
Computation 11 00130 g002
Figure 3. The obtained REFs via GMBPs (left) and QLM-GMBPs (right) methods in Example 1 with ψ , δ = 1 for t [ 0 , 1 ] and different P = 5 , 10 , 15 , 20 .
Figure 3. The obtained REFs via GMBPs (left) and QLM-GMBPs (right) methods in Example 1 with ψ , δ = 1 for t [ 0 , 1 ] and different P = 5 , 10 , 15 , 20 .
Computation 11 00130 g003
Figure 4. Approximate solutions of u ( t ) (left) and v ( t ) (right) obtained via QLM-GMBPs based on domain-splitting ( M = 50 ) in Example 1 with P = 5 , ψ , δ = 1 .
Figure 4. Approximate solutions of u ( t ) (left) and v ( t ) (right) obtained via QLM-GMBPs based on domain-splitting ( M = 50 ) in Example 1 with P = 5 , ψ , δ = 1 .
Computation 11 00130 g004
Figure 5. Approximate solutions of w ( t ) obtained via QLM-GMBPs based domain-splitting ( M = 50 ) in Example 1 with P = 5 , ψ , δ = 1 .
Figure 5. Approximate solutions of w ( t ) obtained via QLM-GMBPs based domain-splitting ( M = 50 ) in Example 1 with P = 5 , ψ , δ = 1 .
Computation 11 00130 g005
Figure 6. 3D uvw-portrait (left) and 2D uv-portrait (right) of approximate solutions obtained via QLM-GMBPs based on domain-splitting ( M = 50 ) in Example 1 with P = 10 , ψ , δ = 1 .
Figure 6. 3D uvw-portrait (left) and 2D uv-portrait (right) of approximate solutions obtained via QLM-GMBPs based on domain-splitting ( M = 50 ) in Example 1 with P = 10 , ψ , δ = 1 .
Computation 11 00130 g006
Figure 7. 2D uw-portrait (left) and 2D vw-portrait (right) of approximate solutions obtained via QLM-GMBPs based on domain-splitting ( M = 50 ) in Example 1 with P = 10 , ψ , δ = 1 .
Figure 7. 2D uw-portrait (left) and 2D vw-portrait (right) of approximate solutions obtained via QLM-GMBPs based on domain-splitting ( M = 50 ) in Example 1 with P = 10 , ψ , δ = 1 .
Computation 11 00130 g007
Figure 8. Approximate solutions obtained via QLM-GMBPs based on domain-splitting ( M = 2 ) in Example 1 with P = 5 , ψ = 0.75 , and δ = 1 , 0.75 .
Figure 8. Approximate solutions obtained via QLM-GMBPs based on domain-splitting ( M = 2 ) in Example 1 with P = 5 , ψ = 0.75 , and δ = 1 , 0.75 .
Computation 11 00130 g008
Figure 9. The obtained REFs via QLM-GMBPs based on domain-splitting ( M = 2 ) in Example 1 with P = 5 , ψ = 0.75 , and δ = 1 , 0.75 .
Figure 9. The obtained REFs via QLM-GMBPs based on domain-splitting ( M = 2 ) in Example 1 with P = 5 , ψ = 0.75 , and δ = 1 , 0.75 .
Computation 11 00130 g009
Figure 10. Approximate solutions of u ( t ) (left) and v ( t ) (right) obtained via QLM-GMBPs based domain-splitting ( M = 3 ) in Example 1 with P = 10 , ψ = 0.25 , 0.5 , 0.75 , 1 , and δ = 1 .
Figure 10. Approximate solutions of u ( t ) (left) and v ( t ) (right) obtained via QLM-GMBPs based domain-splitting ( M = 3 ) in Example 1 with P = 10 , ψ = 0.25 , 0.5 , 0.75 , 1 , and δ = 1 .
Computation 11 00130 g010
Figure 11. Approximate solutions of w ( t ) obtained via QLM-GMBPs based domain-splitting ( M = 3 ) in Example 1 with P = 10 , ψ = 0.25 , 0.5 , 0.75 , 1 , and δ = 1 .
Figure 11. Approximate solutions of w ( t ) obtained via QLM-GMBPs based domain-splitting ( M = 3 ) in Example 1 with P = 10 , ψ = 0.25 , 0.5 , 0.75 , 1 , and δ = 1 .
Computation 11 00130 g011
Table 1. The outcomes of maximum REFs error norms and the related Ord P , q for q = 1 , 2 , 3 achieved via QLM-GMBPs in Example 1 with ψ = 0.5 , γ = 1 , 0.5 , M = 1 , and different P.
Table 1. The outcomes of maximum REFs error norms and the related Ord P , q for q = 1 , 2 , 3 achieved via QLM-GMBPs in Example 1 with ψ = 0.5 , γ = 1 , 0.5 , M = 1 , and different P.
δ = 1 δ = 0.5
P L , 1 Ord P , 1 L , 2 Ord P , 2 L , 3 Ord P , 3 L , 1 Ord P , 1 L , 2 Ord P , 2 L , 3 Ord P , 3
4 4.13 04 1.84 05 4.29 04 5.96 03 2.09 03 1.56 03
8 3.87 07 10.1 1.61 07 6.84 1.11 07 11.9 8.10 05 6.20 9.73 05 4.42 5.82 05 4.75
16 7.09 15 25.7 2.47 14 22.6 1.19 16 29.8 5.73 08 10.5 3.68 08 11.4 3.05 08 10.9
32 9.61 30 49.4 6.73 29 48.4 5.08 31 47.7 1.37 14 22.0 3.90 16 26.5 2.33 15 23.6
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Izadi, M.; Srivastava, H.M. Applications of Modified Bessel Polynomials to Solve a Nonlinear Chaotic Fractional-Order System in the Financial Market: Domain-Splitting Collocation Techniques. Computation 2023, 11, 130. https://doi.org/10.3390/computation11070130

AMA Style

Izadi M, Srivastava HM. Applications of Modified Bessel Polynomials to Solve a Nonlinear Chaotic Fractional-Order System in the Financial Market: Domain-Splitting Collocation Techniques. Computation. 2023; 11(7):130. https://doi.org/10.3390/computation11070130

Chicago/Turabian Style

Izadi, Mohammad, and Hari Mohan Srivastava. 2023. "Applications of Modified Bessel Polynomials to Solve a Nonlinear Chaotic Fractional-Order System in the Financial Market: Domain-Splitting Collocation Techniques" Computation 11, no. 7: 130. https://doi.org/10.3390/computation11070130

APA Style

Izadi, M., & Srivastava, H. M. (2023). Applications of Modified Bessel Polynomials to Solve a Nonlinear Chaotic Fractional-Order System in the Financial Market: Domain-Splitting Collocation Techniques. Computation, 11(7), 130. https://doi.org/10.3390/computation11070130

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop