`stiff3`

is a semi-implicit Runge Kutta solver with adaptive time-stepping originally from the book by Villadsen & Michelsen, *Solutions of Differential Equation Models by Polynomial Approximation*, Prentice-Hall, 1978.

The solver is now available as an `fpm`

package and uses LAPACK for the matrix factorization and back-substitution in the semi-implicit time-stepping.

The solver is designed for autonomous problems of the form

In case of a non-autonomous problem where x appears explicitly, i.e. y' = f(x,y), it suffices to switch to a new integration variable t and introduce the additional ODE

The repository provides an example for the Van der Pol oscillator. Here are some plots of the solutions and phase cycle.

There is one remaining issue with the `stiff3`

solver where I’m seeking help. It has to do with a section of code that sets to zero some elements of the Jacobian matrix. I believe this section is just an outdated practice that has to do with legacy usage patterns. I’ve posted it as a question on SciComp StackExchange, but so far it is still unanswered.