MTK and Jacobian

[BLI]

I have an unbalanced/partial model consisting of 2 differential equations, 11 algebraic equations, and 17 unknowns/states. I want to compute the Jacobian of the model.

I have 3 questions regarding the calculate_jacobian method, see at the bottom.

Here is the model (just to show the structure), named mod_p:

…and here are the 17 unknowns:

I then calculate the Jacobian of the partial model:

J = calculate_jacobian(mod_p)

Questions:

1. When computing J, are the equations preserved in exactly the order shown above (listing of equations for mod_p)?
2. The algebraic equations: I assume all terms are moved to one side of the equation before taking the Jacobian, i.e., that one has algebraic equations 0=g(x)?
3. When computing J, are the unknowns preserved in exactly the order given by states(mod_p)?
4. If the answer to any of questions 1 & 3 is no, is there a way to find the orders used by calculate_jacobian?

[ChrisRackauckas]

@BLI

When computing J, are the equations preserved in exactly the order shown above (listing of equations for mod_p)?

They are computed in the order of states(sys). This ordering can change due to simplifications.

@BLI

The algebraic equations: I assume all terms are moved to one side of the equation before taking the Jacobian, i.e., that one has algebraic equations 0=g(x)0=g(x)0=g(x)?

yes

@BLI

• When computing J, are the unknowns preserved in exactly the order given by states(mod_p)?

yes

[BLI]

@ChrisRackauckas

They are computed in the order of states(sys).

Hm. I need to keep track of the order of the equations.

If I have two expressions (terms?) and an equation expr1 ~ expr2, I can extract the Left Hand Side and the Right Hand Side by

@variables expr1 expr2
(exrp1 ~ expr2).lhs
(expr1 ~ expr2).rhs

Is there a function to do the same, or should I use comprehension if I want to do this for a vector of equations?

[BLI]

OK. I extracted the RHS of equations using field.rhs etc., and created a vector of equations and a vector of variables (states).

Then I used Symbolics.jacobian(rhs,states(sys_p)). The advantage of this is that both rhs and states(sys_p) are vectors, so I know precisely which order the content of each vector has, and therefore I can check each element of the Jacobian.

The result seems to be the same as when using ModelingToolkit.calculate_jacobian(sys_p). I don’t know if this will always be the case, though.

Final question:

• If I build a model sys = ODESystem(eqs) where eqs is a vector of equations, is there a way to extract the equations eqs from sys and get them in exactly the same order as in eqs?
• Will function equations do this, or is there a chance that this function will change the order or change the terms in the equations?

[ChrisRackauckas]

@BLI

• If I build a model sys = ODESystem(eqs) where eqs is a vector of equations, is there a way to extract the equations eqs from sys and get them in exactly the same order as in eqs?

Not necessarily. There won’t even necessarily be the same number of equations, or the same equations, so it’s not a well-defined question.

@BLI

• Will function equations do this, or is there a chance that this function will change the order or change the terms in the equations?

It will give you the current equations within the system, which could be simplified in many ways.

[BLI]

@ChrisRackauckas

Not necessarily. There won’t even necessarily be the same number of equations, or the same equations, so it’s not a well-defined question.

In one case, we put the differential equations at the bottom of eqs_p. Then eqs_p |> ODESystem |> equations had them reordered with the differential equations on top. Also, the order of variables found from get_variables(eqs_p) hand been altered in eqs_p |> ODESystem |> states with the differential variables on top.

Anyways, I can confirm that converting eqs_p into an ODESystem may change order of equations and variables. And as you (CR) say: it can even eliminate variables/equations.