Algorithms for stiff problems

In the following example, the difference between an algorithm for stiff problems and one for nonstiff problems is highlighted. 

Valorani, M., Goussis, D.A. Explicit time-scale splitting algorithm for stiff problems autoignition of gaseous mixtures behind a steady shock. J. Comput. Phys. 169,44—79 (2001) Valorani, M., Najm, FI.N., Goussis, D.A. CSP analysis of a transient flame-vortex interaction time scales and manifolds. Combust. Flame 134, 35-53 (2003)... 

The second approach is a fractional-step method we call asymptotic timestep-splitting. It is developed by consideration of the specific physics of the problem being solved. Stiffness in the governing equations can be handled "asymptotically" as well as implicitly. The individual terms, including those which lead to the stiff behavior, are solved as independently and accurately as possible. Examples of such methods include the Selected Asymptotic Integration Method (4,5) for kinetics problems and the asymptotic slow flow algorithm for hydrodynamic problems where the sound speed is so fast that the pressure is essentially constant (6, 2). ... 

The extension of Gillespie s algorithm to spatially distributed systems is straightforward. A lattice is used to represent binding sites of adsorbates, which correspond to local minima of the potential energy surface. The discrete nature of KMC coupled with possible separation of time scales of various processes could render KMC inefficient. The work of Bortz et al. on the n-fold or continuous time MC CTMC) method can lead to computational speedup of the KMC method, which, however, has been underutilized most probably because of its difficult implementation. This method classifies all atoms in a finite number of classes according to their transition probability. Probabilities are computed a priori and each event is successful, in contrast to the Metropolis method (and other null event algorithms) whose fraction of unsuccessful (null) events increases drastically at low temperatures and for stiff problems. In conjunction with efficient search within a class and dynamic variation of atom coordi-nates, " the CPU time can be practically independent of lattice size. After each event, the time is incremented by a continuous amount. 

The Gear algorithms are stable for stiff problems, whereas the Adams-Moulton are unstable with orders larger than 2. 

The Runge-Kutta algorithm cannot handle so-called stiff problems. Computation times are astronomical and thus the algorithm is useless, for that class of ordinary differential equations, specialised stiff solvers have been developed. In our context, a system of ODEs sometimes becomes stiff if it comprises very fast and also very slow steps and/or very high and very low concentrations. As a typical example we model an oscillating reaction in The Belousov-Zhabotinsky (BZ) Reaction (p.95). 

The main advantage of the implicit algorithms is that they do not become numerically unstable. Very large step sizes can be taken without having to worry about the instability problems that plague the explicit methods. Thus, the implicit methods are very useful for stiff systems. 

It is clear from this expression that the method is stable for all A. and all time steps hn that is, the method is unconditionally stable (for linear problems). A consequence of the strong stability is that the time step can be chosen primarily to maintain accuracy. In the slowly varying regions of stiff problems, the time steps can be very large compared those required to maintain stability for an explicit algorithm. 

A general numerical algorithm of the boundary layer type for stiff systems of differential equations has been proposed by Miranker [173] and applied to a few kinetic problems by Aiken and Lapidus [174,175]. The principle of the method will be briefly described in the case of the following system of differential equations, involving stiff variable x and non-stiff variable y. 

Over the past ten years the numerical simulation of the behavior of complex reaction systems has become a fairly routine procedure, and has been widely used in many areas of chemistry, [l] The most intensive application has been in environmental, atmospheric, and combustion science, where mechanisms often consisting of several hundred reactions are involved. Both deterministic (numerical solution of mass-action differential equations) and stochastic (Monte-Carlo) methods have been used. The former approach is by far the most popular, having been made possible by the development of efficient algorithms for the solution of the "stiff" ODE problem. Edelson has briefly reviewed these developments in a symposium volume which includes several papers on the mathematical techniques and their application. [2]... 

In ref. 144 the author presents the construction of a non-standard explicit algorithm for initial-value problems. The order of the developed method is two and also is A-stable. The new proposed method is proven to be suitable for solving different kind of initial-value problems such as non-singular problems, singular problems, stiff problems and singularly perturbed problems. Some numerical experiments are considered in order to check the behaviour of the method when applied to a variety of initial-value problems. 

For reasons that will be explained in due course, a substitution iterative method is adopted when the system to be solved is nonstiff and, in this case, the algorithm adopted belongs to the Adams-Moulton family. Conversely, in stiff problems, the nonlinear system is solved using the Newton method and the algorithm belongs to the Gear family. 

To implement a multivalue algorithm for solving stiff problems, it is essential to take great care regarding efficiency. The key point with stiff problems is the solution of the nonlinear system (2.232) ... 

For the stiff problems, based on multivalue algorithms of the Gear family ... 

The authors used five algorithms based on the above mentioned techniques. They developed MatLab and Fortran versions of the above formulae and they compared the accuracy and computational efficiency. They used in both cases fixed step size procedure for equality of conditions of all implementations. They applied these methodologies on real problems of sciences and engineering and they expressed the advantages of the proposed algorithms, especially when they are integrating stiff problems. 

---