Matrix-free linear solves

[lassepe]

The Krylov method introduced in #24 lends itself to matrix-free realization via SciMLOperators.jl. This should speed things up by avoiding materialization of the full Jacobian (basically a lazy Jacobian representation)