Innovative Applications of Fractional CalculusView this Special Issue
Qualitative Analysis of Class of Fractional-Order Chaotic System via Bifurcation and Lyapunov Exponents Notions
This paper presents a modified chaotic system under the fractional operator with singularity. The aim of the present subject will be to focus on the influence of the new model’s parameters and its fractional order using the bifurcation diagrams and the Lyapunov exponents. The new fractional model will generate chaotic behaviors. The Lyapunov exponents’ theories in fractional context will be used for the characterization of the chaotic behaviors. In a fractional context, the phase portraits will be obtained with a predictor-corrector numerical scheme method. The details of the numerical scheme will be presented in this paper. The numerical scheme will be used to analyze all the properties addressed in this present paper. The Matignon criterion will also play a fundamental role in the local stability of the presented model’s equilibrium points. We will find a threshold under which the stability will be removed and the chaotic and hyperchaotic behaviors will be generated. An adaptative control will be proposed to correct the instability of the equilibrium points of the model. Sensitive to the initial conditions, we will analyze the influence of the initial conditions on our fractional chaotic system. The coexisting attractors will also be provided for illustrations of the influence of the initial conditions.
In the recent years, modeling chaotic and hyperchaotic systems occupy an important place in the literature and have many applications in physics, biology, electrical circuits, and many other fields [1–4]. The most used fields for the applications of chaos are modeling electrical circuits, and there exist many papers related to the implementation of the chaotic systems in this domain. Many phenomena in the real-world problems are complicated to be predicted and justify the use of chaotic models. Nowadays, there appear many tools for analyzing the chaotic systems as the phase portraits of the system using the numerical discretizations, the bifurcation diagrams to understand the influence of the models’ parameters on the dynamics of the chaotic models, and the Lyapunov exponents used to determine the nature of the chaos. There exists some chaos as chaotic behaviors and hyperchaotic behaviors. As tools, we can also cite the bicoherence and the Poincare map; there is also an algorithm to focus on the initial conditions’ influence. It is known that the chaos systems are sensitive to the variation of the initial conditions. Influencing initial conditions can generate a loss of chaotic behaviors or hyperchaotic behaviors. Fractional calculus has attracted much attention these years, and many fractional operators have been introduced in this new field. As operators in this field, we can cite the Riemann–Liouville derivative . We can cite the Caputo derivative [5, 6], which is the most used operator in fractional calculus due to its physical adequacy with physical problems. Other fractional operators with Mittag–Leffler kernel [7, 8] and exponential kernel exist  and continue to impress all the community in fractional calculus . For the use of the Atangana–Baleanu derivative, see . Many of them have advantages in modeling real-world problems. Modeling chaotic systems and hyperchaotic systems to capture the memories effects have constituted a new direction of research in the recent years; see [8, 12]. For the advancement of fractional calculus and its application, the readers can refer to the following papers: in , the authors address a new numerical scheme for solving fractional differential equations described by Gomez–Atangana–Caputo derivative; in , the authors focus on the characterization of two differential fractional operators without singular kernels; and, in , the authors analyze an epidemic spreading model described by a fractional operator with Mittag–Leffler kernel. For recent works on numerical methods applied to partial differential equations, see [16–19].
Modeling chaotic and hyperchaotic systems are focused on the literature with integer-order derivative and the non-integer-order derivatives. We make a review of the literature in this paragraph. In , the authors focus on the chaotic Chua electrical circuit in the fractional version. In , the authors discuss and prove the algorithm to get the Lyapunov exponents in a fractional version. In , the authors discuss the fractional-order chaotic system and its suppression in some specific order of the fractional derivative. In , the authors discuss the so-called hyperchaotic chameleon system with the aid of fractional-order derivative and propose its electrical implementation. In , the authors propose investigations in the fractional chaotic system used in finance; the findings were interpreted financially and economically. In , finance and chaotic system were also interpreted in this paper. In , the authors investigated on fractional-order exponential jerk system and presented its electrical implementation. In , the authors presented a new chaotic system in integer version with multiple attractors. In , the authors propose a new chaotic system with a self-excited attractor. In , the authors proposed investigations on hidden attractors in the context of dynamical systems. In , the authors proposed a new simple chaotic system but with admitting a line equilibrium. In , the authors propose synchronization investigations using the 4D hyperchaotic jerk system. In , the authors propose a chaotic system with infinite equilibria located on a piecewise linear curve. In , the authors model the hyperchaotic system using fractional derivative with Mittag–Leffler kernel and fractional derivative with exponential derivative. In , the authors introduce chaos in a cancer modeling using fractional derivatives with exponential decay and Mittag–Leffler law. In , Baskonus et al. propose active control to stabilize a fractional-order macroeconomic model using Lyapunov direct method. See more investigations related to fractional modeling of chaotic systems using Caputo derivative, bifurcation, and Lyapunov analysis in , modeling class of fractional-order chaotic or hyperchaotic system with Caputo derivative in , analysis of a four-dimensional hyperchaotic system described by Caputo–Liouville derivative in , modeling Chua’s electrical circuit in the fractional context in , and modeling chaotic processes with Caputo fractional-order derivative in .
In this paper, we model a chaotic system using the Caputo derivative. The main contributions of this paper are mentioned in this paragraph. First, the phase portraits are obtained using the famous predictor-corrector method valid in the fractional differential equations’ discretizations. Second, the numerical scheme as the predictor-corrector method is the main contribution of this present work. Third, the small changes in our introduced fractional chaotic model have been analyzed in terms of the bifurcation diagrams. Different values of the Caputo derivative are considered; at all these values, it will be important to give the chaos’s nature. In other words, we will use the Lyapunov exponents in the fractional context to decide whether we have chaotic behaviors or not. This analysis is fundamental because the Lyapunov exponents’ classical theories are not valid all time in a fractional context. For example, there exist hyperchaotic systems described by the fractional operators with one positive Lyapunov exponent instead of two positive Lyapunov exponents. Fourth, we observed in our investigation that the presence of zero as the Lyapunov exponent is quasi-impossible in a fractional context. The proof of this assumption will be subject to further investigations in the future. Another contribution addressed in this paper is related to the local stability of the fractional chaotic system’s equilibrium points. Due to the chaotic behaviors, all the points are not stable. Alternatively, we propose feedback control to stabilize the fractional error system after combining the slave chaotic system and the master’s chaotic systems. Another contribution of the present paper is that we provide the coexisting attractors for specific values of the model’s parameters at two different initial conditions.
The remainder of this paper is organized as follows. In Section 2, we recall the fractional tools used in the investigations. In Section 3, we introduce the fractional chaotic model by using the Caputo derivative. Section 4 presents the predictor corrector method proposed in the literature to discretize our fractional chaotic model. In Section 5, the phase portraits considering different fractional-order derivative values are proposed for our fractional model. In Section 6, the bifurcation diagrams for the small variation of the model’s parameters are presented. In Section 7, the natures of the chaos are characterized using the Lyapunov exponents’ calculation in the context of fractional calculus. In Section 8, we analyze the influence of the initial condition in the chaotic behaviors. In Section 9, the stability analysis of the fractional chaotic model’s equilibrium point has been proposed, and the feedback control is presented. In Section 10, final remarks for our works are presented.
2. Fractional Operators
In this section, we make a brief recall of the fractional operators which we will use through our investigation. In this section, we will define the Caputo derivative and the Riemann–Liouville derivative.
Definition 2 (see [5, 6]). We define the Riemann–Liouville fractional derivative with the order for the function in the form described bywith representing the Gamma Euler function and we set the order as .
Due to the inconveniences of the Riemann–Liouville operator, we will focus on our paper with the Caputo derivative. The description of this derivative is given in the following definition.
Definition 3 (see [5, 6]). The Caputo fractional derivative operator of order is symbolized as the following form when we consider a function in the form described bywhere the function represents the Gamma Euler function.
For readers’ interest, the Laplace transform is also fundamental for solving the fractional differential equations analytically. The following relation expresses the Laplace transform of the Caputo derivative:
We set the order with the relation . The Caputo derivative’s properties as the composition with the fractional integral, the derivative of a constant function, or the derivative of the Mittag–Leffler function can be found with more pieces of information in [5, 6]. In this paper, the fractional integral in Riemann–Liouville sense will play an important role in the discretization because the proposed method comes from the numerical scheme of this derivative; the discretization details will be found in the next sections.
3. Modeling the Fractional-Order System
In this section, we introduce our model, which we will study in the next section. An integer version of this model is a chaotic system and has been subject to investigation since 2018 in  and is described by the following equations :with initial conditions , , and . The strange attractor is obtained with the following values for the parameters of the previous model , , , , and . Our new advancement here is to study the same chaotic system, but we consider its fractional version described particularly with the Caputo derivative. Therefore, in this paper, we consider the fractional differential system described as follows:
We impose the initial conditions as the following forms:
In the above modeling, the Caputo derivative is used instead of the classical Riemann–Liouville derivative because we want to use physical initial conditions. The Riemann–Liouville derivative does not have physical initial conditions. Furthermore, the Riemann–Liouville derivative of constant value does not give zero. These inconveniences also explain in this section the use of the Caputo derivative. As it will be noticed in the phase portraits, there exist many types of chaos according to the fractional operator’s order. To prove that the fractional order plays an important role in the chaos systems, notably, it allows having new attractors, contrary to the model with the integer-order derivative, where a new type of chaotic system is obtained with the variation of the models’ parameters.
4. Predictor-Corrector Applied on Fractional-Order System
This section describes the numerical method that we will use to obtain the phase portraits of the fractional-order chaotic system (5)–(7). The discretizations are classical in the literature; we try to apply discretization to our model. In fractional calculus, many numerical schemes and analytical methods can be used as the homotopy methods, a domain decomposition method, the Chebyshev method, and many others. But many inconveniences of the cited methods are still to be solved due to the inconveniences in the stability and the convergences of the approximate solutions. The use of the predictor-corrector method in our system has the advantages of having Matlab codes, fundamental in chaotic and hyperchaotic systems. In the rest of this section, we use the predictor-corrector method reported by Garrappa in his review article . The following can describe the solution of the fractional differential system (5)–(7):
We set the following functions from our fractional-order system (5)–(7):as well as point ; then, according to the numerical scheme named the predictor-corrector method, equations (9)–(11) can be rewritten in the following forms:
Furthermore, the predictor has the following form in our fractional system:
In the above formula, denotes the step size, and the parameters of the discretization are defined in the following forms:and when the indices n describe the condition , we set the parameters as the following expressions:
The approximation of the functions in our model step by step is given by the expressions described in the forms
Before ending this section, we give a brief review concerning the method’s stability and convergence; more pieces of information can be found in Garrappa’s paper. We set that , and are the approximate solutions of the fractional system under Caputo derivative (5)–(7) and the exact solution of our model denoted by , and ; then the residual functions as described are given in the following forms:
It is not hard to see that when the step size converges to zero, we get the convergence of the approximate solution to the exact solutions. The predictor-corrector methods’ stability can be obtained from the Lipschitz conditions of the drift functions , , and of our model.
The method of discretization described in this section has many advantages. Firstly, the method is stable and convergent; the Matlab implementation is useful and simple. Note that the convergence and the stability are essential in the numerical methods; comparing our method with the homotopy analysis methods, we can affirm that our method is more useful because, with the homotopy methods, we cannot determine precisely after how many iterations we have the convergence and the stability of the method. The numerical discretization described in this section is also more useful than the Laplace transform method. There are many nonlinear differential equations where the Laplace transform cannot be applied due to the complexity of the equations’ forms, contrary to our described method in this section which is applicable. In the resolution of the diffusion equation, too, the use of the green function is not trivial, and here also the numerical method can be used.
5. Phase Portraits versus Fractional-Order Derivative
This section is devoted to observing via different phase portraits the Caputo derivative’s influence in the dynamics of our fractional system (5)–(7). This section will illustrate the predictor-corrector procedure as well. To arrive at our end, we take different values of the fractional-order derivative. We consider the following values of the fractional derivative , , , , and . The evolution of the chaos will be observed with these different orders. We begin the representation of the phase portraits of the fractional-order system with the order . The graphical representations in different planes are assigned in Figure 1: and planes. In these first graphical representations, the considered order is .
The graphical representations in different planes are assigned in Figure 2: and planes.
We continue this section with the order . We will see the difference existing when the order of the fractional derivative varies. The graphical representations in different planes are assigned in Figure 3: and planes.
The graphical representations in different planes are assigned in Figure 4: and planes.
The first difference in the phase portrait can be observed, and we confirm the existence of new attractors. Therefore, the fractional order can play an interesting role in the behaviors of the solutions. The graphical representations in different planes are assigned in Figure 5: and planes at the order .
The graphical representations in different planes are assigned in Figure 6: and planes.
We finish the phase portrait section by observing the behaviors that happened behind ; in the next figures, we work with the order . We notice in Figures 7 and 8 significant difference from the previous phase portraits, and our remark is that the chaotic behavior is not removed in the dynamics. The graphical representations in different planes are assigned in Figure 7: and planes.
The graphical representations in different planes are assigned in Figure 8: and planes.
The figures in this section illustrate the impact of the order of the Caputo derivative as well; this impact can be seen in the geometries of the attractors which are different in the cases considered in this section. The phase portraits represented in this section will be classified using the bifurcation maps and the Lyapunov exponents. We will notice new characterizations of chaos using the Lyapunov exponents. As will be remarked, the nature of the chaos depends on the Caputo derivative’s order. In other words, the new order of the Caputo derivative generates new types of chaos.
6. Bifurcation Diagrams
In this section, we analyze the sudden qualitative changes in the nature of the solutions due to the variation of the parameters of the fractional-order system equation (5)–(7). In this section, we also illustrate the qualitative changes of the solutions of our model with the phase portraits.
In the first section, we suppose that the first parameter a of our fractional model varies in the small interval precisely in . In Figure 9 the bifurcation diagram according to the variation of parameter is represented.
Figure 9 informs that we notice high chaotic behaviors into the interval . That is, the chaotic behaviors are not removed when the order is maintained to and the parameter has small variation. For more details, we depict phase portraits 10 of the model with parameter . The graphical representations in different planes are assigned in Figure 10: and planes.
The graphical representations in different planes are assigned in Figure 11: and planes.
We continue with the variation of parameter into the interval (9.10). In Figure 12, we represent the bifurcation diagram according to the variation of parameter .
Figure 12 informs that there exist chaotic behaviors into this interval (9.10) when parameter has small variation. To illustrate the chaotic behaviors due to the small changes of , we represent graphically phase portraits 13, with . The graphical representations in different planes are assigned in Figure 13: and planes.
The graphical representations in different planes are assigned in Figure 14: and planes.
The variation of parameter is now considered. We suppose that the parameter varies into the interval . In Figure 15, the bifurcation diagram due to the variation of parameter is represented.
Same conclusion, the chaotic behavior is present in the considered interval and is illustrated in the following phase portraits 16, with the set . The graphical representations in different planes are assigned in Figure 16: and planes.
The graphical representations in different planes are assigned in Figure 17: and planes.
Bifurcation diagram 18 due to the variation of parameter into the interval is represented in Figure 18 and we confirm the changes with phase portraits.
Phase portraits 19 are represented to illustrate the changes in the behaviors of the dynamics of our model when the parameter is into (4.5). The graphical representations in different planes are assigned in Figure 19: and planes.
The graphical representations in different planes are assigned in Figure 20: and planes.
The final observation using bifurcation diagrams 9, 12, 15, and 18 is that our system has chaotic behaviors. But the changes generated by parameters , , , and are approximately the same. This conclusion can be observed in the phase portraits presented in this section, which in general do not have many differences between them.
7. Chaos Detection via Lyapunov Exponents
For chaos detection, we try in this section to characterize the nature of chaos when the order of the fractional derivative varies. We calculate in particular the Lyapunov exponents at the orders , , , , and . In the second part, we will localize the interval where chaotic or hyperchaotic attractors are obtained. We also calculate the dimension of the Lyapunov exponents. According to the values of the Lyapunov exponents at the order previously considered, we will calculate the sum of the Lyapunov exponents to verify whether our fractional model is dissipative or not. Before the calculations of the Lyapunov exponents, we recall the Jacobian matrix necessary for the algorithm to obtain the Lyapunov exponents; we have the following matrix:
The first remark is that the Lyapunov exponents confirm the results in the bifurcation section; that is, our fractional model has chaotic behaviors. It is because there exists one positive Lyapunov exponent. The theory of Lyapunov exponents is very complex in the context of the use of the fractional operators because zero as the value of the Lyapunov exponent seems very difficult to be obtained using the algorithms to get Lyapunov exponents. This result is correct due to the complexity of the numerical scheme of the fractional operators. The second remark is that, for all considered fractional-order derivatives into , the sum of the Lyapunov exponents is negative, which means that the fractional-order chaotic system (5)–(7) considered in this paper is dissipative. We continue our analysis by considering the Lyapunov exponents at the order . The Lyapunov exponents are given as follows:
Their associated Kaplan–Yorke dimension is given as follows:
The second case is the Lyapunov exponents at the order given by the following numbers:
Their associated Kaplan–Yorke dimension is given as follows:
The third case is the Lyapunov exponents at the order given by the following numbers:
Their associated Kaplan–Yorke dimension is given as follows:
The last case is the Lyapunov exponents at the order given by the following numbers:
Their associated Kaplan–Yorke dimension is given as follows:
We can notice that the chaotic attractor is more significant at the order because the positive Lyapunov exponents and the Lyapunov dimension are large. In comparison between the chaotic behaviors at and , we notice by observations of the phase portraits that the chaotic behaviors are more significant when the order converges to . This behavior is explained by the fact that the positive Lyapunov exponent is larger at the order . The same comparison can be made for the order , where the chaotic behaviors are more important than those at the order . These differences can be observed with the Lyapunov exponents’ values and the Lyapunov dimensions, which are larger at .
8. Initial Conditions Influence and Coexistence Attractors
In this section, we analyze the impact of the initial condition. In other words, what the initial conditions give in the nature of the dynamics of our fractional system will be analyzed. The study of the changes of the initial conditions is important because chaotic and hyperchaotic systems are very sensitive to the changes in the initial conditions. We consider many cases in the initial conditions; first, we influence from to . The illustration of this case is represented in Figure 21.
In Figure 21, the initial condition is in blue color and the initial condition is in red color. We notice that significant changes can be generated by the variation of the initial condition related to . We continue by influencing from to . The illustration of this case is represented in Figure 22.
In Figure 22, the initial condition is in blue color and the initial condition is in red color. We notice that significant changes can be generated by the variation of the initial condition related to . We finish by the influence generated at the last variable from to . The illustration of this case is represented in Figure 23.
In Figure 23, the initial condition is in blue color and the initial condition is in red color. We notice that significant changes can be generated by the variation of the initial condition related to .
The general conclusion is that the initial conditions generate many changes in attractors. Thus, due to the fact that the Lyapunov exponents are sensitive to the initial conditions too, the values of the Lyapunov exponents will vary according to the changes in the initial conditions.
Considering the influence of the initial conditions, we end this part by analyzing the coexisting attractors. We have, for example, the presence of two pairs of attractors when parameters , , , , and at the order and two initial conditions which are given by (blue color) and (red color). The figures of the coexisting attractors are represented in Figures 24 and 25. The graphical representations in different planes are assigned in Figure 24: and planes.
The graphical representations in different planes are assigned in Figure 25: and planes.
The order of the fractional operator plays an important role in the existence of pairs of attractors. To observe this influence, we maintain parameters , , , , and and change the order to ; see Figures 26 and 27. The graphical representations in different planes are assigned in Figure 26: and planes.
The graphical representations in different planes are assigned in Figure 27: and planes.
9. Stability Analysis and Feedback Control
In this last section of our investigation, we focus on the local stability of the equilibrium points of the fractional chaotic model (5)–(7). The equilibrium points of our fractional model are given by , , and . At the first point , the Jacobian matrix in the previous section is given as follows:
The eigenvalues are given as follows: , , and . The second and the last eigenvalues have negative real part and thus satisfy the Matignon criterion , but . Thus, the equilibrium point is not stable. At the second equilibrium point , the Jacobian matrix in the previous section is given as follows:
The eigenvalues are given as follows: , , and . The last eigenvalue has negative real part and thus satisfies the Matignon criterion, but the first and the second eigenvalues do not satisfy the Matignon criterion as . Thus, the equilibrium point is not stable when . At the last equilibrium point , the Jacobian matrix in the previous section is given as follows:
We obtain the same eigenvalues as in the previous point. The eigenvalues are given as follows: , , and . The last eigenvalue has negative real part and satisfies the Matignon criterion, but the first and the second eigenvalues do not satisfy the Matignon criterion as . Thus, the equilibrium point is not stable when .
In the last part, we propose a feedback control to stabilize our chaotic system because, as we observe, all the equilibrium points are not stable when the fractional-order derivative exceeds . Let the slave fractional chaotic system be defined by the following equation:and the master system is given by the following equation:where represents the exogenous input, which attracts our attention. Let us define the error terms given by the following equations:
Then, considering the slave system and the master system, we get the following fractional differential error system. Then, considering the slave system and the master system, we get the following fractional differential error system:
Then, here, to stabilize the fractional error equation, we choose feedback control defined by
Using Lyapunov characterization of the global asymptotic stability, we get the global asymptotic stability of the trivial equilibrium point of the fractional system (38)–(40), which in turn implies that
10. Final Remarks
This paper studies the fundamental properties of a class of fractional-order systems in terms of chaotic behaviors, Lyapunov exponents for characterizing the chaotic or hyperchaotic behaviors, the Lyapunov dimensions, and the stability of the equilibrium points of the model in the context of the Matignon criterion. We find that our system admits chaotic behaviors with all fractional-order derivatives into the interval as can be observed in the figures in Section 5. The fractional chaotic system’s equilibrium points are not stable due to the chaotic behaviors, but we find feedback control to stabilize the model’s error term. The different figures of the dynamics of the model represented in this paper were possible with the proposed numerical discretization aid, including the Riemann–Liouville derivative discretization. The numerical method is specially called the predictor-corrector method applied in our system because it is already reported in the literature. The impact of the parameters of the introduced model is analyzed via the bifurcation concept. The main conclusions of this paper are summarized as follows: we find a region where the fractional-order system exhibits chaotic behaviors; the bifurcation diagrams in Section 6 and the Lyapunov exponents inform us that the fractional-order derivative has a significant impact on the dynamics because new attractors are generated when the order of the fractional derivative varies; the effectiveness of all the analysis in the paper is possible with the aid of the numerical scheme. The paper also informs us that the present chaotic system admits coexisting attractors when the initial conditions vary and with specific parameters; see the figures in Section 8. The paper also contributes to proposing adaptative control for global asymptotic stability. For future research directions, the Lyapunov exponents, the bifurcation diagrams, the stability analysis, the synchronization, and the electrical implementation of the model used in this paper can be focused on as regards Caputo derivative terms with different values of the fractional orders. The fractional chaotic systems with the nonsingular derivatives can also be focused on in the future. This paper addresses numerical schemes; it will be interesting in the future research to draw the circuit associated with the present chaotic model and analyze the simulations obtained in the oscilloscopes. The circuits’ schematic can be done with both the integer-order version and the fractional-order version of our present system. The Poincare map of the present chaotic system in fractional version is also the perspective of new research papers.
No data were used to support this study.
Conflicts of Interest
The author declares that there are no conflicts of interest.
A. A. Kilbas, H. M. Srivastava, and J. J. Trujillo, Theory and Applications of Fractional Differential Equations, North-Holland Mathematics Studies, Elsevier, Amsterdam, The Netherlands, 2006.
I. Podlubny, “Fractional differential equations,” Mathematics in Science and Engineering, Academic Press, New York, NY, USA, 1999.View at: Google Scholar
T. Mekkaoui, Z. Hammouch, D. Kumar, and J. Singh, “A new approximation scheme for solving ordinary differential equation with gomez-atangana-caputo fractional derivative,” Methods of Mathematical Modelling, CRC Press, Boca Raton, FL, USA, 2019.View at: Google Scholar
A. Wakif, A. Chamkha, T. Thumma, I. L. Animasaun, and R. Sehaqui, “Thermal radiation and surface roughness effects on the thermo-magneto-hydrodynamic stability of alumina-copper oxide hybrid nanofluids utilizing the generalized buongiorno’s nanofluid model,” Journal of Thermal Analysis and Calorimetry, vol. 143, no. 2, pp. 1201–1220, 2020.View at: Publisher Site | Google Scholar
A. Wakif, Z. Boulahia, F. Ali, M. R. Eid, and R. Sehaqui, “Numerical analysis of the unsteady natural convection MHD couette nanofluid flow in the presence of thermal radiation using single and two-phase nanofluid models for Cu-water nanofluids,” International Journal of Applied and Computational Mathematics, vol. 4, no. 3, 2018.View at: Publisher Site | Google Scholar
D. Matignon, “Stability results on fractional differential equations to control processing,” in Proceedings of the Computational Engineering in Syatems and Application Multiconference IMACS, pp. 963–968, IEEE-SMC, Lille, France, July 1996.View at: Google Scholar