Unbalanced system error when using "compose"

[ctessum]

Hello,

I am trying to compose two systems using the "compose" function. For example, a falling body:

using ModelingToolkit

function newton()
    @variables t
    @variables g(t)
    @variables x(t) v(t)
    D = Differential(t)
    eq = [
        D(v) ~ g,
        D(x) ~ v
    ]

    @named sys = ODESystem(eq)
end

function gravity()
    @variables t
    @variables g(t)
    eq = [g ~ -9.8]
    @named gravity = ODESystem(eq)
end

sys = newton()
grav = gravity()

@named connection = ODESystem([
    sys.g ~ grav.g
])

connected = compose(sys, grav, connection)

simplified_sys = structural_simplify(connected)

In the code above, the dynamics are specified in one system and the amount of gravity is specified in another system, and then I compose them together. The last line of code above gives the error: ERROR: ExtraVariablesSystemException: The system is unbalanced. There are 3 highest order derivative variables and 2 equations.

If I specify everything in one system, however, it works fine:

function newton_combined()
    @variables t
    @variables g(t)
    @variables x(t) v(t)
    D = Differential(t)
    eq = [
        D(v) ~ g,
        D(x) ~ v,
        g ~ -9.8
    ]

    @named sys = ODESystem(eq)
end

sys_combined = newton_combined()
structural_simplify(sys_combined)

I would like to be able to specify things in multiple systems and combine them. In the example above, hypothetically, I might want to use my falling object in different gravity fields or a gravity field that varies with time, and don't want to have to copy and paste the dynamics every time. (Obviously this is a stand-in for a much more complex use case.)

I've moved this over from discourse because based on discussion there it seems that this might be a software bug rather than user error.

Let me know any thoughts anyone might have. Thanks for your help!

[ctessum]

cc @jialinl6

[ChrisRackauckas]

The combined system that you wrote is not the same as composition you created. It would instead be:

function newton_combined()
    @variables t
    @variables g(t) g2(t)
    @variables x(t) v(t)
    D = Differential(t)
    eq = [
        D(v) ~ g,
        D(x) ~ v,
        g ~ -9.8,
        g1 ~ g2
    ]

    @named sys = ODESystem(eq)
end

which you can see is unbalanced and has an extra equation. To handle this, you should probably be using a connection statement instead of building a new equation, or use extend, i.e. inheritance, so that you don't need to define the equality (though standard caveats to inheritance apply). See https://docs.sciml.ai/ModelingToolkit/stable/basics/Composition/#Inheritance-and-Combine for more details.

Let me know if there's anything else. Closing because there's not much more here.

[ctessum]

Thanks! I'm still confused, though, because the example you give above as being unbalanced actually runs just fine (with a minor change from g1 to g which I assume was a typo):

function newton_combined()
    @variables t
    @variables g(t) g2(t)
    @variables x(t) v(t)
    D = Differential(t)
    eq = [
        D(v) ~ g,
        D(x) ~ v,
        g ~ -9.8,
        g ~ g2,
    ]

    @named sys = ODESystem(eq)
end

sys_combined = newton_combined()
structural_simplify(sys_combined)

julia> observed(structural_simplify(sys_combined))
4-element Vector{Equation}:
 xˍtt(t) ~ Differential(t)(v(t))
 g(t) ~ -9.8
 g2(t) ~ g(t)
 xˍtt(t) ~ g(t)

Let me know if I'm missing something here! Thanks.

[ctessum]

In case anyone is looking in the future this is how I figured out to do it with extend:

using ModelingToolkit

@variables t

function newton(t)
    @variables g(t)
    @variables x(t) v(t)
    D = Differential(t)
    eq = [
        D(v) ~ g,
        D(x) ~ v
    ]

    @named sys = ODESystem(eq)
end

function gravity(newton, t; name)
    @unpack g = newton
    eq = [g ~ -9.8]
    extend(ODESystem(eq, t, [], []; name=name), newton)
end

sys = newton(t)
@named grav = gravity(sys, t)

structural_simplify(grav)

However, I'm still not sure why it doesn't work with compose like above.

[ChrisRackauckas]

Interesting, yeah @YingboMa is that not caught by the improved alias elimination?

[YingboMa]

Isn't it just because compose introduces another level of namespacing?

[ctessum]

Yeah, that's it. For posterity, this is the solution:

using ModelingToolkit

function newton()
    @variables t
    @variables g(t)
    @variables x(t) v(t)
    D = Differential(t)
    eq = [
        D(v) ~ g,
        D(x) ~ v
    ]

    @named sys = ODESystem(eq)
end

function gravity()
    @variables t
    @variables g(t)
    eq = [g ~ -9.8]
    @named gravity = ODESystem(eq)
end

sys = newton()
grav = gravity()

@named connection = ODESystem([
    sys.g ~ grav.g
])

connected = extend(sys, compose(grav, connection))

simplified_sys = structural_simplify(connected)

The key here is that it should be connected = extend(sys, compose(grav, connection)) instead of connected = compose(sys, grav, connection). Using extend instead of compose when adding the connector system avoids creating a new namespace for the connector system.

Thanks for your help!