Support time derivative (by supplying a tgrad function) for SciML Rosenbrock solvers

[sjdaines]

SciML Rosenbrock solvers require a time derivative (supplied by tgrad keyword in ODEfunction). Without this, they will attempt to calculate this using autodifferentiation and hence will fail with PALEO models.

Example of what is needed (see SciML/OrdinaryDiffEq.jl#1587) for the case where there is no time dependence PALEOmodel/src/ODE.jl L78:

f = SciMLBase.ODEFunction{true}(m; jac, jac_prototype, mass_matrix=M, tgrad = (dT, u, p, t) -> fill!(dT, 0),)

The more general approach with time dependence would need to calculate this by autodifferentiation (cf how Jacobians are calculated in PALEO).