The definition of stiff ODE system

Question

Consider an IVP for ODE system $y'=f(x,y)$, $y(x_0)=y_0$. Most commonly this problem is considered stiff when Jacobi matrix $\frac{\partial f}{\partial y}(x_0,y_0)$ has both eigenvalues with very large negative real part and eigenvalues with very small negative real part (I consider only the stable case).

On the other hand, in the case of just one equation, for example Prothero-Robinson equation $y'=\lambda y + g'+\lambda g$, it is called stiff when $\lambda\ll -1$.

So there are two questions:

1. Why small eigenvalues are included in the definition of stiffness for ODE systems? I believe that the presence of only very large negative real parts is quite enough for system to be stiff, because this makes us use small timesteps for explicit methods.

2. Yes, I know that the most common stiff problems (e. g. arising from parabolic PDEs) do have both large and small eigenvalues. So the second question: is there a good natural example of large stiff system without very small eigenvalues (or alternatively with mild ratio $\lambda_{\max}/\lambda_{\min}$)?

---

OK, let's modify the question. Consider two two-dimensional linear ODE systems: first with eigenvalues {-1000000,-0.00000001} and second with {-1000000,-999999}. As for me, both of them are stiff. But if we consider stiffness ratio definition, the second system is not. The main question: why stiffness ratio is considered at all?

And the second part of question is still important, lets paraphrase it: I'm looking for a "natural" large ODE system with big negative eigenvalues and mild stiffness ratio (not greater than, say, 100).

Answer 1

Stiffness involves some separation of scales. In general, if you are interested in phase of the fastest mode in the system, then you have to resolve it and the system is not stiff. But frequently, you are interested in the long-term dynamics of a "slow manifold" rather than the precise rate at which a solution off the slow manifold approaches it.

Chemical reactions and reacting flows are common examples of stiff systems. The van der Pol oscillator is a common benchmark problem for ODE integrators that has a tunable stiffness paramater.

An ocean is another example that is perhaps helpful to visualize. Tsunamis (surface gravity waves) travel at the speed of an airplane and produce complex wave structure, but dissipate over long time scales and are mostly inconsequential to the long-term dynamics of the ocean. Eddies, or the other hand, travel about 100 times slower at quite pedestrian speeds, but cause mixing and transport temperature, salinity, and biogeochemical tracers that are relevant. But the same physics that propagates a surface gravity wave also supports an eddy (a quasi-equilibrium structure), so eddy velocity, path under Coriolis, and rate of dissipation are dependent on the gravity wave speed. This presents an opportunity for a time integration scheme designed for stiff systems to step over the time scale of the gravity wave and only resolve the relevant dynamical time scales. See Mousseau, Knoll, and Reisner (2002) for discussion of this problem with a comparison of splitting and fully implicit time integration schemes.

Related: When should implicit methods be used in the integration of hyperbolic PDEs?

Note that diffusive processes are usually considered to be stiff because the fastest time scale in the discrete system is mesh-dependent, scaling with $(\Delta x)^2$, but the time scale of the relevant physics is mesh independent. In fact, the fastest time scales for a given mesh represent spatially local relaxation to the slower manifold on which longer spatial scales evolve, so implicit methods can be very accurate even in strong norms despite not resolving the fastest scales.

Answer 2

Part 1

Small eigenvalues are not included in the definition of stiffness for ODE (initial value problem) systems. There is no satisfying definition of stiffness that I know of, but the best definitions I've come across are:

If a numerical method with a finite region of absolute stability, applied to a system with any initial conditions, is forced to use in a certain interval of integration a step length which is excessively small in relation to the smoothness of the exact solution in that interval, then the system is said to be stiff in that interval. (Lambert, J. D. (1992), Numerical Methods for Ordinary Differential Systems, New York: Wiley.)

An IVP [initial value problem] is stiff in some interval $[0,b]$ if the step size needed to maintain stability of the forward Euler method is much smaller than the step size required to represent the solution accurately. (Ascher, U. M. and Petzold, L. P. (1998), Computer Methods for Ordinary Differential Equations and Differential-Algebraic Equations, Philadelphia:SIAM.)

Stiff equations are equations where certain implicit methods, in particular BDF, perform better, usually tremendously better, than explicit ones. (C. F. Curtiss & J. O. Hirschfelder (1952): Integration of stiff equations. PNAS, vol. 38, p. 235-243)

The Wikipedia article on stiff equations goes on to attribute the following "statements" to Lambert:

1. A linear constant coefficient system is stiff if all of its eigenvalues have negative real part and the stiffness ratio is large.

2. Stiffness occurs when stability requirements, rather than those of accuracy, constrain the step length. [Note that this "observation" is essentially the definition from Ascher and Petzold.]

3. Stiffness occurs when some components of the solution decay much more rapidly than others.

Each of these observations has counterexamples (though admittedly I couldn't produce one off the top of my head).

Part 2

Probably the best example I could come up with would be integrating any sort of large combustion reaction system in chemical kinetics under conditions that result in ignition. The system of equations will be stiff until ignition, and then it will no longer be stiff because the system has passed an initial transient. The ratio of largest to smallest eigenvalue should not be large except around the ignition event, though such systems tend to confound stiff integrators unless you set exceedingly strict integration tolerances.

The book by Hairer and Wanner also gives several other examples in its first section (Part IV, section 1) that illustrate many other examples of stiff equations. (Wanner, G., Hairer, E., Solving Ordinary Differential Equations II: Stiff and Differential-Algebraic Problems (2002), Springer.)

Finally, it's worth pointing out the observation of C. W. Gear:

Although it is common to talk about "stiff differential equations," an equation per se is not stiff, a particular initial value problem for that equation may be stiff, in some regions, but the sizes of these regions depend on the initial values and the error tolerance. (C. W. Gear (1982): Automatic detection and treatment of oscillatory and/or stiff ordinary differential equations. In: Numerical integration of differential equations, Lecture notes in Math., Vol. 968, p. 190-206.)

Answer 3

The absolute magnitude of the eigenvalues (in a linear, autonomous problem) alone has no meaning at all; it's an artifact of the units you choose to express the problem in.

The chain of comments is getting out of control, so I'm making this an answer. I'm not going to answer the full question; as I said, see wikipedia or the other answers here. I'm just answering the bit that says

Consider two two-dimensional linear ODE systems: first with eigenvalues {-1000000,-0.00000001} and second with {-1000000,-999999}. As for me, both of them are stiff. But if we consider stiffness ratio definition, the second system is not. The main question: why stiffness ratio is considered at all?

Okay, let's consider an example of the second case:

$$y_1'(t) = -1000000 y_1(t)$$ $$y_2'(t) = -999999 y_2(t)$$

Now let's consider a time variable with different units: $t^* = 1000000\cdot t$. Then calculus reveals that

$$y_1'(t^*) = - y_1(t^*)$$ $$y_2'(t^*) = -0.999999 y_2(t^*)$$

Note 1: I chose a diagonal system to make it totally obvious, but if you try it with another system with those eigenvalues, you'll see the same effect, since multiplying a matrix by a constant multiplies its eigenvalues by the same constant.

Note 2: I'm not even discussing here whether the system is stiff. I'm just pointing out that your proposed definition of stiffness (i.e., any problem with $|\lambda|\gg1$) makes no sense, since it would mean that stiffness depends on the units in which I choose to express the problem.

Answer 4

In fact Jed Brown has cleared the question for me. What I'm doing now is just putting his words in the context.

1. Both 2d linear ODE systems from above are stiff (i. e. hard to solve with explicit methods) on relative big time intervals (e. g. [0,1]).

2. The linear systems with large stiffness ratio can be considered "more stiff" because most likely one needs to integrate them on large time interval. This is due to slow components corresponding to the smallest eigenvalues: the solution slowly tends to the steady state, and this steady state is usually important to reach.

3. On the other hand, integration of systems with small stiffness ratio on large intervals is not interestng: in this case the steady state is reached very fast and we can just extrapolate it.

Thanks to all for this discussion!