# Put linear ODE system in matrix-vector form
**URL:** https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170
**Category:** Modelling & Simulations
**Tags:** modelingtoolkit
**Created:** [October 11, 2024, 12:59am UTC](https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170 "2024-10-11T00:59:29Z")
**Posts on this page:** 6
**Page:** 1
### Author:  [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [October 11, 2024, 12:59am UTC](https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170/1 "2024-10-11T00:59:30Z")
Hi,
I am trying to put my linear ODE system in state space form.
My model is quite simple:
```julia
function TiThTe(input_funcs::Dict)
@parameters begin
C_i
C_e
C_h
R_ie
R_ih
R_ea
A_w
A_e
T_i_0
T_e_0
T_h_0
end
@variables begin
T_i(t) = T_i_0
T_e(t) = T_e_0
T_h(t) = T_h_0
T_a(t)
P_h(t)
P_s(t)
end
# Define time-varying inputs from the input_funcs dictionary
@named T_a_out = TimeVaryingFunction(input_funcs[:T_a])
@named P_h_out = TimeVaryingFunction(input_funcs[:P_h])
@named P_s_out = TimeVaryingFunction(input_funcs[:P_s])
systems = [T_a_out, P_h_out, P_s_out]
# define the system of differential equations
eqs = [
# dTi equation
C_i * D(T_i) ~ (1 / R_ih) * (T_h - T_i) + (1 / R_ie) * (T_e - T_i) + A_w * P_s,
# dTe equation
C_e * D(T_e) ~ (1 / R_ie) * (T_i - T_e) + (1 / R_ea) * (T_a - T_e) + A_e * P_s,
# dTh equation
C_h * D(T_h) ~ (1 / R_ih) * (T_i - T_h) + P_h
]
# Define the inputs (mapping them to their respective time-varying functions)
inputs = [
T_a ~ T_a_out.output.u,
P_h ~ P_h_out.output.u,
P_s ~ P_s_out.output.u
]
append!(eqs, inputs)
@named sys = ODESystem(eqs, t; systems=systems)
return sys
end
```
My input functions are given by DataInterpolations.jl.
I try to form the matrices in this way, defining my inputs and outputs:
```julia
# inputs
@variables t T_a(t) P_h(t) P_s(t)
inputs = [T_a, P_h, P_s]
# outputs
@variables T_i(t) T_e(t) T_h(t)
outputs= [T_i, T_e, T_h]
# linearize
ss_matrices, ssys = ModelingToolkit.linearize_symbolic(sys, inputs, outputs)
```
Unfortunately this method complains about missing variables which I dont want to see in my linearization:
> ExtraEquationsSystemException: The system is unbalanced. There are 6 highest order derivative variables and 9 equations.
> More equations than variables, here are the potential extra equation(s):
> 0 ~ -T\_a\_out₊output₊u(t) + DataInterpolations.AkimaInterpolation{Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Vector{Float64}, Float64}([8.167, 8.187, 8.077, 7.988, 7.98, 7.814, 7.784, 7.531, 7.638, 8.057
But I tried adding them anyway to see what result I get:
```julia
@variables t T_a(t) P_h(t) P_s(t) T_a_out₊output₊u(t) P_h_out₊output₊u(t) P_s_out₊output₊u(t)
inputs = [T_a, P_h, P_s, T_a_out₊output₊u, P_h_out₊output₊u, P_s_out₊output₊u]
```
This throws a boundserror in AbstractSystem.jl:
> BoundsError: attempt to access 9×3 Matrix{Num} at index [1:3, 4:9]
> Stacktrace:
> [1] throw\_boundserror(A::Matrix{Num}, I::Tuple{UnitRange{Int64}, UnitRange{Int64}})
> @ Base .\abstractarray.jl:737
> [2] checkbounds
> @ .\abstractarray.jl:702 [inlined]
> [3] \_getindex
> @ .\multidimensional.jl:888 [inlined]
> [4] getindex(::Matrix{Num}, ::UnitRange{Int64}, ::UnitRange{Int64})
> @ Base .\abstractarray.jl:1291
> [5] linearize\_symbolic(sys::ODESystem, inputs::Vector{Num}, outputs::Vector{Num}; simplify::Bool, allow\_input\_derivatives::Bool, eval\_expression::Bool, eval\_module::Module, kwargs::@Kwargs{})
> @ ModelingToolkit C:\Users\lange.julia\packages\ModelingToolkit\Vsl3C\src\systems\abstractsystem.jl:2496
> [6] linearize\_symbolic(sys::ODESystem, inputs::Vector{Num}, outputs::Vector{Num})
> @ ModelingToolkit C:\Users\lange.julia\packages\ModelingToolkit\Vsl3C\src\systems\abstractsystem.jl:2472
Could somebody advise me what to do?
---
### Author:  [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [October 12, 2024, 7:15am UTC](https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170/2 "2024-10-12T07:15:34Z")
The first error, extra equations, is expected. This video explains why it appears and how to change your model to make it work
[ ](https://www.youtube.com/watch?v=-XOux-2XDGI)
Tldr is that you cannot have the DataInterpolations connected while linearizing, the inputs must be unconnected.
This other video shared some additional tools for linear analysis you might find interesting
[ ](https://www.youtube.com/watch?v=aThEtmuuazo)
---
### Author:  [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [October 12, 2024, 7:23am UTC](https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170/3 "2024-10-12T07:23:40Z")
Here’s a working example
```julia
using ModelingToolkit, ModelingToolkitStandardLibrary
import ModelingToolkit: t_nounits as t, D_nounits as D
function TiThTe(; name)
@parameters begin
C_i
C_e
C_h
R_ie
R_ih
R_ea
A_w
A_e
T_i_0
T_e_0
T_h_0
end
@variables begin
T_i(t) = T_i_0
T_e(t) = T_e_0
T_h(t) = T_h_0
T_a(t)
P_h(t)
P_s(t)
end
# define the system of differential equations
eqs = [
# dTi equation
C_i * D(T_i) ~ (1 / R_ih) * (T_h - T_i) + (1 / R_ie) * (T_e - T_i) + A_w * P_s,
# dTe equation
C_e * D(T_e) ~ (1 / R_ie) * (T_i - T_e) + (1 / R_ea) * (T_a - T_e) + A_e * P_s,
# dTh equation
C_h * D(T_h) ~ (1 / R_ih) * (T_i - T_h) + P_h
]
ODESystem(eqs, t; name)
end
@named sys = TiThTe()
sys = complete(sys)
# inputs
inputs = [sys.T_a, sys.P_h, sys.P_s]
# outputs
outputs= [sys.T_i, sys.T_e, sys.T_h]
# linearize
ss_matrices, ssys = ModelingToolkit.linearize_symbolic(sys, inputs, outputs)
```
```julia-repl
julia> ss_matrices.A
3×3 Matrix{Num}:
(-1 / R_ie + -1 / R_ih) / C_i 1 / (C_i*R_ie) 1 / (C_i*R_ih)
1 / (C_e*R_ie) (-1 / R_ea + -1 / R_ie) / C_e 0
1 / (C_h*R_ih) 0 -1 / (C_h*R_ih)
```
---
### Author:  [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [October 12, 2024, 12:48pm UTC](https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170/4 "2024-10-12T12:48:42Z")
Thanks a lot, that works great 🙂
Is there perhaps a way to apply this to an already existing model, basically remove the connections? So I don’t have to keep two model definitions.
---
### Author:  [@baggepinnen](https://discourse.julialang.org/u/baggepinnen)
#### Post date: [October 12, 2024, 1:32pm UTC](https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170/5 "2024-10-12T13:32:56Z")
Yes, this is accomplished by the analysis points described in the videos. You can also define the model without connections, and a thin outer model that adds the input connections when they are needed. This is a more extensible approach anyways so I recommend following this approach
---
### Author:  [@langestefan](https://discourse.julialang.org/u/langestefan)
#### Post date: [October 12, 2024, 1:51pm UTC](https://discourse.julialang.org/t/put-linear-ode-system-in-matrix-vector-form/121170/6 "2024-10-12T13:51:51Z")
Yes I also like this approach. I think it can be done using the new MTK v9 syntax for datainterpolations? Like here: [https://github.com/SciML/ModelingToolkit.jl/issues/2823#issuecomment-2366342737](https://github.com/SciML/ModelingToolkit.jl/issues/2823#issuecomment-2366342737)
(last time I checked there weren’t any docs on how to do this exactly and I couldn’t figure out the syntax myself for my model)