Dear Frank, All,
On 26.12.2007, at 23:14, Frank Reininghaus wrote:
I've found some problems with adaptive step size control which can lead to an infinite loop under special circumstances.
Nice write-up of this problem!Although I have not looked through your cases in detail, I can conform that we have had similar infinite-loop problems, for example, on AMD Opteron systems running Linux and operated in 64bit mode (gcc-4.x with -m64).
Because I did not have time to look into it in detail, and because the q&d-hack below works for us, so far we use the following code snippet to work around it:
while(t < tend) {
if(GSL_SUCCESS != gsl_odeiv_evolve_apply(evolve,
control, step,
&system, &t, tend, &h, y)) {
throw Error("Flight::operator(): integration
problem in GSL");
}
if(h < hmin) {
h = hinit;
}
}
with (for example)
hinit = tstep * 0.1
hmin = tstep * 0.01
We use this with gsl_odeiv_step_rkf45, gsl_odeiv_step_rkck, or
gsl_odeiv_step_rk8pd.
We would be quite interested in seeing this fixed in GSL, obviously.I'll try to find some time to test your patch. However, why don't you test what you really describe your comment? This would look like the following (untested):
/* Step size h0 was decreased. Make sure it does not get too small:
It should be positive (i.e., h0 >= GSL_DBL_MIN) and large enough
to make t0 + h0 distinguishable from t0, i.e., at least of size
GSL_DBL_EPSILON*t0 */
h0 = GSL_MAX_DBL(GSL_DBL_MIN, (t0 != t0+h0) ? h0 :
GSL_DBL_EPSILON*t0);
This also avoids the multiplication in most cases. Greetings, Jochen --Einigkeit und Recht und Freiheit http://www.Jochen- Kuepper.de
Liberté, Égalité, Fraternité GnuPG key: CC1B0B4D
Sex, drugs and rock-n-roll
PGP.sig
Description: This is a digitally signed message part
