From 156670de6833b5393ad1ae0791e2140f69445012 Mon Sep 17 00:00:00 2001 From: Uwe Fechner Date: Fri, 27 Sep 2024 11:36:28 +0200 Subject: [PATCH] refactoring --- src/Tether_08.jl | 57 ++++++++++++++++++++++++------------------------ 1 file changed, 28 insertions(+), 29 deletions(-) diff --git a/src/Tether_08.jl b/src/Tether_08.jl index bd97221..03a43c1 100644 --- a/src/Tether_08.jl +++ b/src/Tether_08.jl @@ -82,25 +82,24 @@ function model(se, p1, p2, fix_p1, fix_p2, POS0, VEL0) 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 @parameters rel_compression_stiffness = se.rel_compression_stiffness - @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] - @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 v_apparent(t)[1:3, 1:se.segments] - @variables v_app_perp(t)[1:3, 1:se.segments] - @variables norm_v_app(t)[1:se.segments] - @variables half_drag_force(t)[1:3, 1:se.segments] - @variables total_force(t)[1:3, 1:se.segments+1] + @variables begin + pos(t)[1:3, 1:se.segments+1] = POS0 + vel(t)[1:3, 1:se.segments+1] = VEL0 + acc(t)[1:3, 1:se.segments+1] + segment(t)[1:3, 1:se.segments] + unit_vector(t)[1:3, 1:se.segments] + l_spring(t), c_spring(t), damping(t), m_tether_particle(t) + len(t)[1:se.segments] + rel_vel(t)[1:3, 1:se.segments] + spring_vel(t)[1:se.segments] + c_spr(t)[1:se.segments] + spring_force(t)[1:3, 1:se.segments] + v_apparent(t)[1:3, 1:se.segments] + v_app_perp(t)[1:3, 1:se.segments] + norm_v_app(t)[1:se.segments] + half_drag_force(t)[1:3, 1:se.segments] + total_force(t)[1:3, 1:se.segments+1] + end # basic differential equations eqs1 = vcat(D.(pos) .~ vel, @@ -109,18 +108,18 @@ function model(se, p1, p2, fix_p1, fix_p2, POS0, VEL0) # loop over all segments to calculate the spring and drag forces for i in 1:se.segments eqs = [segment[:, i] ~ pos[:, i+1] - pos[:, i], - norm1[i] ~ norm(segment[:, i]), - unit_vector[:, i] ~ -segment[:, i]/norm1[i], + len[i] ~ norm(segment[:, i]), + unit_vector[:, i] ~ -segment[:, i]/len[i], rel_vel[:, i] ~ vel[:, i+1] - vel[:, i], spring_vel[i] ~ -unit_vector[:, i] ⋅ rel_vel[:, i], - c_spr[i] ~ c_spring/(1+rel_compression_stiffness) - * (rel_compression_stiffness+(norm1[i] > len/se.segments)), - spring_force[:, i] ~ (c_spr[i] * (norm1[i] - (len/se.segments)) + c_spr[i] ~ c_spring / (1+rel_compression_stiffness) + * (rel_compression_stiffness+(len[i] > l_spring)), + spring_force[:, i] ~ (c_spr[i] * (len[i] - l_spring) + damping * spring_vel[i]) * unit_vector[:, i], v_apparent[:, i] ~ se.v_wind_tether .- (vel[:, i] + vel[:, i+1])/2, v_app_perp[:, i] ~ v_apparent[:, i] - (v_apparent[:, i] ⋅ unit_vector[:, i]) .* unit_vector[:, i], norm_v_app[i] ~ norm(v_app_perp[:, i]), - half_drag_force[:, i] ~ 0.25 * se.rho * se.cd_tether * norm_v_app[i] * (norm1[i]*se.d_tether/1000.0) + half_drag_force[:, i] ~ 0.25 * se.rho * se.cd_tether * norm_v_app[i] * (len[i]*se.d_tether/1000.0) * v_app_perp[:, i]] eqs2 = vcat(eqs2, reduce(vcat, eqs)) end @@ -149,10 +148,10 @@ function model(se, p1, p2, fix_p1, fix_p2, POS0, VEL0) eqs2 = vcat(eqs2, reduce(vcat, eqs)) end # scalar equations - eqs = [len ~ se.l0 + se.v_ro*t, - c_spring ~ se.c_spring / (len/se.segments), - m_tether_particle ~ mass_per_meter * (len/se.segments), - damping ~ se.damping / (len/se.segments)] + eqs = [l_spring ~ (se.l0 + se.v_ro*t) / se.segments, + c_spring ~ se.c_spring / l_spring, + m_tether_particle ~ mass_per_meter * l_spring, + damping ~ se.damping / l_spring] eqs2 = vcat(eqs2, reduce(vcat, eqs)) @named sys = ODESystem(reduce(vcat, Symbolics.scalarize.(eqs2)), t)