search

SciML / Sundials.jl

Julia interface to Sundials, including a nonlinear solver (KINSOL), ODE's (CVODE and ARKODE), and DAE's (IDA) in a SciML scientific machine learning enabled manner
https://diffeq.sciml.ai
BSD 2-Clause "Simplified" License
209 stars 78 forks source link

Support for other solvers #51

Closed crbinz closed 6 years ago

crbinz commented 9 years ago

Are there any plans for other solvers (such as ARKode) to be included in this package? Or is the idea to keep Sundials.jl as-is and wait for ODE.jl to mature? Thanks.

stevengj commented 9 years ago

Even if ODE.jl matures, I think it would be nice to have an interface for as much of sundials as possible in this package.

jgoldfar commented 9 years ago

In case any one is interested, I have a very rough wrapper for ARKode (including update to Sundials v2.6.2, provided you can supply your own sources) in the branch at https://github.com/jgoldfar/Sundials.jl/tree/feature-arkode. The ark_brusselator example is completely workable, and added to the test suite.

crbinz commented 9 years ago

I'm definitely still interested, this is great! I'll try to give you some feedback when I get a free moment to play around with it.

ChrisRackauckas commented 6 years ago

This is completed as part of https://github.com/JuliaDiffEq/Sundials.jl/pull/149 .

I do have to warn though that the ARKODE part of Sundials does not live up to the reputation of CVODE. tk;dr: this is done, but don't expect any magic. The CVODE and IDA portions of Sundials are much better tuned than ARKODE.

On the Robertson equation, you actually have to change the defaults:

ARKODE(nonlinear_convergence_coefficient = 1e-7,predictor_method = 2)

to converge. This is actually shown in the actual Sundials C++ example as well (see their example download). For our ARKODE common interface wrapper, we are keeping the exact same defaults as if you used Sundials proper. However, we note that all of the main algorithms from ARKODE have OrdinaryDiffEq.jl equivalents by now, and the benchmarks favor them. Here, ARKODE matches the same tableau as KenCarp4 for Robertson:

arkode_robertson

and here it matches DP5 on a 100x100 linear system (the Dormand-Prince method of each package, and then Tsit5 for comparison):

arkode_linear

The quick explanation is that the ARKODE defaults are weirdly aggressive at expanding the stepsize, far more than any other implementations that I know of, which seems to cause a lot of failed steps. So for everything we are benchmarking, they don't seem to be doing well.

However, there is a catch there. In cases where you know a stepsize limit due to a CFL condition on the (or simply an explicit CFL if IMEX), then these defaults will keep the stepsize right near CFL instability when possible. This probably means they were highly tweaked to do well on reaction-advection-diffusion MOL discretizations, which makes sense because that's what the Kennedy and Carpenter paper refers to. In that case as well, these ARKODE methods are already setup with the same linear solvers as CVODE, so there are some pretty high quality Krylov methods available for large PDEs. Even in that case though, I would check it against an implementation where you swap out the linear solver of KenCarp4 for a Krylov method of IterativeSolvers.jl. It would be an interesting test case if someone wants to setup the benchmark.

As always, the benchmarks backing this up will be uploaded to DiffEqBenchmarks.jl soon. Take care.

crbinz commented 6 years ago

Awesome! My original motivation for asking this question has waned (mostly because of DifferentialEquations.jl!) but I am impressed by the effort, as usual.