From: Di Xiao <dxiao@physics.utexas.edu>
To: Slaven Peles <peles@cns.physics.gatech.edu>
Cc: <gsl-discuss@sources.redhat.com>
Subject: Re: Questions on bsimp
Date: Fri, 02 Aug 2002 08:53:00 -0000 [thread overview]
Message-ID: <Pine.LNX.4.33.0208021040270.18244-100000@linux11.ph.utexas.edu> (raw)
In-Reply-To: <E17aJO1-0002pB-00@kaos.physics.gatech.edu>
Thanks a lot! It works.
However, I still have a question. We use adaptive stepzie control because
we want to deduce the computational effort. The algorithm should be smart
to know at which region there is the smooth uninteresting curve, while in
some region many small steps are needed. My question is why bsimp chose
some unreasonably big steps in my code. It can jump from one period to
another period without showing the details of each one. If you
run the code, you will find from the results givin by bsimp, it is hard
to tell what function it likes. I tried different gsl_odeiv_control, but none
of them gave me good results. I guess the reason is when there are
oscillations, bsimp is not a good one. Can anybody explain it?
Thanks.
Di
On Thu, 1 Aug 2002, Slaven Peles wrote:
> On Thursday 01 August 2002 12:19 pm, you wrote:
> > When I use bsimp as the stepper function, it seems it just jumps too each
> > step. I need bsimp because the speed is an important issue in my project.
>
> When using adaptive step size different integrators will produce different
> number of intermediate values during integration from t to t1. bsimp produces
> fewer intermmediate values so that is why it seems to you like your results
> are messy. Have look at the second example in the documentation at
> http://sources.redhat.com/gsl/ref/gsl-ref_25.html#SEC381
> It's explained there how to "fix" this problem ;-).
>
> Slaven
>
> > --------------------------------------------------------------
> >
> > #include <stdio.h>
> > #include <math.h>
> > #include <gsl/gsl_errno.h>
> > #include <gsl/gsl_matrix.h>
> > #include <gsl/gsl_odeiv.h>
> >
> > int func (double t, const double y[], double f[], void *params)
> > {
> > f[0] = -y[1];
> > f[1] = y[0];
> > return GSL_SUCCESS;
> > }
> >
> >
> > int jac (double t, const double y[], double *dfdy, double dfdt[],
> > void *params)
> > {
> > gsl_matrix_view dfdy_mat
> > = gsl_matrix_view_array (dfdy, 2, 2);
> > gsl_matrix * m = &dfdy_mat.matrix;
> > gsl_matrix_set (m, 0, 0, 0.0);
> > gsl_matrix_set (m, 0, 1, -1.0);
> > gsl_matrix_set (m, 1, 0, 1.0);
> > gsl_matrix_set (m, 1, 1, 0.0);
> > dfdt[0] = 0.0;
> > dfdt[1] = 0.0;
> > return GSL_SUCCESS;
> > }
> >
> > int main (void)
> > {
> > const gsl_odeiv_step_type * T
> > = gsl_odeiv_step_rkf45;
> >
> > gsl_odeiv_step * s
> > = gsl_odeiv_step_alloc (T, 2);
> > gsl_odeiv_control * c
> > = gsl_odeiv_control_standard_new(1e-30, 1e-10,
> > 1.0, 1.0);
> > gsl_odeiv_evolve * e
> > = gsl_odeiv_evolve_alloc (2);
> >
> > gsl_odeiv_system sys = {func, jac, 2, NULL};
> >
> > double t = 0.0, t1 = 10.0;
> > double h = 1e-6;
> > double y[2] = {1.0, 0.0};
> >
> > while (t < t1)
> > {
> > int status = gsl_odeiv_evolve_apply(e, c, s,
> > &sys,
> > &t, t1,
> > &h, y);
> > if (status != GSL_SUCCESS)
> > break;
> > printf ("%.5e %.5e %.5e\n", t, y[0], y[1]);
> > }
> >
> > gsl_odeiv_evolve_free (e);
> > gsl_odeiv_control_free (c);
> > gsl_odeiv_step_free(s);
> > return 0;
> > }
>
next parent reply other threads:[~2002-08-02 15:53 UTC|newest]
Thread overview: 4+ messages / expand[flat|nested] mbox.gz Atom feed top
[not found] <E17aJO1-0002pB-00@kaos.physics.gatech.edu>
2002-08-02 8:53 ` Di Xiao [this message]
2002-08-06 10:42 ` Brian Gough
2002-08-01 10:13 Slaven Peles
-- strict thread matches above, loose matches on Subject: below --
2002-08-01 9:19 Di Xiao
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.33.0208021040270.18244-100000@linux11.ph.utexas.edu \
--to=dxiao@physics.utexas.edu \
--cc=gsl-discuss@sources.redhat.com \
--cc=peles@cns.physics.gatech.edu \
/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).