MSD system with N-DoF
The source files for all examples can be found in /examples.
One of the advantages of using a scripting approach to model dynamic systems is the ability to create and simulate the system algorithmically. To showcase this capability using the Bond Graph Toolkit, let us consider an example of a mass-spring-damper (MSD) system with $n$ degrees of freedom. We will use a for loop to vary the number of masses. The n-degree-of-freedom (n-DoF) MSD system is depicted below:
TODO: Update the image 
The following differential equation governs the motion of this n-DoF harmonic oscillator:
TODO: Update the equation $m\,\frac{d^2 x}{dt^2} + c\,\frac{d^2 x}{dt^2} + k\,x = F(t)$,
where $x_n$ is the position, $m_n$ is the mass, $c_n$ is the damping coefficient, $k_n$ is the spring stiffness, and $F(t)$ is an external force applied to the system.
To model the n-DoF system using the bond graph method, we need to connect inertance (I), compliance (C), and resistance (R) elements to different 1-junctions and 0-junctions, as shown in the bond graph diagram:
TODO: Generate the bond graph for the n-DoF system 
In this representation, the inertance element (I) corresponds to the mass (m) in the mechanical system, the compliance element (C) corresponds to the spring constant (k), and the resistance element (R) represents the damping coefficient (c). For simplicity, we assume they are equal for both degrees of freedom.
Bond Graph Toolkit
We need to import the Bond Graph Toolkit module and the DifferentialEquations.jl package to model the system and solve the resultant ODE. Additionally, we will import the independent variable t from the Bond Graph Toolkit for defining a custom forcing term.
using BondGraph
using BondGraph: t, D
using DifferentialEquations
using ModelingToolkit
using Plots
using Symbolics.LatexifyBuilding the model
Setting the elements
Before defining a system, we need to define the single port elements of the system. Then, we define the mass, spring, and damper elements as follows:
@named m = Mass(m = 1.0) # (kg)
@named s = Spring(k = 1.0) # (N/m)
@named d = Damper(c = 1.0) # (N⋅m⋅s)\[ \begin{align} power_{+}e\left( t \right) =& R power_{+}f\left( t \right) \end{align} \]
Setting the junctions and subsystems
Similarly to the 2-DoF example, in the n-DoF system, we need to define multiple 1-junctions and 0-junctions. To simplify the process, we can define them separately and later connect them accordingly. In the bond graph model of the n-DoF system, there is a section that repeats across the model, a spring, damper, and one-junction. To avoid redundancy, we can create this section as a subsystem and reuse it in other parts of the system where it is needed.
Firstly we define the ground-connected mass ($m_1$) as a subsystem.
@named m1 = Junction1([-1, m], [-1, d], [-1, s])\[ \begin{equation} \left[ \begin{array}{c} \mathrm{connect}\left( power, m_{+}power \right) \\ \mathrm{connect}\left( power, d_{+}power \right) \\ \mathrm{connect}\left( power, s_{+}power \right) \\ \frac{\mathrm{d} m_{+}power_{+}f\left( t \right)}{\mathrm{d}t} = \frac{m_{+}power_{+}e\left( t \right)}{m_{+}I} \\ d_{+}power_{+}e\left( t \right) = d_{+}R d_{+}power_{+}f\left( t \right) \\ s_{+}power_{+}e\left( t \right) = \frac{s_{+}q\left( t \right)}{s_{+}C} \\ \frac{\mathrm{d} s_{+}q\left( t \right)}{\mathrm{d}t} = s_{+}power_{+}f\left( t \right) \\ \end{array} \right] \end{equation} \]
Then, we create the n-DoF spring and damper subsystem.
@named sd = Junction1([-1, s], [-1, d])\[ \begin{equation} \left[ \begin{array}{c} \mathrm{connect}\left( power, s_{+}power \right) \\ \mathrm{connect}\left( power, d_{+}power \right) \\ s_{+}power_{+}e\left( t \right) = \frac{s_{+}q\left( t \right)}{s_{+}C} \\ \frac{\mathrm{d} s_{+}q\left( t \right)}{\mathrm{d}t} = s_{+}power_{+}f\left( t \right) \\ d_{+}power_{+}e\left( t \right) = d_{+}R d_{+}power_{+}f\left( t \right) \\ \end{array} \right] \end{equation} \]
Set the number of degrees of freedom. For this case, it should be greater or equal to $3$.
n = 44Then, we build the model using a for loop. The subsystems are stored inside the Matrix of subsystems (M), while the connections between the subsystems are stored inside the Vector of connections (cons).
M, cons = [], []
for i = 2:n
# Generate the subsystems of one DoF
j0 = Junction0(sd; name = Symbol("m" * string(i) * "_j0"))
j1 = Junction1([-1, m]; name = Symbol("m" * string(i) * "_j1"))
push!(M, [j0, j1])
# Connect the subsystems of one DoF
push!(cons, connect(j0.power, j1.power))
end┌ Warning: power contains 0 flow variables, yet 2 regular (non-flow, non-stream, non-input, non-output) variables. This could lead to imbalanced model that are difficult to debug. Consider marking some of the regular variables as input/output variables.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/GFEA2/src/systems/connectors.jl:40
┌ Warning: power contains 0 flow variables, yet 2 regular (non-flow, non-stream, non-input, non-output) variables. This could lead to imbalanced model that are difficult to debug. Consider marking some of the regular variables as input/output variables.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/GFEA2/src/systems/connectors.jl:40
┌ Warning: power contains 0 flow variables, yet 2 regular (non-flow, non-stream, non-input, non-output) variables. This could lead to imbalanced model that are difficult to debug. Consider marking some of the regular variables as input/output variables.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/GFEA2/src/systems/connectors.jl:40
┌ Warning: power contains 0 flow variables, yet 2 regular (non-flow, non-stream, non-input, non-output) variables. This could lead to imbalanced model that are difficult to debug. Consider marking some of the regular variables as input/output variables.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/GFEA2/src/systems/connectors.jl:40
┌ Warning: power contains 0 flow variables, yet 2 regular (non-flow, non-stream, non-input, non-output) variables. This could lead to imbalanced model that are difficult to debug. Consider marking some of the regular variables as input/output variables.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/GFEA2/src/systems/connectors.jl:40
┌ Warning: power contains 0 flow variables, yet 2 regular (non-flow, non-stream, non-input, non-output) variables. This could lead to imbalanced model that are difficult to debug. Consider marking some of the regular variables as input/output variables.
└ @ ModelingToolkit ~/.julia/packages/ModelingToolkit/GFEA2/src/systems/connectors.jl:40Connect the subsystems
We need now to iterate over the subsystem connecting them.
push!(cons, [connect(M[i][2].power, M[i+1][1].power) for i = 1:(length(M)-1)]...)5-element Vector{Any}:
connect(m2_j0.power, m2_j1.power)
connect(m3_j0.power, m3_j1.power)
connect(m4_j0.power, m4_j1.power)
connect(m2_j1.power, m3_j0.power)
connect(m3_j1.power, m4_j0.power)Finally, we need to connect the ground-connected mass ($m_1$) to the second subsystem.
push!(cons, connect(m1.power, M[1][1].power))6-element Vector{Any}:
connect(m2_j0.power, m2_j1.power)
connect(m3_j0.power, m3_j1.power)
connect(m4_j0.power, m4_j1.power)
connect(m2_j1.power, m3_j0.power)
connect(m3_j1.power, m4_j0.power)
connect(m1.power, m2_j0.power)Build the system
We flatten the subsystem Matrix into a Vector for simplicity.
M_flat = vcat(M...);We build the system by creating an ODESystem with the connections and setting the independent variable as t.
@named mdl = ODESystem(cons, t)\[ \begin{equation} \left[ \begin{array}{c} \mathrm{connect}\left( m2_{j0_{+}power}, m2_{j1_{+}power} \right) \\ \mathrm{connect}\left( m3_{j0_{+}power}, m3_{j1_{+}power} \right) \\ \mathrm{connect}\left( m4_{j0_{+}power}, m4_{j1_{+}power} \right) \\ \mathrm{connect}\left( m2_{j1_{+}power}, m3_{j0_{+}power} \right) \\ \mathrm{connect}\left( m3_{j1_{+}power}, m4_{j0_{+}power} \right) \\ \mathrm{connect}\left( m1_{+}power, m2_{j0_{+}power} \right) \\ \end{array} \right] \end{equation} \]
Then, we need to add the subsystems to the model.
mdl = compose(mdl, m1, M_flat...)\[ \begin{equation} \left[ \begin{array}{c} \mathrm{connect}\left( m2_{j0_{+}power}, m2_{j1_{+}power} \right) \\ \mathrm{connect}\left( m3_{j0_{+}power}, m3_{j1_{+}power} \right) \\ \mathrm{connect}\left( m4_{j0_{+}power}, m4_{j1_{+}power} \right) \\ \mathrm{connect}\left( m2_{j1_{+}power}, m3_{j0_{+}power} \right) \\ \mathrm{connect}\left( m3_{j1_{+}power}, m4_{j0_{+}power} \right) \\ \mathrm{connect}\left( m1_{+}power, m2_{j0_{+}power} \right) \\ \mathrm{connect}\left( power, m_{+}power \right) \\ \mathrm{connect}\left( power, d_{+}power \right) \\ \mathrm{connect}\left( power, s_{+}power \right) \\ \frac{\mathrm{d} m1_{+}m_{+}power_{+}f\left( t \right)}{\mathrm{d}t} = \frac{m1_{+}m_{+}power_{+}e\left( t \right)}{m1_{+}m_{+}I} \\ m1_{+}d_{+}power_{+}e\left( t \right) = m1_{+}d_{+}R m1_{+}d_{+}power_{+}f\left( t \right) \\ m1_{+}s_{+}power_{+}e\left( t \right) = \frac{m1_{+}s_{+}q\left( t \right)}{m1_{+}s_{+}C} \\ \frac{\mathrm{d} m1_{+}s_{+}q\left( t \right)}{\mathrm{d}t} = m1_{+}s_{+}power_{+}f\left( t \right) \\ \mathrm{connect}\left( sd_{+}power, power \right) \\ \mathrm{connect}\left( power, s_{+}power \right) \\ \mathrm{connect}\left( power, d_{+}power \right) \\ \mathrm{m2}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) = \frac{\mathrm{m2}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m2_{j0_{+}sd_{+}s_{+}C}} \\ \frac{\mathrm{d} \mathrm{m2}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} = \mathrm{m2}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ \mathrm{m2}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) = m2_{j0_{+}sd_{+}d_{+}R} \mathrm{m2}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ \mathrm{connect}\left( power, m_{+}power \right) \\ \frac{\mathrm{d} \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} = \frac{\mathrm{m2}_{j1_{+}m_{+}power_{+}e}\left( t \right)}{m2_{j1_{+}m_{+}I}} \\ \mathrm{connect}\left( sd_{+}power, power \right) \\ \mathrm{connect}\left( power, s_{+}power \right) \\ \mathrm{connect}\left( power, d_{+}power \right) \\ \mathrm{m3}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) = \frac{\mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m3_{j0_{+}sd_{+}s_{+}C}} \\ \frac{\mathrm{d} \mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} = \mathrm{m3}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ \mathrm{m3}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) = m3_{j0_{+}sd_{+}d_{+}R} \mathrm{m3}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ \mathrm{connect}\left( power, m_{+}power \right) \\ \frac{\mathrm{d} \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} = \frac{\mathrm{m3}_{j1_{+}m_{+}power_{+}e}\left( t \right)}{m3_{j1_{+}m_{+}I}} \\ \mathrm{connect}\left( sd_{+}power, power \right) \\ \mathrm{connect}\left( power, s_{+}power \right) \\ \mathrm{connect}\left( power, d_{+}power \right) \\ \mathrm{m4}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) = \frac{\mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}} \\ \frac{\mathrm{d} \mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} = \mathrm{m4}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ \mathrm{m4}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) = m4_{j0_{+}sd_{+}d_{+}R} \mathrm{m4}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ \mathrm{connect}\left( power, m_{+}power \right) \\ \frac{\mathrm{d} \mathrm{m4}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} = \frac{\mathrm{m4}_{j1_{+}m_{+}power_{+}e}\left( t \right)}{m4_{j1_{+}m_{+}I}} \\ \end{array} \right] \end{equation} \]
Analysing the model
We can visualize the bond-graph connections by generating a graph plot of the model.
generate_graph(mdl)We can obtain the bond-graph equations using the code below.
equations(expand_connections(mdl))\[ \begin{align} \frac{\mathrm{d} m1_{+}m_{+}power_{+}f\left( t \right)}{\mathrm{d}t} =& \frac{m1_{+}m_{+}power_{+}e\left( t \right)}{m1_{+}m_{+}I} \\ m1_{+}d_{+}power_{+}e\left( t \right) =& m1_{+}d_{+}R m1_{+}d_{+}power_{+}f\left( t \right) \\ m1_{+}s_{+}power_{+}e\left( t \right) =& \frac{m1_{+}s_{+}q\left( t \right)}{m1_{+}s_{+}C} \\ \frac{\mathrm{d} m1_{+}s_{+}q\left( t \right)}{\mathrm{d}t} =& m1_{+}s_{+}power_{+}f\left( t \right) \\ \mathrm{m2}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) =& \frac{\mathrm{m2}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m2_{j0_{+}sd_{+}s_{+}C}} \\ \frac{\mathrm{d} \mathrm{m2}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} =& \mathrm{m2}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ \mathrm{m2}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) =& m2_{j0_{+}sd_{+}d_{+}R} \mathrm{m2}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ \frac{\mathrm{d} \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} =& \frac{\mathrm{m2}_{j1_{+}m_{+}power_{+}e}\left( t \right)}{m2_{j1_{+}m_{+}I}} \\ \mathrm{m3}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) =& \frac{\mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m3_{j0_{+}sd_{+}s_{+}C}} \\ \frac{\mathrm{d} \mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} =& \mathrm{m3}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ \mathrm{m3}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) =& m3_{j0_{+}sd_{+}d_{+}R} \mathrm{m3}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ \frac{\mathrm{d} \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} =& \frac{\mathrm{m3}_{j1_{+}m_{+}power_{+}e}\left( t \right)}{m3_{j1_{+}m_{+}I}} \\ \mathrm{m4}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) =& \frac{\mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}} \\ \frac{\mathrm{d} \mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} =& \mathrm{m4}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ \mathrm{m4}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) =& m4_{j0_{+}sd_{+}d_{+}R} \mathrm{m4}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ \frac{\mathrm{d} \mathrm{m4}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} =& \frac{\mathrm{m4}_{j1_{+}m_{+}power_{+}e}\left( t \right)}{m4_{j1_{+}m_{+}I}} \\ - \mathrm{m2}_{j0_{+}power_{+}e}\left( t \right) =& m1_{+}power_{+}e\left( t \right) \\ - \mathrm{m2}_{j0_{+}power_{+}e}\left( t \right) =& \mathrm{m2}_{j0_{+}sd_{+}power_{+}e}\left( t \right) \\ - \mathrm{m3}_{j0_{+}power_{+}e}\left( t \right) =& \mathrm{m2}_{j1_{+}power_{+}e}\left( t \right) \\ - \mathrm{m3}_{j0_{+}power_{+}e}\left( t \right) =& \mathrm{m3}_{j0_{+}sd_{+}power_{+}e}\left( t \right) \\ 0 =& - \mathrm{m2}_{j1_{+}m_{+}power_{+}e}\left( t \right) - \mathrm{m2}_{j1_{+}power_{+}e}\left( t \right) + \mathrm{m2}_{j0_{+}power_{+}e}\left( t \right) \\ 0 =& - \mathrm{m3}_{j0_{+}power_{+}f}\left( t \right) + \mathrm{m2}_{j1_{+}power_{+}f}\left( t \right) + \mathrm{m3}_{j0_{+}sd_{+}power_{+}f}\left( t \right) \\ 0 =& - \mathrm{m4}_{j1_{+}m_{+}power_{+}e}\left( t \right) + \mathrm{m4}_{j0_{+}power_{+}e}\left( t \right) \\ 0 =& - \mathrm{m2}_{j0_{+}power_{+}f}\left( t \right) + m1_{+}power_{+}f\left( t \right) + \mathrm{m2}_{j0_{+}sd_{+}power_{+}f}\left( t \right) \\ 0 =& - \mathrm{m3}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) - \mathrm{m3}_{j0_{+}sd_{+}power_{+}e}\left( t \right) - \mathrm{m3}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) \\ - \mathrm{m4}_{j0_{+}power_{+}e}\left( t \right) =& \mathrm{m3}_{j1_{+}power_{+}e}\left( t \right) \\ - \mathrm{m4}_{j0_{+}power_{+}e}\left( t \right) =& \mathrm{m4}_{j0_{+}sd_{+}power_{+}e}\left( t \right) \\ - \mathrm{m4}_{j0_{+}sd_{+}power_{+}f}\left( t \right) =& - \mathrm{m4}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ - \mathrm{m4}_{j0_{+}sd_{+}power_{+}f}\left( t \right) =& - \mathrm{m4}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ - \mathrm{m4}_{j1_{+}m_{+}power_{+}f}\left( t \right) =& \mathrm{m4}_{j0_{+}power_{+}f}\left( t \right) \\ 0 =& - m1_{+}d_{+}power_{+}e\left( t \right) - m1_{+}m_{+}power_{+}e\left( t \right) - m1_{+}power_{+}e\left( t \right) - m1_{+}s_{+}power_{+}e\left( t \right) \\ - \mathrm{m2}_{j0_{+}sd_{+}power_{+}f}\left( t \right) =& - \mathrm{m2}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ - \mathrm{m2}_{j0_{+}sd_{+}power_{+}f}\left( t \right) =& - \mathrm{m2}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ 0 =& - \mathrm{m3}_{j1_{+}m_{+}power_{+}e}\left( t \right) - \mathrm{m3}_{j1_{+}power_{+}e}\left( t \right) + \mathrm{m3}_{j0_{+}power_{+}e}\left( t \right) \\ - m1_{+}power_{+}f\left( t \right) =& - m1_{+}s_{+}power_{+}f\left( t \right) \\ - m1_{+}power_{+}f\left( t \right) =& - m1_{+}m_{+}power_{+}f\left( t \right) \\ - m1_{+}power_{+}f\left( t \right) =& - m1_{+}d_{+}power_{+}f\left( t \right) \\ - \mathrm{m3}_{j1_{+}power_{+}f}\left( t \right) =& \mathrm{m3}_{j0_{+}power_{+}f}\left( t \right) \\ - \mathrm{m3}_{j1_{+}power_{+}f}\left( t \right) =& - \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right) \\ 0 =& - \mathrm{m4}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) - \mathrm{m4}_{j0_{+}sd_{+}power_{+}e}\left( t \right) - \mathrm{m4}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) \\ - \mathrm{m2}_{j1_{+}power_{+}f}\left( t \right) =& \mathrm{m2}_{j0_{+}power_{+}f}\left( t \right) \\ - \mathrm{m2}_{j1_{+}power_{+}f}\left( t \right) =& - \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right) \\ - \mathrm{m3}_{j0_{+}sd_{+}power_{+}f}\left( t \right) =& - \mathrm{m3}_{j0_{+}sd_{+}s_{+}power_{+}f}\left( t \right) \\ - \mathrm{m3}_{j0_{+}sd_{+}power_{+}f}\left( t \right) =& - \mathrm{m3}_{j0_{+}sd_{+}d_{+}power_{+}f}\left( t \right) \\ 0 =& - \mathrm{m4}_{j0_{+}power_{+}f}\left( t \right) + \mathrm{m3}_{j1_{+}power_{+}f}\left( t \right) + \mathrm{m4}_{j0_{+}sd_{+}power_{+}f}\left( t \right) \\ 0 =& - \mathrm{m2}_{j0_{+}sd_{+}d_{+}power_{+}e}\left( t \right) - \mathrm{m2}_{j0_{+}sd_{+}power_{+}e}\left( t \right) - \mathrm{m2}_{j0_{+}sd_{+}s_{+}power_{+}e}\left( t \right) \end{align} \]
The model's equations above were not simplified and represented the bond-graphs equations directly derived from the connections. Then, we simplify the DAE above and obtain the ODE. The simplified set of equations is required for using solve from DifferentialEquations.jl.
@named sys = simplifysys(mdl);Print system states.
states(sys)8-element Vector{SymbolicUtils.BasicSymbolic{Real}}:
m1₊m₊power₊f(t)
m1₊s₊q(t)
m2_j0₊sd₊s₊q(t)
m2_j1₊m₊power₊f(t)
m3_j0₊sd₊s₊q(t)
m3_j1₊m₊power₊f(t)
m4_j0₊sd₊s₊q(t)
m4_j1₊m₊power₊f(t)Print system parameters.
parameters(sys)12-element Vector{SymbolicUtils.BasicSymbolic{Real}}:
m4_j0₊sd₊s₊C
m3_j0₊sd₊s₊C
m1₊m₊I
m2_j0₊sd₊s₊C
m1₊s₊C
m1₊d₊R
m2_j0₊sd₊d₊R
m3_j0₊sd₊d₊R
m2_j1₊m₊I
m4_j0₊sd₊d₊R
m3_j1₊m₊I
m4_j1₊m₊IPrint the system's simplified equations.
equations(sys)\[ \begin{align} \frac{\mathrm{d} m1_{+}m_{+}power_{+}f\left( t \right)}{\mathrm{d}t} =& \frac{m2_{j0_{+}sd_{+}d_{+}R} \left( - m1_{+}m_{+}power_{+}f\left( t \right) - \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right) \right) + \frac{ - m1_{+}s_{+}q\left( t \right)}{m1_{+}s_{+}C} + \frac{ - \mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m3_{j0_{+}sd_{+}s_{+}C}} + \frac{ - \mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}} + \frac{\mathrm{m2}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m2_{j0_{+}sd_{+}s_{+}C}} + \frac{\mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m3_{j0_{+}sd_{+}s_{+}C}} + \frac{\mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}} - m1_{+}d_{+}R m1_{+}m_{+}power_{+}f\left( t \right)}{m1_{+}m_{+}I} \\ \frac{\mathrm{d} m1_{+}s_{+}q\left( t \right)}{\mathrm{d}t} =& m1_{+}m_{+}power_{+}f\left( t \right) \\ \frac{\mathrm{d} \mathrm{m2}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} =& - m1_{+}m_{+}power_{+}f\left( t \right) - \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right) \\ \frac{\mathrm{d} \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} =& \frac{m2_{j0_{+}sd_{+}d_{+}R} \left( - m1_{+}m_{+}power_{+}f\left( t \right) - \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right) \right) + m3_{j0_{+}sd_{+}d_{+}R} \left( - \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right) - \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right) \right) + \frac{ - \mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}} + \frac{\mathrm{m2}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m2_{j0_{+}sd_{+}s_{+}C}} + \frac{\mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m3_{j0_{+}sd_{+}s_{+}C}} + \frac{\mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}}}{m2_{j1_{+}m_{+}I}} \\ \frac{\mathrm{d} \mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} =& - \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right) - \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right) \\ \frac{\mathrm{d} \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} =& \frac{m3_{j0_{+}sd_{+}d_{+}R} \left( - \mathrm{m2}_{j1_{+}m_{+}power_{+}f}\left( t \right) - \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right) \right) + m4_{j0_{+}sd_{+}d_{+}R} \left( - \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right) - \mathrm{m4}_{j1_{+}m_{+}power_{+}f}\left( t \right) \right) + \frac{\mathrm{m3}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m3_{j0_{+}sd_{+}s_{+}C}} + \frac{\mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}}}{m3_{j1_{+}m_{+}I}} \\ \frac{\mathrm{d} \mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{\mathrm{d}t} =& - \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right) - \mathrm{m4}_{j1_{+}m_{+}power_{+}f}\left( t \right) \\ \frac{\mathrm{d} \mathrm{m4}_{j1_{+}m_{+}power_{+}f}\left( t \right)}{\mathrm{d}t} =& \frac{m4_{j0_{+}sd_{+}d_{+}R} \left( - \mathrm{m3}_{j1_{+}m_{+}power_{+}f}\left( t \right) - \mathrm{m4}_{j1_{+}m_{+}power_{+}f}\left( t \right) \right) + \frac{\mathrm{m4}_{j0_{+}sd_{+}s_{+}q}\left( t \right)}{m4_{j0_{+}sd_{+}s_{+}C}}}{m4_{j1_{+}m_{+}I}} \end{align} \]
Simulate the system
Unforced case
In this case, we will consider an initial displacement of $1.0 (m)$ on the second mass of the system.
Define the simulation time.
tspan = (0, 10)(0, 10)When defining the ODEProblem as in the ModelingToolkit, we can define the initial value of the states.
prob = ODEProblem(sys, [M[1][2].m.power.f => 0.0, M[1][1].sd.s.q => 1.0], tspan)ODEProblem with uType Vector{Float64} and tType Int64. In-place: true
timespan: (0, 10)
u0: 8-element Vector{Float64}:
0.0
0.0
1.0
0.0
0.0
0.0
0.0
0.0Use the solve to simulate the system.
sol = solve(prob)
plot(sol, xlabel="Time", ylabel="Amplitude")This page was generated using Literate.jl.