Models Methods Software

Chaos and ODEs Part 1a: The Literature Sources

I have way too much material for a single post. I have spent days trying to force a good fit for all the material into a single document. I have put that aside for a while. So these discussions will be broken into several parts. At some future time I might try to tie all the pieces together by use of HTML/PDF.

This post gives the literature that has been the subject of the present investigations. Some of the material listed below is background and fundamental material. Other of the citations contain the actual equation systems that I have investigated. A range of such systems have been the focus of these investigations including the original Lorenz systems, Lorenz-like systems, predator-prey systems, porous media, systems involving membrane impulses, and systems derived to exhibit the characteristics of chaotic response. The numerical solution methods used included numerous standard ODE solution methods, finite-difference approximations, and the expansion method developed by Adomian [1986, 1988, and a very expensive book available at Amazon] and extended by others. I will eventually post the results of the numerical solution methods, and the specific calculations done for the selected equation systems.

The equation systems investigated are:
(1) the original system by Saltzman [1962] that seems to at the foundation of many Lorenz-like systems,
(2) the original Lorenz system of 1963 which has been the subject of this paper, this post, and a comment by Young [1966],
(3) additional systems by Lorenz [1984, 1990] that have been also solved by Pielke and Zeng [1994]; these systems are the subject of recent publications by Lorenz [2006a, 2006b],
(4) the Rossler system as reported in Stuart and Humphries, [1996],
(5) the Chen and Ueta system [1999]
(6) the general predator-prey Lotka-Volterra system reported by Vano et al. [2006],
(7) a porous-media system investigated by Vadasz and Olek [2000b, 1999, 2000a],
(8) the bursting membrane system reported by Terman [1991], and
(9) the generalized Lorenz systems developed by Roy and Musielak [I don’t have actual publication dates for these drafts].

I don’t think I’ve forgotten any, but I might have.

I have not been successful in demonstrating that numerical solution methods will converge for any of these systems. The Saltzman, Lotka-Volterra, and membrane systems are special cases relative to this statement. The Saltzman system because it does not and never has exhibited chaotic response, the Lotka-Volterra system because the lack of convergence is marginal in a practical sense, but nonetheless not converged, and the latter because I suspect the system does not in fact exhibit chaotic response. I have not yet calculated the generalized Lorenz systems of (9), but by extrapolation I suspect convergence cannot be demonstrated for those. I will calculate these real soon now.

Here’s a rough estimate of the other parts of the posts on this subject. There will be a part that gives the actual equations systems that were solved including parameter values and initial conditions, a part that gives the specific numerical solution methods used, and a part that gives the results of the calculations (these might be separate posts for each system). I plan to give sufficient information that all my results can be replicated by others. It has been very interesting to me that for a subject that is based on numerical calculations some of the papers listed above do not provide sufficient information to allow replication by others. The fourth-order Runge-Kutta method for example comes in several flavors.

I will send my computer software for any specific calculations to anyone who would like to look it over and/or use it. Be aware that I make and use research-grade software and do not claim to produce production-grade plug-n-play software.

If you have calculated any of these systems and either agree or disagree with my conclusions I would like to hear from you.

REFERENCES

Books

Uri M. Ascher and Linda R. Petzold, Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations, Society for Industrial and Applied Mathematics, SIAM, New York, 1998.

Arieh Iserles, A First Course in the Numerical Analysis of Differential Equations, Cambridge Texts in Applied Mathematics, Cambridge University Press, Cambridge, 1996.

Ronald E. Mickens, Nonstandard Finite Difference Models of Differential Equations, World Scientific, River Edge, NJ, 1994.

Ronald E. Mickens, Advances in Nonstandard Finite Difference Models of Differential Equations, World Scientific, River Edge, NJ, 2001.

A. M. Stuart and A. R. Humphries, Dynamical Systems and Numerical Analysis, Cambridge university Press, Cambridge, 1996.

Charles R. Doering and J. D. Gibbon, Applied Analysis of the Navier-Stokes Equations, Cambridge University Press, 1995.

D. J. Tritton, Physical Fluid Dynamics, Oxford University Press, Oxford, 1988.

Uriel Frisch, Turbulence, Cambridge University Press, Cambridge, 1995.

J. M. Thompson and H. B. Stewart, Nonlinear Dynamics and Chaos, John Wiley & Sons, Ltd., West Sussex, 2002.

E. Adams and U. Kulisch, Scientific Computing with Automatic Result Verification, Academic Press, San Diego, 1993.

Papers

The Classic Lorenz Systems

Barry Saltzman, “Finite Amplitude Free Convection as an Initial Value Problem”, Journal of the Atmospheric Sciences, Vol. 19, pp. 329-341, 1962.
Abstract
The Oberbeck-Boussinesq equations are reduced to a two-dimensional form governing “roll� convection between two free surfaces maintained at a constant temperature difference. These equations are then transformed to a set of ordinary differential equations governing the time variations of the double-Fourier coefficients for the motion and temperature fields. Non-linear transfer processes are retained and appear as quadratic interactions between Fourier coefficients. Energy and heat transfer relations appropriate to this Fourier resolution, and a numerical method for solution from arbitrary initial conditions are given. These solutions, which are for a fixed Prandtl number and a variable Rayleigh number, show the transient growth of convection from small perturbations, and in all cases studied approach steady states. The steady-states obtained agree favorably with steady state solutions obtained by previous investigators.

Edward N. Lorenz, “Deterministic Nonperiodic Flow”, Journal of the Atmospheric Sciences, Vol. 20, pp. 130-141, 1963.
Abstract
Finite systems of deterministic ordinary differential equations may be designed to represent forced dissipative hydrodynamic flow. Solutions of these equations can be identified with trajectories in phase space. For those systems with bounded solutions, it is found that nonperiodic solutions are ordinarily unstable with respect to small modifications, so that slightly differing initial states can evolve into considerably different states. Systems with bounded solutions are shown to process bounded numerical solutions.
A simple system representing cellular convection is solved numerically. All of the numerical solutions are found to be stable, and almost all of them are nonperiodic.
The feasibility of very-long-range weather prediction is examined in the light of these results.

Charles Young, “Comment on “Deterministic Nonperiodic Flow””, Journal of Atmospheric the Sciences, Vol. 23, pp. 628-629, 1966.

Edward N. Lorenz, “Climate Change as a Mathematical Problem”, Journal of Applied Meteorology”, Vol. 9. No. 3, pp. 325-329, 1970.
Abstract
Formulating reasonable hypotheses regarding climate change requires physical insight and ingenuity, but subsequently testing these hypotheses demands quantitative computation. Many features of today’s climate have been reproduced by mathematical models (equations arranged for numerical solution by digital computers), similar to those used in weather prediction. Models currently in use generally predict only the atmosphere, and pre-specify the state of the environment (oceans, land surfaces, sun, etc.). Newer models, where certain environmental conditions enter as additional dependent variables, should be suitable for testing climatic-change hypotheses. Aspects of the atmosphere which play no role in these hypotheses may be highly simplified. A super-model where virtually all not-strictly-constant features of the atmosphere and its environment enter as variables may ultimately lead to an acceptable theory of climate change.

E. N. Lorenz, “Irregularity: A Fundamental Property of the Atmosphere”, Tellus, Vol. 36A, pp. 98-110, 1984.

8. Edward N. Lorenz, “Computational Chaos – A Prelude to Computational Instability”, Physica D, Vol. 35, pp. 299-317, 1989.
Abstract
Chaotic behavior sometimes occurs when difference equations used as approximations to ordinary differential equations are solved numerically with an excessively large time increment r. In two simple examples we find that, as r increases, chaos first sets in when an attractor A acquires two distinct points that map to the same point. This happens when A acquires slopes of the same sign, in a rectifying coordinate system, at two consecutive intersections with the critical curve. Chaotic and quasi-periodic behavior may then alternate within a range of rc before computational instability finally prevails. Bifurcations to and from chaos and transitions to computational instability are highly scheme-dependent, even among differencing schemes of the same order. Systems exhibiting computational chaos can serve as illustrative examples in more general studies of noninvertible mappings.

E. N. Lorenz, “Can Chaos and Intransitivity Lead to Interannual Variability?”, Tellus, Vol. 42A, pp. 378-389, 1990.
Abstract
We suggest that the atmosphere-ocean-earth system is unlikely to be intransitive, i.e., to admit two or more possible climates, any one of which, once established, will persist forever. Our reasoning is that even if the system would be intransitive if the external heating could be held fixed, say as in summer, the new heating patterns that actually accompany the advance of the seasons will break up any established summer circulation, and an alternative circulation may develop during the following summer, particularly if chaos has prevailed during the intervening winter. We introduce a very-low-order geostrophic baroclinic “general circulation” model, which may be run with or without seasonal variations of heating. Under perpetual summer conditions the model is intransitive, admitting either weakly oscillating or strongly oscillating westerly flow, while under perpetual winter conditions it is chaotic. When seasonal variations of heating are introduced, weak oscillations prevail through some summers and strong oscillations prevail through others, thus lending support to our original suggestion. We develop some additional properties of the model as a dynamical system, and we speculate as to whether its behavior has a counterpart in the real world.

Edward N. Lorenz, “An Attractor Embedded in the Atmosphere”, Tellus, Vol. 58A, pp. 425-429, 2006.
Abstract
A procedure for deriving a three-variable, one-level baroclinic model by truncating the familiar two-level quasigeostrophic model is described. The longitude, latitude and isobaric height of the low centre are introduced into the new model as alternative dependent variables. State space then becomes equivalent to geographical space, and the attractor becomes a structure within the atmosphere.

Edward N. Lorenz, “Computational Periodicity as Observed in a Simple System”, Tellus, Vol. 58A, pp. 549-557, 2006.
Abstract
When the exact time-dependent solutions of a system of ordinary differential equations are chaotic, numerical solutions obtained by using particular schemes for approximating time derivatives by finite differences, with particular values of the time increment _ , are sometimes stably periodic. It is suggested that this phenomenon be called computational periodicity.
A particular system of three equations with a chaotic exact solution is solved numerically with an Nth-order Taylor series scheme, with various values of N, and with values of _ ranging from near zero to just below the critical value for computational instability. When N = 1, the value of _ below which computational periodicity never appears is extremely small, and frequent alternations between chaos and periodicity occur within the range of _ . Computational periodicity occupies most of the range when N = 2 or 3, and about half when N = 4.
These solutions are compared with those produced by fourth-order Runge-Kutta and Adams-Bashforth schemes, and with numerical solutions of two other simple systems. There is some evidence that computational periodicity will more likely occur when the chaos in the exact solutions is not very robust, that is, if relatively small changes in the values of the constants can replace the chaos by periodicity.

Roger A. Pielke and Xubin Zeng, “Long-Term Variability of Climate”, Journal of the Atmospheric Sciences, Vol. 51, No. 1, pp. 155-159, 1994.
Abstract
In this research note, we address the following general question: In a nonlinear dynamical system (such as the climate system), can a known short-periodic variation lead to significant long-term variability? It is known from chaos studies (e.g., Lorenz 1991) that any perturbations in chaotic dynamic systems can lead to a red-noise spectrum; however, whether a significant long-term variability can be induced is unknown. To perform this study, an idealized nonlinear model developed by Lorenz (1984,1990) is used. The model and the results are presented in sections 2 and 3, respectively. Finally, the implications of our research to the understanding of the natural variability of the climate system due to internal dynamics will be discussed in section 4.

Roger A. Pielke, “Climate Prediction as an Initial Value Problem”, Bulletin of the American Meteorology Society, Vol. 79, pp. 2743-2746, 1998.

Generalized Lorenz Equation Systems

D. Roy and Z. E. Musielak, “Generalized Lorenz Models and their Routes to Chaos, II: Energy-Conserving Horizontal Mode Truncations”.
Abstract
Many attempts have been made to extend 3D Lorenz model in higher dimension to
describe 2D Rayleigh-Benard convection. 3D Lorenz model was originally derived by taking into account only three Fourier modes. We make a novel attempt here to generalize this 3D model to higher dimensions by picking up the fourier modes in consistent way. In this paper we have derived a generalized Lorenz by selecting horizontal modes that conserve energy in the dissipationless limit and lead to the systems that have bounded solutions. An interesting result is that the lowest-order generalized Lorenz model, which satisfies these criteria, is a 5D model and that its route to chaos is the same as that observed in the original 3D Lorenz model. A rather interesting result is that the onset of fully developed chaos is at low Rayleigh ratio.

D. Roy and Z. E. Musielak, “Generalized Lorenz Models and their Routes to Chaos, I: Energy-Conserving Vertical Mode Truncations”.
Abstract
A two-dimensional and dissipative Rayleigh-Benard convection can be approximated by Lorenz model, which was originally derived by taking into account only three Fourier modes. Numerous attempts have been made to generalize this 3D model to higher dimensions and several different methods of selecting Fourier modes have been proposed. In this paper, generalized Lorenz models with dimensions ranging from four to nine are constructed by selecting vertical modes that conserve energy in the dissipationless limit and lead to the systems that have bounded solutions. An interesting result is that the lowest-order generalized Lorenz model, which satisfies these criteria, is a 9D model and that its route to chaos is the same as that observed in the original 3D Lorenz model. The latter is in contradiction to some previous results that imply different routes to chaos in several generalized Lorenz systems. The discrepancy is explained by the fact that previously obtained generalized Lorenz models were derived by using different methods of selecting Fourier modes and, as result, most of them did not obey the principle of the conservation of energy in the dissipationless limit.

D. Roy and Z.E. Musielak, “Generalized Lorenz Models and their Routes to Chaos, III: Energy-conserving horizontal and vertical mode truncations”.
Abstract
To construct generalized Lorenz systems, higher-order modes in doubled Fourier expansions of a stream function and temperature variations must be considered. Selection of these modes is guided by the requirements that they conserve energy in the dissipationless limit and lead to systems that have bounded solutions. The previous study showed how to select the modes by using either vertical or horizontal mode truncations. In this paper, the most general method of horizontal and vertical mode truncations is presented and it is shown that the lowest-order generalized Lorenz system derived by this method is an eight dimensional system. An interesting result is that a route to chaos in this system is different than that observed in the original Lorenz model. Possible physical consequences of this result are discussed.

Effects of Numerical Methods

A. B. Gumel, “Removal of Contrived Chaos in Finite Difference Methods”, Intern. J. Computer Math., Vol. 79(9), pp. 1033-1041, 2002.
Abstract
A chaos-free numerical method will be developed for the solution of a system of non-linear initial-value problems (IVP’s) associated with the transmission dynamics of two HIV subtypes. It will be shown that integration of this four dimensional system with some standard numerical integrators like the fourth-order Runge–Kutta algorithm leads to scheme-dependent numerical instabilities.

R.M. Corless, C. Essex and M.A.H. Nerenberg, “Numerical Methods can Suppress Chaos”, Physics Letters A, Vol. 157, No. 1, pp. 27-36, 1991.
Abstract
Numerical methods for the solution of ordinary differential equations are one of the main tools used in the theoretical investigation of nonlinear continuous dynamical systems. These replace the continuous dynamical system under study by a discrete dynamical system that is then usually simulated on a digital computer. It is well known that such discrete dynamical systems may be chaotic even when the underlying continuous dynamical system is not chaotic. We here show that some numerical methods may produce discrete dynamical systems that are not chaotic, even when the underlying continuous dynamical system is thought to be chaotic. We find in this case that the transition to chaos from false stability mimics the transition to chaos that has been previously observed as parameters were changed in the Rössler system.

H.C. Yee and P.K. Sweby, “Nonlinear Dynamics and Numerical Uncertainties in CFD”, Invited review paper for Journal of Computational Physics (invitation from the Editor-in-Chief of Journal of Computational Physics Dr. Jerry Brackbill). See also NASA Technical Memorandum 110398, 1996.
Abstract
The application of nonlinear dynamics to improve the understanding of numerical uncertainties in computational fluid dynamics (CFD) is reviewed. Elementary examples in the use of dynamics to explain the nonlinear phenomena and spurious behavior that occur in numerics are given. The role of dynamics in the understanding of long time behavior of numerical integrations and the nonlinear stability, convergence, and reliability of using time-marching approaches for obtaining steady-state numerical solutions in CFD is explained. The study is complemented with examples of spurious behavior observed in CFD computations.

H. C. Yee, J. R. Torczynski, S. A. Morton, M. R. Visbal, and P. K. Sweby, “On Spurious Behavior of CFD Simulations”, Invited paper AIAA-97-1869 for the 13th AIAA CFD Conference, June 29 – July 2, 1997, Snowmass, Colorado. CONF-970666-2, also Sandia National Laboratory Report SAND97-1062C, 1997.
Abstract
Spurious behavior in underresolved grids and/or semi-implicit temporal discretizations for four computational fluid dynamics (CFD) simulations are studied. The numerical simulations consist of (a) a 1-D chemically relaxed nonequilibrium model, (b) the direct numerical simulation (DNS) of 2-D incompressible flow over a backward facing step, (c) a loosely-coupled approach doe a 2-D fluid-structure interaction, and (d) a 3-D compressible simulations of vortex breakdown in delta wings. Using knowledge from dynamical systems theory, various types of spurious behaviors that are numerical artifacts were systematically identified. These studies reveal the various possible dangers of misinterpreting numerical simulation of realistic complex flows that are constrained by the available computing power. In large scale computations underresolved grids, semi-implicit procedures, loosely-coupled implicit procedures, and insufficiently long time integration in DNS are most often unavoidable. Consequently, care must be taken in both computation and in interpretation of the numerical data. The results presented confirm the important role that dynamical system theory can play in understanding of the nonlinear behavior of numerical algorithms and in aiding the identification of the sources of numerical uncertainties in CFD.

H. C, Yee, P. K. Sweby, and D. F. Griffths, “Dynamical Approach Study of Spurious Steady-State Numerical Solutions of Nonlinear Differential Equations Part I: The ODE Connection and Its Implications for Algorithm Development in Computational Fluid Dynamics”, NASA Ames Research Center Report, Technical Memorandum 10820, 1990.
Abstract
The goal of this paper is to attempt to give some insight and guidelines on the application of nonlinear dynamic theory to the better understanding of steady-state numerical solutions and nonlinear instability in algorithm development for nonlinear differential equations that display genuinely nonlinear behavior in computational sciences and, in particular, computational fluid dynamics (CFD). This stems from the fact that, although the study of nonlinear dynamics and chaotic dynamics for nonlinear differential equations and for discrete maps have independently flourished rapidly for the last decade, there are very few investigators addressing the issue on the connection between the nonlinear dynamical behavior of the continuous systems and the corresponding discrete map resulting from finite-difference discretizations. This issue is especially vital for computational sciences since nonlinear differential equations in applied sciences can rarely be solved in closed form and it is often necessary to replace them by finite dimensional nonlinear discrete maps. In addition, it is also important to realize that these nonlinear discrete maps can exhibit a much richer range of dynamical behavior than their continuum counterparts. Furthermore, it is also very important to identify some of the implications of what happens when linear stability breaks down for problems with genuinely nonlinear behavior. Studies indicate that for relatively simple nonlinear ordinary differential (ODEs) and well-known time-discretization with modest step-sizes some schemes can converge to a spurious (false) steady-state solution in a deceptively “smooth” manner. In some instances, spurious steady states may appear below the linearized stability limit of the schemes, and consequently computation may lead to erroneous results. Our preliminary studies on partial differential equations (PDEs) also show that much of nonlinear dynamic (e.g. chaotic) phenomena have a direct relation for problems containing nonlinear source terms such as the reaction-diffusion, the reaction-convection or the reaction-convection-diffusion equations. Here our object is neither to provide theory nor to illustrate with realistic examples the connection of the dynamical behavior of practical PDEs with their discretized counterparts, but rather to give insight into the nonlinear features unconventional to this type of study and to concentrate on the fundamental ideas. Thus, in order to bring out the special properties, the illustrations center on simple scalar differential equation (DE) examples in which the exact solutions of the DEs are known.

H. C, Yee, P. K. Sweby, and D. F. Griffths, “On Spurious Steady-State Solutions of Explicit Runge-Kutta Schemes”, NASA Ames Research Center Report, Technical Memorandum 10819, 1990.
Abstract
The bifurcation diagram associated with the logistic equation v n+l = avn(1 – v n) is by now well known, as is its equivalence to solving the ordinary differential equation (ODE) u’ = au(1 – u) by the explicit Euler difference scheme. It has also been noted by Iserles that other popular difference schemes may not only exhibit period doubling and chaotic phenomena but also possess spurious fixed points. We investigate computationally and analytically Runge-Kutta schemes applied to both the equation u’ = au(1 – u) and the cubic equation u’ = au(1 -u)(b-u), contrasting their behaviour with the explicit Euler scheme. We note their spurious fixed points and periodic orbits. In particular we observe that these may appear below the linearized stability limit of the scheme, and, consequently computation may lead to erroneous results.

H. C. Yee and P. K. Sweby, “Dynamical Approach Study of Spurious Steady-State Numerical Solutions of Nonlinear Differential Equations II. The Dynamics of Numerics of Systems of ODEs and Its Connection to Finite Discretizations of Nonlinear PDEs”, NASA Applied Research Technical Report NASA Ames Research Center, RNR-92-008, March 1992.
Abstract
The objective of this paper is to study the dynamics of numerics of time discretizations for various systems of first-order autonomous nonlinear ordinary differential equations (ODEs) with known analytic solutions. The main emphasis is to gain a basic understanding of the difference in the dynamics of numerics between the scalars and systems of nonlinear ODEs. It is found that in addition to the phenomenon of stable spurious steady-state numerical solutions occurring below and above the linearized stability limit of the exact steady states, more complex phenomena such as stable spurious limit cycles, stable spurious higher dimensional tori and the changing type and stability of numerical asymptotes were observed. With the aid of a parallel Connection Machine (CM2), the complex behavior and sometimes fractal like structure of the associated basins of attraction of the various widely used time discretizations in computational fluid dynamics (CFD) are revealed and compared. The underlying purpose of this study is to set the baseline global asymptotic solution behavior of the schemes so that one can use them more wisely in other more complicated settings such as when nonlinear systems of partial differential equations (PDEs) for which the exact solutions are not known are encountered in nonlinear sciences and in particular in CFD. This is especially important when there are no experimental data for comparison and/or when the numerical solution indicates a new flow structure not easily understood. The results of this investigation can be used as an explanation for possible causes of error, and slow convergence and nonconvergence of steady-state numerical solution when using the time-dependent approach for nonlinear hyperbolic or parabolic PDEs. The knowledge gained can also aid the construction of appropriate iteration methods, relaxation procedures, or preconditioners for convergence acceleration strategies in numerically solved boundary-value problems of nonlinear PDEs, since most of these procedures can be viewed as approximations of time-dependent PDEs. It can also enhance the understanding of flow patterns in 2-D and 3-D flow visualizations of numerical data.

Professor J. C. Sprott’s Web Site shows computational chaos and computational periodicity as artifacts of applications of the explicit Euler method to simple 2-dimensional systems of ODEs.

L. D. Cloutman, “A Note on the Stability and Accuracy of Finite Difference Approximations to Differential Equations”, Lawrence Livermore National Laboratory, UCRL-ID-125549, 1996.
Abstract
There are many finite difference approximations to ordinary and partial differential equations, and these vary in their accuracy and stability properties. We examine selected commonly used methods and illustrate their stability and accuracy using both stability analysis and numerical examples. We find that the formal order of accuracy alone gives an incomplete picture of the accuracy of the method. Specifically, the Adams-Bashforth and Crank-Nicholson methods are shown to have some undesirable features for both ordinary and partial differential equations.

L. D. Cloutman, “Chaos and Instabilities in Finite Difference Approximations to Nonlinear Differential Equations”, Lawrence Livermore National Laboratory, UCRL- ID-131333, 1998.
Abstract
The numerical solution of time-dependent ordinary and partial differential equations
by finite difference techniques is a common task in computational physics and engineering. The rate equations for chemical kinetics in combustion modeling are an important example. They not only are nonlinear, but they tend to be stiff, which makes their solution a challenge for transient problems. We show that one must be very careful how such equations are solved In addition to the danger of large time-marching errors, there can be unphysical chaotic solutions that remain numerically stable for a range of time steps that depends on the particular finite difference method used. We point out that the solutions of the finite difference equations converge to those of the differential equations only in the limit as the time step approaches zero for stable and consistent finite difference approximations. The chaotic behavior observed for finite time steps in some nonlinear difference equations is unrelated to solutions of the differential equations, but is connected with the onset of numerical instabilities of the finite difference equations. This behavior suggests that the use of the theory of chaos in nonlinear iterated maps may be useful in stability analysis of finite difference approximations to nonlinear differential equations, providing more stringent time step limits than the formal linear stability analysis that tests only for unbounded solutions. This observation implies that apparently stable numerical solutions of nonlinear differential equations by finite difference techniques may in fact be contaminated (if not dominated) by nonphysical chaotic parasitic solutions that degrade the accuracy of the numerical solution. We demonstrate this phenomenon with some solutions of the logistic equation and a simple two-dimensional computational fluid dynamics example.

Predator-Prey Equation Systems

J. A. Vano, J. C. Wildenberg, M. B. Anderson, J. K. Noel and J. C. Sprott, “Chaos in low-dimensional Lotka-Volterra models of Competition”, Nonlinearity, Vol. 19, pp. 2391-2404, 2006. doi:10.1088/0951-7715/19/10/006
Abstract
The occurrence of chaos in basic Lotka-Volterra models of four competing species is studied. A brute-force numerical search conditioned on the largest Lyapunov exponent (LE) indicates that chaos occurs in a narrow region of parameter space but is robust to perturbations. The dynamics of the attractor for a maximally chaotic case are studied using symbolic dynamics, and the question of self-organized critical behaviour (scale-invariance) of the solution is considered.

Khaled Batiha, “Numerical Solutions of the Multispecies Predator-Prey Model by Variational Iteration Method”, Journal of Computer Science, Vol. 3, No. 7, pp. 523-527, 2007.
Abstract
The main objective of the current work was to solve the multispecies predator-prey model. The techniques used here were called the variational iteration method (VIM) and the Adomian decomposition method (ADM). The advantage of this work is twofold. Firstly, the VIM reduces the computational work. Secondly, in comparison with existing techniques, the VIM is an improvement with regard to its accuracy and rapid convergence. The VIM has the advantage of being more concise for analytical and numerical purposes. Comparisons with the exact solution and the fourth-order Runge-Kutta method (RK4) show that the VIM is a powerful method for the solution of nonlinear equations.

Adomian Decomposition Method (ADM) Applied to Solving ODEs

G. Adomian, “Decomposition Solution For Duffing And Van Der Pol Oscillators”, Internat. J. Math. & Math. Sci. Vol. 9, No. 4, pp. 731-732, 1986.
Abstract
The decomposition method is applied to solve the Duffing and Van der Pol oscillators without customary restrictive assumptions [I-4] and without resort to perturbation methods.

G. Adomian, “A Review of the Decomposition Method in Applied Mathematics”, Journal Of Mathematical Analysis And Applications , Vol. 135, pp. 501-544 , 1988.
Abstract
The decomposition method can be an effective procedure for analytical solution of a wide class of dynamical systems without linearization or weak nonlinearity assumptions, closure approximations, perturbation theory, or restrictive assumptions on stochasticitiy.

Yonggui Zhu, Qianshun Chang, and Shengchang Wu, “A New Algorithm For Calculating Adomian Polynomials”, Applied Mathematics and Computation, Vol. 169, pp. 402-416, 2005. doi:10.1016/j.amc.2004.09.082.
Abstract
In this paper, a new algorithm for calculating Adomian polynomials for nonlinear operators will be established by parameterization. The algorithm requires less formula than the previous method developed by Adomian [Nonlinear Stochastic Operator Equations, Academic Press, 1986, G. Adomian, R. Rach, On composite nonlinearities and decomposition method. J. Math. Anal. Appl. 113 (1986) 504-509, G. Adomian, Applications of Nonlinear Stochastic Systems Theory to Physics, Kluwer, 1988].Many forms of nonlinearity will be studied to illustrate the new algorithm. The new algorithm will be extended to calculate Adomian polynomials for nonlinearity of several variables.

J. Biazar , M. Ilie, and A. Khoshkenar, “An Improvement To An Alternate Algorithm For Computing Adomian Polynomials In Special Cases”, Applied Mathematics and Computation, Vol. 173, pp. 582-592, 2006. doi:10.1016/j.amc.2005.04.052
Abstract
In a published article entitled “An alternate algorithm for computing Adomian polynomials, in special cases” by the first author and some of his collaborators, the algorithm were applied in special cases which restricts the application of the algorithm. In this article new ideas are suggested for applying the algorithm in any cases without any more assumptions. Some examples are presented to illustrate the applications of the algorithm in more comprehensive situations, and also to show the ability and simplicity of the algorithm.

Applications of the Adomian Decomposition Method (ADM)

Jafar Biazar and Leila Farrokhi, “Solving The Lorenz System By Adomian Decomposition Method And Comparing The Results With The Results Of Runge-Kutta Method”, Far East Journal of Applied Mathematics, Vol. 24, No. 1, pp. 39-44, 2006.
Abstract
In this article the Lorenz system is solved by two different methods. The first one is Adomian decomposition method which is a well-known method for solving functional equations nowadays. The second one is Runge-Kutta method which is a numerical method for solving differential equations and systems of differential equations. Theoretical considerations are discussed.

S. Guellal, P. Crimalt, and Y. Cherruault, “Numerical Study of Lorenz’s Equation by the Adomian Method”, Computers Math. Applic. Vol. 33, No. 3, pp. 25-29, 1997
Abstract
In this paper, we use the decomposition method (of Adomian) for solving differential systems coming from physics (meteorology). We also give a comparison between the Runge-Kutta method and the decomposition technique. Furthermore, we reconfirm the famous “Butterfly effect.”

I. Hashim, M.S.M. Noorani, R. Ahmad, S.A. Bakar, E.S. Ismail, and A.M. Zakaria, “Accuracy Of The Adomian Decomposition Method Applied To The Lorenz System”, Chaos, Solitons and Fractals, Vol. 28, pp. 1149-1158, 2006.
Abstract
In this paper, the Adomian decomposition method (ADM) is applied to the famous Lorenz system. The ADM yields an analytical solution in terms of a rapidly convergent infinite power series with easily computable terms. Comparisons between the decomposition solutions and the fourth-order Runge-Kutta (RK4) numerical solutions are made for various time steps. In particular we look at the accuracy of the ADM as the Lorenz system changes from a non-chaotic system to a chaotic one.

Porous Media

Peter Vadasza and Shmuel Olek, “Convergence and Accuracy of Adomian’s Decomposition Method for the Solution of Lorenz Equations”, International Journal of Heat and Mass Transfer, Vol. 43, pp. 1715-1734, 2000a.
Abstract
The convergence and accuracy of Adomian’s decomposition method of solution is analyzed in the context of its application to the solution of Lorenz equations which govern at lower order the convection in a porous layer (or respectively in a pure fluid layer) heated from below. Adomian’s decomposition method provides an analytical solution in terms of an infinite power series and is applicable to a much wider range of heat transfer problems. The practical need to evaluate the solution and obtain numerical values from the infinite power series, the consequent series truncation, and the practical procedure to accomplish this task, transform the analytical results into a computational solution evaluated up to a finite accuracy. The analysis indicates that the series converges within a sufficiently small time domain, a result that proves to be significant in the derivation of the practical procedure to compute the infinite power series. Comparison of the results obtained by using Adomian’s decomposition method with corresponding results obtained by using a numerical Runge-Kutta-Verner method show that both solutions agree up to 12-13 significant digits at subcritical conditions, and up to 8-9 significant digits at certain supercritical conditions, the critical conditions being associated with the loss of linear stability of the steady convection solution.
The difference between the two solutions is presented as projections of trajectories in the state space, producing similar shapes that preserve under scale reduction or magnification, and are presumed to be of a fractal form.

Peter Vadasz and Shmuel Olek, “Weak Turbulence and Chaos for Low Prandtl Number Gravity Driven Convection in Porous Media”, Transport in Porous Media, Vol. 37, pp. 69-91, 1999.
Abstract
Low Prandtl number convection in porous media is relevant to modern applications of transport phenomena in porous media such as the process of solidification of binary alloys. The transition from steady convection to chaos is analysed by using Adomian’s decomposition method to obtain an analytical solution in terms of infinite power series. The practical need to evaluate the solution and obtain numerical values from the infinite power series, the consequent series truncation, and the practical procedure to accomplish this task, transform the analytical results into a computational solution evaluated up to a finite accuracy. The solution shows a transition from steady convection to chaos via a Hopf bifurcation producing a ‘solitary limit cycle’ which may be associated with an homoclinic explosion. This occurs at a slightly subcritical value of Rayleigh number, the critical value being associated with the loss of linear stability of the steady convection solution. Periodic windows within the broad band of parameter regime where the chaotic solution persists are identified and analysed. It is evident that the further transition from chaos to a high Rayleigh number periodic convection occurs via a period halving sequence of bifurcations.

Peter Vadasz, “Local and Global Transitions to Chaos and Hysteresis in a Porous Layer Heated from Below”, Transport in Porous Media, Vol. 37, pp. 213-245, 1999.
Abstract
The routes to chaos in a fluid saturated porous layer heated from below are investigated by using the weak nonlinear theory as well as Adomian’s decomposition method to solve a system of ordinary differential equations which result from a truncated Galerkin representation of the governing equations. This representation is equivalent to the familiar Lorenz equations with different coefficients which correspond to the porous media convection. While the weak nonlinear method of solution provides significant insight to the problem, to its solution and corresponding bifurcations and other transitions, it is limited because of its local domain of validity, which in the present case is in the neighbourhood of any one of the two steady state convective solutions.
On the other hand, the Adomian’s decomposition method provides an analytical solution to the problem in terms of infinite power series. The practical need to evaluate numerical values from the infinite power series, the consequent series truncation, and the practical procedure to accomplish this task transform the otherwise analytical results into a computational solution achieved up to a finite accuracy. The transition from the steady solution to chaos is analysed by using both methods and their results are compared, showing a very good agreement in the neighbourhood of the convective steady solutions.
The analysis explains previously obtained computational results for low Prandtl number convection in porous media suggesting a transition from steady convection to chaos via a Hopf bifurcation, represented by a solitary limit cycle at a sub-critical value of Rayleigh number. A simple explanation of the well known experimental phenomenon of Hysteresis in the transition from steady convection to chaos and backwards from chaos to steady state is provided in terms of the present analysis results.

Peter Vadasz, “Heat transfer Regimes and Hysteresis in Porous Media Convection”, ASME Journal of Heat Transfer, Vol. 123, pp. 145-156, 2001. Doi 10.1115/1.1336505.
Abstract
Results of an investigation of different heat transfer rgimes in porous media convection are presented using a truncated Galerkin representation of the governing equations that yields the familiar Lorenz equations ofr the variation of the amplitude in the time domain. The solution to this system is obtained analytically by using a weal non-linear analysis and computationally by using Adomian’s decomposition method. Expressions for the average Nusselt number are derived for steady, periodic, as well as weak-turbulent (temporal-chaotic) convection. The phenomena of Hysteresis in the transition from steady to weak-turbulent convection, and backwards, is particularly investigated, identifying analytically its mechanism, which is confirmed by the computational results. While the post-transient chaotic solution in terms of the dependent variables is very sensitive to the initial conditions, the affinity of the averaged values of these variables to initial conditions is very weak. Therefore, long-term predictability of these averaged variables, and in particular the Nusselt number, becomes possible, a result of substantial practical significance. Actually, the only impact that the transition to chaos causes on the predicted results in terms of the average heat flux is a minor loss of accuracy. Therefore, the predictability of the results in the sense of the averaged heat flux is not significantly affected by the transition from steady to weak-turbulent convection. The transition point is shown to be very sensitive to a particular scaling of the equations, which leads the solution to an invariant value of steady-state for sub-transitional conditions, a result that affects the transition point in some cases.

Peter Vadasz and Shmuel Olek, “Route to Chaos for Moderate Prandtl Number Convection in a Porous Layer Heated from Below”, Transport in Porous Media, Vol. 41, pp. 211-239, 2000.
Abstract
The route to chaos for moderate Prandtl number gravity driven convection in porous media is analysed by using Adomian’s decomposition method which provides an accurate analytical solution in terms of infinite power series. The practical need to evaluate numerical values from the infinite power series, the consequent series truncation, and the practical procedure to accomplish this task, transform the otherwise analytical results into a computational solution achieved up to a desired but finite accuracy. The solution shows a transition to chaos via a period doubling sequence of bifurcatons at a Rayleigh number value far beyond the critical value associated with the loss of stability of the convection steady solution. This result is extremely distinct from the sequence of events leading to chaos in low Prandtl number convection in porous media, where a sudden transition from steady convection to chaos associated with an homoclinic explosion occurs in the neighbourhood of the critical Rayleigh number (unless mentioned otherwise by ‘the critical Rayleigh number’ we mean the value associated with the loss of stability of the convection steady solution). In the present case of moderate Prandtl number convection the homoclinic explosion leads to a transition from steady convection to a period-2 periodic solution in the neighbourhood of the critical Rayleigh number. This occurs at a slightly sub-critical value of Rayleigh number via a transition associated with a period-1 limit cycle which seem to belong to the sub-critical Hopf bifurcation around the point where the convection steady solution looses its stability. The different regimes are analysed and periodic windows within the chaotic regime are identified. The significance of including a time derivative term in Darcy’s equation when wave phenomena are being investigated becomes evident from the results.

Lotka-Volterra Species System

Shmuel Olek, “An Accurate Solution To The Multispecies Lotka-Volterra Equations”, SIAM Review, Vol. 36, No. 3, pp. 480-488, September 1994.
Abstract
The decomposition method is applied to the solution of the Lotka-Volterra equations modelling the dynamic behaviour of an arbitrary number of species. The analytical solution derived is an infinite power series for each species, where the n term is given by a recurrence relation. As particular examples, the cases of one, two, and three species are considered. For these cases, comparisons between the present semi-analytical solution and a fully numerical solution (or an exact one for one species) show excellent agreement.

The Chen ODE System

G. Chen and T. Ueta, “Yet another chaotic attractor”, International Journal Bifurcation and Chaos, Vol. 9, No. 7, pp. 1465-1466, 1999.
Abstract

M.S.M. Noorani, I. Hashim, R. Ahmad, S.A. Bakar, E.S. Ismail, and A.M. Zakaria, “Comparing Numerical Methods For The Solutions Of The Chen System”, Chaos, Solitons and Fractals, Vol. 32, pp. 1296-1304, 2007.
Abstract
In this paper, the Adomian decomposition method (ADM) is applied to the Chen system which is a three-dimensional system of ODEs with quadratic nonlinearities. The ADM yields an analytical solution in terms of a rapidly convergent infinite power series with easily computable terms. Comparisons between the decomposition solutions and the classical fourth-order Runge–Kutta (RK4) numerical solutions are made. In particular we look at the accuracy of the ADM as the Chen system changes from a non-chaotic system to a chaotic one. To highlight some computational difficulties due to a high Lyapunov exponent, a comparison with the Lorenz system is given.

Non-Standard Finite Difference Methods

Christophe Letellier, Saber Elaydi, Luis A. Aguirre, and Aziz Alaoui, “Difference Equations Versus Differential Equations, A Possible Equivalence For The Rossler System?”, Physica D 195 (2004) 29-49
Abstract
When a set of nonlinear differential equations is investigated, most of time there is no analytical solution and only numerical integration techniques can provide accurate numerical solutions. In a general way the process of numerical integration is the replacement of a set of differential equations with a continuous dependence on the time by a model for which the time variable is discrete. In numerical investigations a fourth-order Runge–Kutta integration scheme is usually sufficient. Nevertheless, sometimes a set of difference equations may be required and, in this case, standard schemes like the forward Euler, backward Euler or central difference schemes are used. The major problem encountered with these schemes is that they offer numerical solutions equivalent to those of the set of differential equations only for sufficiently small integration time steps. In some cases, it may be of interest to obtain difference equations with the same type of solutions as for the differential equations but with significantly large time steps. Nonstandard schemes as introduced by Mickens [Nonstandard Finite Difference Models of Differential Equations, World Scientific, 1994] allow to obtain more robust difference equations. In this paper, using such nonstandard scheme, we propose some difference equations as discrete analogues of the Rössler system for which it is shown that the dynamics is less dependent on the time step size than when a nonstandard scheme is used. In particular, it has been observed that the solutions to the discrete models are topologically equivalent to the solutions to the Rössler system as long as the time step is less than the threshold value associated with the Nyquist criterion.

Membranes

David Terman, “Chaotic Spikes arising from a Model of Bursting in Membranes”, SIAM J. Appl. Math., Vol. 51, No. 5, pp. 1418-1450, October 1991.
Abstract
A class of differential equations that model electrical activity in pancreatic beta cells is considered. It is demonstrated that these equations must give rise to both bursting solutions and, for different values of the parameters, continuous spiking. We also consider how the number of spikes per burst changes as parameters in the equations are varied. This transition may be continuous, in which case the period of the bursting solution increases significantly and then decreases. Hence, small perturbations may cause macroscopic changes in the bursting solution. This transition may also give rise to chaotic dynamics due to the existence of a Smale horseshoe.

November 15, 2007 - Posted by Dan Hughes | Chaos and Lorenz | Chaos, Chaos and Lorenz, Numerical methods Verification, ODEs