Choice of adjoint sensitivity

[jiweiqi]

Currently, RMS is using InterpolatingAdjoint(autojacvec=ReverseDiffVJP(false)) for adjoint sensitivity.

ReactionMechanismSimulator.jl/src/Simulation.jl

Line 141 in 53fe904

function getadjointsensitivities(bsol::Q,target::String,solver::W;sensalg::W2=InterpolatingAdjoint(autojacvec=ReverseDiffVJP(false)),abstol::Float64=1e-6,reltol::Float64=1e-3,kwargs...) where {Q,W,W2}

For stiff systems like combustion chemistry, call ForwardDiff.gradient(fsol, u0) might be more efficient, something like.
https://github.com/DENG-MIT/Arrhenius.jl/blob/99dfd44f78c0c26743087e423e763f54d93d631d/example/pyrolysis.jl#L102

But this is case dependent, and unfortunately, Arrhenius.jl doesn't support reverse-mode AD currently, due to mutable array. Therefore, I haven't really benchmarked the forward and reverse mode AD + adjoint yet.

@ChrisRackauckas might have better explanations for this point and can suggest what to try.

[ChrisRackauckas]

The recent manuscript better explains some of https://arxiv.org/abs/1812.01892 . Basically, ForwardDiff.gradient is going to scale better with stiff ODE solvers unless the adjoint is QuadratureAdjoint (which is what's in this preprint). InterpolatingAdjoint with chunking could be an option, but I think QuadratureAdjoint is then going to win whenever you have the memory for it. The cutoff of 100 shown there on the stiff Brusselator PDE semi-discretization shows about where to expect the cutoff, but when you allow ForwardDiff.jl to do chunking it's probably a little higher. Especially if you parallelize the chunks. So somewhere in the 100-1000 zone is probably the cutoff zone on modern computers, going to the higher end if you multithread.

(This is all in the stiff case of course, which acts very different from the non-stiff case)

[jiweiqi]

Therefore, what can be explored are

• QuadratureAdjoint
• ForwardDiff.jl

[ChrisRackauckas]

Right now yes.

[ChrisRackauckas]

with the right VJP choice.