Skip to content

AutoFiniteDiff JVP recomputes f(x) on every Krylov matvec; nf accounting double-counts it #4162

Description

@ChrisRackauckas-Claude

Summary

On the matrix-free Krylov path with autodiff = AutoFiniteDiff(), the base evaluation f(x) is recomputed on every Krylov matvec, even though x is fixed for the entire linear solve. FiniteDiff.finite_difference_jvp! accepts an f_in argument precisely so the caller can supply a cached f(x), but DifferentiationInterface's pushforward! (DifferentiationInterfaceFiniteDiffExt/twoarg.jl:113-131) never passes it.

Measured cost

2D Brusselator, N=32 (2048 ODEs), FBDF(linsolve = KrylovJL_GMRES()), reltol = abstol = 1e-6, Julia 1.12.4, BLAS.set_num_threads(1):

time/matvec raw f calls per matvec
plain f 9.79 µs
AutoForwardDiff JVP 27.67 µs 1 (dual)
AutoFiniteDiff JVP 24.96 µs 2

One of those two f calls is the redundant base evaluation. In an end-to-end solve with a correctly scaled JVP step, caching it moved wall time from 0.97× to 0.83× relative to AutoForwardDiff — about a 1.2× improvement — with bit-identical step statistics.

It also halves the reported nf. dolinsolve accounts 2 * linres.iters for the finite-difference path versus 1 * linres.iters for AD, so nf currently overstates the real cost difference by 2× and makes AutoFiniteDiff look worse in solver statistics than it is.

Context

Found while root-causing a separate and much larger problem on the same path: FiniteDiff's JVP step size is dimensionally wrong (sqrt(|x⋅v|) has units √([x][v]), not [x]/[v]), which pins the step at the absstep floor regardless of state scale and drives Newton convergence failures worth several hundred× in wall time. That is being fixed separately in FiniteDiff.jl. This issue is the smaller, independent inefficiency — worth ~1.2× on its own, and worth fixing regardless of the step-size work.

Where the fix belongs

Primarily in DifferentiationInterface's FiniteDiff extension: thread a cached base evaluation through the prep object and pass it as f_in. That is a different repo; filing here because this is where the cost lands and where the nf accounting lives, and because the nf double-count in dolinsolve is fixable here independently.

For reference, the pre-DifferentiationInterface implementation handled this: SparseDiffTools.num_jacvec! carried a compute_f0 flag to skip the redundant evaluation. Both that and the correct step scaling were lost in the migration (commit 13c036b2f4).

Environment

Julia 1.12.4
OrdinaryDiffEq v7.3.0
OrdinaryDiffEqDifferentiation 3.7.0
DifferentiationInterface 0.7.20
FiniteDiff 2.32.1

prepare_jvp and dolinsolve on master are byte-identical to the released versions, so this applies to master unchanged.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions