From: Timo Laitinen <timolai@green.srl.utu.fi>
To: gsl-discuss@sources.redhat.com
Subject: Stepping to negative direction in ODE
Date: Wed, 19 Dec 2001 13:20:00 -0000 [thread overview]
Message-ID: <Pine.LNX.4.10.10110031239010.15403-100000@green.srl.utu.fi> (raw)
Hi,
I found a problem with the step size evolution mechanism in the functions
gsl_odeiv_evolve_apply
and
std_control_hadjust
The problem is that they cannot handle negative step sizes, ie they cannot
be used for solving an equation backwards in time (or in space!!!). My
solution for the problem is to change in
evolve.c, line 147
if (h0 > dt)
to
if (GSL_SIGN(h0)*h0 > GSL_SIGN(dt)*dt)
and in cstd.c, lines 109 and 118,
*h = GSL_MAX_DBL(0.2 * h_old, h_new);
to
*h = (double) GSL_SIGN(*h)*GSL_MAX_DBL(0.2 * h_old, h_new);
and
*h = GSL_MIN_DBL(5.0 * h_old, h_new);
*h = (double) GSL_SIGN(*h)*GSL_MIN_DBL(5.0 * h_old, h_new);
, respectively. There may be a better/faster way of doing this so think
first (offhand, i can't remember why I chose not to use fabs in the first
one...)
There is also another problem in function gsl_odeiv_evolve_apply.
If the initial h0 is larger than the remaining time dt, the final_step
flag is set, on line 150. BUT, if the step is now decreased by
gsl_odeiv_control_hadjust, line 182, and the new step h0<dt, the
final_step flag is NOT RESET!. Thus, you end up taking a small step *h,
and updating the time *t to t1, which would give strange results.
Of course one should take care of not starting with a step size larger
than the integration period, but accidents happen (as in my case :-) ),
and tracking the problem is a real pain...
Otherwise it seems to work fine, thank you for a great work, and for
saving my time considerably :-)
----------------
Timo FABRICATI DIEM, PVNK.
email:timo.laitinen@utu.fi
next reply other threads:[~2001-12-19 13:20 UTC|newest]
Thread overview: 2+ messages / expand[flat|nested] mbox.gz Atom feed top
2001-12-19 13:20 Timo Laitinen [this message]
2001-12-19 13:20 ` Brian Gough
Reply instructions:
You may reply publicly to this message via plain-text email
using any one of the following methods:
* Save the following mbox file, import it into your mail client,
and reply-to-all from there: mbox
Avoid top-posting and favor interleaved quoting:
https://en.wikipedia.org/wiki/Posting_style#Interleaved_style
* Reply using the --to, --cc, and --in-reply-to
switches of git-send-email(1):
git send-email \
--in-reply-to=Pine.LNX.4.10.10110031239010.15403-100000@green.srl.utu.fi \
--to=timolai@green.srl.utu.fi \
--cc=gsl-discuss@sources.redhat.com \
/path/to/YOUR_REPLY
https://kernel.org/pub/software/scm/git/docs/git-send-email.html
* If your mail client supports setting the In-Reply-To header
via mailto: links, try the mailto: link
Be sure your reply has a Subject: header at the top and a blank line
before the message body.
This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).