Stiffness ratio evolution for odes

Question

I had a general query related to calculation of stiffness ratio evolution for a set of coupled odes over a certain time interval. My question is, while calculating the eigen values from the jacobian matrix do we need to know the exact solution at the different time points or can we use solutions approximated by an integration method?

Answer 1

I made a number of comments on stiffness in response to another question that are relevant here; I won't repeat them here in their entirety. The most important point is that stiffness ratio is only correlated with stiffness, and not a definition of stiffness. Consequently, it does not seem critical to calculate eigenvalues of a Jacobian matrix based on an exact solution; calculating approximate eigenvalues (for example, via a Krylov subspace method that provides eigenvalue estimates, like GMRES) of a reasonably accurate approximate solution, then calculating the stiffness ratio is probably fine for diagnostic purposes.

Another criterion that could work, if you can calculate the necessary data, is to compare the time steps your integrator is taking at the end of each iterate to the maximum stable time step it could take. If stability is limiting time step size rather than accuracy, your system is probably stiff. (Again, this sort of behavior tends to occur in stiff systems, but it is not a definition.)

Answer 2

As Geoff says, stiffness is more pragmatic in its definition and use than mathematical. A good practical way to define it is that a stiff problem is one where simple explicit (i.e. Runge-Kutta) methods need a really small time step to be stable, so your simulation never seems to finish. One way to see this is if the eigenvalues of the Jacobian are large, but if you can't tell by eye, just run the simulation with an adaptive Runge-Kutta method (i.e. ode45). If it's taking for ever and seems to get stuck, that means you have to move to a stiff solver (usually more work).

[If you don't have an easy adaptive method to implement, just try a Runge-Kutta method with some reasonable step sizes. If the solution keeps blowing up, it's unstable, which means go to implicit methods].