Conditional exact integration in ODE solving

Say in my ODE problem I have an equation like this:

\frac{\text{d}u_i}{\text{d}t} = r(\mathbf{u})A(t)

Here A(t) is a known non-negative interpolation function and r(\mathbf{u}) is what we call a ‘reduction factor’; it’s a function whose value is in the interval [0,1] used to maintain physical feasiblity, and generally depends on a small number of other states.

There can be large time periods in the simulation where simply r(\mathbf{u}) =1, in which case the above equation can be integrated exactly. My question is this: in OrdinaryDiffEq.jl, is there a nice way to apply the exact integration where applicable? In theory it is possible to rig something together with VectorContinuousCallback, but that is too costly because our simulations can contain hundreds of equations of the above form.

Now that I think of it, it might be nicer to have the equation

\frac{\text{d}u_i}{\text{d}t} = (1-r(\mathbf{u}))A(t)

so that this state represents the deviation from the exact integration case.

How about rigging a composite algorithm for this purpose? ODE Solvers · DifferentialEquations.jl

That looks interesting. I don’t think it’s applicable though, because this decision has to be made on a per-state basis.