From 98f475ae6eccec59de715b67bdd2fe5adde457f0 Mon Sep 17 00:00:00 2001
From: Mauco03 <marco.gaibotti@skywarder.eu>
Date: Thu, 16 May 2024 09:59:21 +0200
Subject: [PATCH] [refactoring-ode][functions] Added uncertanty to ode, moved
old files to OLD
---
classes/DataWrapper.m | 6 +-
functions/odeFunctions/ballistic.m | 1 +
functions/odeFunctions/descentParachute.m | 1 +
functions/odeFunctions/descentParafoil.m | 1 +
.../simulations/{ => OLD}/quickApogeeOnly.m | 0
.../simulations/{ => OLD}/stochAscentLoop.m | 0
.../{ => OLD}/stochDescentBalLoop.m | 0
.../simulations/{ => OLD}/stochDescentLoop.m | 0
.../simulations/{ => OLD}/stochPostProcess.m | 0
.../{ => OLD}/stochStabilityLoop.m | 0
.../simulations/{ => OLD}/stochWindData.m | 0
functions/simulations/stdStability.m | 128 ++++++++++++++++++
settings/odeConfig.m | 3 +-
13 files changed, 135 insertions(+), 5 deletions(-)
rename functions/simulations/{ => OLD}/quickApogeeOnly.m (100%)
rename functions/simulations/{ => OLD}/stochAscentLoop.m (100%)
rename functions/simulations/{ => OLD}/stochDescentBalLoop.m (100%)
rename functions/simulations/{ => OLD}/stochDescentLoop.m (100%)
rename functions/simulations/{ => OLD}/stochPostProcess.m (100%)
rename functions/simulations/{ => OLD}/stochStabilityLoop.m (100%)
rename functions/simulations/{ => OLD}/stochWindData.m (100%)
create mode 100644 functions/simulations/stdStability.m
diff --git a/classes/DataWrapper.m b/classes/DataWrapper.m
index 62981eb..3fa54c6 100644
--- a/classes/DataWrapper.m
+++ b/classes/DataWrapper.m
@@ -15,12 +15,12 @@ classdef DataWrapper < handle
% improve code maintainability
properties(SetAccess = private)
- data struct % Cached data
- dummy struct % Empty struct for field reference
+ data struct % Cached data
+ dummy struct % Empty struct for field reference
+ counter int32 % Current size of data stored
end
properties(Access = private)
- counter int32 % Current size of data stored
shift int32 % Number of fields to skip in cache, useful to preallocate initial values / starting conditions
chunkSize int32 % Size of preallocated memory
dataLen int32 % Current preallocated size of cache
diff --git a/functions/odeFunctions/ballistic.m b/functions/odeFunctions/ballistic.m
index 5992de8..6c4a1f8 100644
--- a/functions/odeFunctions/ballistic.m
+++ b/functions/odeFunctions/ballistic.m
@@ -350,6 +350,7 @@ if nargout == 2 || ~isempty(wrapper)
parout.coeff.XCPlat = XCPlat;
parout.coeff.XCPtot = XCPtot;
+ parout.uncertanty = [];
if ~isempty(wrapper)
wrapper.setCache(parout);
end
diff --git a/functions/odeFunctions/descentParachute.m b/functions/odeFunctions/descentParachute.m
index 6ec9537..556fbbd 100644
--- a/functions/odeFunctions/descentParachute.m
+++ b/functions/odeFunctions/descentParachute.m
@@ -150,6 +150,7 @@ if nargout == 2 || ~isempty(wrapper)
parout.accelerations.body = [du; dv; dw];
parout.coeff = [];
+ parout.uncertanty = [];
if ~isempty(wrapper)
wrapper.setCache(parout);
diff --git a/functions/odeFunctions/descentParafoil.m b/functions/odeFunctions/descentParafoil.m
index 66ceca0..279a785 100644
--- a/functions/odeFunctions/descentParafoil.m
+++ b/functions/odeFunctions/descentParafoil.m
@@ -227,6 +227,7 @@ if nargout == 2 || ~isempty(wrapper)
parout.accelerations.body = bodyAcc;
parout.coeff = [];
+ parout.uncertanty = [];
if ~isempty(wrapper)
wrapper.setCache(parout);
diff --git a/functions/simulations/quickApogeeOnly.m b/functions/simulations/OLD/quickApogeeOnly.m
similarity index 100%
rename from functions/simulations/quickApogeeOnly.m
rename to functions/simulations/OLD/quickApogeeOnly.m
diff --git a/functions/simulations/stochAscentLoop.m b/functions/simulations/OLD/stochAscentLoop.m
similarity index 100%
rename from functions/simulations/stochAscentLoop.m
rename to functions/simulations/OLD/stochAscentLoop.m
diff --git a/functions/simulations/stochDescentBalLoop.m b/functions/simulations/OLD/stochDescentBalLoop.m
similarity index 100%
rename from functions/simulations/stochDescentBalLoop.m
rename to functions/simulations/OLD/stochDescentBalLoop.m
diff --git a/functions/simulations/stochDescentLoop.m b/functions/simulations/OLD/stochDescentLoop.m
similarity index 100%
rename from functions/simulations/stochDescentLoop.m
rename to functions/simulations/OLD/stochDescentLoop.m
diff --git a/functions/simulations/stochPostProcess.m b/functions/simulations/OLD/stochPostProcess.m
similarity index 100%
rename from functions/simulations/stochPostProcess.m
rename to functions/simulations/OLD/stochPostProcess.m
diff --git a/functions/simulations/stochStabilityLoop.m b/functions/simulations/OLD/stochStabilityLoop.m
similarity index 100%
rename from functions/simulations/stochStabilityLoop.m
rename to functions/simulations/OLD/stochStabilityLoop.m
diff --git a/functions/simulations/stochWindData.m b/functions/simulations/OLD/stochWindData.m
similarity index 100%
rename from functions/simulations/stochWindData.m
rename to functions/simulations/OLD/stochWindData.m
diff --git a/functions/simulations/stdStability.m b/functions/simulations/stdStability.m
new file mode 100644
index 0000000..7ff4dfb
--- /dev/null
+++ b/functions/simulations/stdStability.m
@@ -0,0 +1,128 @@
+function stability = stdStability(Y0, rocket, env, wind, settings)
+% stochStabilityLoop - This function runs a simulation to compute the rocket
+% stability inside the parfor loop
+%
+% INPUTS:
+% - settings, struct: stores data of the rocket and of the simulation.
+% - Y0, double [13, 1]: initial state vector of the integration
+% [x, y, z, u, v, w, p, q, r, q0, q1, q2, q3]
+%
+% OUTPUTS:
+% - data_stability,cell [1, k]: cell containg the struct with the rocket
+% stability data at the launchpad exit
+%
+% CALLED FUNCTIONS:
+%
+% REVISIONS:
+% - 0 14/06/2023, Release, Riccardo Cadamuro
+%
+% Copyright © 2023, Skyward Experimental Rocketry, AFD department
+% All rights reserved
+%
+% SPDX-License-Identifier: GPL-3.0-or-later
+
+if rocket.airbrakes.enabled && length(rocket.airbrakes.extension) > 1 && ...
+ length(rocket.airbrakes.deltaTime) < length(rocket.airbrakes.extension) - 1
+ error('In airbrakes smooth opening simulations, AB configuration usage time vector must be at least of length length(airbrakes.extension)-1, check rocketConfig.m')
+end
+
+if any(rocket.airbrakes.extension > 3) || any(rocket.airbrakes.extension < 1)
+ error('At least a value in airbrakes.extension is set incorrectly in rocketConfig, they must be between 1 and 3')
+end
+
+%% ASCENT
+t0 = 0;
+tf = 10;
+solution = ode113(@ballistic, [t0, tf], Y0, options, ...
+ rocket, env, wind);
+
+if nargout == 2
+ t = solution.x;
+ Y = solution.y;
+
+ [~, stability] = ballistic(t(end), Y(:, end), rocket, env, wind);
+
+ ascent = getOdeFcn(solution);
+ ascent.state.Y = Y;
+ ascent.state.time = t;
+end
+
+%%
+
+%% comupute launchpad run
+[Tpad, Ypad] = ode113(@ascent, [0, 10], Y0, settings.ode.optionspad, settings);
+
+z = -Ypad(end, 3);
+[~, a, ~, ~] = atmosisa(settings.z0 + z);
+
+%% retrive the correct wind data
+if settings.wind.model
+ Hour = settings.stoch.Hour;
+ Day = settings.stoch.Day;
+ t = Tpad(end);
+
+ if settings.stoch.N > 1
+ [uw, vw, ww] = windMatlabGenerator(settings, z, t, Hour, Day);
+ else
+ [uw, vw, ww] = windMatlabGenerator(settings, z, t);
+ end
+
+elseif settings.wind.input
+ [uw, vw, ww] = windInputGenerator(settings, z, uncert);
+elseif settings.wind.variable
+
+ [uw, vw, ww] = windVariableGenerator(z, settings.stoch);
+else
+ uw = settings.stoch.uw; vw = settings.stoch.vw; ww = settings.stoch.ww;
+end
+inertialWind = [uw, vw, ww];
+
+%% compute the exit conditions in the aerodynamic frame
+vExit = norm(Ypad(end, 4));
+machExit = vExit/a;
+
+bodyWind = quatrotate(Y0(10:13)', inertialWind);
+bodyVelocity = [Ypad(end, 4), 0, 0];
+Vr = bodyVelocity - bodyWind;
+ur = Vr(1); vr = Vr(2); wr = Vr(3);
+
+if ur > 1e-9
+ alpha = atan(wr/ur)*180/pi;
+ beta = atan(vr/ur)*180/pi;
+else
+ alpha = 0;
+ beta = 0;
+end
+
+%% create dissileMatcom input data struct
+dissileVars.Alpha = sort([-alpha, alpha, 0, -1, 1]);
+dissileVars.Beta = beta;
+dissileVars.Alt = settings.z0 + z;
+dissileVars.Mach = machExit;
+dissileVars.xcg = interpLinear(settings.motor.expTime, settings.xcg, Tpad(end));
+dissileVars.hprot = [];
+input = createDissileInput(settings, dissileVars);
+
+input.flagXCP = true;
+coeff = dissileMatcom(input);
+
+vars.useAlphaTot = false;
+vars.checkOpenRocket = false;
+vars.modHaack = true;
+
+if dissileVars.Alpha(1) == alpha
+ iAlpha = 1;
+else
+ iAlpha = 5;
+end
+
+%% compute the stability margin and coresponding quality indices
+[XCPvect, qualityIndex, stabilityCoeffs] = getXCP(coeff, iAlpha, 1, input, vars);
+stability.time = Tpad;
+stability.state = Ypad;
+stability.wind = [uw, vw, ww];
+stability.XCPvect = XCPvect;
+stability.qualityIndex = qualityIndex;
+stability.stabilityCoeffs = stabilityCoeffs;
+
+end
\ No newline at end of file
diff --git a/settings/odeConfig.m b/settings/odeConfig.m
index 031e351..8e1fb78 100644
--- a/settings/odeConfig.m
+++ b/settings/odeConfig.m
@@ -3,10 +3,9 @@
% ODE settings
ode.finalTime = 2000; % [s] Final integration time
ode.optAscent = odeset('Events', @eventApogee, 'InitialStep', 1); % ODE options for ascend
-ode.optAscentDelayPara= odeset('InitialStep', 1); % ODE options for due to the opening delay of the parachute
ode.optParachute = odeset('Events', @eventParaCut, 'InitialStep', 1); % ODE options for the parachutes
ode.optDescent = odeset('Events', @eventLanding,...
-'RelTol', 1e-3, 'AbsTol', 1e-3, 'OutputFcn', @storeData); % ODE options for ballistic descent
+'RelTol', 1e-3, 'AbsTol', 1e-3); % ODE options for ballistic descent
ode.optLaunchpad = odeset('Events', @eventPad); % ODE options for the launchpad phase
% Settings for descent 6dof simulation
--
GitLab