Draft: Performance fixes for DualLinearProblems

[jClugstor]

Checklist

• Appropriate tests were added
• Any code changes were done in a way that does not break public API
• All documentation related to code changes were updated
• The new code follows the
  contributor guidelines, in particular the SciML Style Guide and
  COLPRAC.
• Any new documentation only uses public API

Additional context

[jClugstor]

@oscardssmith this caches a lot more stuff during init, so should help performance.

The biggest source of allocations and something that takes a large chunk of time is the creation of the solution here:

LinearSolve.jl/ext/LinearSolveForwardDiffExt.jl

Line 115 in 547db8a

function linearsolve_dual_solution(u::AbstractArray, partials,

which can definitely be improved.

[jClugstor]

Now I'm down to solve!(::DualLinearCache) only taking one alloc.

The getfields everywhere are because I overrided getproperty for DualLinearCache to go to getproperty of LinearCache if appropriate, but of course that's not type stable. Maybe using Val would help if we really wanted to change back to getproperty, not sure.

So with this branch:

using LinearSolve
using ForwardDiff
using Test
using Chairmarks

function h(p)
    (A=[p[1] p[2]+1 p[2]^3;
            3*p[1] p[1]+5 p[2]*p[1]-4;
            p[2]^2 9*p[1] p[2]],
        b=[p[1] + 1, p[2] * 2, p[1]^2])
end

A, b = h([ForwardDiff.Dual(5.0, 1.0, 0.0), ForwardDiff.Dual(5.0, 0.0, 1.0)])

prob = LinearProblem(A, b)
@b solve(prob, LUFactorization())
cache = init(prob, LUFactorization())
@b solve!(cache)
371.835 ns (1 allocs: 48 bytes)

Then LinearSolve 3.17:

@b solve!(cache)
74.021 ns (1 allocs: 48 bytes)

So this is a lot better than current main but still quite slower than before.

[oscardssmith]

delete this and just return sol?

[oscardssmith]

these functions seem to be just shufling data back and forth. Can you not do this without the shuffling?

[jClugstor]

The issue is that we need to take a matrix of ForwardDiff.Partials and convert them in to a list of matrices where each matrix holds corresponding entries of the Partials. We need to do this in order to have matrices that mul! and LinearSolve can actually act on, because arithmetic isn't defined for ForwardDiff.Partials. So I'm not sure if there's a way to accomplish that without this.

[jClugstor]

I guess one thing we could do is precompute the list of the necessary RHS matrices during init, since that will only need to be computed once, but would need to be recomputed if A or b change

[jClugstor]

FYI these changes also broke the nested dual number handling so I'm trying to figure that out

[jClugstor]

Nested Duals are back to working

[ChrisRackauckas]

what's left here?