Varying length of residuals in NonlinearLeastSquaresProblem

[Ickaser]

Is your feature request related to a problem? Please describe.

In my application, I fit parameters to measurements from a drying process. This means that the model has a varying end time (i.e. when drying is complete), and so the number of data points with a meaningful comparison (and therefor number of residuals) may vary as the parameters vary.

One way of handling this might be to pad out a cached array with the maximum possible number of data points, and strategically fill it, which I could do, but it would be nice to just let the algorithm handle it. In person @ChrisRackauckas suggested that the varying number of residuals should be fine with the out-of-place form of the algorithms, so I tried that out, but it doesn't seem to work.

MRE

using Random
x = range(0, 10, length=20)
y0 = rand(20) .+ x.^2
function varying_length(u, p)
    len = rand(19:29)
    pred_y = @. u[1] + u[2]*x + u[3]*x^2
    if len > length(x)
        return pred_y .- y0
    else
        return pred_y[1:len] .- y0[1:len]
    end
end

Random.seed!(1234)
nlvary = NonlinearLeastSquaresProblem{false}(varying_length, [0.0, 0.0, 0.0])
solvary = solve(nlvary, LevenbergMarquardt(), show_trace=Val(true))

Results of MRE:

Algorithm: LevenbergMarquardt(
    trustregion = LevenbergMarquardtTrustRegion(
        β_uphill = 1.0
    ),
    descent = GeodesicAcceleration(
        descent = DampedNewtonDescent(
            initial_damping = 1.0,
            damping_fn = LevenbergMarquardtDampingFunction(
                increase_factor = 2.0,
                decrease_factor = 3.0,
                min_damping = 1.0e-8
            )
        ),
        finite_diff_step_geodesic = 0.1,
        α = 0.75
    ),
    autodiff = AutoForwardDiff(),
    vjp_autodiff = AutoFiniteDiff(
        fdtype = Val{:forward}(),
        fdjtype = Val{:forward}(),
        fdhtype = Val{:hcentral}(),
        dir = true
    ),
    jvp_autodiff = AutoForwardDiff(),
    concrete_jac = Val{true}()
)

----            -------------           -----------
Iter            f(u) 2-norm             Step 2-norm
----            -------------           -----------
0               2.09308241e+02          0.00000000e+00
1               7.32421822e+01          6.78876231e+00
2               2.94288416e+01          3.74900805e+00
3               1.32998631e+01          5.37657850e+00
ERROR: DimensionMismatch: new dimensions (19,) must be consistent with array size 20
Stacktrace:
  [1] (::Base.var"#throw_dmrsa#365")(dims::Tuple{Int64}, len::Int64)
    @ Base .\reshapedarray.jl:41
  [2] reshape
    @ .\reshapedarray.jl:45 [inlined]
  [3] reshape
    @ .\reshapedarray.jl:127 [inlined]
  [4] restructure
    @ C:\Users\iwheeler\.julia\packages\ArrayInterface\9um5f\src\ArrayInterface.jl:642 [inlined]
  [5] restructure
    @ C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\utils.jl:100 [inlined]
  [6] solve!(cache::NonlinearSolveBase.GeodesicAccelerationCache{…}, J::Matrix{…}, fu::Vector{…}, u::Vector{…}, idx::Val{…}; skip_solve::Bool, kwargs::@Kwargs{…})
    @ NonlinearSolveBase C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\descent\geodesic_acceleration.jl:116
  [7] solve! (repeats 2 times)
    @ C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\descent\geodesic_acceleration.jl:98 [inlined]
  [8] step!(cache::NonlinearSolveFirstOrder.GeneralizedFirstOrderAlgorithmCache{…}; recompute_jacobian::Nothing)
    @ NonlinearSolveFirstOrder C:\Users\iwheeler\.julia\packages\NonlinearSolveFirstOrder\mQbg1\src\solve.jl:246
  [9] step!
    @ C:\Users\iwheeler\.julia\packages\NonlinearSolveFirstOrder\mQbg1\src\solve.jl:224 [inlined]
 [10] #step!#155
    @ C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\solve.jl:253 [inlined]
 [11] step!
    @ C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\solve.jl:247 [inlined]
 [12] solve!(cache::NonlinearSolveFirstOrder.GeneralizedFirstOrderAlgorithmCache{…})
    @ NonlinearSolveBase C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\solve.jl:18
 [13] __solve(::NonlinearLeastSquaresProblem{…}, ::GeneralizedFirstOrderAlgorithm{…}; kwargs::@Kwargs{…})
    @ NonlinearSolveBase C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\solve.jl:6
 [14] __solve
    @ C:\Users\iwheeler\.julia\packages\NonlinearSolveBase\NIRyZ\src\solve.jl:1 [inlined]
 [15] #solve_call#36
    @ C:\Users\iwheeler\.julia\packages\DiffEqBase\cu9Gd\src\solve.jl:670 [inlined]
 [16] solve_call
    @ C:\Users\iwheeler\.julia\packages\DiffEqBase\cu9Gd\src\solve.jl:627 [inlined]
 [17] #solve_up#45
    @ C:\Users\iwheeler\.julia\packages\DiffEqBase\cu9Gd\src\solve.jl:1208 [inlined]
 [18] solve_up
    @ C:\Users\iwheeler\.julia\packages\DiffEqBase\cu9Gd\src\solve.jl:1201 [inlined]
 [19] #solve#43
    @ C:\Users\iwheeler\.julia\packages\DiffEqBase\cu9Gd\src\solve.jl:1096 [inlined]
 [20] top-level scope
    @ c:\Users\iwheeler\.julia\dev\LyoPronto\src\mwe.jl:38
Some type information was truncated. Use `show(err)` to see complete types.

(In my application the change in length is deterministic and the result of a DAE solution, but I think this gets to the heart of it more easily.)