I’m struggling to figure out how to loop through numbered variables in ModelingToolkit [MTK].
As a simple, yet concrete example, consider a model of a telegraph line, which consists in a telegraph line module component
which is cascaded + an input voltage and a terminal load model:
# Importing packages
using ModelingToolkit
using ModelingToolkit: t_nounits as t, D_nounits as Dt
using DifferentialEquations
using Plots, LaTeXStrings
#
# Some plot constants
LW1 = 2.5
#
LS1 = :solid
LS2 = :dash
LS3 = :dot
#
LC1 = :blue
LC2 = :red
LC3 = :green
#
# Telegraph Line module
@mtkmodel Telegraph_Module begin
# Model parameters
@parameters begin
# Default parameters X_, X ∈ {R,L,G,C} taken from https://en.wikipedia.org/wiki/Telegrapher%27s_equations
R_=172.28, [description = "Unit resistance, Ω/km"]
L_=612.5e-6, [description = "Unit inductance, H/km"]
G_=0.072e-6, [description = "Unit conductance, S/km"]
C_=51.57e-9, [description = "Unit capacitance, F/km"]
Δx=1, [description = "Module length, km"]
R=R_*Δx, [description = "Resistance, Ω"]
L=L_*Δx, [description = "Inductance, H"]
G=G_*Δx, [description = "Conductance, S"]
C=C_*Δx, [description = "Capacitance, F"]
end
#
# Model variables, with initial values needed
@variables begin
# Inputs
u(t), [description = "Voltage at input, V"]
# Variables
u₊(t)=0.5, [description = "Voltage at output, V"]
i(t)=1e-3, [description = "Current at input, A"]
i₊(t), [description = "Current at output, A"]
end
#
# Telegraph line module equations
#
@equations begin
u ~ R*i + L*Dt(i) + u₊
i ~ G*u₊ + C*Dt(u₊) + i₊
end
end
# Telegraph Line system: cascaded modules with input voltage + terminal load
@mtkmodel Telegraph_Line begin
#
@parameters begin
R_L=600, [description = "Load resistance, Ω"]
end
#
@structural_parameters begin
ℓ = 10 # km
N = 2
end
# Components used
@components begin
tm = [Telegraph_Module(; Δx=ℓ/N) for i in 1:N]
end
# Equations for connecting components
@equations begin
# Distribution manifold splitting
tm[1].u ~ sin(2π*1_000*t) # Input at 1 kHz
[tm[j].u ~ tm[j-1].u₊ for j in 2:N]...
[tm[j].i ~ tm[j-1].i₊ for j in 2:N]...
tm[N].i₊ ~ tm[N].u₊/R_L # Terminal load
end
end
One would then instantiate the model, and solve the model, e.g.,
ℓ = 18 # length in km
N = 3 # number of cascaded modules
# Instantiation
@mtkbuild tL = Telegraph_Line(; ℓ=ℓ, N=N)
#
# Simulation set-up
tspan = (0, 1e-2)
#
prob = ODEProblem(tL, [], tspan)
sol = solve(prob)
One can then plot the solution, e.g., voltages:
plot(sol, idxs=[tL.tm_1.u₊, tL.tm_2.u₊, tL.tm_3.u₊], lw=LW1, lc=[LC1 LC2 LC3], ls=[LS1 LS2 LS3],
xlabel=L"t \mathrm{\ [s]}", ylabel=L"u \mathrm{\ [V]}",
title="Telegraph line voltage",
frame=:box,
ylims=(-0.5,0.75))
leading to (example):
Suppose I now want to use more cascaded modules, say N=10 or something? And that I want to plot all voltages? Then it becomes very impractical to explicitly write every voltage variable as:
plot(sol, idxs=[tL.tm_1.u₊, tL.tm_2.u₊, tL.tm_3.u₊, tL.tm_4.u₊, ....]
Question: How can I construct the vector of indices [tL.tm_1.u₊, tL.tm_2.u₊, tL.tm_3.u₊, tL.tm_4.u₊, ....] using comprehension, or similar techniques?
