error tolerance for integrate_ode?

[Bill Gillespie]

Hi Bob,

I don't have posting privileges on the developers group, so I am responding to a post from Sept 22 here. I am a regular user of numerical ODE solvers in the context of PKPD modeling.

Right now, the error tolerance in the ODE solvers are set to 
1e-6 absolute and 1e-6 relative. 
Question 1:  Why would we ever want an absolute tolerance? 
  Can I just make it the max double value? 

That will be fine in cases where relative tolerance dominates, but it may lead to failures when absolute tolerance dominates, e.g., when the solution is near 0. By failure I mean inability to meet the tolerance specification with the specified limits on step size or number of steps.

Question 2:  Can we just stick with 1e-6 relative tolerance, which is 
   what all the Boost examples use and pretty close to what everyone 
   seems to use for all these integrators across packages.  It seems 
   like 1e-6 is plenty of tolerance for stats, but I'm not an ODE 
   solving expert. 

Relative tolerance = 1e-6 is probably a good default, but there will be occasions where a user will want to change it either to prevent integration failures or to increase precision.

The error control is now for both the system itself and the 
sensitivities (derivatives of solutions w.r.t. parameters). 
Currently we have a single function w/o any user-specified tolerances: 
   integrate_ode(...); 
To allow control of both tolerances we'ed need three 
more functions since there's ambiguity in dropping 
arguments.  For example, 
   integrate_ode_tol(..., abs_tol, rel_tol) 
   integrate_ode_abstol(..., abs_tol) 
   integrate_ode_reltol(..., rel_tol) 
Question 3:  Do we need these? 

I think it would be sufficient to have two functions: one with only the defaults and one that requires specification of all options. If an expert user really wants other choices, they can create them via user-defined functions.

Question 4:  If so, can we put off adding the tolerances until later 
   or are these dealbreakers? 

I'll take what I can get for now, but eventually want control of tolerances and possibly limits on step size or number of steps.

If the answer to (3) is yes, we could use suggestions on how to test. 
Do we just set up some equations, set different tolerances, and check the 
ODEs determinstically? 

For formal testing, deterministic tests using ODE's with known solutions seems reasonable. On the other hand when ODE solvers are exercised by sometimes extreme sets of parameter values within an MCMC run we learn more about what will break them. I encourage you to put prototypes in the hands of users who often require numerical ODE solvers and see what happens. I would be glad to participate.

- Bob 

[Bob Carpenter]

Thanks for the clarification on absolute/relative tolerances.
I think someone else wrote in recently and said the same thing.

There's no way to limit number of steps with the integrator
we're currently using (Boost's odeint).

We can eventually add in something to change tolerances, and it
makes sense to only have the fully specified version if both are
always relevant.

We're going to roll out our first version soon, which I think
will be more of a "prototype" than a finished product. Maybe we'll
stick a "beta" sticker on it or something.

We're really bogged down in testing at the moment with four
different run conditions of known/unknown initial state and
needing to test all the gradients and error conditions (about four
different ways for each to raise errors).

We hope to have it done in another couple weeks, but who knows
what'll pop up next. Right now, Daniel's wrestling with a bug in
make that's blocking all our testing.

- Bob

> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

[Avraham Adler]

Yes, Bob, I sent you a direct e-mail with similar content as I don't have posting privileges on stan-devel either.

Avi

[Bob Carpenter]

Thanks to both of you. I hadn't thought about relative
error being unstable with values near 0.

And I can see how you'd need to tune if you had diff eqs
on radically different scales, because if your variables scale
is 1e9, then I'm not likely to get 1e-6 absolute error tolerances,
whereas if the variable scale is 1e-9, then 1e-6 is too large.
out.

- Bob