Hi! I’m trying to use ModelingToolkit with matrix variables/parameters. I’ve got it working, but I need to be very careful about how I define things, and a couple behaviors feel non-intuitive. I’m wondering if I’m missing an idiomatic way to do this, or if what I’m seeing is expected.
Minimal example 1 (works, but with quirks)
using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D
using OrdinaryDiffEq
@mtkmodel FOL begin
@parameters begin
A1[1:2,1:2] = [-1 0;
0 -2]
v[1:1,1:2] = [1 2]
k[1:2,1:1] = [1; 1]
end
@variables begin
x(t)[1:2,1:1]
x1(t) = 1
x2(t) = 1
u(t)
end
@equations begin
x ~ [x1 x2]'
D(x) ~ A1 * x * u
u ~ (v * k)[1,1] # (v*k) is 1×1, so I index [1,1] to get a scalar for u
end
end
@mtkcompile fol = FOL()
t2 = 3.0
prob = ODEProblem(fol, Dict(fol.x => [defaults(fol)[fol.x1], defaults(fol)[fol.x2]]), (0.0, t2))
sol = solve(prob, saveat=0.01)
Observations:
- Defaults don’t seem to propagate from
x1(t),x2(t)tox(t)even thoughx ~ [x1 x2]'. With scalar variables I’m used to defaults flowing, but with matrix variables I have to manually supplyfol.xin theODEProbleminit. - Because
vis 1×2 andkis 2×1,v*kis 1×1 (a matrix), not a scalar. I have to do(v*k)[1,1]to feedu(t).
Minimal example 2 (also works, but more plumbing)
If I instead define u(t) as a 1×1 matrix, then I don’t need [1,1] on the RHS—but I have to provide u’s initial condition explicitly and in matrix shape:
using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as D
using OrdinaryDiffEq
@mtkmodel FOL begin
@parameters begin
A1[1:2,1:2] = [-1 0;
0 -2]
v[1:1,1:2] = [1 2]
k[1:2,1:1] = [1; 1]
end
@variables begin
x(t)[1:2,1:1]
x1(t) = 1
x2(t) = 1
u(t)[1:1,1:1]
end
@equations begin
x ~ [x1 x2]'
D(x) ~ A1 * x * u
u ~ (v * k)
end
end
@mtkcompile fol = FOL()
t2 = 3.0
prob = ODEProblem(
fol,
Dict(
fol.x => [defaults(fol)[fol.x1], defaults(fol)[fol.x2]],
fol.u => 3 * ones(1,1) # must supply as 1×1 matrix
),
(0.0, t2)
)
sol = solve(prob, saveat=0.01)
This works, but feels a bit clunky compared to the scalar case. Also now it would fail if I had
D(x) ~ u * A1 * x
since the matrix dimension would be wrong (even though it works fine math wise since u is a scalar)
Questions
-
Default propagation for array variables:
Is it expected that defaults on the scalar components (x1,x2) do not propagate to the array variable (x) when using an equation likex ~ [x1 x2]'? Is there a recommended way to have defaults “flow” intox(t)[…]without manually building the vector forODEProblem? -
Scalar vs 1×1 matrix ergonomics:
When expressions likev*knaturally produce a 1×1 matrix, what’s the idiomatic way to feed a scalar variable on the LHS? Is(v*k)[1,1]the intended approach? -
Defining array variables/parameters idiomatically:
Is there a best practice for declaring vectors vs 2D “column vectors” in@variables/@parameters? For example, should “state vectors” generally bex(t)[1:2](1D) rather thanx(t)[1:2,1:1](2D)? Any pitfalls to be aware of when mixing 1D and 2D shapes in equations?
Side note on Julia array syntax (source of confusion)
Just to show where part of my confusion came from:
[1 1] # 1×2 Matrix
[1,1] # 2-element Vector
[1;1] # 2-element Vector
ones(2,1) # 2×1 Matrix
I get why [1 1] vs [1,1] differ, but I still find myself tripping over it when switching between vector and matrix forms in symbolic equations.
Thanks a lot! Any guidance, docs pointers, or examples would be super appreciated.
