Parameters of IDA DAE solver

I have the following code:

using DifferentialEquations, Sundials, StaticArrays
using Revise

if ! @isdefined KPS3
    using .KPS3

my_state = KPS3.get_state()
y0, yd0 = KPS3.init(my_state)

tspan = (0.0, 0.05)         # time span; fails when changed to (0.0, 0.06)

differential_vars =  ones(Bool, 36)
prob = DAEProblem(residual!, yd0, y0, tspan, differential_vars=differential_vars)

sol = solve(prob, IDA(), saveat=0.001, abstol=0.01, reltol=0.001)

time = sol.t

pos_x = sol[3*5+1, :]
pos_z = sol[3*5+3, :]

It works fine. You can check it out here: GitHub - ufechner7/KiteViewer at sim
(careful, use branch sim)

But if I increase the simulation time to 0.06 seconds it fails with the following error

julia> include("src/RTSim.jl")

  At t = 0.0599924 and h = 2.40789e-162, the corrector convergence failed repeatedly or with |h| = hmin.

  3.462572 seconds (10.40 M allocations: 599.441 MiB, 15.08% gc time)

I would like to set the following parameters of the solver to make it work:

        inith = 0.002
        maxord = 3

        pos_tol = 2.0/ 100.0
        vel_tol = pos_tol / 60.0
        atol = vel_tol * ones(36)
        for i in 1:18
            atol[i] = pos_tol
        rtol = 1e-3

The names I am using here are coming from the Python/ Assimulo code that I am converting to Julia,
so I might have to use slightly different names. (See: IDA — Assimulo 3.0 documentation)

Special question: Can I use arrays for the absolute end the relative tolerance?

Any ideas how to set inith, maxord, rtol and atol when using IDA from DifferentialEquations?


I created an issue: Document IDA options · Issue #306 · SciML/Sundials.jl · GitHub

OK, I found the functions:
function IDASVtolerances(ida_mem, reltol::realtype, abstol::N_Vector)
function IDASetMaxOrd(ida_mem, maxord::Cint)
function IDASetInitStep(ida_mem, hin::realtype)

In the file:

But how can I use them?

Where do I get ida_mem from?

Most of those are common arguments documented in the same place they always are:

The order though is specific to IDA, so it’s documented in the same place all of the algorithm-specific arguments are, i.e. in the algorithm’s documentation.

It all follows a pattern.

Ok, but how can I specify a vector for abstol?

Just make it a vector. solve(prob,alg,abstol=[1e-3,1e-4])

The scaled error is guaranteed to be <1 for a given local error estimate (note: error estimates are local unless the method specifies otherwise). abstol controls the non-scaling error and thus can be though of as the error around zero. reltol scales with the size of the dependent variables and so one can interpret reltol=1e-3 as roughly being (locally) correct to 3 digits. Note tolerances can be specified element-wise by passing a vector whose size matches u0 .

OK, lets see if that works.

It should, for reference the source code shows how the common interface lowers to using the functions you’re talking about:

1 Like