Once again we stayed true to form and didn’t solve the problems in the development list but adding a ton of new features anyways. Now that Google Summer of Code (GSoC) is in full force, a lot of these updates are due to our very awesome and productive students. Here’s what we got.
Approximate Bayesian Computation (ABC)
Marc Williams (@marcjwilliams1) contributed the
abc_inference function to
DiffEqBayes.jl which utilizes approximate Bayesian computation (also known as
ABC) from ApproxBayes.jl
to perform the estimation of parameter posteriors. Compared to standard
Bayesian methods, ABC is computationally cheaper and faster at the cost of
making some approximations (the name is quite appropriate). For models with lots
of parameters where you have a good prior guess for the posterior parameter
point estimates, ABC can be a great way to get posterior distributions.
Maximum A Posteriori Estimation
GSoC student Vaibhav Dixit (@Vaibhavdixit02) added maximum a posteriori estimation to our existing parameter estimation routines in DiffEqParamEstim.jl. Now these optimization-based methods can take into account prior distributions which can help global optimizers stay in the right parameter ranges and improve fitting. All it takes is passing an optional prior distribution to any of the existing methods, so it’s very easy to add it your current work!
Multiple shooting objective
In other parameter estimation news, Vaibhav Dixit (@Vaibhavdixit02) added the
multiple_shooting_objective which, like the
allows for fitting ODE results to any loss function. However, multiple shooting
methods are naturally more robust by solving simultaneously from many different
SteadyStateProblem types, a new solver
DynamicSS has been added that
utilizes the ODE solvers to find steady states. It builds in a callback that
will halt when the derivative is sufficiently small, allowing it to be a
robust automated steady-state finding machine.
Unified retcode handling
Takafumi Arakaki (@tkf) performed an awesome refactoring in the DiffEq core which allows all of the DiffEq solvers to utilize the same retcode handling code. This means that we have unified exit warning messages and the resulting retcodes. Each package should act the same for the common retcodes.
Fixed order variable time step Adams methods
Shubham Maddhashiya (@sipah00) contributed a variety of new adaptive Adams-Bashforth and Adams-Bashforth-Moulton methods. These methods are fixed order 3-5 which minimize function evaluations and are designed for large non-stiff ODE discretizations. They utilize Runge-Kutta methods to hotstart the Adams methods, making them more efficient than variable order versions when only lower orders are required or when there are a lot of events.
Yingbo Ma (@YingboMa) contributed
ABDF2, an adaptive implementation of the
Backwards Differentiation Formula (BDF) order 2 method. BDF schemes are the
go-to choice for PDE discretizations since they only have one step (one
nonlinear equation to solve) for each time step, minimizing the number of
function calculations. Unlike adaptive order BDF discretizations, this adaptive
BDF2 is A-B-L-stable, meaning that it should be stable for any stiff ODE that is
thrown at it. This makes it a nice testing or anchor method: one that can always
be relied on to do well-enough.
A lot of the next developments will come from our GSoC students. Here’s a list of things we are aiming for:
SABDF2, which is a strong order 0.5 adaptive BDF2 implementation for stochastic differential equations which is order 2 for small noise SDEs. This will be the first implicit adaptive integrator for small noise SDEs and will be a great choice for SPDEs.
Adaptive order Adams-Bashforth-Moulton. This will be JuliaDiffEq’s native Julia implementation of an AOAT Adams method. It will utilize variable coefficient formulas like Shampine’s
ddeabmwhich has been shown in sources like Hairer’s Solving Ordinary Differential Equations to be more efficient than the more standard
CVODE_Adamssince it’s more easily able to utilize higher order steps.
Adaptive order Nordsieck methods. Nordsieck methods are the special implementation multistep methods which Sundials’ CVODE uses. This implementation utilizes a fixed leading coefficient which makes it more optimized for solving large PDE discretizations since it can re-use the Jacobian between steps while adapting time steps. Other problem solving environments such as MATLAB and SciPy utilize quasi-fixed timesteps to get a similar benefit, but this has the added cost of fixing time steps for many steps which reduces stepping efficiency. The reason this is done is simple: Nordsieck implementations are hard which is why the only existing ones (that I know of) are EPISODE and VODE/CVODE (and VODE is an adaption of EPISODE, so it’s really one code!). However, our GSoC student has been hammering away at this for months and we’ve made great progress, with the fixed order Adams Nordsieck method having already merged. The last step is to figure out order adaptivity and then the translation to BDF coefficients is trivial. While this methodology isn’t as optimized for non-stiff ODEs in its Adams form, this BDF form will likely be the new goto method in DiffEq for large stiff PDE discretizations.
Until exponential integrators! Xingjian Guo (@MSeeker1340) has been doing an extensive study of the numerical implementations of
phimvfor building efficient exponential integrators. If you’re curious, you can read up on some of the discussions and here. He is quite close to having both efficient and numerically stable methods for lazy adaptive Krylov
phimvwhich will allow for fast c alculations of exponential integrators without computing full matrix exponentials. This will then be used in the current exponential integrator implementations, along with the adaptive high order EPIRK methods. This class of integrators will be one to keep aware of if you’re interested in time-dependent PDEs. It will be great to compare and contrast between these methods and BDF integrators in these problems.
Yiannis Simillides (@ysimillides) keeps making improvements to FEniCS.jl. At this part a large portion (a majority?) of the tutorial works from Julia.
Vaibhav Dixit (@Vaibhavdixit02) finished most of the parameter estimation methods on our wish list, so he’s onto implementing global sensitivity analysis methods, starting with the Morris method and then going on to the eFAST and Sobol method.
Mikhail Vaganov (@Mikhail-Vaganov) is making good progress on his N-body modeling language. This will make it easy to utilize DiffEq as a backend for molecular dynamics simulation. Follow the progress in DiffEqPhysics.jl
And here’s a quick view of the rest of our “in development” list:
- Preconditioner choices for Sundials methods
- Adaptivity in the MIRK BVP solvers
- More general Banded and sparse Jacobian support outside of Sundials
- IMEX methods
- Function input for initial conditions and time span (
- LSODA integrator interface
Are you a student who is interested in working on differential equations software and modeling? If so, please get in touch with us since we may have some funding after August for some student developers to contribute towards some related goals. It’s not guaranteed yet, but getting in touch never hurts!