Is there a stiff differential equation that cannot be solved by the Runge-Kutta 4th order method, but which has an analytical solution for testing?
$\begingroup$
$\endgroup$
Add a comment
|
2 Answers
$\begingroup$
$\endgroup$
1
How about "impractical"? Here's an example from Numerical Recipes in C, Second Ed., Sec. 16.6:
$$u'=998 u+1998 v$$
$$v'=-999 u -1999 v$$
where $u(0) = 1$ and $v(0) = 0$.
Then
$$u(x) = 2 e^{-x}- e^{-1000 x}$$
$$v(x) = -e^{-x} + e^{-1000 x}$$
Runge-Kutta 4th order would require a stepsize of less than $0.001$ for any accuracy, which is obviously aggravating for the one of the solutions which could be acquired with a far coarser step.
answered Apr 6, 2013 at 14:42
-
$\begingroup$ Thanks, I like your equations and the book that was referenced, but if you calculate eigenvalues of the Jacobian, you will find $ \lambda_1 = 2560 $ and $\lambda_2 = -1559$. A system of differential equations is stiff, when \begin{equation} S = \dfrac {max_{1 \leqslant i \leqslant n} \vert{Re(-\lambda_i)}\vert} {min_{1 \leqslant i \leqslant n} \vert{Re(-\lambda_i)}\vert} >> 1. \end{equation} For this equations $S = 1.64$, please correct me if I wrong. $\endgroup$– WarlockApr 7, 2013 at 8:29
$\begingroup$
$\endgroup$
5
Cleve Moler, in this note, gives an innocuous-looking DE he attributes to Larry Shampine that models flame propagation. The differential equation is
$$y^\prime=y^2-y^3$$
with initial condition $y(0)=\frac1{h}$, and integrated over the interval $[0,2h]$. The exact solution of this differential equation is
$$y(t)=\frac1{1+W((h-1)\exp(h-t-1))}$$
where $W(t)$ is the Lambert function, which is the inverse of the function $t\exp\,t$. (Whether the Lambert function is considered an analytical solution might well be up for debate, but I'm in the camp that considers it a closed form.)
For instance, in Mathematica, the following code throws a warning about possible stiffness:
With[{h = 40}, y /.
First @ NDSolve[{y'[t] == y[t]^2 - y[t]^3, y[0] == 1/h}, y,
{t, 0, 2 h}, Method -> {"ExplicitRungeKutta",
"DifferenceOrder" -> 4}]]
answered Apr 6, 2013 at 16:17
J. M. ain't a mathematicianJ. M. ain't a mathematician
75.1k8 gold badges210 silver badges340 bronze badges
-
$\begingroup$ But why is it a stiff problem, beside the fact that Matlab says so? $\endgroup$ Apr 6, 2013 at 19:54
-
$\begingroup$ @Ron: the working definition of "stiff" that I use is "...equations where certain implicit methods, in particular BDF, perform better, usually tremendously better, than explicit ones". Having also tried this out with my own implementations of explicit and implicit RK with stiffness detection, I agree with the assessment of Shampine and Moler. $\endgroup$ Apr 6, 2013 at 20:02
-
$\begingroup$ I don't disagree with that assessment, it is just very nonobvious to me, especially given that I have a very dim grasp of how the Lambert W function behaves. Perhaps it is also that my numerical education came from, for better or worse, Numerical Recipes. There, the definition of "stiff" - and one that agrees with that given in the reference you provide - is when there are two dimension scales in a problem that are vastly different. I of course may be missing something, but this is why the stiffness here isn't clear to me. $\endgroup$ Apr 6, 2013 at 20:08
-
1$\begingroup$ @Ron, maybe the asymptotics might make things clearer for you: $$W((h-1)\exp(h-t-1))=\exp(h-t-1)(h-1)-\exp(2 h-2 t-2) (h-1)^2+\frac32 \exp(3 h-3 t-3) (h-1)^3+\cdots$$ $\endgroup$ Apr 7, 2013 at 2:20
-
$\begingroup$ Thank you for the pretty cool equation, but it is required an additional research. $\endgroup$– WarlockApr 8, 2013 at 15:36