Update new fields in SysState

Project.toml

@@ -42,7 +42,7 @@ DocStringExtensions = "0.8, 0.9"
 Documenter = "1.0"
 Interpolations = "0.15"
 KitePodModels = "0.3.6"
-KiteUtils = "0.8.4"
+KiteUtils = "0.9.0"
 LaTeXStrings = "1.4.0"
 ModelingToolkit = "~9.42.0"
 NLsolve = "4.5"
@@ -51,7 +51,7 @@ OrdinaryDiffEqCore = "1.4.0"
 OrdinaryDiffEqSDIRK = "1.1.0"
 PackageCompiler = "2.0"
 Parameters = "0.12"
-Pkg = "1.10"
+Pkg = "1.10, 1.11"
 PrecompileTools = "1.2"
 Reexport = "1.1, 1.2"
 Rotations = "1.7"
@@ -61,7 +61,7 @@ SteadyStateDiffEq = "2.3.0"
 Sundials = "4.24"
 SymbolicIndexingInterface = "0.3.27"
 Timers = "0.1.5"
-WinchModels = "0.3.4"
+WinchModels = "0.3.5"
 Xfoil = "0.5.0"
 julia = "1.10, 1.11"
 

examples/test_init_1p.jl

@@ -62,6 +62,13 @@ println("azimuth: ", round(rad2deg(azimuth), digits = 2), "°")
 # println("alpha2, alpha3, alpha4: ", kps3.alpha_2, " ", kps3.alpha_3, " ", kps3.alpha_4)
 println("heading: ", round(rad2deg(calc_heading(kps3)), digits = 2), "°")
 
+ss = SysState(kps3)
+println("AoA: ", rad2deg(ss.AoA), "°")
+println("CL2: ", ss.CL2, " CD2: ", ss.CD2)
+println("v_wind_gnd: ", ss.v_wind_gnd)
+println("v_wind_200m: ", ss.v_wind_200m)
+println("v_wind_kite: ", ss.v_wind_kite)
+
 # output on main branch
 # Lift: 1047.1795339611076, Drag: 281.39765463928745, elev: 72.77014, Iterations: 1552
 # roll: -0.0 pitch: 5.755691707832273 yaw: 90.00000250447816

examples/test_init_4p.jl

@@ -64,6 +64,13 @@ pos_C, pos_D = kps4.pos[kps4.set.segments+4], kps4.pos[kps4.set.segments+5]
 println("alpha2, alpha3, alpha4: ", kps4.alpha_2, " ", kps4.alpha_3, " ", kps4.alpha_4)
 println("heading: ", round(rad2deg(calc_heading(kps4)), digits = 2), "°")
 
+ss = SysState(kps4)
+println("AoA: ", rad2deg(ss.AoA), "°")
+println("CL2: ", ss.CL2, " CD2: ", ss.CD2)
+println("v_wind_gnd: ", ss.v_wind_gnd)
+println("v_wind_200m: ", ss.v_wind_200m)
+println("v_wind_kite: ", ss.v_wind_kite)
+
 # output on main branch
 # Lift: 1047.1795339611076, Drag: 281.39765463928745, elev: 72.77014, Iterations: 1552
 # roll: -0.0 pitch: 5.755691707832273 yaw: 90.00000250447816

src/KPS3.jl

@@ -70,6 +70,8 @@ $(TYPEDFIELDS)
     v_apparent::T =       zeros(S, 3)
     "vector, perpendicular to v_apparent; output of calc_drag"
     v_app_perp::T =       zeros(S, 3)
+    "angle of attack of the kite; output of set_cl_cd!"
+    alpha_2::S =          0.0
     "drag force of kite and bridle; output of calc_aero_forces"
     drag_force::T =       zeros(S, 3)
     "lift force of the kite; output of calc_aero_forces"
@@ -345,6 +347,14 @@ function calc_set_cl_cd!(s::KPS3, vec_c, v_app)
     alpha = calc_alpha(v_app, s.vec_z) - s.alpha_depower
     set_cl_cd!(s, alpha)
 end
+"""
+    cl_cd(s::KPS3)
+
+Calculate the lift and drag coefficients of the kite, based on the current angles of attack.
+"""
+function cl_cd(s::KPS3)
+    CL2, CD2 = s.calc_cl(s.alpha_2), s.calc_cd(s.alpha_2)
+end
 
 """
     residual!(res, yd, y::MVector{S, SimFloat}, s::KPS3, time) where S

src/KiteModels.jl

@@ -120,8 +120,9 @@ include("KPS3.jl") # include code, specific for the one point kite model
 include("init.jl") # functions to calculate the inital state vector, the inital masses and initial springs
 
 # Calculate the lift and drag coefficient as a function of the angle of attack alpha.
-function set_cl_cd!(s::AKM, alpha)   
-    angle =  alpha * 180.0 / π
+function set_cl_cd!(s::AKM, alpha)
+    angle =  rad2deg(alpha)
+    s.alpha_2 = angle
     if angle > 180.0
         angle -= 360.0
     end
@@ -211,7 +212,7 @@ function set_v_wind_ground!(s::AKM, height, v_wind_gnd=s.set.v_wind; upwind_dir=
     wind_dir = -upwind_dir - pi/2
     s.v_wind .= v_wind_gnd * calc_wind_factor(s.am, height) .* [cos(wind_dir), sin(wind_dir), 0]
     s.v_wind_gnd .= [v_wind_gnd * cos(wind_dir), v_wind_gnd * sin(wind_dir), 0.0]
-    s.v_wind_tether .= s.v_wind_gnd * calc_wind_factor(s.am, height / 2.0) # .* [cos(wind_dir), sin(wind_dir), 0]
+    s.v_wind_tether .= s.v_wind_gnd * calc_wind_factor(s.am, height / 2.0)
     s.rho = calc_rho(s.am, height)
     nothing
 end
@@ -443,6 +444,14 @@ function update_sys_state!(ss::SysState, s::AKM, zoom=1.0)
     ss.depower = s.depower
     ss.steering = s.steering/s.set.cs_4p
     ss.vel_kite .= s.vel_kite
+    ss.t_sim = 0.0
+    ss.AoA = deg2rad(s.alpha_2)
+    cl, cd = cl_cd(s)
+    ss.CL2 = cl
+    ss.CD2 = cd
+    ss.v_wind_gnd  .= s.v_wind_gnd
+    ss.v_wind_200m .= s.v_wind_gnd * calc_wind_factor(s.am, 200.0)
+    ss.v_wind_kite .= s.v_wind
     nothing
 end
 
@@ -455,28 +464,9 @@ system state in a viewer. Optionally the position arrays can be zoomed
 according to the requirements of the viewer.
 """
 function SysState(s::AKM, zoom=1.0)
-    pos = s.pos
-    P = length(pos)
-    X = zeros(MVector{P, MyFloat})
-    Y = zeros(MVector{P, MyFloat})
-    Z = zeros(MVector{P, MyFloat})
-    for i in 1:P
-        X[i] = pos[i][1] * zoom
-        Y[i] = pos[i][2] * zoom
-        Z[i] = pos[i][3] * zoom
-    end
-
-    orient = calc_orient_quat(s)
-
-    elevation = calc_elevation(s)
-    azimuth = calc_azimuth(s)
-    force = winch_force(s)
-    heading = calc_heading(s)
-    course = calc_course(s)
-    v_app_norm = norm(s.v_apparent)
-    t_sim = 0
-    KiteUtils.SysState{P}(s.t_0, t_sim, 0, 0, orient, elevation, azimuth, s.l_tether, s.v_reel_out, force, s.depower, s.steering/s.set.cs_4p, 
-                          heading, course, v_app_norm, s.vel_kite, X, Y, Z, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0)
+    ss = SysState{length(s.pos)}()
+    update_sys_state!(ss, s, zoom)
+    ss
 end
 
 function calc_pre_tension(s::AKM)