Determined complete system raises an overdetermined warn when creating ODEProblem

[pmc4]

Describe the bug 🐞

I have created a system with structural_simplify that has 2 equations, 2 unknowns and 8 observed variables. The system is marked as complete and it's a determined system of equations. When I create the problem with ODEProblem, a warning stating that the initialization system is overdetermined appears, saying that there are 7 equations for 0 unknowns.

I can still solve the problem and get a success retcode, though.

Expected behavior

I should be able to pass a complete determined system to the ODEProblem constructor and create the problem without warnings.

Minimal Reproducible Example 👇

I have simplified a lot the structure of my models, but in essence is what I've written below. I know I'm using @independent_variables x instead of t_nounits.

using ModelingToolkit, DifferentialEquations, Integrals

@independent_variables x
D = Differential(x)

rhoA(z) = z^4
rhoB(z) = z^3
PB(z) = z^4

# Define models and solve
@mtkmodel ModelA begin
    @variables begin
        z(x)
        ρ(x)
        P(x)
        ρ̇(x)
    end
    @equations begin
        ρ ~ rhoA(z)
        P ~ ρ / 3
        ρ̇ ~ 0
    end
end


@mtkmodel ModelB begin
    @variables begin
        z(x)
        ρ(x)
        P(x)
        ρ̇(x)
    end
    @equations begin
        ρ ~ rhoB(z)
        P ~ PB(z)
        ρ̇ ~ (ρ - 3*P) / x
    end
end


function Model()
    z0 = 1.0

    # Model A
    ρ0A = rhoA(z0)
    P0A = ρ0A / 3
    z0A = z0
    @named A = ModelA(ρ=ρ0A, P=P0A, z=z0A)

    # Model B
    ρ0B = rhoB(z0)
    P0B = PB(z0)
    z0B = z0
    @named B = ModelB(ρ=ρ0B, P=P0B, z=z0B)

    system = [A, B]
    ρtot0 = ρ0A + ρ0B

    @variables z(x) = z0 ρtot(x)

    eqs = [
        ρtot ~ A.ρ + B.ρ
        D(ρtot) ~ A.ρ̇ + B.ρ̇
        A.z ~ z
        B.z ~ z
    ]

    @named _model = ODESystem(eqs, x)
    @named model = compose(_model, system)

    return model, ρtot0
end


model, ρtot0 = Model()
sys = structural_simplify(model)


xspan = (1e-2, 1e2)
prob = ODEProblem(sys, [sys.ρtot => ρtot0], xspan, []; warn_initialize_determined = true)
sol = solve(prob)

Error & Stacktrace ⚠️

┌ Warning: Initialization system is overdetermined. 7 equations for 0 unknowns. Initialization will default to using least squares. `SCCNonlinearProblem` can only be used for initialization of fully determined systems and hence will not be used here. To suppress this warning pass warn_initialize_determined = false. To make this warning into an error, pass fully_determined = true
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/41wGH/src/systems/diffeqs/abstractodesystem.jl:1358

Environment (please complete the following information):

• Output of using Pkg; Pkg.status()

Status `~/software/test/Project.toml`
  [0c46a032] DifferentialEquations v7.15.0
  [de52edbc] Integrals v4.5.0
  [961ee093] ModelingToolkit v9.61.0

• Output of versioninfo()

Julia Version 1.11.2
Commit 5e9a32e7af2 (2024-12-01 20:02 UTC)
Build Info:
  Official https://julialang.org/ release
Platform Info:
  OS: Linux (x86_64-linux-gnu)
  CPU: 20 × 12th Gen Intel(R) Core(TM) i7-12700H
  WORD_SIZE: 64
  LLVM: libLLVM-16.0.6 (ORCJIT, alderlake)
Threads: 1 default, 0 interactive, 1 GC (on 20 virtual cores)
Environment:
  JULIA_EDITOR = code
  JULIA_NUM_THREADS =