`MonteCarloMeasurements` for a `NonlinearProblem`

[bradcarman]

Not able to solve a NonlinearProblem that contains MonteCarloMeasurements Particles type. Also, split=false is required (see #2746)

MWE

using ModelingToolkit
using DifferentialEquations
using MonteCarloMeasurements

pars = @parameters a, b
vars = @variables x, y

eqs = [
    x ~ a
    y ~ x + b
]

@mtkbuild model = NonlinearSystem(eqs, vars, pars) split=false

u0 = Dict(
    model.x => 1, 
    model.y => 2
    )

p = Dict(
    model.a => 1 ± 0.1,
    model.b => 1 ± 0.1
    )

prob = NonlinearProblem(model, u0, p)
sol = solve(prob)

This returns the error

ERROR: MethodError: promote_u0(::Nothing, ::Vector{Particles{Float64, 2000}}, ::Nothing) is ambiguous.

Environment:

pkg> st
  [0c46a032] DifferentialEquations v7.13.0
⌅ [06fc5a27] DynamicQuantities v0.13.2
  [961ee093] ModelingToolkit v9.15.0
  [0987c9cc] MonteCarloMeasurements v1.1.6

[baggepinnen]

I wouldn't consider this a bug, it's not expected to work. To make this work, an explicit overload would be required that breaks the problem up into N individual problems and then recombine the result in the end.

[baggepinnen]

This is a temporary workaround

using ModelingToolkit
using MonteCarloMeasurements

pars = @parameters a, b
vars = @variables x, y

eqs = [
    x ~ a
    y ~ x + b
]

@mtkbuild model = NonlinearSystem(eqs, vars, pars) #split=false

u0 = Dict(
    x => 1, 
    y => 2
)


function solve_wrapper(a0, b0)
    p = Dict(
        model.a => a0,
        model.b => b0
        )

    prob = NonlinearProblem(model, u0, p)
    sol = solve(prob)
    [sol[model.x], sol[model.y]]
end


a0 = 1
b0 = 1
solve_wrapper(a0, b0)

a = 1 ± 0.1
b = 1 ± 0.1
xy = @bymap solve_wrapper(a, b)

julia> xy = @bymap solve_wrapper(a, b)
2-element Vector{Particles{Float64, 2000}}:
 1.0 ± 0.1
 2.0 ± 0.14