Merge pull request #14 from ufechner7/broken

Update MTK to v9.40.1 and remove redundant initialization

Project.toml

@@ -25,7 +25,7 @@ Conda = "1.10.2"
 CondaPkg = "0.2.23"
 ControlPlots = "0.2.1"
 LiveServer = "1.3.1"
-ModelingToolkit = "~9.38.0"
+ModelingToolkit = "~9.40.1"
 OrdinaryDiffEq = "6.87.0"
 PackageCompiler = "2.1.17"
 Parameters = "0.12.3"

src/Tether_01.jl

@@ -1,13 +1,13 @@
 # Example one: Falling mass.
-using ModelingToolkit, OrdinaryDiffEq
+using ModelingToolkit, OrdinaryDiffEq, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
 G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81]    # gravitational acceleration     [m/s²]
 
 # definiting the model
 @variables   pos(t)[1:3] = [0.0, 0.0,  0.0]
 @variables   vel(t)[1:3] = [0.0, 0.0, 50.0] 
-@variables   acc(t)[1:3] = [0.0, 0.0, -9.81] 
+@variables   acc(t)[1:3]
 
 eqs = vcat(D(pos) ~ vel,
            D(vel) ~ acc,

src/Tether_02.jl

@@ -1,5 +1,5 @@
 # Example two: Falling mass, attached to linear spring
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
 G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81]    # gravitational acceleration     [m/s²]
@@ -9,10 +9,10 @@ L0::Float64 = -10.0                             # initial spring length      [m]
 @parameters mass=1.0 c_spring=50.0 damping=0.5 l0=L0
 @variables   pos(t)[1:3] = [0.0, 0.0,  L0]
 @variables   vel(t)[1:3] = [0.0, 0.0,  0.0] 
-@variables   acc(t)[1:3] = G_EARTH
-@variables unit_vector(t)[1:3]  = [0.0, 0.0, -sign(L0)]
-@variables spring_force(t)[1:3] = [0.0, 0.0, 0.0]
-@variables norm1(t) = abs(l0) spring_vel(t) = 0.0
+@variables   acc(t)[1:3]
+@variables unit_vector(t)[1:3]
+@variables spring_force(t)[1:3]
+@variables norm1(t) spring_vel(t)
 
 eqs = vcat(D(pos)       ~ vel,
            D(vel)       ~ acc,

src/Tether_03.jl

@@ -1,6 +1,6 @@
 # Example three: Falling mass, attached to non-linear spring without compression stiffness,
 # initially moving upwards with 4 m/s.
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
 G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81]    # gravitational acceleration     [m/s²]
@@ -11,11 +11,11 @@ V0::Float64 = 4                                 # initial velocity           [m/
 @parameters mass=1.0 c_spring0=50.0 damping=0.5 l0=L0
 @variables   pos(t)[1:3] = [0.0, 0.0,  L0]
 @variables   vel(t)[1:3] = [0.0, 0.0,  V0] 
-@variables   acc(t)[1:3] = G_EARTH
-@variables unit_vector(t)[1:3]  = [0.0, 0.0, -sign(L0)]
-@variables c_spring(t) = c_spring0
-@variables spring_force(t)[1:3] = [0.0, 0.0, 0.0]
-@variables force(t) = 0.0 norm1(t) = abs(l0) spring_vel(t) = 0.0
+@variables   acc(t)[1:3]
+@variables unit_vector(t)[1:3]
+@variables c_spring(t)
+@variables spring_force(t)[1:3]
+@variables force(t) norm1(t) spring_vel(t)
 
 eqs = vcat(D(pos)       ~ vel,
            D(vel)       ~ acc,

src/Tether_03b.jl

@@ -1,6 +1,6 @@
 # Example three: Falling mass, attached to non-linear spring without compression stiffness
 # with callback, initially moving upwards with 4 m/s.
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
 G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81]    # gravitational acceleration     [m/s²]
@@ -11,11 +11,11 @@ V0::Float64 = 4                                 # initial velocity           [m/
 @parameters mass=1.0 c_spring0=50.0 damping=0.5 l0=L0
 @variables   pos(t)[1:3] = [0.0, 0.0,  L0]
 @variables   vel(t)[1:3] = [0.0, 0.0,  V0] 
-@variables   acc(t)[1:3] = G_EARTH
-@variables unit_vector(t)[1:3]  = [0.0, 0.0, -sign(L0)]
-@variables c_spring(t) = c_spring0
-@variables spring_force(t)[1:3] = [0.0, 0.0, 0.0]
-@variables force(t) = 0.0 norm1(t) = abs(l0) spring_vel(t) = 0.0
+@variables   acc(t)[1:3]
+@variables unit_vector(t)[1:3]
+@variables c_spring(t)
+@variables spring_force(t)[1:3]
+@variables force(t) norm1(t) spring_vel(t)
 
 eqs = vcat(D.(pos)      ~ vel,
            D.(vel)      ~ acc,

src/Tether_03c.jl

@@ -1,6 +1,6 @@
 # Example three: Falling mass, attached to non-linear spring without compression stiffness,
 # initially moving upwards with 4 m/s. Comparing results with and without callbacks.
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
 G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81]    # gravitational acceleration     [m/s²]
@@ -12,11 +12,11 @@ function model3(G_EARTH, L0, V0)
     @parameters mass=1.0 c_spring0=50.0 damping=0.5 l0=L0
     @variables   pos(t)[1:3] = [0.0, 0.0,  L0]
     @variables   vel(t)[1:3] = [0.0, 0.0,  V0] 
-    @variables   acc(t)[1:3] = G_EARTH
-    @variables unit_vector(t)[1:3]  = [0.0, 0.0, -sign(L0)]
-    @variables c_spring(t) = c_spring0
-    @variables spring_force(t)[1:3] = [0.0, 0.0, 0.0]
-    @variables force(t) = 0.0 norm1(t) = abs(l0) spring_vel(t) = 0.0
+    @variables   acc(t)[1:3]
+    @variables unit_vector(t)[1:3]
+    @variables c_spring(t)
+    @variables spring_force(t)[1:3]
+    @variables force(t) norm1(t) spring_vel(t)
 
     eqs = vcat(D.(pos)   ~ vel,
             D.(vel)      ~ acc,

src/Tether_04.jl

@@ -5,7 +5,7 @@ for l < l_0) and n tether segments.
 # TODO: Distribute force correctly
 # TODO: Add 2D plot
 
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
 G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81]    # gravitational acceleration     [m/s²]
@@ -14,33 +14,23 @@ V0::Float64 = 4                                 # initial velocity of lowest mas
 segments::Int64 = 2                             # number of tether segments         [-]
 POS0 = zeros(3, segments+1)
 VEL0 = zeros(3, segments+1)
-ACC0 = zeros(3, segments+1)
-SEGMENTS0 = zeros(3, segments) 
-UNIT_VECTORS0 = zeros(3, segments)
 for i in 1:segments+1
     POS0[:, i] .= [0.0, 0, -(i-1)*L0]
     VEL0[:, i] .= [0.0, 0, (i-1)*V0/segments]
 end
-for i in 2:segments+1
-    ACC0[:, i] .= G_EARTH
-end
-for i in 1:segments
-    UNIT_VECTORS0[:, i] .= [0, 0, 1.0]
-    SEGMENTS0[:, i] .= POS0[:, i+1] - POS0[:, i]
-end
 
 # defining the model, Z component upwards
 @parameters mass=1.0 c_spring0=50.0 damping=0.5 l_seg=L0
 @variables pos(t)[1:3, 1:segments+1]  = POS0
 @variables vel(t)[1:3, 1:segments+1]  = VEL0
-@variables acc(t)[1:3, 1:segments+1]  = ACC0
-@variables segment(t)[1:3, 1:segments]  = SEGMENTS0
-@variables unit_vector(t)[1:3, 1:segments]  = UNIT_VECTORS0
-@variables norm1(t)[1:segments] = l_seg * ones(segments)
-@variables rel_vel(t)[1:3, 1:segments]  = zeros(3, segments)
-@variables spring_vel(t)[1:segments] = zeros(segments)
-@variables c_spring(t)[1:segments] = c_spring0 * ones(segments)
-@variables spring_force(t)[1:3, 1:segments] = zeros(3, segments)
+@variables acc(t)[1:3, 1:segments+1]
+@variables segment(t)[1:3, 1:segments]
+@variables unit_vector(t)[1:3, 1:segments]
+@variables norm1(t)[1:segments]
+@variables rel_vel(t)[1:3, 1:segments]
+@variables spring_vel(t)[1:segments]
+@variables c_spring(t)[1:segments]
+@variables spring_force(t)[1:3, 1:segments]
 
 # basic differential equations
 eqs1 = vcat(D.(pos) .~ vel,

src/Tether_05.jl

@@ -1,6 +1,6 @@
 # Tutorial example simulating a 3D mass-spring system with a nonlinear spring (no spring forces
 # for l < l_0) and n tether segments. 
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, Timers
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, Timers, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 
 G_EARTH::Vector{Float64} = [0.0, 0.0, -9.81]    # gravitational acceleration     [m/s²]
@@ -13,37 +13,27 @@ segments::Int64 = 5                             # number of tether segments
 duration = 10.0                                 # duration of the simulation        [s]
 POS0 = zeros(3, segments+1)
 VEL0 = zeros(3, segments+1)
-ACC0 = zeros(3, segments+1)
-SEGMENTS0 = zeros(3, segments) 
-UNIT_VECTORS0 = zeros(3, segments)
 for i in 1:segments+1
     local l0
     l0 = -(i-1) * L0
     v0 = (i-1) * V0/segments
     POS0[:, i] .= [sin(α0) * l0, 0, cos(α0) * l0]
     VEL0[:, i] .= [sin(α0) * v0, 0, cos(α0) * v0]
 end
-for i in 2:segments+1
-    ACC0[:, i] .= G_EARTH
-end
-for i in 1:segments
-    UNIT_VECTORS0[:, i] .= [0, 0, 1.0]
-    SEGMENTS0[:, i] .= POS0[:, i+1] - POS0[:, i]
-end
 
 # defining the model, Z component upwards
 @parameters mass=M0 c_spring0=C_SPRING damping=0.5 l_seg=L0
 @variables pos(t)[1:3, 1:segments+1]  = POS0
 @variables vel(t)[1:3, 1:segments+1]  = VEL0
-@variables acc(t)[1:3, 1:segments+1]  = ACC0
-@variables segment(t)[1:3, 1:segments]  = SEGMENTS0
-@variables unit_vector(t)[1:3, 1:segments]  = UNIT_VECTORS0
-@variables norm1(t)[1:segments] = l_seg * ones(segments)
-@variables rel_vel(t)[1:3, 1:segments]  = zeros(3, segments)
-@variables spring_vel(t)[1:segments] = zeros(segments)
-@variables c_spring(t)[1:segments] = c_spring0 * ones(segments)
-@variables spring_force(t)[1:3, 1:segments] = zeros(3, segments)
-@variables total_force(t)[1:3, 1:segments+1] = zeros(3, segments+1)
+@variables acc(t)[1:3, 1:segments+1]
+@variables segment(t)[1:3, 1:segments]
+@variables unit_vector(t)[1:3, 1:segments]
+@variables norm1(t)[1:segments]
+@variables rel_vel(t)[1:3, 1:segments]
+@variables spring_vel(t)[1:segments]
+@variables c_spring(t)[1:segments]
+@variables spring_force(t)[1:3, 1:segments]
+@variables total_force(t)[1:3, 1:segments+1]
 
 # basic differential equations
 eqs1 = vcat(D.(pos) .~ vel,

src/Tether_06.jl

@@ -1,7 +1,7 @@
 # Tutorial example simulating a 3D mass-spring system with a nonlinear spring (1% spring forces
 # for l < l_0), n tether segments and reel-in and reel-out. The compression stiffness is hardcoded
 # to 1% of the spring stiffness.
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, Timers
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, Timers, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 using ControlPlots
 
@@ -21,37 +21,29 @@ mass_per_meter::Float64 = RHO_TETHER * π * (D_TETHER/2000.0)^2
 # calculating consistant initial conditions
 POS0 = zeros(3, SEGMENTS+1)
 VEL0 = zeros(3, SEGMENTS+1)
-ACC0 = zeros(3, SEGMENTS+1)
-SEGMENTS0 = zeros(3, SEGMENTS) 
-UNIT_VECTORS0 = zeros(3, SEGMENTS)
 for i in 1:SEGMENTS+1
     l0_ = -(i-1)*L0/SEGMENTS
     POS0[:, i] .= [sin(α0) * l0_, 0, cos(α0) * l0_]
     VEL0[:, i] .= [0.0, 0, 0]
 end
-for i in 1:SEGMENTS
-    ACC0[:, i+1] .= G_EARTH
-    UNIT_VECTORS0[:, i] .= [0, 0, 1.0]
-    SEGMENTS0[:, i] .= POS0[:, i+1] - POS0[:, i]
-end
 
 # defining the model, Z component upwards
 @parameters c_spring0=C_SPRING/(L0/SEGMENTS) l_seg=L0/SEGMENTS
 @variables pos(t)[1:3, 1:SEGMENTS+1]  = POS0
 @variables vel(t)[1:3, 1:SEGMENTS+1]  = VEL0
-@variables acc(t)[1:3, 1:SEGMENTS+1]  = ACC0
-@variables segment(t)[1:3, 1:SEGMENTS]  = SEGMENTS0
-@variables unit_vector(t)[1:3, 1:SEGMENTS]  = UNIT_VECTORS0
-@variables len(t) = L0
-@variables c_spring(t) = c_spring0
-@variables damping(t) = DAMPING  / l_seg
-@variables m_tether_particle(t) = mass_per_meter * l_seg
-@variables norm1(t)[1:SEGMENTS] = l_seg * ones(SEGMENTS)
-@variables rel_vel(t)[1:3, 1:SEGMENTS]  = zeros(3, SEGMENTS)
-@variables spring_vel(t)[1:SEGMENTS] = zeros(SEGMENTS)
-@variables c_spr(t)[1:SEGMENTS] = c_spring0 * ones(SEGMENTS)
-@variables spring_force(t)[1:3, 1:SEGMENTS] = zeros(3, SEGMENTS)
-@variables total_force(t)[1:3, 1:SEGMENTS+1] = zeros(3, SEGMENTS+1)
+@variables acc(t)[1:3, 1:SEGMENTS+1]
+@variables segment(t)[1:3, 1:SEGMENTS]
+@variables unit_vector(t)[1:3, 1:SEGMENTS]
+@variables len(t)
+@variables c_spring(t)
+@variables damping(t)
+@variables m_tether_particle(t)
+@variables norm1(t)[1:SEGMENTS]
+@variables rel_vel(t)[1:3, 1:SEGMENTS]
+@variables spring_vel(t)[1:SEGMENTS]
+@variables c_spr(t)[1:SEGMENTS]
+@variables spring_force(t)[1:3, 1:SEGMENTS]
+@variables total_force(t)[1:3, 1:SEGMENTS+1]
 
 # basic differential equations
 eqs1 = vcat(D.(pos) .~ vel,

src/Tether_06b.jl

@@ -1,7 +1,7 @@
 # Tutorial example simulating a 3D mass-spring system with a nonlinear spring (no spring forces
 # for l < l_0), n tether segments and reel-in and reel-out. Can create a video of the simulation.
 # For creating the video, set save=true in the Settings struct.
-using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, Timers, Parameters
+using ModelingToolkit, OrdinaryDiffEq, LinearAlgebra, Timers, Parameters, ControlPlots
 using ModelingToolkit: t_nounits as t, D_nounits as D
 using ControlPlots
 
@@ -22,41 +22,33 @@ end
 function calc_initial_state(se)
     POS0 = zeros(3, se.segments+1)
     VEL0 = zeros(3, se.segments+1)
-    ACC0 = zeros(3, se.segments+1)
-    SEGMENTS0 = zeros(3, se.segments) 
-    UNIT_VECTORS0 = zeros(3, se.segments)
     for i in 1:se.segments+1
         l0 = -(i-1)*se.l0/se.segments
         POS0[:, i] .= [sin(se.α0) * l0, 0, cos(se.α0) * l0]
         VEL0[:, i] .= [0, 0, 0]
     end
-    for i in 1:se.segments
-        ACC0[:, i+1] .= se.g_earth
-        UNIT_VECTORS0[:, i] .= [0, 0, 1.0]
-        SEGMENTS0[:, i] .= POS0[:, i+1] - POS0[:, i]
-    end
-    POS0, VEL0, ACC0, SEGMENTS0, UNIT_VECTORS0
+    POS0, VEL0
 end
 
 function model(se)
-    POS0, VEL0, ACC0, SEGMENTS0, UNIT_VECTORS0 = calc_initial_state(se)
+    POS0, VEL0 = calc_initial_state(se)
     mass_per_meter = se.rho_tether * π * (se.d_tether/2000.0)^2
     @parameters c_spring0=se.c_spring/(se.l0/se.segments) l_seg=se.l0/se.segments
     @variables pos(t)[1:3, 1:se.segments+1]  = POS0
     @variables vel(t)[1:3, 1:se.segments+1]  = VEL0
-    @variables acc(t)[1:3, 1:se.segments+1]  = ACC0
-    @variables segment(t)[1:3, 1:se.segments]  = SEGMENTS0
-    @variables unit_vector(t)[1:3, 1:se.segments]  = UNIT_VECTORS0
-    @variables len(t) = se.l0
-    @variables c_spring(t) = c_spring0
-    @variables damping(t) = se.damping  / l_seg
-    @variables m_tether_particle(t) = mass_per_meter * l_seg
-    @variables norm1(t)[1:se.segments] = l_seg * ones(se.segments)
-    @variables rel_vel(t)[1:3, 1:se.segments]  = zeros(3, se.segments)
-    @variables spring_vel(t)[1:se.segments] = zeros(se.segments)
-    @variables c_spr(t)[1:se.segments] = c_spring0 * ones(se.segments)
-    @variables spring_force(t)[1:3, 1:se.segments] = zeros(3, se.segments)
-    @variables total_force(t)[1:3, 1:se.segments+1] = zeros(3, se.segments+1)
+    @variables acc(t)[1:3, 1:se.segments+1]
+    @variables segment(t)[1:3, 1:se.segments]
+    @variables unit_vector(t)[1:3, 1:se.segments]
+    @variables len(t)
+    @variables c_spring(t)
+    @variables damping(t)
+    @variables m_tether_particle(t)
+    @variables norm1(t)[1:se.segments]
+    @variables rel_vel(t)[1:3, 1:se.segments]
+    @variables spring_vel(t)[1:se.segments]
+    @variables c_spr(t)[1:se.segments]
+    @variables spring_force(t)[1:3, 1:se.segments]
+    @variables total_force(t)[1:3, 1:se.segments+1]
 
     # basic differential equations
     eqs1 = vcat(D.(pos) .~ vel,