ODE solver stats

New to Julia
1

I would really appreciate if someone can briefly explain stats 3, 7-9. I did find any documentation on this:

sol.destats
DiffEqBase.DEStats
Number of function 1 evaluations:                  261
Number of function 2 evaluations:                  0
Number of W matrix evaluations:                    37
Number of linear solves:                           222
Number of Jacobians created:                       37
Number of nonlinear solver iterations:             0
Number of nonlinear solver convergence failures:   0
Number of rootfind condition calls:                0
Number of accepted steps:                          37
Number of rejected steps:                          0

Thanks!

1 Like
2

In the computation of an implicit solve, W = I - gamma*J is what’s actually used in the updates. In the Newton method, then this matrix is used and reused. (7) is then the number of times the Newton iterations diverged, the accepted steps are the number of steps which were accepted because the error was below the tolerance. The rootfind calls refers to the number of condition calls in the ContinuousCallbacks

3

Thanks!
So the function 1 refers to the objective function containing the ode problem? What does function 2 refer to?
So the implicit solver used here used the Newton’s method which is a non-linear. What do the ‘linear solves imply’? Because aren’t the jacobians used for the non-linear iterations?( Which are zero here)?

4

Function 2 is the second function when using a split ODE.

You are correct on that second part, but maybe the ODE solver you are using doesn’t give us that information (or we haven’t hooked into it)

5

Okay, so taking a concrete example (Lorentz system):

using DifferentialEquations, Sundials

function lorenz!(du,u,p,t)
 du[1] = 10.0*(u[2]-u[1])
 du[2] = u[1]*(28.0-u[3]) - u[2]
 du[3] = u[1]*u[2] - (8/3)*u[3]
end

u0 = [1.0;0.0;0.0]
tspan = (0.0,100.0)
prob = ODEProblem(lorenz!,u0,tspan)
sol = solve(prob, CVODE_BDF(), reltol = 1e-4, abstol = 1e-7);
println(sol.destats);

This gives an output:

Number of function 1 evaluations:                  8302
Number of function 2 evaluations:                  0
Number of W matrix evaluations:                    828
Number of linear solves:                           0
Number of Jacobians created:                       113
Number of nonlinear solver iterations:             8298
Number of nonlinear solver convergence failures:   0
Number of rootfind condition calls:                0
Number of accepted steps:                          5889
Number of rejected steps:                          416

So since CVODE_BDF uses Newton iteration the “number of linear solves” is zero. Then,

  1. “number of non-linear solver iterations” refers to the total number of iterations done in all the steps taken by the solver in the given tspan?
  2. Why don’t the (accepted + rejected steps) add up to the total number of function evaluations?
  3. So the W matrix is obtained from the Jacobian as W = I - gamma*J? What do I and gamma refer to here?

If I use Rodas4 for the same problem, I get:

println(sol.destats);
DiffEqBase.DEStats
Number of function 1 evaluations:                  35818
Number of function 2 evaluations:                  0
Number of W matrix evaluations:                    4477
Number of linear solves:                           26862
Number of Jacobians created:                       4477
Number of nonlinear solver iterations:             0
Number of nonlinear solver convergence failures:   0
Number of rootfind condition calls:                0
Number of accepted steps:                          4051
Number of rejected steps:                          426

So is Rodas4 an example where the solver doesn’t return the number of “non-linear iterations” as you referred to? Why does it show “linear solves” then?

Thanks

6

Calculation the Jacobian will also give some function evaluations, if no sparsity is involved it will cost size(u)^2 calls per jacobian (so 3*3*113=1017 in your example ). You see that this is the point where large systems become expensive to compute. The remaining 152 function calls… I don’t know…

7

Okay so the Jacobians account for the additional function calls per ODE function call. Then I guess the rest of it has to do with the W matrix evaluation. You considered 1 function call per evaluation to arrive at this number right?

8

I used call and evaluation synonymously but good point. The number of calls may depend on the way the algorithm takes the derivatives. And since autodiff is the default I assumed one call per derivative (how ever finite difference approximation should only give one additional call). I am not an expert in such details and not shure if this true at all.

9

No it’s length(u) since it’s computed column by column

10

Let me start by noting that the stats we get out of the Julia solvers are known to be correct and varified, while ones from Sundials and ODEInterface are just the stats we know how to get out from their interface and some of them seem a little odd. This is a known issue and if there are ways to improve there stats I’d be interested to see.

That is correct.

Because every step of BDF is solving an implicit system, where the implicit system has to call f every step. 8298 nonlinear solver iterations where each has an f call in it, and a few for the automatic dt selection algorithm would make 8302. These nonlinear solves include the accepted and rejected steps.

That f then cannot be a total f because the Jacobians would add 113*3 f calls. This differs from how the *DiffEq libraries calculate it, and it would be worth double checking the exact number of f calls and adjusting it the number given from Sundials.

I is the mass matrix of an ODE without a mass matrix defined Mu'=f(u,p,t). gamma is a constant that is proportional to dt.

11

Thanks,
Then the “function calls” only refer to the calls to the objective function. So by the same logic then in case of Rodas4 for example (as above) the “number of linear solves” + function calls for optimizing dt = “number of function 1 evaluations”? But then that is almost 10000 additional calls. and why are the Jacobians evaluated if the nonlinear solver iterations are zero?

12

No, with Rodas4 and the rest of the *DiffEq methods it’s counting the calls in the Jacobian construction.

Rosenbrock methods are not implicit, but they still use Jacobians and linear solves.