From 83609755b45a78afc02da5997e023f8c751eb781 Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Sun, 13 Nov 2022 15:01:22 -0300 Subject: [PATCH 01/13] MAINT: remove unnecessary private attributes --- rocketpy/Flight.py | 81 +--------------------------------------------- 1 file changed, 1 insertion(+), 80 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index 8be8a2342..79b13ffb3 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -628,10 +628,6 @@ def __init__( self.impactState = np.array([0]) self.parachuteEvents = [] self.postProcessed = False - self._drift = Function(0) - self._bearing = Function(0) - self._latitude = Function(0) - self._longitude = Function(0) # Initialize solver monitors self.functionEvaluations = [] self.functionEvaluationsPerTimeStep = [] @@ -1112,83 +1108,8 @@ def __init__( def __init_post_process_variables(self): """Initialize post-process variables.""" - # Initialize all variables created during Flight.postProcess() + # Initialize all variables calculated after initialization. # Important to do so that MATLAB® can access them - self._windVelocityX = Function(0) - self._windVelocityY = Function(0) - self._density = Function(0) - self._pressure = Function(0) - self._dynamicViscosity = Function(0) - self._speedOfSound = Function(0) - self._ax = Function(0) - self._ay = Function(0) - self._az = Function(0) - self._alpha1 = Function(0) - self._alpha2 = Function(0) - self._alpha3 = Function(0) - self._speed = Function(0) - self._maxSpeed = 0 - self._maxSpeedTime = 0 - self._horizontalSpeed = Function(0) - self._Acceleration = Function(0) - self._maxAcceleration = 0 - self._maxAccelerationTime = 0 - self._pathAngle = Function(0) - self._attitudeVectorX = Function(0) - self._attitudeVectorY = Function(0) - self._attitudeVectorZ = Function(0) - self._attitudeAngle = Function(0) - self._lateralAttitudeAngle = Function(0) - self._phi = Function(0) - self._theta = Function(0) - self._psi = Function(0) - self._R1 = Function(0) - self._R2 = Function(0) - self._R3 = Function(0) - self._M1 = Function(0) - self._M2 = Function(0) - self._M3 = Function(0) - self._aerodynamicLift = Function(0) - self._aerodynamicDrag = Function(0) - self._aerodynamicBendingMoment = Function(0) - self._aerodynamicSpinMoment = Function(0) - self._railButton1NormalForce = Function(0) - self._maxRailButton1NormalForce = 0 - self._railButton1ShearForce = Function(0) - self._maxRailButton1ShearForce = 0 - self._railButton2NormalForce = Function(0) - self._maxRailButton2NormalForce = 0 - self._railButton2ShearForce = Function(0) - self._maxRailButton2ShearForce = 0 - self._rotationalEnergy = Function(0) - self._translationalEnergy = Function(0) - self._kineticEnergy = Function(0) - self._potentialEnergy = Function(0) - self._totalEnergy = Function(0) - self._thrustPower = Function(0) - self._dragPower = Function(0) - self._attitudeFrequencyResponse = Function(0) - self._omega1FrequencyResponse = Function(0) - self._omega2FrequencyResponse = Function(0) - self._omega3FrequencyResponse = Function(0) - self._streamVelocityX = Function(0) - self._streamVelocityY = Function(0) - self._streamVelocityZ = Function(0) - self._freestreamSpeed = Function(0) - self._apogeeFreestreamSpeed = 0 - self._MachNumber = Function(0) - self._maxMachNumber = 0 - self._maxMachNumberTime = 0 - self._ReynoldsNumber = Function(0) - self._maxReynoldsNumber = 0 - self._maxReynoldsNumberTime = 0 - self._dynamicPressure = Function(0) - self._maxDynamicPressure = 0 - self._maxDynamicPressureTime = 0 - self._totalPressure = Function(0) - self._maxTotalPressure = 0 - self._maxTotalPressureTime = 0 - self._angleOfAttack = Function(0) self.flutterMachNumber = Function(0) self.difference = Function(0) self.safetyFactor = Function(0) From 9ad1423427f874429e3a8b60e60c33a636525e7e Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Sun, 13 Nov 2022 15:56:24 -0300 Subject: [PATCH 02/13] ENH: add funcify_method decorator, missing docs and create helper attributes outOfRailTimeIndex, solutionArray and time --- rocketpy/Flight.py | 1503 ++++++++++++-------------------------------- 1 file changed, 419 insertions(+), 1084 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index 79b13ffb3..e32a58698 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -9,6 +9,7 @@ import math import time import warnings +from copy import deepcopy from functools import cached_property import matplotlib.pyplot as plt @@ -16,7 +17,7 @@ import simplekml from scipy import integrate -from .Function import Function +from .Function import Function, funcify_method class Flight: @@ -615,6 +616,7 @@ def __init__( self.__init_post_process_variables() # Initialize solution monitors self.outOfRailTime = 0 + self.outOfRailTimeIndex = 0 self.outOfRailState = np.array([0]) self.outOfRailVelocity = 0 self.apogeeState = np.array([0]) @@ -671,6 +673,7 @@ def __init__( # TODO: Check if rocket is actually out of rail. Otherwise, start at rail self.outOfRailState = self.initialSolution[1:] self.outOfRailTime = self.initialSolution[0] + self.outOfRailTimeIndex = 0 self.initialDerivative = self.uDot self.tInitial = self.initialSolution[0] @@ -893,12 +896,17 @@ def __init__( raise ValueError( "Multiple roots found when solving for rail exit time." ) + elif len(valid_t_root) == 0: + raise ValueError( + "No valid roots found when solving for rail exit time." + ) # Determine final state when upper button is going out of rail self.t = valid_t_root[0] + self.solution[-2][0] interpolator = phase.solver.dense_output() self.ySol = interpolator(self.t) self.solution[-1] = [self.t, *self.ySol] self.outOfRailTime = self.t + self.outOfRailTimeIndex = len(self.solution) - 1 self.outOfRailState = self.ySol self.outOfRailVelocity = ( self.ySol[3] ** 2 + self.ySol[4] ** 2 + self.ySol[5] ** 2 @@ -1526,709 +1534,267 @@ def uDotParachute(self, t, u, postProcessing=False): return [vx, vy, vz, ax, ay, az, 0, 0, 0, 0, 0, 0, 0] - # Process first type of outputs - state vector - # Transform solution array into Functions - @cached_property - def x(self, interpolation="spline", extrapolation="natural"): - """Rocket x position as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - x: Function - Rocket x position as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 1]], - "Time (s)", - "X (m)", - interpolation, - extrapolation, - ) - - @cached_property - def y(self, interpolation="spline", extrapolation="natural"): - """ocket y position as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - y: Function - Rocket y position as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 2]], - "Time (s)", - "Y (m)", - interpolation, - extrapolation, - ) - @cached_property - def z(self, interpolation="spline", extrapolation="natural"): - """Rocket z position as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - z: Function - Rocket z position as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 3]], - "Time (s)", - "Z (m)", - interpolation, - extrapolation, - ) + def solutionArray(self): + """Returns solution array of the rocket flight.""" + return np.array(self.solution) @cached_property - def vx(self, interpolation="spline", extrapolation="natural"): - """Rocket x velocity as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - vx: Function - Rocket x velocity as a function of time. - """ - - return Function( - np.array(self.solution)[:, [0, 4]], - "Time (s)", - "Vx (m/s)", - interpolation, - extrapolation, - ) - - @cached_property - def vy(self, interpolation="spline", extrapolation="natural"): - """Rocket y velocity as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. + def time(self): + """Returns time array from solution.""" + return self.solutionArray[:, 0] - Returns - ------- - vy: Function - Rocket y velocity as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 5]], - "Time (s)", - "Vy (m/s)", - interpolation, - extrapolation, - ) - - @cached_property - def vz(self, interpolation="spline", extrapolation="natural"): - """Rocket z velocity as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - vz: Function - Rocket z velocity as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 6]], - "Time (s)", - "Vz (m/s)", - interpolation, - extrapolation, - ) - - @cached_property - def e0(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 7]], - "Time (s)", - "e0", - interpolation, - extrapolation, - ) - - @cached_property - def e1(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 8]], - "Time (s)", - "e1", - interpolation, - extrapolation, - ) - - @cached_property - def e2(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 9]], - "Time (s)", - "e2", - interpolation, - extrapolation, - ) - - @cached_property - def e3(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 10]], - "Time (s)", - "e3", - interpolation, - extrapolation, - ) - - @cached_property - def w1(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 11]], - "Time (s)", - "ω1 (rad/s)", - interpolation, - extrapolation, - ) - - @cached_property - def w2(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 12]], - "Time (s)", - "ω2 (rad/s)", - interpolation, - extrapolation, - ) - - @cached_property - def w3(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 13]], - "Time (s)", - "ω3 (rad/s)", - interpolation, - extrapolation, - ) + # Process first type of outputs - state vector + # Transform solution array into Functions + @funcify_method("Time (s)", "X (m)", "spline", "constant") + def x(self): + """Rocket x position as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 1]] + + @funcify_method("Time (s)", "Y (m)", "spline", "constant") + def y(self): + """Rocket y position as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 2]] + + @funcify_method("Time (s)", "Z (m)", "spline", "constant") + def z(self): + """Rocket z position as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 3]] + + @funcify_method("Time (s)", "Vx (m/s)", "spline", "constant") + def vx(self): + """Rocket x velocity as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 4]] + + @funcify_method("Time (s)", "Vy (m/s)", "spline", "constant") + def vy(self): + """Rocket y velocity as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 5]] + + @funcify_method("Time (s)", "Vz (m/s)", "spline", "constant") + def vz(self): + """Rocket z velocity as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 6]] + + @funcify_method("Time (s)", "e0", "spline", "constant") + def e0(self): + """Rocket quaternion e0 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 7]] + + @funcify_method("Time (s)", "e1", "spline", "constant") + def e1(self): + """Rocket quaternion e1 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 8]] + + @funcify_method("Time (s)", "e2", "spline", "constant") + def e2(self): + """Rocket quaternion e2 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 9]] + + @funcify_method("Time (s)", "e3", "spline", "constant") + def e3(self): + """Rocket quaternion e3 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 10]] + + @funcify_method("Time (s)", "ω1 (rad/s)", "spline", "constant") + def w1(self): + """Rocket angular velocity ω1 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 11]] + + @funcify_method("Time (s)", "ω2 (rad/s)", "spline", "constant") + def w2(self): + """Rocket angular velocity ω2 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 12]] + + @funcify_method("Time (s)", "ω3 (rad/s)", "spline", "constant") + def w3(self): + """Rocket angular velocity ω3 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 13]] # Process second type of outputs - accelerations components - @cached_property - def ax(self, interpolation="spline", extrapolation="natural"): - ax = self.retrieve_acceleration_arrays[0] - # Convert accelerations to functions - ax = Function(ax, "Time (s)", "Ax (m/s2)", interpolation, extrapolation) - return ax - - @cached_property - def ay(self, interpolation="spline", extrapolation="natural"): - ay = self.retrieve_acceleration_arrays[1] - # Convert accelerations to functions - ay = Function(ay, "Time (s)", "Ay (m/s2)", interpolation, extrapolation) - return ay - - @cached_property - def az(self, interpolation="spline", extrapolation="natural"): - az = self.retrieve_acceleration_arrays[2] - # Convert accelerations to functions - az = Function(az, "Time (s)", "Az (m/s2)", interpolation, extrapolation) - return az - - @cached_property - def alpha1(self, interpolation="spline", extrapolation="natural"): - alpha1 = self.retrieve_acceleration_arrays[3] - # Convert accelerations to functions - alpha1 = Function( - alpha1, "Time (s)", "α1 (rad/s2)", interpolation, extrapolation - ) - return alpha1 - - @cached_property - def alpha2(self, interpolation="spline", extrapolation="natural"): - alpha2 = self.retrieve_acceleration_arrays[4] - # Convert accelerations to functions - alpha2 = Function( - alpha2, "Time (s)", "α2 (rad/s2)", interpolation, extrapolation - ) - return alpha2 - - @cached_property - def alpha3(self, interpolation="spline", extrapolation="natural"): - alpha3 = self.retrieve_acceleration_arrays[5] - # Convert accelerations to functions - alpha3 = Function( - alpha3, "Time (s)", "α3 (rad/s2)", interpolation, extrapolation - ) - return alpha3 + @funcify_method("Time (s)", "Ax (m/s²)", "spline", "constant") + def ax(self): + """Rocket x acceleration as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[0] + + @funcify_method("Time (s)", "Ay (m/s²)", "spline", "constant") + def ay(self): + """Rocket y acceleration as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[1] + + @funcify_method("Time (s)", "Az (m/s²)", "spline", "constant") + def az(self): + """Rocket z acceleration as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[2] + + @funcify_method("Time (s)", "α1 (rad/s²)", "spline", "constant") + def alpha1(self): + """Rocket angular acceleration α1 as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[3] + + @funcify_method("Time (s)", "α2 (rad/s²)", "spline", "constant") + def alpha2(self): + """Rocket angular acceleration α2 as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[4] + + @funcify_method("Time (s)", "α3 (rad/s²)", "spline", "constant") + def alpha3(self): + """Rocket angular acceleration α3 as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[5] # Process third type of outputs - Temporary values - @cached_property - def R1(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "R1 (N)", "spline", "constant") + def R1(self): """Aerodynamic force along the first axis that is perpendicular to the - rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. + rocket's axis of symmetry as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[0] - Returns - ------- - R1: Function - Aero force along the first axis that is perpendicular to the - rocket's axis of symmetry. - """ - R1 = self.retrieve_temporary_values_arrays[0] - R1 = Function(R1, "Time (s)", "R1 (m)", interpolation, extrapolation) - - return R1 - - @cached_property - def R2(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "R2 (N)", "spline", "constant") + def R2(self): """Aerodynamic force along the second axis that is perpendicular to the - rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - R2: Function - Aero force along the second axis that is perpendicular to the - rocket's axis of symmetry. - """ - R2 = self.retrieve_temporary_values_arrays[1] - R2 = Function(R2, "Time (s)", "R2 (m)", interpolation, extrapolation) + rocket's axis of symmetry as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[1] - return R2 + @funcify_method("Time (s)", "R3 (N)", "spline", "constant") + def R3(self): + """Aerodynamic force along the rocket's axis of symmetry as a rocketpy.Function + of time.""" + return self.retrieve_temporary_values_arrays[2] - @cached_property - def R3(self, interpolation="spline", extrapolation="natural"): - """Aerodynamic force along the rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - R3: Function - Aerodynamic force along the rocket's axis of symmetry. - """ - R3 = self.retrieve_temporary_values_arrays[2] - R3 = Function(R3, "Time (s)", "R3 (m)", interpolation, extrapolation) - - return R3 - - @cached_property - def M1(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "M1 (Nm)", "spline", "constant") + def M1(self): """Aerodynamic bending moment in the same direction as the axis that is - perpendicular to the rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - M1: Function - Aero moment along the first axis that is perpendicular to the - rocket's axis of symmetry. - """ - M1 = self.retrieve_temporary_values_arrays[3] - M1 = Function(M1, "Time (s)", "M1 (Nm)", interpolation, extrapolation) - - return M1 - - @cached_property - def M2(self, interpolation="spline", extrapolation="natural"): - """Aerodynamic moment in the same direction as the second axis that is - perpendicular to the rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - M2: Function - Aero moment along the second axis that is perpendicular to the - rocket's axis of symmetry. - """ - M2 = self.retrieve_temporary_values_arrays[4] - M2 = Function(M2, "Time (s)", "M2 (Nm)", interpolation, extrapolation) - - return M2 - - @cached_property - def M3(self, interpolation="spline", extrapolation="natural"): - """Aerodynamic bending in the rocket's axis of symmetry direction. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - M3: Function - Aero moment in the same direction as the rocket's axis of symmetry. - """ - M3 = self.retrieve_temporary_values_arrays[5] - M3 = Function(M3, "Time (s)", "M3 (Nm)", interpolation, extrapolation) - - return M3 - - @cached_property - def pressure(self, interpolation="spline", extrapolation="natural"): - """Air pressure at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - pressure: Function - Air pressure at each time step. - """ - pressure = self.retrieve_temporary_values_arrays[6] - pressure = Function( - pressure, "Time (s)", "Pressure (Pa)", interpolation, extrapolation - ) - - return pressure - - @cached_property - def density(self, interpolation="spline", extrapolation="natural"): - """Air density at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - density: Function - Air density at each time step. - """ - density = self.retrieve_temporary_values_arrays[7] - density = Function( - density, "Time (s)", "Density (kg/m³)", interpolation, extrapolation - ) - - return density - - @cached_property - def dynamicViscosity(self, interpolation="spline", extrapolation="natural"): - """Air dynamic viscosity at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - dynamicViscosity: Function - Air dynamic viscosity at each time step. - """ - dynamicViscosity = self.retrieve_temporary_values_arrays[8] - dynamicViscosity = Function( - dynamicViscosity, - "Time (s)", - "Dynamic Viscosity (Pa s)", - interpolation, - extrapolation, - ) - - return dynamicViscosity - - @cached_property - def speedOfSound(self, interpolation="spline", extrapolation="natural"): - """Speed of sound at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - speedOfSound: Function - Speed of sound at each time step. - """ - speedOfSound = self.retrieve_temporary_values_arrays[9] - speedOfSound = Function( - speedOfSound, - "Time (s)", - "Speed of Sound (m/s)", - interpolation, - extrapolation, - ) - - return speedOfSound - - @cached_property - def windVelocityX(self, interpolation="spline", extrapolation="natural"): - """Wind velocity in the X direction at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - windVelocityX: Function - Wind velocity in the X direction at each time step. + perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time. """ - windVelocityX = self.retrieve_temporary_values_arrays[10] - windVelocityX = Function( - windVelocityX, - "Time (s)", - "Wind Velocity X (East) (m/s)", - interpolation, - extrapolation, - ) - - return windVelocityX - - @cached_property - def windVelocityY(self, interpolation="spline", extrapolation="natural"): - """Wind velocity in the Y direction at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. + return self.retrieve_temporary_values_arrays[3] - Returns - ------- - windVelocityY: Function - Wind velocity in the Y direction at each time step. + @funcify_method("Time (s)", "M2 (Nm)", "spline", "constant") + def M2(self): + """Aerodynamic bending moment in the same direction as the axis that is + perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time. """ - windVelocityY = self.retrieve_temporary_values_arrays[10] - windVelocityY = Function( - windVelocityY, - "Time (s)", - "Wind Velocity Y (North) (m/s)", - interpolation, - extrapolation, - ) - - return windVelocityY + return self.retrieve_temporary_values_arrays[4] + + @funcify_method("Time (s)", "M3 (Nm)", "spline", "constant") + def M3(self): + """Aerodynamic bending moment in the same direction as the rocket's axis of + symmetry as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[5] + + @funcify_method("Time (s)", "Pressure (Pa)", "spline", "constant") + def pressure(self): + """Air pressure felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[6] + + @funcify_method("Time (s)", "Density (kg/m³)", "spline", "constant") + def density(self): + """Air density felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[7] + + @funcify_method("Time (s)", "Dynamic Viscosity (Pa s)", "spline", "constant") + def dynamicViscosity(self): + """Air dynamic viscosity felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[7] + + @funcify_method("Time (s)", "Speed of Sound (m/s)", "spline", "constant") + def speedOfSound(self): + """Speed of sound in the air felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[9] + + @funcify_method("Time (s)", "Wind Velocity X (East) (m/s)", "spline", "constant") + def windVelocityX(self): + """Wind velocity in the X direction (east) as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[10] + + @funcify_method("Time (s)", "Wind Velocity Y (North) (m/s)", "spline", "constant") + def windVelocityY(self): + """Wind velocity in the y direction (north) as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[11] # Process fourth type of output - values calculated from previous outputs # Kinematics functions and values # Velocity Magnitude - @cached_property + @funcify_method("Time (s)", "Speed - Velocity Magnitude (m/s)") def speed(self): - speed = (self.vx**2 + self.vy**2 + self.vz**2) ** 0.5 - speed.setOutputs("Speed - Velocity Magnitude (m/s)") - return speed + """Rocket speed, or velocity magnitude, as a rocketpy.Function of time.""" + return (self.vx**2 + self.vy**2 + self.vz**2) ** 0.5 @cached_property def maxSpeedTime(self): + """Time at which the rocket reaches its maximum speed.""" maxSpeedTimeIndex = np.argmax(self.speed[:, 1]) return self.speed[maxSpeedTimeIndex, 0] @cached_property def maxSpeed(self): + """Maximum speed reached by the rocket.""" return self.speed(self.maxSpeedTime) # Accelerations - @cached_property + @funcify_method("Time (s)", "Acceleration Magnitude (m/s²)") def acceleration(self): - acceleration = (self.ax**2 + self.ay**2 + self.az**2) ** 0.5 - acceleration.setOutputs("Acceleration Magnitude (m/s²)") - return acceleration + """Rocket acceleration magnitude as a rocketpy.Function of time.""" + return (self.ax**2 + self.ay**2 + self.az**2) ** 0.5 @cached_property def maxAcceleration(self): + """Maximum acceleration reached by the rocket.""" maxAccelerationTimeIndex = np.argmax(self.acceleration[:, 1]) return self.acceleration[maxAccelerationTimeIndex, 1] @cached_property def maxAccelerationTime(self): + """Time at which the rocket reaches its maximum acceleration.""" maxAccelerationTimeIndex = np.argmax(self.acceleration[:, 1]) return self.acceleration[maxAccelerationTimeIndex, 0] - # Path Angle - @cached_property + @funcify_method("Time (s)", "Horizontal Speed (m/s)") def horizontalSpeed(self): + """Rocket horizontal speed as a rocketpy.Function of time.""" return (self.vx**2 + self.vy**2) ** 0.5 - @cached_property + # Path Angle + @funcify_method("Time (s)", "Path Angle (°)", "spline", "constant") def pathAngle(self): + """Rocket path angle as a rocketpy.Function of time.""" pathAngle = (180 / np.pi) * np.arctan2( self.vz[:, 1], self.horizontalSpeed[:, 1] ) - pathAngle = np.column_stack([self.vz[:, 0], pathAngle]) - return Function(pathAngle, "Time (s)", "Path Angle (°)") + return np.column_stack([self.time, pathAngle]) # Attitude Angle - @cached_property + @funcify_method("Time (s)", "Attitude Vector X Component") def attitudeVectorX(self): + """Rocket attitude vector X component as a rocketpy.Function of time.""" return 2 * (self.e1 * self.e3 + self.e0 * self.e2) # a13 - @cached_property + @funcify_method("Time (s)", "Attitude Vector Y Component") def attitudeVectorY(self): + """Rocket attitude vector Y component as a rocketpy.Function of time.""" return 2 * (self.e2 * self.e3 - self.e0 * self.e1) # a23 - @cached_property + @funcify_method("Time (s)", "Attitude Vector Z Component") def attitudeVectorZ(self): + """Rocket attitude vector Z component as a rocketpy.Function of time.""" return 1 - 2 * (self.e1**2 + self.e2**2) # a33 - @cached_property + @funcify_method("Time (s)", "Attitude Angle (°)") def attitudeAngle(self): + """Rocket attitude angle as a rocketpy.Function of time.""" horizontalAttitudeProj = ( self.attitudeVectorX**2 + self.attitudeVectorY**2 ) ** 0.5 attitudeAngle = (180 / np.pi) * np.arctan2( self.attitudeVectorZ[:, 1], horizontalAttitudeProj[:, 1] ) - attitudeAngle = np.column_stack([self.attitudeVectorZ[:, 0], attitudeAngle]) - return Function(attitudeAngle, "Time (s)", "Attitude Angle (°)") + attitudeAngle = np.column_stack([self.time, attitudeAngle]) + return attitudeAngle # Lateral Attitude Angle - @cached_property + @funcify_method("Time (s)", "Lateral Attitude Angle (°)") def lateralAttitudeAngle(self): + """Rocket lateral attitude angle as a rocketpy.Function of time.""" lateralVectorAngle = (np.pi / 180) * (self.heading - 90) lateralVectorX = np.sin(lateralVectorAngle) lateralVectorY = np.cos(lateralVectorAngle) @@ -2249,204 +1815,178 @@ def lateralAttitudeAngle(self): lateralAttitudeAngle = (180 / np.pi) * np.arctan2( attitudeLateralProj, attitudeLateralPlaneProj ) - lateralAttitudeAngle = np.column_stack( - [self.attitudeVectorZ[:, 0], lateralAttitudeAngle] - ) - return Function(lateralAttitudeAngle, "Time (s)", "Lateral Attitude Angle (°)") + lateralAttitudeAngle = np.column_stack([self.time, lateralAttitudeAngle]) + return lateralAttitudeAngle # Euler Angles - @cached_property + @funcify_method("Time (s)", "Precession Angle - ψ (°)", "spline", "constant") def psi(self): + """Precession angle as a rocketpy.Function of time.""" psi = (180 / np.pi) * ( np.arctan2(self.e3[:, 1], self.e0[:, 1]) + np.arctan2(-self.e2[:, 1], -self.e1[:, 1]) ) # Precession angle - psi = np.column_stack([self.e1[:, 0], psi]) # Precession angle - return Function(psi, "Time (s)", "Precession Angle - ψ (°)") + psi = np.column_stack([self.time, psi]) # Precession angle + return psi - @cached_property + @funcify_method("Time (s)", "Spin Angle - φ (°)", "spline", "constant") def phi(self): + """Spin angle as a rocketpy.Function of time.""" phi = (180 / np.pi) * ( np.arctan2(self.e3[:, 1], self.e0[:, 1]) - np.arctan2(-self.e2[:, 1], -self.e1[:, 1]) ) # Spin angle - phi = np.column_stack([self.e1[:, 0], phi]) # Spin angle - return Function(phi, "Time (s)", "Spin Angle - φ (°)") + phi = np.column_stack([self.time, phi]) # Spin angle + return phi - @cached_property + @funcify_method("Time (s)", "Nutation Angle - θ (°)", "spline", "constant") def theta(self): + """Nutation angle as a rocketpy.Function of time.""" theta = ( (180 / np.pi) * 2 * np.arcsin(-((self.e1[:, 1] ** 2 + self.e2[:, 1] ** 2) ** 0.5)) ) # Nutation angle - theta = np.column_stack([self.e1[:, 0], theta]) # Nutation angle - return Function(theta, "Time (s)", "Nutation Angle - θ (°)") + theta = np.column_stack([self.time, theta]) # Nutation angle + return theta # Fluid Mechanics variables # Freestream Velocity - @cached_property - def streamVelocityX(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "Freestream Velocity X (m/s)", "spline", "constant") + def streamVelocityX(self): + """Freestream velocity X component as a rocketpy.Function of time.""" streamVelocityX = np.column_stack( - (self.vx[:, 0], self.windVelocityX[:, 1] - self.vx[:, 1]) - ) - streamVelocityX = Function( - streamVelocityX, - "Time (s)", - "Freestream Velocity X (m/s)", - interpolation, - extrapolation, + (self.time, self.windVelocityX[:, 1] - self.vx[:, 1]) ) return streamVelocityX - @cached_property - def streamVelocityY(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "Freestream Velocity Y (m/s)", "spline", "constant") + def streamVelocityY(self): + """Freestream velocity Y component as a rocketpy.Function of time.""" streamVelocityY = np.column_stack( - (self.vy[:, 0], self.windVelocityY[:, 1] - self.vy[:, 1]) - ) - streamVelocityY = Function( - streamVelocityY, - "Time (s)", - "Freestream Velocity Y (m/s)", - interpolation, - extrapolation, + (self.time, self.windVelocityY[:, 1] - self.vy[:, 1]) ) return streamVelocityY - @cached_property + @funcify_method("Time (s)", "Freestream Velocity Z (m/s)", "spline", "constant") def streamVelocityZ(self, interpolation="spline", extrapolation="natural"): - streamVelocityZ = np.column_stack((self.vz[:, 0], -self.vz[:, 1])) - streamVelocityZ = Function( - streamVelocityZ, - "Time (s)", - "Freestream Velocity Z (m/s)", - interpolation, - extrapolation, - ) + """Freestream velocity Z component as a rocketpy.Function of time.""" + streamVelocityZ = np.column_stack((self.time, -self.vz[:, 1])) return streamVelocityZ - @cached_property + @funcify_method("Time (s)", "Freestream Speed (m/s)", "spline", "constant") def freestreamSpeed(self): + """Freestream speed as a rocketpy.Function of time.""" freestreamSpeed = ( self.streamVelocityX**2 + self.streamVelocityY**2 + self.streamVelocityZ**2 ) ** 0.5 - freestreamSpeed.setInputs("Time (s)") - freestreamSpeed.setOutputs("Freestream Speed (m/s)") - return freestreamSpeed + return freestreamSpeed.source # Apogee Freestream speed @cached_property def apogeeFreestreamSpeed(self): + """Freestream speed at apogee in m/s.""" return self.freestreamSpeed(self.apogeeTime) # Mach Number - @cached_property + @funcify_method("Time (s)", "Mach Number", "spline", "constant") def MachNumber(self): - MachNumber = self.freestreamSpeed / self.speedOfSound - MachNumber.setInputs("Time (s)") - MachNumber.setOutputs("Mach Number") - return MachNumber + """Mach number as a rocketpy.Function of time.""" + return self.freestreamSpeed / self.speedOfSound @cached_property def maxMachNumberTime(self): + """Time of maximum Mach number.""" maxMachNumberTimeIndex = np.argmax(self.MachNumber[:, 1]) return self.MachNumber[maxMachNumberTimeIndex, 0] @cached_property def maxMachNumber(self): + """Maximum Mach number.""" return self.MachNumber(self.maxMachNumberTime) # Reynolds Number - @cached_property + @funcify_method("Time (s)", "Reynolds Number", "spline", "constant") def ReynoldsNumber(self): - ReynoldsNumber = ( - self.density * self.freestreamSpeed / self.dynamicViscosity - ) * (2 * self.rocket.radius) - ReynoldsNumber.setInputs("Time (s)") - ReynoldsNumber.setOutputs("Reynolds Number") - return ReynoldsNumber + """Reynolds number as a rocketpy.Function of time.""" + return (self.density * self.freestreamSpeed / self.dynamicViscosity) * ( + 2 * self.rocket.radius + ) @cached_property def maxReynoldsNumberTime(self): + """Time of maximum Reynolds number.""" maxReynoldsNumberTimeIndex = np.argmax(self.ReynoldsNumber[:, 1]) return self.ReynoldsNumber[maxReynoldsNumberTimeIndex, 0] @cached_property def maxReynoldsNumber(self): + """Maximum Reynolds number.""" return self.ReynoldsNumber(self.maxReynoldsNumberTime) # Dynamic Pressure - @cached_property + @funcify_method("Time (s)", "Dynamic Pressure (Pa)", "spline", "constant") def dynamicPressure(self): - dynamicPressure = 0.5 * self.density * self.freestreamSpeed**2 - dynamicPressure.setInputs("Time (s)") - dynamicPressure.setOutputs("Dynamic Pressure (Pa)") - return dynamicPressure + """Dynamic pressure as a rocketpy.Function of time.""" + return 0.5 * self.density * self.freestreamSpeed**2 @cached_property def maxDynamicPressureTime(self): + """Time of maximum dynamic pressure.""" maxDynamicPressureTimeIndex = np.argmax(self.dynamicPressure[:, 1]) return self.dynamicPressure[maxDynamicPressureTimeIndex, 0] @cached_property def maxDynamicPressure(self): + """Maximum dynamic pressure.""" return self.dynamicPressure(self.maxDynamicPressureTime) # Total Pressure - @cached_property + @funcify_method("Time (s)", "Total Pressure (Pa)", "spline", "constant") def totalPressure(self): - totalPressure = self.pressure * (1 + 0.2 * self.MachNumber**2) ** (3.5) - totalPressure.setInputs("Time (s)") - totalPressure.setOutputs("Total Pressure (Pa)") - return totalPressure + return self.pressure * (1 + 0.2 * self.MachNumber**2) ** (3.5) @cached_property - def maxtotalPressureTime(self): - maxtotalPressureTimeIndex = np.argmax(self.totalPressure[:, 1]) - return self.totalPressure[maxtotalPressureTimeIndex, 0] + def maxTotalPressureTime(self): + """Time of maximum total pressure.""" + maxTotalPressureTimeIndex = np.argmax(self.totalPressure[:, 1]) + return self.totalPressure[maxTotalPressureTimeIndex, 0] @cached_property - def maxtotalPressure(self): - return self.totalPressure(self.maxtotalPressureTime) + def maxTotalPressure(self): + """Maximum total pressure.""" + return self.totalPressure(self.maxTotalPressureTime) # Dynamics functions and variables # Aerodynamic Lift and Drag - @cached_property + @funcify_method("Time (s)", "Aerodynamic Lift Force (N)", "spline", "constant") def aerodynamicLift(self): - aerodynamicLift = (self.R1**2 + self.R2**2) ** 0.5 - aerodynamicLift.setInputs("Time (s)") - aerodynamicLift.setOutputs("Aerodynamic Lift Force (N)") - return aerodynamicLift + """Aerodynamic lift force as a rocketpy.Function of time.""" + return (self.R1**2 + self.R2**2) ** 0.5 - @cached_property + @funcify_method("Time (s)", "Aerodynamic Drag Force (N)", "spline", "constant") def aerodynamicDrag(self): - aerodynamicDrag = -1 * self.R3 - aerodynamicDrag.setInputs("Time (s)") - aerodynamicDrag.setOutputs("Aerodynamic Drag Force (N)") - return aerodynamicDrag + """Aerodynamic drag force as a rocketpy.Function of time.""" + return -1 * self.R3 - @cached_property + @funcify_method("Time (s)", "Aerodynamic Bending Moment (Nm)", "spline", "constant") def aerodynamicBendingMoment(self): + """Aerodynamic bending moment as a rocketpy.Function of time.""" + return (self.M1**2 + self.M2**2) ** 0.5 - aerodynamicBendingMoment = (self.M1**2 + self.M2**2) ** 0.5 - aerodynamicBendingMoment.setInputs("Time (s)") - aerodynamicBendingMoment.setOutputs("Aerodynamic Bending Moment (N m)") - return aerodynamicBendingMoment - - @cached_property + @funcify_method("Time (s)", "Aerodynamic Spin Moment (Nm)", "spline", "constant") def aerodynamicSpinMoment(self): - aerodynamicSpinMoment = self.M3 - aerodynamicSpinMoment.setInputs("Time (s)") - aerodynamicSpinMoment.setOutputs("Aerodynamic Spin Moment (N m)") - return aerodynamicSpinMoment + """Aerodynamic spin moment as a rocketpy.Function of time.""" + return self.M3 # Energy # Kinetic Energy - @cached_property + @funcify_method("Time (s)", "Rotational Kinetic Energy (J)", "spline", "constant") def rotationalEnergy(self): + """Rotational kinetic energy as a rocketpy.Function of time.""" b = -self.rocket.distanceRocketPropellant mu = self.rocket.reducedMass Rz = self.rocket.inertiaZ @@ -2454,103 +1994,80 @@ def rotationalEnergy(self): Tz = self.rocket.motor.inertiaZ Ti = self.rocket.motor.inertiaI I1, I2, I3 = (Ri + Ti + mu * b**2), (Ri + Ti + mu * b**2), (Rz + Tz) - # Redefine I1, I2 and I3 grid - grid = self.vx[:, 0] - I1 = Function(np.column_stack([grid, I1(grid)]), "Time (s)") - I2 = Function(np.column_stack([grid, I2(grid)]), "Time (s)") - I3 = Function(np.column_stack([grid, I3(grid)]), "Time (s)") + # Redefine I1, I2 and I3 time grid to allow for efficient Function algebra + I1.setDiscreteBasedOnModel(self.w1) + I2.setDiscreteBasedOnModel(self.w1) + I3.setDiscreteBasedOnModel(self.w1) rotationalEnergy = 0.5 * ( I1 * self.w1**2 + I2 * self.w2**2 + I3 * self.w3**2 ) - rotationalEnergy.setInputs("Time (s)") - rotationalEnergy.setOutputs("Rotational Kinetic Energy (J)") return rotationalEnergy - @cached_property + @funcify_method( + "Time (s)", "Translational Kinetic Energy (J)", "spline", "constant" + ) def translationalEnergy(self): - # Redefine total mass grid - totalMass = self.rocket.totalMass - grid = self.vx[:, 0] - totalMass = Function(np.column_stack([grid, totalMass(grid)]), "Time (s)") + """Translational kinetic energy as a rocketpy.Function of time.""" + # Redefine totalMass time grid to allow for efficient Function algebra + totalMass = deepcopy(self.rocket.totalMass) + totalMass.setDiscreteBasedOnModel(self.vz) translationalEnergy = ( 0.5 * totalMass * (self.vx**2 + self.vy**2 + self.vz**2) ) - translationalEnergy.setInputs("Time (s)") - translationalEnergy.setOutputs("Translational Kinetic Energy (J)") return translationalEnergy - @cached_property + @funcify_method("Time (s)", "Kinetic Energy (J)", "spline", "constant") def kineticEnergy(self): - kineticEnergy = self.rotationalEnergy + self.translationalEnergy - kineticEnergy.setInputs("Time (s)") - kineticEnergy.setOutputs("Kinetic Energy (J)") - return kineticEnergy + """Total kinetic energy as a rocketpy.Function of time.""" + return self.rotationalEnergy + self.translationalEnergy # Potential Energy - @cached_property + @funcify_method("Time (s)", "Potential Energy (J)", "spline", "constant") def potentialEnergy(self): - # Redefine total mass grid - totalMass = self.rocket.totalMass - grid = self.vx[:, 0] - totalMass = Function(np.column_stack([grid, totalMass(grid)]), "Time (s)") + """Potential energy as a rocketpy.Function of time.""" + # Redefine totalMass time grid to allow for efficient Function algebra + totalMass = deepcopy(self.rocket.totalMass) + totalMass.setDiscreteBasedOnModel(self.z) + # TODO: change calculation method to account for variable gravity potentialEnergy = totalMass * self.env.g * self.z - potentialEnergy.setInputs("Time (s)") return potentialEnergy # Total Mechanical Energy - @cached_property + @funcify_method("Time (s)", "Mechanical Energy (J)", "spline", "constant") def totalEnergy(self): - totalEnergy = self.kineticEnergy + self.potentialEnergy - totalEnergy.setInputs("Time (s)") - totalEnergy.setOutputs("Total Mechanical Energy (J)") - return totalEnergy + """Total mechanical energy as a rocketpy.Function of time.""" + return self.kineticEnergy + self.potentialEnergy # Thrust Power - @cached_property + @funcify_method("Time (s)", "Thrust Power (W)", "spline", "constant") def thrustPower(self): - grid = self.vx[:, 0] - # Redefine thrust grid - thrust = Function( - np.column_stack([grid, self.rocket.motor.thrust(grid)]), "Time (s)" - ) + """Thrust power as a rocketpy.Function of time.""" + thrust = deepcopy(self.rocket.motor.thrust) + thrust = thrust.setDiscreteBasedOnModel(self.speed) thrustPower = thrust * self.speed - thrustPower.setOutputs("Thrust Power (W)") return thrustPower # Drag Power - @cached_property + @funcify_method("Time (s)", "Drag Power (W)", "spline", "constant") def dragPower(self): + """Drag power as a rocketpy.Function of time.""" dragPower = self.R3 * self.speed dragPower.setOutputs("Drag Power (W)") return dragPower # Angle of Attack - @cached_property - def angleOfAttack(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "Angle of Attack (°)", "spline", "constant") + def angleOfAttack(self): """Angle of attack of the rocket with respect to the freestream - velocity vector. - - Parameters - ---------- - interpolation : str, optional - Interpolation method, by default "spline" - extrapolation : str, optional - Extrapolation method, by default "natural" - - Returns - ------- - angleOfAttack: Function - Angle of attack of the rocket with respect to the freestream - velocity vector. - """ + velocity vector.""" dotProduct = [ -self.attitudeVectorX.getValueOpt(i) * self.streamVelocityX.getValueOpt(i) - self.attitudeVectorY.getValueOpt(i) * self.streamVelocityY.getValueOpt(i) - self.attitudeVectorZ.getValueOpt(i) * self.streamVelocityZ.getValueOpt(i) - for i in self.x[:, 0] + for i in self.time ] # Define freestream speed list - freestreamSpeed = [self.freestreamSpeed(i) for i in self.x[:, 0]] + freestreamSpeed = [self.freestreamSpeed(i) for i in self.time] freestreamSpeed = np.nan_to_num(freestreamSpeed) # Normalize dot product @@ -2562,146 +2079,85 @@ def angleOfAttack(self, interpolation="spline", extrapolation="natural"): # Calculate angle of attack and convert to degrees angleOfAttack = np.rad2deg(np.arccos(dotProductNormalized)) - angleOfAttack = np.column_stack([self.x[:, 0], angleOfAttack]) - - # Convert to Function - angleOfAttack = Function( - angleOfAttack, - "Time (s)", - "Angle Of Attack (°)", - interpolation, - extrapolation, - ) + angleOfAttack = np.column_stack([self.time, angleOfAttack]) return angleOfAttack # Stability and Control variables + def _frequency_analysis(self, func, samplingFrequency=1000): + """Frequency analysis of any function of time between out of rail time + and final simulation time. + + Parameters + ---------- + func : rocketpy.Function, callable + Function of time to be analyzed. + + Returns + ------- + FourierFrequencies : numpy.ndarray + Array of frequencies. + FourierAmplitude : numpy.ndarray + Array of amplitude values. + """ + samplingTimeStep = 1.0 / samplingFrequency + samplingRange = np.arange(self.outOfRailTime, self.tFinal, samplingTimeStep) + sampledPoints = func(samplingRange) + sampledPoints -= np.mean(sampledPoints) + numberOfSamples = len(sampledPoints) + FourierAmplitude = np.abs(np.fft.fft(sampledPoints) / numberOfSamples) + FourierFrequencies = np.fft.fftfreq(numberOfSamples, samplingTimeStep) + return FourierFrequencies, FourierAmplitude + # Angular velocities frequency response - Fourier Analysis - @cached_property + @funcify_method( + "Frequency (Hz)", "ω1 Angle Fourier Amplitude", "spline", "constant" + ) def omega1FrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.w1(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - omega1FrequencyResponse = np.column_stack([frq, abs(Y)]) - omega1FrequencyResponse = Function( - omega1FrequencyResponse, "Frequency (Hz)", "Omega 1 Angle Fourier Amplitude" - ) - return omega1FrequencyResponse + """Angular velocity 1 frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.w1)) - @cached_property + @funcify_method( + "Frequency (Hz)", "ω2 Angle Fourier Amplitude", "spline", "constant" + ) def omega2FrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.w2(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - omega2FrequencyResponse = np.column_stack([frq, abs(Y)]) - omega2FrequencyResponse = Function( - omega2FrequencyResponse, "Frequency (Hz)", "Omega 2 Angle Fourier Amplitude" - ) - return omega2FrequencyResponse + """Angular velocity 2 frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.w2)) - @cached_property + @funcify_method( + "Frequency (Hz)", "ω3 Angle Fourier Amplitude", "spline", "constant" + ) def omega3FrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.w3(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - omega3FrequencyResponse = np.column_stack([frq, abs(Y)]) - omega3FrequencyResponse = Function( - omega3FrequencyResponse, "Frequency (Hz)", "Omega 3 Angle Fourier Amplitude" - ) - return omega3FrequencyResponse + """Angular velocity 3 frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.w3)) - @cached_property + @funcify_method( + "Frequency (Hz)", "Attitude Angle Fourier Amplitude", "spline", "constant" + ) def attitudeFrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.attitudeAngle(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - attitudeFrequencyResponse = np.column_stack([frq, abs(Y)]) - attitudeFrequencyResponse = Function( - attitudeFrequencyResponse, - "Frequency (Hz)", - "Attitude Angle Fourier Amplitude", - ) - return attitudeFrequencyResponse + """Attitude frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.attitudeAngle)) # Static Margin @cached_property def staticMargin(self): + """Static margin of the rocket.""" return self.rocket.staticMargin # Rail Button Forces - @cached_property + @funcify_method("Time (s)", "Upper Rail Button Normal Force (N)", "spline", "zero") def railButton1NormalForce(self): - """Upper rail button normal force as a function of time - - Returns - ------- - railButton1NormalForce: Function - Upper rail button normal force as a function of time - """ - + """Upper rail button normal force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F11, F12 = self.calculate_rail_button_forces[0:2] else: F11, F12 = self.calculate_rail_button_forces()[0:2] alpha = self.rocket.railButtons.angularPosition * (np.pi / 180) - railButton1NormalForce = F11 * np.cos(alpha) + F12 * np.sin(alpha) - railButton1NormalForce.setOutputs("Upper Rail Button Normal Force (N)") - - return railButton1NormalForce + return F11 * np.cos(alpha) + F12 * np.sin(alpha) - @cached_property + @funcify_method("Time (s)", "Upper Rail Button Shear Force (N)", "spline", "zero") def railButton1ShearForce(self): - """Upper rail button shear force as a function of time - - Returns - ------- - _type_ - _description_ - """ + """Upper rail button shear force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F11, F12 = self.calculate_rail_button_forces[0:2] else: @@ -2709,159 +2165,91 @@ def railButton1ShearForce(self): alpha = self.rocket.railButtons.angularPosition * ( np.pi / 180 ) # Rail buttons angular position - railButton1ShearForce = F11 * -np.sin(alpha) + F12 * np.cos(alpha) - railButton1ShearForce.setOutputs("Upper Rail Button Shear Force (N)") - - return railButton1ShearForce + return F11 * -np.sin(alpha) + F12 * np.cos(alpha) - @cached_property + @funcify_method("Time (s)", "Lower Rail Button Normal Force (N)", "spline", "zero") def railButton2NormalForce(self): - """Lower rail button normal force as a function of time - - Returns - ------- - railButton2NormalForce: Function - Lower rail button normal force as a function of time - """ + """Lower rail button normal force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F21, F22 = self.calculate_rail_button_forces[2:4] else: F21, F22 = self.calculate_rail_button_forces()[2:4] alpha = self.rocket.railButtons.angularPosition * (np.pi / 180) - railButton2NormalForce = F21 * np.cos(alpha) + F22 * np.sin(alpha) - railButton2NormalForce.setOutputs("Lower Rail Button Normal Force (N)") - - return railButton2NormalForce + return F21 * np.cos(alpha) + F22 * np.sin(alpha) - @cached_property + @funcify_method("Time (s)", "Lower Rail Button Shear Force (N)", "spline", "zero") def railButton2ShearForce(self): - """Lower rail button shear force as a function of time - - Returns - ------- - railButton2ShearForce: Function - Lower rail button shear force as a function of time - """ + """Lower rail button shear force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F21, F22 = self.calculate_rail_button_forces[2:4] else: F21, F22 = self.calculate_rail_button_forces()[2:4] alpha = self.rocket.railButtons.angularPosition * (np.pi / 180) - railButton2ShearForce = F21 * -np.sin(alpha) + F22 * np.cos(alpha) - railButton2ShearForce.setOutputs("Lower Rail Button Shear Force (N)") - - return railButton2ShearForce + return F21 * -np.sin(alpha) + F22 * np.cos(alpha) @cached_property def maxRailButton1NormalForce(self): - """Maximum upper rail button normal force, in Newtons - - Returns - ------- - maxRailButton1NormalForce: float - Maximum upper rail button normal force, in Newtons - """ + """Maximum upper rail button normal force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1NormalForce[:outOfRailTimeIndex]) + return np.max(self.railButton1NormalForce[: self.outOfRailTimeIndex]) @cached_property def maxRailButton1ShearForce(self): - """Maximum upper rail button shear force, in Newtons - - Returns - ------- - maxRailButton1ShearForce: float - Maximum upper rail button shear force, in Newtons - """ + """Maximum upper rail button shear force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1ShearForce[:outOfRailTimeIndex]) + return np.max(self.railButton1ShearForce[: self.outOfRailTimeIndex]) @cached_property def maxRailButton2NormalForce(self): - """Maximum lower rail button normal force, in Newtons - - Returns - ------- - maxRailButton2NormalForce: float - Maximum lower rail button normal force, in Newtons - """ + """Maximum lower rail button normal force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2NormalForce[:outOfRailTimeIndex]) + return np.max(self.railButton2NormalForce[: self.outOfRailTimeIndex]) @cached_property def maxRailButton2ShearForce(self): - """Maximum lower rail button shear force, in Newtons - - Returns - ------- - maxRailButton2ShearForce: float - Maximum lower rail button shear force, in Newtons - """ + """Maximum lower rail button shear force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2ShearForce[:outOfRailTimeIndex]) + return np.max(self.railButton2ShearForce[: self.outOfRailTimeIndex]) - @cached_property + @funcify_method( + "Time (s)", "Horizontal Distance to Launch Point (m)", "spline", "constant" + ) def drift(self): - """Rocket horizontal distance to tha launch point, in meters, at each - time step. - - Returns - ------- - drift: Function - Rocket horizontal distance to tha launch point, in meters, at each - time step. - """ - drift = np.column_stack( - (self.x[:, 0], (self.x[:, 1] ** 2 + self.y[:, 1] ** 2) ** 0.5) - ) - drift = Function( - drift, - "Time (s)", - "Horizontal Distance to Launch Point (m)", + """Rocket horizontal distance to tha launch point, in meters, as a + rocketpy.Function of time.""" + return np.column_stack( + (self.time, (self.x[:, 1] ** 2 + self.y[:, 1] ** 2) ** 0.5) ) - return drift - - @cached_property - def bearing(self, interpolation="spline", extrapolation="natural"): - """Rocket bearing compass, in degrees, at each time step. - - Returns - ------- - bearing: Function - Rocket bearing, in degrees, at each time step. - """ - + @funcify_method("Time (s)", "Bearing (°)", "spline", "constant") + def bearing(self): + """Rocket bearing compass, in degrees, as a rocketpy.Function of time.""" # Get some nicknames - t = self.x[:, 0] + t = self.time x, y = self.x[:, 1], self.y[:, 1] bearing = [] for i in range(len(t)): @@ -2891,26 +2279,11 @@ def bearing(self, interpolation="spline", extrapolation="natural"): bearing = np.rad2deg(bearing) bearing = np.column_stack((t, bearing)) - print(bearing) - bearing = Function( - bearing, - "Time (s)", - "Bearing (deg)", - interpolation, - extrapolation, - ) - return bearing - @cached_property + @funcify_method("Time (s)", "Latitude (°)", "linear") def latitude(self): - """Rocket latitude coordinate, in degrees, at each time step. - - Returns - ------- - latitude: Function - Rocket latitude coordinate, in degrees, at each time step. - """ + """Rocket latitude coordinate, in degrees, as a rocketpy.Function of time.""" lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians # Applies the haversine equation to find final lat/lon coordinates @@ -2922,20 +2295,11 @@ def latitude(self): * np.cos(np.deg2rad(self.bearing[:, 1])) ) ) - latitude = np.column_stack((self.x[:, 0], latitude)) - latitude = Function(latitude, "Time (s)", "Latitude (°)", "linear") + return np.column_stack((self.time, latitude)) - return latitude - - @cached_property + @funcify_method("Time (s)", "Longitude (°)", "linear") def longitude(self): - """Rocket longitude coordinate, in degrees, at each time step. - - Returns - ------- - longitude: Function - Rocket longitude coordinate, in degrees, at each time step. - """ + """Rocket longitude coordinate, in degrees, as a rocketpy.Function of time.""" lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians lon1 = np.deg2rad(self.env.lon) # Launch lon point converted to radians @@ -2951,10 +2315,7 @@ def longitude(self): ) ) - longitude = np.column_stack((self.x[:, 0], longitude)) - longitude = Function(longitude, "Time (s)", "Longitude (°)", "linear") - - return longitude + return np.column_stack((self.time, longitude)) @cached_property def retrieve_acceleration_arrays(self): @@ -3104,7 +2465,7 @@ def retrieve_temporary_values_arrays(self): @cached_property def calculate_rail_button_forces(self): """Calculate the forces applied to the rail buttons while rocket is still - on the launch rail. It will return 0 if none rail buttons are defined. + on the launch rail. It will return 0 if no rail buttons are defined. Returns ------- @@ -3117,7 +2478,6 @@ def calculate_rail_button_forces(self): F22: Function Rail Button 2 force in the 2 direction """ - if self.rocket.railButtons is None: warnings.warn( "Trying to calculate rail button forces without rail buttons defined." @@ -3137,10 +2497,11 @@ def calculate_rail_button_forces(self): F22 = (self.R2 * D1 - self.M1) / (D1 + D2) F22.setOutputs("Lower button force direction 2 (m)") - # F11 = F11[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail - # F12 = F12[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail - # F21 = F21[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail - # F22 = F22[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail + # Limit force calculation to when rocket is in rail + F11 = F11[: self.outOfRailTimeIndex, :] + F12 = F12[: self.outOfRailTimeIndex, :] + F21 = F21[: self.outOfRailTimeIndex, :] + F22 = F22[: self.outOfRailTimeIndex, :] return F11, F12, F21, F22 @@ -3180,6 +2541,9 @@ def postProcess(self, interpolation="spline", extrapolation="natural"): calculation of secondary values such as energy and conversion of lists to Function objects to facilitate plotting. + * This method is deprecated and is only kept here for backwards compatibility. + * All attributes that need to be post processed are computed just in time. + Parameters ---------- None @@ -3188,7 +2552,6 @@ def postProcess(self, interpolation="spline", extrapolation="natural"): ------ None """ - # Register post processing self.postProcessed = True @@ -3206,16 +2569,10 @@ def info(self): None """ - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3645,7 +3002,7 @@ def plotAttitudeData(self): # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3710,7 +3067,7 @@ def plotFlightPathAngleData(self): # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3762,7 +3119,7 @@ def plotAngularKinematicsData(self): # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3840,16 +3197,10 @@ def plotTrajectoryForceData(self): None """ - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3859,13 +3210,13 @@ def plotTrajectoryForceData(self): ax1 = plt.subplot(211) ax1.plot( - self.railButton1NormalForce[:outOfRailTimeIndex, 0], - self.railButton1NormalForce[:outOfRailTimeIndex, 1], + self.railButton1NormalForce[: self.outOfRailTimeIndex, 0], + self.railButton1NormalForce[: self.outOfRailTimeIndex, 1], label="Upper Rail Button", ) ax1.plot( - self.railButton2NormalForce[:outOfRailTimeIndex, 0], - self.railButton2NormalForce[:outOfRailTimeIndex, 1], + self.railButton2NormalForce[: self.outOfRailTimeIndex, 0], + self.railButton2NormalForce[: self.outOfRailTimeIndex, 1], label="Lower Rail Button", ) ax1.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) @@ -3877,13 +3228,13 @@ def plotTrajectoryForceData(self): ax2 = plt.subplot(212) ax2.plot( - self.railButton1ShearForce[:outOfRailTimeIndex, 0], - self.railButton1ShearForce[:outOfRailTimeIndex, 1], + self.railButton1ShearForce[: self.outOfRailTimeIndex, 0], + self.railButton1ShearForce[: self.outOfRailTimeIndex, 1], label="Upper Rail Button", ) ax2.plot( - self.railButton2ShearForce[:outOfRailTimeIndex, 0], - self.railButton2ShearForce[:outOfRailTimeIndex, 1], + self.railButton2ShearForce[: self.outOfRailTimeIndex, 0], + self.railButton2ShearForce[: self.outOfRailTimeIndex, 1], label="Lower Rail Button", ) ax2.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) @@ -3963,17 +3314,10 @@ def plotEnergyData(self): ------- None """ - - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -4059,12 +3403,6 @@ def plotFluidMechanicsData(self): None """ - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Trajectory Fluid Mechanics Plots fig10 = plt.figure(figsize=(9, 12)) @@ -4126,8 +3464,8 @@ def calculateFinFlutterAnalysis(self, finThickness, shearModulus): NACA Technical Paper 4197. Be careful, these results are only estimates of a real problem and may not be useful for fins made from non-isotropic materials. These results - should not be used as a way to fully prove the safety of any rocket’s fins. - IMPORTANT: This function works if only a single set of fins is added + should not be used as a way to fully prove the safety of any rocket's fins. + IMPORTANT: This function works if only a single set of fins is added. Parameters ---------- @@ -4442,8 +3780,9 @@ def exportData(self, fileName, *variables, timeStep=None): '<, >, :, ", /, \\, | ?, *' in Windows. variables : strings, optional Names of the data variables which shall be exported. Must be Flight - classes attribute which are an instance of the Function class. Usage - example: TestFlight.exportData('test.csv', 'z', 'angleOfAttack', 'machNumber'). + class attributes which are instances of the Function class. Usage + example: TestFlight.exportData('test.csv', 'z', 'angleOfAttack', + 'machNumber'). timeStep : float, optional Time step desired for the data. If None, all integration time steps will be exported. Otherwise, linear interpolation is carried out to @@ -4491,8 +3830,7 @@ def exportData(self, fileName, *variables, timeStep=None): ] if timeStep is None: - # Get time from any Function, should all be the same - timePoints = self.x[:, 0] + timePoints = self.time else: timePoints = np.arange(self.tInitial, self.tFinal, timeStep) @@ -4503,13 +3841,11 @@ def exportData(self, fileName, *variables, timeStep=None): for variable in variables: if variable in self.__dict__.keys(): variableFunction = self.__dict__[variable] - # Deal with variables that are not Flight attributes (e.g. 'angleOfAttack') - # Usually, these are properties of the Flight class + # Deal with decorated Flight methods else: try: obj = getattr(self.__class__, variable) - if isinstance(obj, property) or isinstance(obj, cached_property): - variableFunction = obj.__get__(self, self.__class__) + variableFunction = obj.__get__(self, self.__class__) except AttributeError: raise AttributeError( "Variable '{}' not found in Flight class".format(variable) @@ -4570,8 +3906,7 @@ def exportKML( """ # Define time points vector if timeStep is None: - # Get time from any Function, should all be the same - timePoints = self.z[:, 0] + timePoints = self.time else: timePoints = np.arange(self.tInitial, self.tFinal + timeStep, timeStep) # Open kml file with simplekml library From 9ea7d9c07cc833f0d7eaceb7bb335463d6cc35bc Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Sun, 13 Nov 2022 16:47:29 -0300 Subject: [PATCH 03/13] Fix: solve rail buttons force bug introduced in last commit --- rocketpy/Flight.py | 41 +++++++++++++++++++++++++---------------- 1 file changed, 25 insertions(+), 16 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index e32a58698..55f8e8711 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -2197,7 +2197,7 @@ def maxRailButton1NormalForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1NormalForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton1NormalForce) @cached_property def maxRailButton1ShearForce(self): @@ -2209,7 +2209,7 @@ def maxRailButton1ShearForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1ShearForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton1ShearForce) @cached_property def maxRailButton2NormalForce(self): @@ -2221,7 +2221,7 @@ def maxRailButton2NormalForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2NormalForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton2NormalForce) @cached_property def maxRailButton2ShearForce(self): @@ -2233,7 +2233,7 @@ def maxRailButton2ShearForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2ShearForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton2ShearForce) @funcify_method( "Time (s)", "Horizontal Distance to Launch Point (m)", "spline", "constant" @@ -2483,6 +2483,10 @@ def calculate_rail_button_forces(self): "Trying to calculate rail button forces without rail buttons defined." ) return 0, 0, 0, 0 + if self.outOfRailTimeIndex == 0: + # No rail phase, no rail button forces + nullForce = 0 * self.R1 + return nullForce, nullForce, nullForce, nullForce # Distance from Rail Button 1 (upper) to CM D1 = self.rocket.railButtons.distanceToCM[0] @@ -2497,11 +2501,16 @@ def calculate_rail_button_forces(self): F22 = (self.R2 * D1 - self.M1) / (D1 + D2) F22.setOutputs("Lower button force direction 2 (m)") + model = Function( + F11.source[: self.outOfRailTimeIndex + 1, :], + interpolation=F11.__interpolation__, + ) + # Limit force calculation to when rocket is in rail - F11 = F11[: self.outOfRailTimeIndex, :] - F12 = F12[: self.outOfRailTimeIndex, :] - F21 = F21[: self.outOfRailTimeIndex, :] - F22 = F22[: self.outOfRailTimeIndex, :] + F11.setDiscreteBasedOnModel(model) + F12.setDiscreteBasedOnModel(model) + F21.setDiscreteBasedOnModel(model) + F22.setDiscreteBasedOnModel(model) return F11, F12, F21, F22 @@ -3210,13 +3219,13 @@ def plotTrajectoryForceData(self): ax1 = plt.subplot(211) ax1.plot( - self.railButton1NormalForce[: self.outOfRailTimeIndex, 0], - self.railButton1NormalForce[: self.outOfRailTimeIndex, 1], + self.railButton1NormalForce[:, 0], + self.railButton1NormalForce[:, 1], label="Upper Rail Button", ) ax1.plot( - self.railButton2NormalForce[: self.outOfRailTimeIndex, 0], - self.railButton2NormalForce[: self.outOfRailTimeIndex, 1], + self.railButton2NormalForce[:, 0], + self.railButton2NormalForce[:, 1], label="Lower Rail Button", ) ax1.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) @@ -3228,13 +3237,13 @@ def plotTrajectoryForceData(self): ax2 = plt.subplot(212) ax2.plot( - self.railButton1ShearForce[: self.outOfRailTimeIndex, 0], - self.railButton1ShearForce[: self.outOfRailTimeIndex, 1], + self.railButton1ShearForce[:, 0], + self.railButton1ShearForce[:, 1], label="Upper Rail Button", ) ax2.plot( - self.railButton2ShearForce[: self.outOfRailTimeIndex, 0], - self.railButton2ShearForce[: self.outOfRailTimeIndex, 1], + self.railButton2ShearForce[:, 0], + self.railButton2ShearForce[:, 1], label="Lower Rail Button", ) ax2.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) From 8689d6cf4d3d0f58e8d12c6bf7aac96c73d7b654 Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Sun, 13 Nov 2022 15:01:22 -0300 Subject: [PATCH 04/13] MAINT: remove unnecessary private attributes --- rocketpy/Flight.py | 81 +--------------------------------------------- 1 file changed, 1 insertion(+), 80 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index cb908dcab..f5cd12d6b 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -632,10 +632,6 @@ def __init__( self.impactState = np.array([0]) self.parachuteEvents = [] self.postProcessed = False - self._drift = Function(0) - self._bearing = Function(0) - self._latitude = Function(0) - self._longitude = Function(0) # Initialize solver monitors self.functionEvaluations = [] self.functionEvaluationsPerTimeStep = [] @@ -1116,83 +1112,8 @@ def __init__( def __init_post_process_variables(self): """Initialize post-process variables.""" - # Initialize all variables created during Flight.postProcess() + # Initialize all variables calculated after initialization. # Important to do so that MATLAB® can access them - self._windVelocityX = Function(0) - self._windVelocityY = Function(0) - self._density = Function(0) - self._pressure = Function(0) - self._dynamicViscosity = Function(0) - self._speedOfSound = Function(0) - self._ax = Function(0) - self._ay = Function(0) - self._az = Function(0) - self._alpha1 = Function(0) - self._alpha2 = Function(0) - self._alpha3 = Function(0) - self._speed = Function(0) - self._maxSpeed = 0 - self._maxSpeedTime = 0 - self._horizontalSpeed = Function(0) - self._Acceleration = Function(0) - self._maxAcceleration = 0 - self._maxAccelerationTime = 0 - self._pathAngle = Function(0) - self._attitudeVectorX = Function(0) - self._attitudeVectorY = Function(0) - self._attitudeVectorZ = Function(0) - self._attitudeAngle = Function(0) - self._lateralAttitudeAngle = Function(0) - self._phi = Function(0) - self._theta = Function(0) - self._psi = Function(0) - self._R1 = Function(0) - self._R2 = Function(0) - self._R3 = Function(0) - self._M1 = Function(0) - self._M2 = Function(0) - self._M3 = Function(0) - self._aerodynamicLift = Function(0) - self._aerodynamicDrag = Function(0) - self._aerodynamicBendingMoment = Function(0) - self._aerodynamicSpinMoment = Function(0) - self._railButton1NormalForce = Function(0) - self._maxRailButton1NormalForce = 0 - self._railButton1ShearForce = Function(0) - self._maxRailButton1ShearForce = 0 - self._railButton2NormalForce = Function(0) - self._maxRailButton2NormalForce = 0 - self._railButton2ShearForce = Function(0) - self._maxRailButton2ShearForce = 0 - self._rotationalEnergy = Function(0) - self._translationalEnergy = Function(0) - self._kineticEnergy = Function(0) - self._potentialEnergy = Function(0) - self._totalEnergy = Function(0) - self._thrustPower = Function(0) - self._dragPower = Function(0) - self._attitudeFrequencyResponse = Function(0) - self._omega1FrequencyResponse = Function(0) - self._omega2FrequencyResponse = Function(0) - self._omega3FrequencyResponse = Function(0) - self._streamVelocityX = Function(0) - self._streamVelocityY = Function(0) - self._streamVelocityZ = Function(0) - self._freestreamSpeed = Function(0) - self._apogeeFreestreamSpeed = 0 - self._MachNumber = Function(0) - self._maxMachNumber = 0 - self._maxMachNumberTime = 0 - self._ReynoldsNumber = Function(0) - self._maxReynoldsNumber = 0 - self._maxReynoldsNumberTime = 0 - self._dynamicPressure = Function(0) - self._maxDynamicPressure = 0 - self._maxDynamicPressureTime = 0 - self._totalPressure = Function(0) - self._maxTotalPressure = 0 - self._maxTotalPressureTime = 0 - self._angleOfAttack = Function(0) self.flutterMachNumber = Function(0) self.difference = Function(0) self.safetyFactor = Function(0) From b332e8a1042c881d30111c05fe93619130062cae Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Sun, 13 Nov 2022 15:56:24 -0300 Subject: [PATCH 05/13] ENH: add funcify_method decorator, missing docs and create helper attributes outOfRailTimeIndex, solutionArray and time --- rocketpy/Flight.py | 1504 +++++++++++++------------------------------- 1 file changed, 420 insertions(+), 1084 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index f5cd12d6b..711002778 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -9,13 +9,15 @@ import math import time import warnings +from copy import deepcopy +from functools import cached_property import matplotlib.pyplot as plt import numpy as np import simplekml from scipy import integrate -from .Function import Function +from .Function import Function, funcify_method try: from functools import cached_property @@ -619,6 +621,7 @@ def __init__( self.__init_post_process_variables() # Initialize solution monitors self.outOfRailTime = 0 + self.outOfRailTimeIndex = 0 self.outOfRailState = np.array([0]) self.outOfRailVelocity = 0 self.apogeeState = np.array([0]) @@ -675,6 +678,7 @@ def __init__( # TODO: Check if rocket is actually out of rail. Otherwise, start at rail self.outOfRailState = self.initialSolution[1:] self.outOfRailTime = self.initialSolution[0] + self.outOfRailTimeIndex = 0 self.initialDerivative = self.uDot self.tInitial = self.initialSolution[0] @@ -897,12 +901,17 @@ def __init__( raise ValueError( "Multiple roots found when solving for rail exit time." ) + elif len(valid_t_root) == 0: + raise ValueError( + "No valid roots found when solving for rail exit time." + ) # Determine final state when upper button is going out of rail self.t = valid_t_root[0] + self.solution[-2][0] interpolator = phase.solver.dense_output() self.ySol = interpolator(self.t) self.solution[-1] = [self.t, *self.ySol] self.outOfRailTime = self.t + self.outOfRailTimeIndex = len(self.solution) - 1 self.outOfRailState = self.ySol self.outOfRailVelocity = ( self.ySol[3] ** 2 + self.ySol[4] ** 2 + self.ySol[5] ** 2 @@ -1530,709 +1539,267 @@ def uDotParachute(self, t, u, postProcessing=False): return [vx, vy, vz, ax, ay, az, 0, 0, 0, 0, 0, 0, 0] - # Process first type of outputs - state vector - # Transform solution array into Functions - @cached_property - def x(self, interpolation="spline", extrapolation="natural"): - """Rocket x position as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - x: Function - Rocket x position as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 1]], - "Time (s)", - "X (m)", - interpolation, - extrapolation, - ) - - @cached_property - def y(self, interpolation="spline", extrapolation="natural"): - """ocket y position as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - y: Function - Rocket y position as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 2]], - "Time (s)", - "Y (m)", - interpolation, - extrapolation, - ) - @cached_property - def z(self, interpolation="spline", extrapolation="natural"): - """Rocket z position as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - z: Function - Rocket z position as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 3]], - "Time (s)", - "Z (m)", - interpolation, - extrapolation, - ) + def solutionArray(self): + """Returns solution array of the rocket flight.""" + return np.array(self.solution) @cached_property - def vx(self, interpolation="spline", extrapolation="natural"): - """Rocket x velocity as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - vx: Function - Rocket x velocity as a function of time. - """ - - return Function( - np.array(self.solution)[:, [0, 4]], - "Time (s)", - "Vx (m/s)", - interpolation, - extrapolation, - ) - - @cached_property - def vy(self, interpolation="spline", extrapolation="natural"): - """Rocket y velocity as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. + def time(self): + """Returns time array from solution.""" + return self.solutionArray[:, 0] - Returns - ------- - vy: Function - Rocket y velocity as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 5]], - "Time (s)", - "Vy (m/s)", - interpolation, - extrapolation, - ) - - @cached_property - def vz(self, interpolation="spline", extrapolation="natural"): - """Rocket z velocity as a function of time. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps - interpolation, 'constant', which returns the value of the function - at the edge of the interval, and 'zero', which returns zero for all - points outside of source range. Default is 'natural'. - - Returns - ------- - vz: Function - Rocket z velocity as a function of time. - """ - return Function( - np.array(self.solution)[:, [0, 6]], - "Time (s)", - "Vz (m/s)", - interpolation, - extrapolation, - ) - - @cached_property - def e0(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 7]], - "Time (s)", - "e0", - interpolation, - extrapolation, - ) - - @cached_property - def e1(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 8]], - "Time (s)", - "e1", - interpolation, - extrapolation, - ) - - @cached_property - def e2(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 9]], - "Time (s)", - "e2", - interpolation, - extrapolation, - ) - - @cached_property - def e3(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 10]], - "Time (s)", - "e3", - interpolation, - extrapolation, - ) - - @cached_property - def w1(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 11]], - "Time (s)", - "ω1 (rad/s)", - interpolation, - extrapolation, - ) - - @cached_property - def w2(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 12]], - "Time (s)", - "ω2 (rad/s)", - interpolation, - extrapolation, - ) - - @cached_property - def w3(self, interpolation="spline", extrapolation="natural"): - return Function( - np.array(self.solution)[:, [0, 13]], - "Time (s)", - "ω3 (rad/s)", - interpolation, - extrapolation, - ) + # Process first type of outputs - state vector + # Transform solution array into Functions + @funcify_method("Time (s)", "X (m)", "spline", "constant") + def x(self): + """Rocket x position as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 1]] + + @funcify_method("Time (s)", "Y (m)", "spline", "constant") + def y(self): + """Rocket y position as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 2]] + + @funcify_method("Time (s)", "Z (m)", "spline", "constant") + def z(self): + """Rocket z position as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 3]] + + @funcify_method("Time (s)", "Vx (m/s)", "spline", "constant") + def vx(self): + """Rocket x velocity as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 4]] + + @funcify_method("Time (s)", "Vy (m/s)", "spline", "constant") + def vy(self): + """Rocket y velocity as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 5]] + + @funcify_method("Time (s)", "Vz (m/s)", "spline", "constant") + def vz(self): + """Rocket z velocity as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 6]] + + @funcify_method("Time (s)", "e0", "spline", "constant") + def e0(self): + """Rocket quaternion e0 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 7]] + + @funcify_method("Time (s)", "e1", "spline", "constant") + def e1(self): + """Rocket quaternion e1 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 8]] + + @funcify_method("Time (s)", "e2", "spline", "constant") + def e2(self): + """Rocket quaternion e2 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 9]] + + @funcify_method("Time (s)", "e3", "spline", "constant") + def e3(self): + """Rocket quaternion e3 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 10]] + + @funcify_method("Time (s)", "ω1 (rad/s)", "spline", "constant") + def w1(self): + """Rocket angular velocity ω1 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 11]] + + @funcify_method("Time (s)", "ω2 (rad/s)", "spline", "constant") + def w2(self): + """Rocket angular velocity ω2 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 12]] + + @funcify_method("Time (s)", "ω3 (rad/s)", "spline", "constant") + def w3(self): + """Rocket angular velocity ω3 as a rocketpy.Function of time.""" + return self.solutionArray[:, [0, 13]] # Process second type of outputs - accelerations components - @cached_property - def ax(self, interpolation="spline", extrapolation="natural"): - ax = self.retrieve_acceleration_arrays[0] - # Convert accelerations to functions - ax = Function(ax, "Time (s)", "Ax (m/s2)", interpolation, extrapolation) - return ax - - @cached_property - def ay(self, interpolation="spline", extrapolation="natural"): - ay = self.retrieve_acceleration_arrays[1] - # Convert accelerations to functions - ay = Function(ay, "Time (s)", "Ay (m/s2)", interpolation, extrapolation) - return ay - - @cached_property - def az(self, interpolation="spline", extrapolation="natural"): - az = self.retrieve_acceleration_arrays[2] - # Convert accelerations to functions - az = Function(az, "Time (s)", "Az (m/s2)", interpolation, extrapolation) - return az - - @cached_property - def alpha1(self, interpolation="spline", extrapolation="natural"): - alpha1 = self.retrieve_acceleration_arrays[3] - # Convert accelerations to functions - alpha1 = Function( - alpha1, "Time (s)", "α1 (rad/s2)", interpolation, extrapolation - ) - return alpha1 - - @cached_property - def alpha2(self, interpolation="spline", extrapolation="natural"): - alpha2 = self.retrieve_acceleration_arrays[4] - # Convert accelerations to functions - alpha2 = Function( - alpha2, "Time (s)", "α2 (rad/s2)", interpolation, extrapolation - ) - return alpha2 - - @cached_property - def alpha3(self, interpolation="spline", extrapolation="natural"): - alpha3 = self.retrieve_acceleration_arrays[5] - # Convert accelerations to functions - alpha3 = Function( - alpha3, "Time (s)", "α3 (rad/s2)", interpolation, extrapolation - ) - return alpha3 + @funcify_method("Time (s)", "Ax (m/s²)", "spline", "constant") + def ax(self): + """Rocket x acceleration as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[0] + + @funcify_method("Time (s)", "Ay (m/s²)", "spline", "constant") + def ay(self): + """Rocket y acceleration as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[1] + + @funcify_method("Time (s)", "Az (m/s²)", "spline", "constant") + def az(self): + """Rocket z acceleration as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[2] + + @funcify_method("Time (s)", "α1 (rad/s²)", "spline", "constant") + def alpha1(self): + """Rocket angular acceleration α1 as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[3] + + @funcify_method("Time (s)", "α2 (rad/s²)", "spline", "constant") + def alpha2(self): + """Rocket angular acceleration α2 as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[4] + + @funcify_method("Time (s)", "α3 (rad/s²)", "spline", "constant") + def alpha3(self): + """Rocket angular acceleration α3 as a rocketpy.Function of time.""" + return self.retrieve_acceleration_arrays[5] # Process third type of outputs - Temporary values - @cached_property - def R1(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "R1 (N)", "spline", "constant") + def R1(self): """Aerodynamic force along the first axis that is perpendicular to the - rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. + rocket's axis of symmetry as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[0] - Returns - ------- - R1: Function - Aero force along the first axis that is perpendicular to the - rocket's axis of symmetry. - """ - R1 = self.retrieve_temporary_values_arrays[0] - R1 = Function(R1, "Time (s)", "R1 (m)", interpolation, extrapolation) - - return R1 - - @cached_property - def R2(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "R2 (N)", "spline", "constant") + def R2(self): """Aerodynamic force along the second axis that is perpendicular to the - rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - R2: Function - Aero force along the second axis that is perpendicular to the - rocket's axis of symmetry. - """ - R2 = self.retrieve_temporary_values_arrays[1] - R2 = Function(R2, "Time (s)", "R2 (m)", interpolation, extrapolation) + rocket's axis of symmetry as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[1] - return R2 + @funcify_method("Time (s)", "R3 (N)", "spline", "constant") + def R3(self): + """Aerodynamic force along the rocket's axis of symmetry as a rocketpy.Function + of time.""" + return self.retrieve_temporary_values_arrays[2] - @cached_property - def R3(self, interpolation="spline", extrapolation="natural"): - """Aerodynamic force along the rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - R3: Function - Aerodynamic force along the rocket's axis of symmetry. - """ - R3 = self.retrieve_temporary_values_arrays[2] - R3 = Function(R3, "Time (s)", "R3 (m)", interpolation, extrapolation) - - return R3 - - @cached_property - def M1(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "M1 (Nm)", "spline", "constant") + def M1(self): """Aerodynamic bending moment in the same direction as the axis that is - perpendicular to the rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - M1: Function - Aero moment along the first axis that is perpendicular to the - rocket's axis of symmetry. - """ - M1 = self.retrieve_temporary_values_arrays[3] - M1 = Function(M1, "Time (s)", "M1 (Nm)", interpolation, extrapolation) - - return M1 - - @cached_property - def M2(self, interpolation="spline", extrapolation="natural"): - """Aerodynamic moment in the same direction as the second axis that is - perpendicular to the rocket's axis of symmetry. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - M2: Function - Aero moment along the second axis that is perpendicular to the - rocket's axis of symmetry. - """ - M2 = self.retrieve_temporary_values_arrays[4] - M2 = Function(M2, "Time (s)", "M2 (Nm)", interpolation, extrapolation) - - return M2 - - @cached_property - def M3(self, interpolation="spline", extrapolation="natural"): - """Aerodynamic bending in the rocket's axis of symmetry direction. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - M3: Function - Aero moment in the same direction as the rocket's axis of symmetry. - """ - M3 = self.retrieve_temporary_values_arrays[5] - M3 = Function(M3, "Time (s)", "M3 (Nm)", interpolation, extrapolation) - - return M3 - - @cached_property - def pressure(self, interpolation="spline", extrapolation="natural"): - """Air pressure at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - pressure: Function - Air pressure at each time step. - """ - pressure = self.retrieve_temporary_values_arrays[6] - pressure = Function( - pressure, "Time (s)", "Pressure (Pa)", interpolation, extrapolation - ) - - return pressure - - @cached_property - def density(self, interpolation="spline", extrapolation="natural"): - """Air density at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - density: Function - Air density at each time step. - """ - density = self.retrieve_temporary_values_arrays[7] - density = Function( - density, "Time (s)", "Density (kg/m³)", interpolation, extrapolation - ) - - return density - - @cached_property - def dynamicViscosity(self, interpolation="spline", extrapolation="natural"): - """Air dynamic viscosity at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - dynamicViscosity: Function - Air dynamic viscosity at each time step. - """ - dynamicViscosity = self.retrieve_temporary_values_arrays[8] - dynamicViscosity = Function( - dynamicViscosity, - "Time (s)", - "Dynamic Viscosity (Pa s)", - interpolation, - extrapolation, - ) - - return dynamicViscosity - - @cached_property - def speedOfSound(self, interpolation="spline", extrapolation="natural"): - """Speed of sound at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - speedOfSound: Function - Speed of sound at each time step. - """ - speedOfSound = self.retrieve_temporary_values_arrays[9] - speedOfSound = Function( - speedOfSound, - "Time (s)", - "Speed of Sound (m/s)", - interpolation, - extrapolation, - ) - - return speedOfSound - - @cached_property - def windVelocityX(self, interpolation="spline", extrapolation="natural"): - """Wind velocity in the X direction at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. - - Returns - ------- - windVelocityX: Function - Wind velocity in the X direction at each time step. + perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time. """ - windVelocityX = self.retrieve_temporary_values_arrays[10] - windVelocityX = Function( - windVelocityX, - "Time (s)", - "Wind Velocity X (East) (m/s)", - interpolation, - extrapolation, - ) - - return windVelocityX - - @cached_property - def windVelocityY(self, interpolation="spline", extrapolation="natural"): - """Wind velocity in the Y direction at each time step. - - Parameters - ---------- - extrapolation : str, optional - Function extrapolation mode. Options are 'linear', 'polynomial', - 'akima' and 'spline'. Default is 'spline'. - interpolation : str, optional - Function extrapolation mode. Options are 'natural', which keeps interpolation, 'constant', - which returns the value of the function at the edge of the interval, - and 'zero', which returns zero for all points outside of source - range. Default is 'natural'. + return self.retrieve_temporary_values_arrays[3] - Returns - ------- - windVelocityY: Function - Wind velocity in the Y direction at each time step. + @funcify_method("Time (s)", "M2 (Nm)", "spline", "constant") + def M2(self): + """Aerodynamic bending moment in the same direction as the axis that is + perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time. """ - windVelocityY = self.retrieve_temporary_values_arrays[10] - windVelocityY = Function( - windVelocityY, - "Time (s)", - "Wind Velocity Y (North) (m/s)", - interpolation, - extrapolation, - ) - - return windVelocityY + return self.retrieve_temporary_values_arrays[4] + + @funcify_method("Time (s)", "M3 (Nm)", "spline", "constant") + def M3(self): + """Aerodynamic bending moment in the same direction as the rocket's axis of + symmetry as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[5] + + @funcify_method("Time (s)", "Pressure (Pa)", "spline", "constant") + def pressure(self): + """Air pressure felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[6] + + @funcify_method("Time (s)", "Density (kg/m³)", "spline", "constant") + def density(self): + """Air density felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[7] + + @funcify_method("Time (s)", "Dynamic Viscosity (Pa s)", "spline", "constant") + def dynamicViscosity(self): + """Air dynamic viscosity felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[7] + + @funcify_method("Time (s)", "Speed of Sound (m/s)", "spline", "constant") + def speedOfSound(self): + """Speed of sound in the air felt by the rocket as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[9] + + @funcify_method("Time (s)", "Wind Velocity X (East) (m/s)", "spline", "constant") + def windVelocityX(self): + """Wind velocity in the X direction (east) as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[10] + + @funcify_method("Time (s)", "Wind Velocity Y (North) (m/s)", "spline", "constant") + def windVelocityY(self): + """Wind velocity in the y direction (north) as a rocketpy.Function of time.""" + return self.retrieve_temporary_values_arrays[11] # Process fourth type of output - values calculated from previous outputs # Kinematics functions and values # Velocity Magnitude - @cached_property + @funcify_method("Time (s)", "Speed - Velocity Magnitude (m/s)") def speed(self): - speed = (self.vx**2 + self.vy**2 + self.vz**2) ** 0.5 - speed.setOutputs("Speed - Velocity Magnitude (m/s)") - return speed + """Rocket speed, or velocity magnitude, as a rocketpy.Function of time.""" + return (self.vx**2 + self.vy**2 + self.vz**2) ** 0.5 @cached_property def maxSpeedTime(self): + """Time at which the rocket reaches its maximum speed.""" maxSpeedTimeIndex = np.argmax(self.speed[:, 1]) return self.speed[maxSpeedTimeIndex, 0] @cached_property def maxSpeed(self): + """Maximum speed reached by the rocket.""" return self.speed(self.maxSpeedTime) # Accelerations - @cached_property + @funcify_method("Time (s)", "Acceleration Magnitude (m/s²)") def acceleration(self): - acceleration = (self.ax**2 + self.ay**2 + self.az**2) ** 0.5 - acceleration.setOutputs("Acceleration Magnitude (m/s²)") - return acceleration + """Rocket acceleration magnitude as a rocketpy.Function of time.""" + return (self.ax**2 + self.ay**2 + self.az**2) ** 0.5 @cached_property def maxAcceleration(self): + """Maximum acceleration reached by the rocket.""" maxAccelerationTimeIndex = np.argmax(self.acceleration[:, 1]) return self.acceleration[maxAccelerationTimeIndex, 1] @cached_property def maxAccelerationTime(self): + """Time at which the rocket reaches its maximum acceleration.""" maxAccelerationTimeIndex = np.argmax(self.acceleration[:, 1]) return self.acceleration[maxAccelerationTimeIndex, 0] - # Path Angle - @cached_property + @funcify_method("Time (s)", "Horizontal Speed (m/s)") def horizontalSpeed(self): + """Rocket horizontal speed as a rocketpy.Function of time.""" return (self.vx**2 + self.vy**2) ** 0.5 - @cached_property + # Path Angle + @funcify_method("Time (s)", "Path Angle (°)", "spline", "constant") def pathAngle(self): + """Rocket path angle as a rocketpy.Function of time.""" pathAngle = (180 / np.pi) * np.arctan2( self.vz[:, 1], self.horizontalSpeed[:, 1] ) - pathAngle = np.column_stack([self.vz[:, 0], pathAngle]) - return Function(pathAngle, "Time (s)", "Path Angle (°)") + return np.column_stack([self.time, pathAngle]) # Attitude Angle - @cached_property + @funcify_method("Time (s)", "Attitude Vector X Component") def attitudeVectorX(self): + """Rocket attitude vector X component as a rocketpy.Function of time.""" return 2 * (self.e1 * self.e3 + self.e0 * self.e2) # a13 - @cached_property + @funcify_method("Time (s)", "Attitude Vector Y Component") def attitudeVectorY(self): + """Rocket attitude vector Y component as a rocketpy.Function of time.""" return 2 * (self.e2 * self.e3 - self.e0 * self.e1) # a23 - @cached_property + @funcify_method("Time (s)", "Attitude Vector Z Component") def attitudeVectorZ(self): + """Rocket attitude vector Z component as a rocketpy.Function of time.""" return 1 - 2 * (self.e1**2 + self.e2**2) # a33 - @cached_property + @funcify_method("Time (s)", "Attitude Angle (°)") def attitudeAngle(self): + """Rocket attitude angle as a rocketpy.Function of time.""" horizontalAttitudeProj = ( self.attitudeVectorX**2 + self.attitudeVectorY**2 ) ** 0.5 attitudeAngle = (180 / np.pi) * np.arctan2( self.attitudeVectorZ[:, 1], horizontalAttitudeProj[:, 1] ) - attitudeAngle = np.column_stack([self.attitudeVectorZ[:, 0], attitudeAngle]) - return Function(attitudeAngle, "Time (s)", "Attitude Angle (°)") + attitudeAngle = np.column_stack([self.time, attitudeAngle]) + return attitudeAngle # Lateral Attitude Angle - @cached_property + @funcify_method("Time (s)", "Lateral Attitude Angle (°)") def lateralAttitudeAngle(self): + """Rocket lateral attitude angle as a rocketpy.Function of time.""" lateralVectorAngle = (np.pi / 180) * (self.heading - 90) lateralVectorX = np.sin(lateralVectorAngle) lateralVectorY = np.cos(lateralVectorAngle) @@ -2253,204 +1820,178 @@ def lateralAttitudeAngle(self): lateralAttitudeAngle = (180 / np.pi) * np.arctan2( attitudeLateralProj, attitudeLateralPlaneProj ) - lateralAttitudeAngle = np.column_stack( - [self.attitudeVectorZ[:, 0], lateralAttitudeAngle] - ) - return Function(lateralAttitudeAngle, "Time (s)", "Lateral Attitude Angle (°)") + lateralAttitudeAngle = np.column_stack([self.time, lateralAttitudeAngle]) + return lateralAttitudeAngle # Euler Angles - @cached_property + @funcify_method("Time (s)", "Precession Angle - ψ (°)", "spline", "constant") def psi(self): + """Precession angle as a rocketpy.Function of time.""" psi = (180 / np.pi) * ( np.arctan2(self.e3[:, 1], self.e0[:, 1]) + np.arctan2(-self.e2[:, 1], -self.e1[:, 1]) ) # Precession angle - psi = np.column_stack([self.e1[:, 0], psi]) # Precession angle - return Function(psi, "Time (s)", "Precession Angle - ψ (°)") + psi = np.column_stack([self.time, psi]) # Precession angle + return psi - @cached_property + @funcify_method("Time (s)", "Spin Angle - φ (°)", "spline", "constant") def phi(self): + """Spin angle as a rocketpy.Function of time.""" phi = (180 / np.pi) * ( np.arctan2(self.e3[:, 1], self.e0[:, 1]) - np.arctan2(-self.e2[:, 1], -self.e1[:, 1]) ) # Spin angle - phi = np.column_stack([self.e1[:, 0], phi]) # Spin angle - return Function(phi, "Time (s)", "Spin Angle - φ (°)") + phi = np.column_stack([self.time, phi]) # Spin angle + return phi - @cached_property + @funcify_method("Time (s)", "Nutation Angle - θ (°)", "spline", "constant") def theta(self): + """Nutation angle as a rocketpy.Function of time.""" theta = ( (180 / np.pi) * 2 * np.arcsin(-((self.e1[:, 1] ** 2 + self.e2[:, 1] ** 2) ** 0.5)) ) # Nutation angle - theta = np.column_stack([self.e1[:, 0], theta]) # Nutation angle - return Function(theta, "Time (s)", "Nutation Angle - θ (°)") + theta = np.column_stack([self.time, theta]) # Nutation angle + return theta # Fluid Mechanics variables # Freestream Velocity - @cached_property - def streamVelocityX(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "Freestream Velocity X (m/s)", "spline", "constant") + def streamVelocityX(self): + """Freestream velocity X component as a rocketpy.Function of time.""" streamVelocityX = np.column_stack( - (self.vx[:, 0], self.windVelocityX[:, 1] - self.vx[:, 1]) - ) - streamVelocityX = Function( - streamVelocityX, - "Time (s)", - "Freestream Velocity X (m/s)", - interpolation, - extrapolation, + (self.time, self.windVelocityX[:, 1] - self.vx[:, 1]) ) return streamVelocityX - @cached_property - def streamVelocityY(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "Freestream Velocity Y (m/s)", "spline", "constant") + def streamVelocityY(self): + """Freestream velocity Y component as a rocketpy.Function of time.""" streamVelocityY = np.column_stack( - (self.vy[:, 0], self.windVelocityY[:, 1] - self.vy[:, 1]) - ) - streamVelocityY = Function( - streamVelocityY, - "Time (s)", - "Freestream Velocity Y (m/s)", - interpolation, - extrapolation, + (self.time, self.windVelocityY[:, 1] - self.vy[:, 1]) ) return streamVelocityY - @cached_property + @funcify_method("Time (s)", "Freestream Velocity Z (m/s)", "spline", "constant") def streamVelocityZ(self, interpolation="spline", extrapolation="natural"): - streamVelocityZ = np.column_stack((self.vz[:, 0], -self.vz[:, 1])) - streamVelocityZ = Function( - streamVelocityZ, - "Time (s)", - "Freestream Velocity Z (m/s)", - interpolation, - extrapolation, - ) + """Freestream velocity Z component as a rocketpy.Function of time.""" + streamVelocityZ = np.column_stack((self.time, -self.vz[:, 1])) return streamVelocityZ - @cached_property + @funcify_method("Time (s)", "Freestream Speed (m/s)", "spline", "constant") def freestreamSpeed(self): + """Freestream speed as a rocketpy.Function of time.""" freestreamSpeed = ( self.streamVelocityX**2 + self.streamVelocityY**2 + self.streamVelocityZ**2 ) ** 0.5 - freestreamSpeed.setInputs("Time (s)") - freestreamSpeed.setOutputs("Freestream Speed (m/s)") - return freestreamSpeed + return freestreamSpeed.source # Apogee Freestream speed @cached_property def apogeeFreestreamSpeed(self): + """Freestream speed at apogee in m/s.""" return self.freestreamSpeed(self.apogeeTime) # Mach Number - @cached_property + @funcify_method("Time (s)", "Mach Number", "spline", "constant") def MachNumber(self): - MachNumber = self.freestreamSpeed / self.speedOfSound - MachNumber.setInputs("Time (s)") - MachNumber.setOutputs("Mach Number") - return MachNumber + """Mach number as a rocketpy.Function of time.""" + return self.freestreamSpeed / self.speedOfSound @cached_property def maxMachNumberTime(self): + """Time of maximum Mach number.""" maxMachNumberTimeIndex = np.argmax(self.MachNumber[:, 1]) return self.MachNumber[maxMachNumberTimeIndex, 0] @cached_property def maxMachNumber(self): + """Maximum Mach number.""" return self.MachNumber(self.maxMachNumberTime) # Reynolds Number - @cached_property + @funcify_method("Time (s)", "Reynolds Number", "spline", "constant") def ReynoldsNumber(self): - ReynoldsNumber = ( - self.density * self.freestreamSpeed / self.dynamicViscosity - ) * (2 * self.rocket.radius) - ReynoldsNumber.setInputs("Time (s)") - ReynoldsNumber.setOutputs("Reynolds Number") - return ReynoldsNumber + """Reynolds number as a rocketpy.Function of time.""" + return (self.density * self.freestreamSpeed / self.dynamicViscosity) * ( + 2 * self.rocket.radius + ) @cached_property def maxReynoldsNumberTime(self): + """Time of maximum Reynolds number.""" maxReynoldsNumberTimeIndex = np.argmax(self.ReynoldsNumber[:, 1]) return self.ReynoldsNumber[maxReynoldsNumberTimeIndex, 0] @cached_property def maxReynoldsNumber(self): + """Maximum Reynolds number.""" return self.ReynoldsNumber(self.maxReynoldsNumberTime) # Dynamic Pressure - @cached_property + @funcify_method("Time (s)", "Dynamic Pressure (Pa)", "spline", "constant") def dynamicPressure(self): - dynamicPressure = 0.5 * self.density * self.freestreamSpeed**2 - dynamicPressure.setInputs("Time (s)") - dynamicPressure.setOutputs("Dynamic Pressure (Pa)") - return dynamicPressure + """Dynamic pressure as a rocketpy.Function of time.""" + return 0.5 * self.density * self.freestreamSpeed**2 @cached_property def maxDynamicPressureTime(self): + """Time of maximum dynamic pressure.""" maxDynamicPressureTimeIndex = np.argmax(self.dynamicPressure[:, 1]) return self.dynamicPressure[maxDynamicPressureTimeIndex, 0] @cached_property def maxDynamicPressure(self): + """Maximum dynamic pressure.""" return self.dynamicPressure(self.maxDynamicPressureTime) # Total Pressure - @cached_property + @funcify_method("Time (s)", "Total Pressure (Pa)", "spline", "constant") def totalPressure(self): - totalPressure = self.pressure * (1 + 0.2 * self.MachNumber**2) ** (3.5) - totalPressure.setInputs("Time (s)") - totalPressure.setOutputs("Total Pressure (Pa)") - return totalPressure + return self.pressure * (1 + 0.2 * self.MachNumber**2) ** (3.5) @cached_property - def maxtotalPressureTime(self): - maxtotalPressureTimeIndex = np.argmax(self.totalPressure[:, 1]) - return self.totalPressure[maxtotalPressureTimeIndex, 0] + def maxTotalPressureTime(self): + """Time of maximum total pressure.""" + maxTotalPressureTimeIndex = np.argmax(self.totalPressure[:, 1]) + return self.totalPressure[maxTotalPressureTimeIndex, 0] @cached_property - def maxtotalPressure(self): - return self.totalPressure(self.maxtotalPressureTime) + def maxTotalPressure(self): + """Maximum total pressure.""" + return self.totalPressure(self.maxTotalPressureTime) # Dynamics functions and variables # Aerodynamic Lift and Drag - @cached_property + @funcify_method("Time (s)", "Aerodynamic Lift Force (N)", "spline", "constant") def aerodynamicLift(self): - aerodynamicLift = (self.R1**2 + self.R2**2) ** 0.5 - aerodynamicLift.setInputs("Time (s)") - aerodynamicLift.setOutputs("Aerodynamic Lift Force (N)") - return aerodynamicLift + """Aerodynamic lift force as a rocketpy.Function of time.""" + return (self.R1**2 + self.R2**2) ** 0.5 - @cached_property + @funcify_method("Time (s)", "Aerodynamic Drag Force (N)", "spline", "constant") def aerodynamicDrag(self): - aerodynamicDrag = -1 * self.R3 - aerodynamicDrag.setInputs("Time (s)") - aerodynamicDrag.setOutputs("Aerodynamic Drag Force (N)") - return aerodynamicDrag + """Aerodynamic drag force as a rocketpy.Function of time.""" + return -1 * self.R3 - @cached_property + @funcify_method("Time (s)", "Aerodynamic Bending Moment (Nm)", "spline", "constant") def aerodynamicBendingMoment(self): + """Aerodynamic bending moment as a rocketpy.Function of time.""" + return (self.M1**2 + self.M2**2) ** 0.5 - aerodynamicBendingMoment = (self.M1**2 + self.M2**2) ** 0.5 - aerodynamicBendingMoment.setInputs("Time (s)") - aerodynamicBendingMoment.setOutputs("Aerodynamic Bending Moment (N m)") - return aerodynamicBendingMoment - - @cached_property + @funcify_method("Time (s)", "Aerodynamic Spin Moment (Nm)", "spline", "constant") def aerodynamicSpinMoment(self): - aerodynamicSpinMoment = self.M3 - aerodynamicSpinMoment.setInputs("Time (s)") - aerodynamicSpinMoment.setOutputs("Aerodynamic Spin Moment (N m)") - return aerodynamicSpinMoment + """Aerodynamic spin moment as a rocketpy.Function of time.""" + return self.M3 # Energy # Kinetic Energy - @cached_property + @funcify_method("Time (s)", "Rotational Kinetic Energy (J)", "spline", "constant") def rotationalEnergy(self): + """Rotational kinetic energy as a rocketpy.Function of time.""" b = -self.rocket.distanceRocketPropellant mu = self.rocket.reducedMass Rz = self.rocket.inertiaZ @@ -2458,103 +1999,80 @@ def rotationalEnergy(self): Tz = self.rocket.motor.inertiaZ Ti = self.rocket.motor.inertiaI I1, I2, I3 = (Ri + Ti + mu * b**2), (Ri + Ti + mu * b**2), (Rz + Tz) - # Redefine I1, I2 and I3 grid - grid = self.vx[:, 0] - I1 = Function(np.column_stack([grid, I1(grid)]), "Time (s)") - I2 = Function(np.column_stack([grid, I2(grid)]), "Time (s)") - I3 = Function(np.column_stack([grid, I3(grid)]), "Time (s)") + # Redefine I1, I2 and I3 time grid to allow for efficient Function algebra + I1.setDiscreteBasedOnModel(self.w1) + I2.setDiscreteBasedOnModel(self.w1) + I3.setDiscreteBasedOnModel(self.w1) rotationalEnergy = 0.5 * ( I1 * self.w1**2 + I2 * self.w2**2 + I3 * self.w3**2 ) - rotationalEnergy.setInputs("Time (s)") - rotationalEnergy.setOutputs("Rotational Kinetic Energy (J)") return rotationalEnergy - @cached_property + @funcify_method( + "Time (s)", "Translational Kinetic Energy (J)", "spline", "constant" + ) def translationalEnergy(self): - # Redefine total mass grid - totalMass = self.rocket.totalMass - grid = self.vx[:, 0] - totalMass = Function(np.column_stack([grid, totalMass(grid)]), "Time (s)") + """Translational kinetic energy as a rocketpy.Function of time.""" + # Redefine totalMass time grid to allow for efficient Function algebra + totalMass = deepcopy(self.rocket.totalMass) + totalMass.setDiscreteBasedOnModel(self.vz) translationalEnergy = ( 0.5 * totalMass * (self.vx**2 + self.vy**2 + self.vz**2) ) - translationalEnergy.setInputs("Time (s)") - translationalEnergy.setOutputs("Translational Kinetic Energy (J)") return translationalEnergy - @cached_property + @funcify_method("Time (s)", "Kinetic Energy (J)", "spline", "constant") def kineticEnergy(self): - kineticEnergy = self.rotationalEnergy + self.translationalEnergy - kineticEnergy.setInputs("Time (s)") - kineticEnergy.setOutputs("Kinetic Energy (J)") - return kineticEnergy + """Total kinetic energy as a rocketpy.Function of time.""" + return self.rotationalEnergy + self.translationalEnergy # Potential Energy - @cached_property + @funcify_method("Time (s)", "Potential Energy (J)", "spline", "constant") def potentialEnergy(self): - # Redefine total mass grid - totalMass = self.rocket.totalMass - grid = self.vx[:, 0] - totalMass = Function(np.column_stack([grid, totalMass(grid)]), "Time (s)") + """Potential energy as a rocketpy.Function of time.""" + # Redefine totalMass time grid to allow for efficient Function algebra + totalMass = deepcopy(self.rocket.totalMass) + totalMass.setDiscreteBasedOnModel(self.z) + # TODO: change calculation method to account for variable gravity potentialEnergy = totalMass * self.env.g * self.z - potentialEnergy.setInputs("Time (s)") return potentialEnergy # Total Mechanical Energy - @cached_property + @funcify_method("Time (s)", "Mechanical Energy (J)", "spline", "constant") def totalEnergy(self): - totalEnergy = self.kineticEnergy + self.potentialEnergy - totalEnergy.setInputs("Time (s)") - totalEnergy.setOutputs("Total Mechanical Energy (J)") - return totalEnergy + """Total mechanical energy as a rocketpy.Function of time.""" + return self.kineticEnergy + self.potentialEnergy # Thrust Power - @cached_property + @funcify_method("Time (s)", "Thrust Power (W)", "spline", "constant") def thrustPower(self): - grid = self.vx[:, 0] - # Redefine thrust grid - thrust = Function( - np.column_stack([grid, self.rocket.motor.thrust(grid)]), "Time (s)" - ) + """Thrust power as a rocketpy.Function of time.""" + thrust = deepcopy(self.rocket.motor.thrust) + thrust = thrust.setDiscreteBasedOnModel(self.speed) thrustPower = thrust * self.speed - thrustPower.setOutputs("Thrust Power (W)") return thrustPower # Drag Power - @cached_property + @funcify_method("Time (s)", "Drag Power (W)", "spline", "constant") def dragPower(self): + """Drag power as a rocketpy.Function of time.""" dragPower = self.R3 * self.speed dragPower.setOutputs("Drag Power (W)") return dragPower # Angle of Attack - @cached_property - def angleOfAttack(self, interpolation="spline", extrapolation="natural"): + @funcify_method("Time (s)", "Angle of Attack (°)", "spline", "constant") + def angleOfAttack(self): """Angle of attack of the rocket with respect to the freestream - velocity vector. - - Parameters - ---------- - interpolation : str, optional - Interpolation method, by default "spline" - extrapolation : str, optional - Extrapolation method, by default "natural" - - Returns - ------- - angleOfAttack: Function - Angle of attack of the rocket with respect to the freestream - velocity vector. - """ + velocity vector.""" dotProduct = [ -self.attitudeVectorX.getValueOpt(i) * self.streamVelocityX.getValueOpt(i) - self.attitudeVectorY.getValueOpt(i) * self.streamVelocityY.getValueOpt(i) - self.attitudeVectorZ.getValueOpt(i) * self.streamVelocityZ.getValueOpt(i) - for i in self.x[:, 0] + for i in self.time ] # Define freestream speed list - freestreamSpeed = [self.freestreamSpeed(i) for i in self.x[:, 0]] + freestreamSpeed = [self.freestreamSpeed(i) for i in self.time] freestreamSpeed = np.nan_to_num(freestreamSpeed) # Normalize dot product @@ -2566,146 +2084,85 @@ def angleOfAttack(self, interpolation="spline", extrapolation="natural"): # Calculate angle of attack and convert to degrees angleOfAttack = np.rad2deg(np.arccos(dotProductNormalized)) - angleOfAttack = np.column_stack([self.x[:, 0], angleOfAttack]) - - # Convert to Function - angleOfAttack = Function( - angleOfAttack, - "Time (s)", - "Angle Of Attack (°)", - interpolation, - extrapolation, - ) + angleOfAttack = np.column_stack([self.time, angleOfAttack]) return angleOfAttack # Stability and Control variables + def _frequency_analysis(self, func, samplingFrequency=1000): + """Frequency analysis of any function of time between out of rail time + and final simulation time. + + Parameters + ---------- + func : rocketpy.Function, callable + Function of time to be analyzed. + + Returns + ------- + FourierFrequencies : numpy.ndarray + Array of frequencies. + FourierAmplitude : numpy.ndarray + Array of amplitude values. + """ + samplingTimeStep = 1.0 / samplingFrequency + samplingRange = np.arange(self.outOfRailTime, self.tFinal, samplingTimeStep) + sampledPoints = func(samplingRange) + sampledPoints -= np.mean(sampledPoints) + numberOfSamples = len(sampledPoints) + FourierAmplitude = np.abs(np.fft.fft(sampledPoints) / numberOfSamples) + FourierFrequencies = np.fft.fftfreq(numberOfSamples, samplingTimeStep) + return FourierFrequencies, FourierAmplitude + # Angular velocities frequency response - Fourier Analysis - @cached_property + @funcify_method( + "Frequency (Hz)", "ω1 Angle Fourier Amplitude", "spline", "constant" + ) def omega1FrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.w1(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - omega1FrequencyResponse = np.column_stack([frq, abs(Y)]) - omega1FrequencyResponse = Function( - omega1FrequencyResponse, "Frequency (Hz)", "Omega 1 Angle Fourier Amplitude" - ) - return omega1FrequencyResponse + """Angular velocity 1 frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.w1)) - @cached_property + @funcify_method( + "Frequency (Hz)", "ω2 Angle Fourier Amplitude", "spline", "constant" + ) def omega2FrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.w2(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - omega2FrequencyResponse = np.column_stack([frq, abs(Y)]) - omega2FrequencyResponse = Function( - omega2FrequencyResponse, "Frequency (Hz)", "Omega 2 Angle Fourier Amplitude" - ) - return omega2FrequencyResponse + """Angular velocity 2 frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.w2)) - @cached_property + @funcify_method( + "Frequency (Hz)", "ω3 Angle Fourier Amplitude", "spline", "constant" + ) def omega3FrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.w3(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - omega3FrequencyResponse = np.column_stack([frq, abs(Y)]) - omega3FrequencyResponse = Function( - omega3FrequencyResponse, "Frequency (Hz)", "Omega 3 Angle Fourier Amplitude" - ) - return omega3FrequencyResponse + """Angular velocity 3 frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.w3)) - @cached_property + @funcify_method( + "Frequency (Hz)", "Attitude Angle Fourier Amplitude", "spline", "constant" + ) def attitudeFrequencyResponse(self): - Fs = 100.0 - # sampling rate - Ts = 1.0 / Fs - # sampling interval - t = np.arange(1, self.tFinal, Ts) # time vector - y = self.attitudeAngle(t) - y -= np.mean(y) - n = len(y) # length of the signal - k = np.arange(n) - T = n / Fs - frq = k / T # two sides frequency range - frq = frq[range(n // 2)] # one side frequency range - Y = np.fft.fft(y) / n # fft computing and normalization - Y = Y[range(n // 2)] - attitudeFrequencyResponse = np.column_stack([frq, abs(Y)]) - attitudeFrequencyResponse = Function( - attitudeFrequencyResponse, - "Frequency (Hz)", - "Attitude Angle Fourier Amplitude", - ) - return attitudeFrequencyResponse + """Attitude frequency response as a rocketpy.Function of frequency.""" + return np.column_stack(self._frequency_analysis(self.attitudeAngle)) # Static Margin @cached_property def staticMargin(self): + """Static margin of the rocket.""" return self.rocket.staticMargin # Rail Button Forces - @cached_property + @funcify_method("Time (s)", "Upper Rail Button Normal Force (N)", "spline", "zero") def railButton1NormalForce(self): - """Upper rail button normal force as a function of time - - Returns - ------- - railButton1NormalForce: Function - Upper rail button normal force as a function of time - """ - + """Upper rail button normal force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F11, F12 = self.calculate_rail_button_forces[0:2] else: F11, F12 = self.calculate_rail_button_forces()[0:2] alpha = self.rocket.railButtons.angularPosition * (np.pi / 180) - railButton1NormalForce = F11 * np.cos(alpha) + F12 * np.sin(alpha) - railButton1NormalForce.setOutputs("Upper Rail Button Normal Force (N)") - - return railButton1NormalForce + return F11 * np.cos(alpha) + F12 * np.sin(alpha) - @cached_property + @funcify_method("Time (s)", "Upper Rail Button Shear Force (N)", "spline", "zero") def railButton1ShearForce(self): - """Upper rail button shear force as a function of time - - Returns - ------- - _type_ - _description_ - """ + """Upper rail button shear force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F11, F12 = self.calculate_rail_button_forces[0:2] else: @@ -2713,159 +2170,91 @@ def railButton1ShearForce(self): alpha = self.rocket.railButtons.angularPosition * ( np.pi / 180 ) # Rail buttons angular position - railButton1ShearForce = F11 * -np.sin(alpha) + F12 * np.cos(alpha) - railButton1ShearForce.setOutputs("Upper Rail Button Shear Force (N)") - - return railButton1ShearForce + return F11 * -np.sin(alpha) + F12 * np.cos(alpha) - @cached_property + @funcify_method("Time (s)", "Lower Rail Button Normal Force (N)", "spline", "zero") def railButton2NormalForce(self): - """Lower rail button normal force as a function of time - - Returns - ------- - railButton2NormalForce: Function - Lower rail button normal force as a function of time - """ + """Lower rail button normal force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F21, F22 = self.calculate_rail_button_forces[2:4] else: F21, F22 = self.calculate_rail_button_forces()[2:4] alpha = self.rocket.railButtons.angularPosition * (np.pi / 180) - railButton2NormalForce = F21 * np.cos(alpha) + F22 * np.sin(alpha) - railButton2NormalForce.setOutputs("Lower Rail Button Normal Force (N)") - - return railButton2NormalForce + return F21 * np.cos(alpha) + F22 * np.sin(alpha) - @cached_property + @funcify_method("Time (s)", "Lower Rail Button Shear Force (N)", "spline", "zero") def railButton2ShearForce(self): - """Lower rail button shear force as a function of time - - Returns - ------- - railButton2ShearForce: Function - Lower rail button shear force as a function of time - """ + """Lower rail button shear force as a rocketpy.Function of time.""" if isinstance(self.calculate_rail_button_forces, tuple): F21, F22 = self.calculate_rail_button_forces[2:4] else: F21, F22 = self.calculate_rail_button_forces()[2:4] alpha = self.rocket.railButtons.angularPosition * (np.pi / 180) - railButton2ShearForce = F21 * -np.sin(alpha) + F22 * np.cos(alpha) - railButton2ShearForce.setOutputs("Lower Rail Button Shear Force (N)") - - return railButton2ShearForce + return F21 * -np.sin(alpha) + F22 * np.cos(alpha) @cached_property def maxRailButton1NormalForce(self): - """Maximum upper rail button normal force, in Newtons - - Returns - ------- - maxRailButton1NormalForce: float - Maximum upper rail button normal force, in Newtons - """ + """Maximum upper rail button normal force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1NormalForce[:outOfRailTimeIndex]) + return np.max(self.railButton1NormalForce[: self.outOfRailTimeIndex]) @cached_property def maxRailButton1ShearForce(self): - """Maximum upper rail button shear force, in Newtons - - Returns - ------- - maxRailButton1ShearForce: float - Maximum upper rail button shear force, in Newtons - """ + """Maximum upper rail button shear force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1ShearForce[:outOfRailTimeIndex]) + return np.max(self.railButton1ShearForce[: self.outOfRailTimeIndex]) @cached_property def maxRailButton2NormalForce(self): - """Maximum lower rail button normal force, in Newtons - - Returns - ------- - maxRailButton2NormalForce: float - Maximum lower rail button normal force, in Newtons - """ + """Maximum lower rail button normal force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2NormalForce[:outOfRailTimeIndex]) + return np.max(self.railButton2NormalForce[: self.outOfRailTimeIndex]) @cached_property def maxRailButton2ShearForce(self): - """Maximum lower rail button shear force, in Newtons - - Returns - ------- - maxRailButton2ShearForce: float - Maximum lower rail button shear force, in Newtons - """ + """Maximum lower rail button shear force, in Newtons.""" if isinstance(self.calculate_rail_button_forces, tuple): F11 = self.calculate_rail_button_forces[0] else: F11 = self.calculate_rail_button_forces()[0] - outOfRailTimeIndex = np.searchsorted(F11[:, 0], self.outOfRailTime) - if outOfRailTimeIndex == 0: + if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2ShearForce[:outOfRailTimeIndex]) + return np.max(self.railButton2ShearForce[: self.outOfRailTimeIndex]) - @cached_property + @funcify_method( + "Time (s)", "Horizontal Distance to Launch Point (m)", "spline", "constant" + ) def drift(self): - """Rocket horizontal distance to tha launch point, in meters, at each - time step. - - Returns - ------- - drift: Function - Rocket horizontal distance to tha launch point, in meters, at each - time step. - """ - drift = np.column_stack( - (self.x[:, 0], (self.x[:, 1] ** 2 + self.y[:, 1] ** 2) ** 0.5) - ) - drift = Function( - drift, - "Time (s)", - "Horizontal Distance to Launch Point (m)", + """Rocket horizontal distance to tha launch point, in meters, as a + rocketpy.Function of time.""" + return np.column_stack( + (self.time, (self.x[:, 1] ** 2 + self.y[:, 1] ** 2) ** 0.5) ) - return drift - - @cached_property - def bearing(self, interpolation="spline", extrapolation="natural"): - """Rocket bearing compass, in degrees, at each time step. - - Returns - ------- - bearing: Function - Rocket bearing, in degrees, at each time step. - """ - + @funcify_method("Time (s)", "Bearing (°)", "spline", "constant") + def bearing(self): + """Rocket bearing compass, in degrees, as a rocketpy.Function of time.""" # Get some nicknames - t = self.x[:, 0] + t = self.time x, y = self.x[:, 1], self.y[:, 1] bearing = [] for i in range(len(t)): @@ -2895,26 +2284,11 @@ def bearing(self, interpolation="spline", extrapolation="natural"): bearing = np.rad2deg(bearing) bearing = np.column_stack((t, bearing)) - print(bearing) - bearing = Function( - bearing, - "Time (s)", - "Bearing (deg)", - interpolation, - extrapolation, - ) - return bearing - @cached_property + @funcify_method("Time (s)", "Latitude (°)", "linear") def latitude(self): - """Rocket latitude coordinate, in degrees, at each time step. - - Returns - ------- - latitude: Function - Rocket latitude coordinate, in degrees, at each time step. - """ + """Rocket latitude coordinate, in degrees, as a rocketpy.Function of time.""" lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians # Applies the haversine equation to find final lat/lon coordinates @@ -2926,20 +2300,11 @@ def latitude(self): * np.cos(np.deg2rad(self.bearing[:, 1])) ) ) - latitude = np.column_stack((self.x[:, 0], latitude)) - latitude = Function(latitude, "Time (s)", "Latitude (°)", "linear") + return np.column_stack((self.time, latitude)) - return latitude - - @cached_property + @funcify_method("Time (s)", "Longitude (°)", "linear") def longitude(self): - """Rocket longitude coordinate, in degrees, at each time step. - - Returns - ------- - longitude: Function - Rocket longitude coordinate, in degrees, at each time step. - """ + """Rocket longitude coordinate, in degrees, as a rocketpy.Function of time.""" lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians lon1 = np.deg2rad(self.env.lon) # Launch lon point converted to radians @@ -2955,10 +2320,7 @@ def longitude(self): ) ) - longitude = np.column_stack((self.x[:, 0], longitude)) - longitude = Function(longitude, "Time (s)", "Longitude (°)", "linear") - - return longitude + return np.column_stack((self.time, longitude)) @cached_property def retrieve_acceleration_arrays(self): @@ -3108,7 +2470,7 @@ def retrieve_temporary_values_arrays(self): @cached_property def calculate_rail_button_forces(self): """Calculate the forces applied to the rail buttons while rocket is still - on the launch rail. It will return 0 if none rail buttons are defined. + on the launch rail. It will return 0 if no rail buttons are defined. Returns ------- @@ -3121,7 +2483,6 @@ def calculate_rail_button_forces(self): F22: Function Rail Button 2 force in the 2 direction """ - if self.rocket.railButtons is None: warnings.warn( "Trying to calculate rail button forces without rail buttons defined." @@ -3141,10 +2502,11 @@ def calculate_rail_button_forces(self): F22 = (self.R2 * D1 - self.M1) / (D1 + D2) F22.setOutputs("Lower button force direction 2 (m)") - # F11 = F11[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail - # F12 = F12[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail - # F21 = F21[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail - # F22 = F22[:outOfRailTimeIndex + 1, :] # Limit force calculation to when rocket is in rail + # Limit force calculation to when rocket is in rail + F11 = F11[: self.outOfRailTimeIndex, :] + F12 = F12[: self.outOfRailTimeIndex, :] + F21 = F21[: self.outOfRailTimeIndex, :] + F22 = F22[: self.outOfRailTimeIndex, :] return F11, F12, F21, F22 @@ -3184,6 +2546,9 @@ def postProcess(self, interpolation="spline", extrapolation="natural"): calculation of secondary values such as energy and conversion of lists to Function objects to facilitate plotting. + * This method is deprecated and is only kept here for backwards compatibility. + * All attributes that need to be post processed are computed just in time. + Parameters ---------- None @@ -3192,7 +2557,6 @@ def postProcess(self, interpolation="spline", extrapolation="natural"): ------ None """ - # Register post processing self.postProcessed = True @@ -3210,16 +2574,10 @@ def info(self): None """ - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3649,7 +3007,7 @@ def plotAttitudeData(self): # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3714,7 +3072,7 @@ def plotFlightPathAngleData(self): # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3766,7 +3124,7 @@ def plotAngularKinematicsData(self): # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3844,16 +3202,10 @@ def plotTrajectoryForceData(self): None """ - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -3863,13 +3215,13 @@ def plotTrajectoryForceData(self): ax1 = plt.subplot(211) ax1.plot( - self.railButton1NormalForce[:outOfRailTimeIndex, 0], - self.railButton1NormalForce[:outOfRailTimeIndex, 1], + self.railButton1NormalForce[: self.outOfRailTimeIndex, 0], + self.railButton1NormalForce[: self.outOfRailTimeIndex, 1], label="Upper Rail Button", ) ax1.plot( - self.railButton2NormalForce[:outOfRailTimeIndex, 0], - self.railButton2NormalForce[:outOfRailTimeIndex, 1], + self.railButton2NormalForce[: self.outOfRailTimeIndex, 0], + self.railButton2NormalForce[: self.outOfRailTimeIndex, 1], label="Lower Rail Button", ) ax1.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) @@ -3881,13 +3233,13 @@ def plotTrajectoryForceData(self): ax2 = plt.subplot(212) ax2.plot( - self.railButton1ShearForce[:outOfRailTimeIndex, 0], - self.railButton1ShearForce[:outOfRailTimeIndex, 1], + self.railButton1ShearForce[: self.outOfRailTimeIndex, 0], + self.railButton1ShearForce[: self.outOfRailTimeIndex, 1], label="Upper Rail Button", ) ax2.plot( - self.railButton2ShearForce[:outOfRailTimeIndex, 0], - self.railButton2ShearForce[:outOfRailTimeIndex, 1], + self.railButton2ShearForce[: self.outOfRailTimeIndex, 0], + self.railButton2ShearForce[: self.outOfRailTimeIndex, 1], label="Lower Rail Button", ) ax2.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) @@ -3967,17 +3319,10 @@ def plotEnergyData(self): ------- None """ - - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Get index of time before parachute event if len(self.parachuteEvents) > 0: eventTime = self.parachuteEvents[0][0] + self.parachuteEvents[0][1].lag - eventTimeIndex = np.nonzero(self.x[:, 0] == eventTime)[0][0] + eventTimeIndex = np.nonzero(self.time == eventTime)[0][0] else: eventTime = self.tFinal eventTimeIndex = -1 @@ -4063,12 +3408,6 @@ def plotFluidMechanicsData(self): None """ - # Get index of out of rail time - outOfRailTimeIndexes = np.nonzero(self.x[:, 0] == self.outOfRailTime) - outOfRailTimeIndex = ( - -1 if len(outOfRailTimeIndexes) == 0 else outOfRailTimeIndexes[0][0] - ) - # Trajectory Fluid Mechanics Plots fig10 = plt.figure(figsize=(9, 12)) @@ -4130,8 +3469,8 @@ def calculateFinFlutterAnalysis(self, finThickness, shearModulus): NACA Technical Paper 4197. Be careful, these results are only estimates of a real problem and may not be useful for fins made from non-isotropic materials. These results - should not be used as a way to fully prove the safety of any rocket’s fins. - IMPORTANT: This function works if only a single set of fins is added + should not be used as a way to fully prove the safety of any rocket's fins. + IMPORTANT: This function works if only a single set of fins is added. Parameters ---------- @@ -4446,8 +3785,9 @@ def exportData(self, fileName, *variables, timeStep=None): '<, >, :, ", /, \\, | ?, *' in Windows. variables : strings, optional Names of the data variables which shall be exported. Must be Flight - classes attribute which are an instance of the Function class. Usage - example: TestFlight.exportData('test.csv', 'z', 'angleOfAttack', 'machNumber'). + class attributes which are instances of the Function class. Usage + example: TestFlight.exportData('test.csv', 'z', 'angleOfAttack', + 'machNumber'). timeStep : float, optional Time step desired for the data. If None, all integration time steps will be exported. Otherwise, linear interpolation is carried out to @@ -4495,8 +3835,7 @@ def exportData(self, fileName, *variables, timeStep=None): ] if timeStep is None: - # Get time from any Function, should all be the same - timePoints = self.x[:, 0] + timePoints = self.time else: timePoints = np.arange(self.tInitial, self.tFinal, timeStep) @@ -4507,13 +3846,11 @@ def exportData(self, fileName, *variables, timeStep=None): for variable in variables: if variable in self.__dict__.keys(): variableFunction = self.__dict__[variable] - # Deal with variables that are not Flight attributes (e.g. 'angleOfAttack') - # Usually, these are properties of the Flight class + # Deal with decorated Flight methods else: try: obj = getattr(self.__class__, variable) - if isinstance(obj, property) or isinstance(obj, cached_property): - variableFunction = obj.__get__(self, self.__class__) + variableFunction = obj.__get__(self, self.__class__) except AttributeError: raise AttributeError( "Variable '{}' not found in Flight class".format(variable) @@ -4574,8 +3911,7 @@ def exportKML( """ # Define time points vector if timeStep is None: - # Get time from any Function, should all be the same - timePoints = self.z[:, 0] + timePoints = self.time else: timePoints = np.arange(self.tInitial, self.tFinal + timeStep, timeStep) # Open kml file with simplekml library From b527cce748942dddc0da9aecd9defa2be9c00889 Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Sun, 13 Nov 2022 16:47:29 -0300 Subject: [PATCH 06/13] Fix: solve rail buttons force bug introduced in last commit --- rocketpy/Flight.py | 41 +++++++++++++++++++++++++---------------- 1 file changed, 25 insertions(+), 16 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index 711002778..b8773a044 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -2202,7 +2202,7 @@ def maxRailButton1NormalForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1NormalForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton1NormalForce) @cached_property def maxRailButton1ShearForce(self): @@ -2214,7 +2214,7 @@ def maxRailButton1ShearForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton1ShearForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton1ShearForce) @cached_property def maxRailButton2NormalForce(self): @@ -2226,7 +2226,7 @@ def maxRailButton2NormalForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2NormalForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton2NormalForce) @cached_property def maxRailButton2ShearForce(self): @@ -2238,7 +2238,7 @@ def maxRailButton2ShearForce(self): if self.outOfRailTimeIndex == 0: return 0 else: - return np.max(self.railButton2ShearForce[: self.outOfRailTimeIndex]) + return np.max(self.railButton2ShearForce) @funcify_method( "Time (s)", "Horizontal Distance to Launch Point (m)", "spline", "constant" @@ -2488,6 +2488,10 @@ def calculate_rail_button_forces(self): "Trying to calculate rail button forces without rail buttons defined." ) return 0, 0, 0, 0 + if self.outOfRailTimeIndex == 0: + # No rail phase, no rail button forces + nullForce = 0 * self.R1 + return nullForce, nullForce, nullForce, nullForce # Distance from Rail Button 1 (upper) to CM D1 = self.rocket.railButtons.distanceToCM[0] @@ -2502,11 +2506,16 @@ def calculate_rail_button_forces(self): F22 = (self.R2 * D1 - self.M1) / (D1 + D2) F22.setOutputs("Lower button force direction 2 (m)") + model = Function( + F11.source[: self.outOfRailTimeIndex + 1, :], + interpolation=F11.__interpolation__, + ) + # Limit force calculation to when rocket is in rail - F11 = F11[: self.outOfRailTimeIndex, :] - F12 = F12[: self.outOfRailTimeIndex, :] - F21 = F21[: self.outOfRailTimeIndex, :] - F22 = F22[: self.outOfRailTimeIndex, :] + F11.setDiscreteBasedOnModel(model) + F12.setDiscreteBasedOnModel(model) + F21.setDiscreteBasedOnModel(model) + F22.setDiscreteBasedOnModel(model) return F11, F12, F21, F22 @@ -3215,13 +3224,13 @@ def plotTrajectoryForceData(self): ax1 = plt.subplot(211) ax1.plot( - self.railButton1NormalForce[: self.outOfRailTimeIndex, 0], - self.railButton1NormalForce[: self.outOfRailTimeIndex, 1], + self.railButton1NormalForce[:, 0], + self.railButton1NormalForce[:, 1], label="Upper Rail Button", ) ax1.plot( - self.railButton2NormalForce[: self.outOfRailTimeIndex, 0], - self.railButton2NormalForce[: self.outOfRailTimeIndex, 1], + self.railButton2NormalForce[:, 0], + self.railButton2NormalForce[:, 1], label="Lower Rail Button", ) ax1.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) @@ -3233,13 +3242,13 @@ def plotTrajectoryForceData(self): ax2 = plt.subplot(212) ax2.plot( - self.railButton1ShearForce[: self.outOfRailTimeIndex, 0], - self.railButton1ShearForce[: self.outOfRailTimeIndex, 1], + self.railButton1ShearForce[:, 0], + self.railButton1ShearForce[:, 1], label="Upper Rail Button", ) ax2.plot( - self.railButton2ShearForce[: self.outOfRailTimeIndex, 0], - self.railButton2ShearForce[: self.outOfRailTimeIndex, 1], + self.railButton2ShearForce[:, 0], + self.railButton2ShearForce[:, 1], label="Lower Rail Button", ) ax2.set_xlim(0, self.outOfRailTime if self.outOfRailTime > 0 else self.tFinal) From 3543ec806650ce773e5f4014f9e6a063fe4514d8 Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Wed, 16 Nov 2022 20:33:57 -0300 Subject: [PATCH 07/13] ENH: improve extrapolation option for each funcified method --- rocketpy/Flight.py | 82 +++++++++++++++++++++------------------------- 1 file changed, 37 insertions(+), 45 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index 653bea134..649665a9a 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -1565,17 +1565,17 @@ def z(self): """Rocket z position as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 3]] - @funcify_method("Time (s)", "Vx (m/s)", "spline", "constant") + @funcify_method("Time (s)", "Vx (m/s)", "spline", "zero") def vx(self): """Rocket x velocity as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 4]] - @funcify_method("Time (s)", "Vy (m/s)", "spline", "constant") + @funcify_method("Time (s)", "Vy (m/s)", "spline", "zero") def vy(self): """Rocket y velocity as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 5]] - @funcify_method("Time (s)", "Vz (m/s)", "spline", "constant") + @funcify_method("Time (s)", "Vz (m/s)", "spline", "zero") def vz(self): """Rocket z velocity as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 6]] @@ -1600,86 +1600,86 @@ def e3(self): """Rocket quaternion e3 as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 10]] - @funcify_method("Time (s)", "ω1 (rad/s)", "spline", "constant") + @funcify_method("Time (s)", "ω1 (rad/s)", "spline", "zero") def w1(self): """Rocket angular velocity ω1 as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 11]] - @funcify_method("Time (s)", "ω2 (rad/s)", "spline", "constant") + @funcify_method("Time (s)", "ω2 (rad/s)", "spline", "zero") def w2(self): """Rocket angular velocity ω2 as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 12]] - @funcify_method("Time (s)", "ω3 (rad/s)", "spline", "constant") + @funcify_method("Time (s)", "ω3 (rad/s)", "spline", "zero") def w3(self): """Rocket angular velocity ω3 as a rocketpy.Function of time.""" return self.solutionArray[:, [0, 13]] # Process second type of outputs - accelerations components - @funcify_method("Time (s)", "Ax (m/s²)", "spline", "constant") + @funcify_method("Time (s)", "Ax (m/s²)", "spline", "zero") def ax(self): """Rocket x acceleration as a rocketpy.Function of time.""" return self.retrieve_acceleration_arrays[0] - @funcify_method("Time (s)", "Ay (m/s²)", "spline", "constant") + @funcify_method("Time (s)", "Ay (m/s²)", "spline", "zero") def ay(self): """Rocket y acceleration as a rocketpy.Function of time.""" return self.retrieve_acceleration_arrays[1] - @funcify_method("Time (s)", "Az (m/s²)", "spline", "constant") + @funcify_method("Time (s)", "Az (m/s²)", "spline", "zero") def az(self): """Rocket z acceleration as a rocketpy.Function of time.""" return self.retrieve_acceleration_arrays[2] - @funcify_method("Time (s)", "α1 (rad/s²)", "spline", "constant") + @funcify_method("Time (s)", "α1 (rad/s²)", "spline", "zero") def alpha1(self): """Rocket angular acceleration α1 as a rocketpy.Function of time.""" return self.retrieve_acceleration_arrays[3] - @funcify_method("Time (s)", "α2 (rad/s²)", "spline", "constant") + @funcify_method("Time (s)", "α2 (rad/s²)", "spline", "zero") def alpha2(self): """Rocket angular acceleration α2 as a rocketpy.Function of time.""" return self.retrieve_acceleration_arrays[4] - @funcify_method("Time (s)", "α3 (rad/s²)", "spline", "constant") + @funcify_method("Time (s)", "α3 (rad/s²)", "spline", "zero") def alpha3(self): """Rocket angular acceleration α3 as a rocketpy.Function of time.""" return self.retrieve_acceleration_arrays[5] # Process third type of outputs - Temporary values - @funcify_method("Time (s)", "R1 (N)", "spline", "constant") + @funcify_method("Time (s)", "R1 (N)", "spline", "zero") def R1(self): """Aerodynamic force along the first axis that is perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time.""" return self.retrieve_temporary_values_arrays[0] - @funcify_method("Time (s)", "R2 (N)", "spline", "constant") + @funcify_method("Time (s)", "R2 (N)", "spline", "zero") def R2(self): """Aerodynamic force along the second axis that is perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time.""" return self.retrieve_temporary_values_arrays[1] - @funcify_method("Time (s)", "R3 (N)", "spline", "constant") + @funcify_method("Time (s)", "R3 (N)", "spline", "zero") def R3(self): """Aerodynamic force along the rocket's axis of symmetry as a rocketpy.Function of time.""" return self.retrieve_temporary_values_arrays[2] - @funcify_method("Time (s)", "M1 (Nm)", "spline", "constant") + @funcify_method("Time (s)", "M1 (Nm)", "spline", "zero") def M1(self): """Aerodynamic bending moment in the same direction as the axis that is perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time. """ return self.retrieve_temporary_values_arrays[3] - @funcify_method("Time (s)", "M2 (Nm)", "spline", "constant") + @funcify_method("Time (s)", "M2 (Nm)", "spline", "zero") def M2(self): """Aerodynamic bending moment in the same direction as the axis that is perpendicular to the rocket's axis of symmetry as a rocketpy.Function of time. """ return self.retrieve_temporary_values_arrays[4] - @funcify_method("Time (s)", "M3 (Nm)", "spline", "constant") + @funcify_method("Time (s)", "M3 (Nm)", "spline", "zero") def M3(self): """Aerodynamic bending moment in the same direction as the rocket's axis of symmetry as a rocketpy.Function of time.""" @@ -1895,7 +1895,7 @@ def apogeeFreestreamSpeed(self): return self.freestreamSpeed(self.apogeeTime) # Mach Number - @funcify_method("Time (s)", "Mach Number", "spline", "constant") + @funcify_method("Time (s)", "Mach Number", "spline", "zero") def MachNumber(self): """Mach number as a rocketpy.Function of time.""" return self.freestreamSpeed / self.speedOfSound @@ -1912,7 +1912,7 @@ def maxMachNumber(self): return self.MachNumber(self.maxMachNumberTime) # Reynolds Number - @funcify_method("Time (s)", "Reynolds Number", "spline", "constant") + @funcify_method("Time (s)", "Reynolds Number", "spline", "zero") def ReynoldsNumber(self): """Reynolds number as a rocketpy.Function of time.""" return (self.density * self.freestreamSpeed / self.dynamicViscosity) * ( @@ -1931,7 +1931,7 @@ def maxReynoldsNumber(self): return self.ReynoldsNumber(self.maxReynoldsNumberTime) # Dynamic Pressure - @funcify_method("Time (s)", "Dynamic Pressure (Pa)", "spline", "constant") + @funcify_method("Time (s)", "Dynamic Pressure (Pa)", "spline", "zero") def dynamicPressure(self): """Dynamic pressure as a rocketpy.Function of time.""" return 0.5 * self.density * self.freestreamSpeed**2 @@ -1948,7 +1948,7 @@ def maxDynamicPressure(self): return self.dynamicPressure(self.maxDynamicPressureTime) # Total Pressure - @funcify_method("Time (s)", "Total Pressure (Pa)", "spline", "constant") + @funcify_method("Time (s)", "Total Pressure (Pa)", "spline", "zero") def totalPressure(self): return self.pressure * (1 + 0.2 * self.MachNumber**2) ** (3.5) @@ -1966,29 +1966,29 @@ def maxTotalPressure(self): # Dynamics functions and variables # Aerodynamic Lift and Drag - @funcify_method("Time (s)", "Aerodynamic Lift Force (N)", "spline", "constant") + @funcify_method("Time (s)", "Aerodynamic Lift Force (N)", "spline", "zero") def aerodynamicLift(self): """Aerodynamic lift force as a rocketpy.Function of time.""" return (self.R1**2 + self.R2**2) ** 0.5 - @funcify_method("Time (s)", "Aerodynamic Drag Force (N)", "spline", "constant") + @funcify_method("Time (s)", "Aerodynamic Drag Force (N)", "spline", "zero") def aerodynamicDrag(self): """Aerodynamic drag force as a rocketpy.Function of time.""" return -1 * self.R3 - @funcify_method("Time (s)", "Aerodynamic Bending Moment (Nm)", "spline", "constant") + @funcify_method("Time (s)", "Aerodynamic Bending Moment (Nm)", "spline", "zero") def aerodynamicBendingMoment(self): """Aerodynamic bending moment as a rocketpy.Function of time.""" return (self.M1**2 + self.M2**2) ** 0.5 - @funcify_method("Time (s)", "Aerodynamic Spin Moment (Nm)", "spline", "constant") + @funcify_method("Time (s)", "Aerodynamic Spin Moment (Nm)", "spline", "zero") def aerodynamicSpinMoment(self): """Aerodynamic spin moment as a rocketpy.Function of time.""" return self.M3 # Energy # Kinetic Energy - @funcify_method("Time (s)", "Rotational Kinetic Energy (J)", "spline", "constant") + @funcify_method("Time (s)", "Rotational Kinetic Energy (J)", "spline", "zero") def rotationalEnergy(self): """Rotational kinetic energy as a rocketpy.Function of time.""" b = -self.rocket.distanceRocketPropellant @@ -2007,9 +2007,7 @@ def rotationalEnergy(self): ) return rotationalEnergy - @funcify_method( - "Time (s)", "Translational Kinetic Energy (J)", "spline", "constant" - ) + @funcify_method("Time (s)", "Translational Kinetic Energy (J)", "spline", "zero") def translationalEnergy(self): """Translational kinetic energy as a rocketpy.Function of time.""" # Redefine totalMass time grid to allow for efficient Function algebra @@ -2020,7 +2018,7 @@ def translationalEnergy(self): ) return translationalEnergy - @funcify_method("Time (s)", "Kinetic Energy (J)", "spline", "constant") + @funcify_method("Time (s)", "Kinetic Energy (J)", "spline", "zero") def kineticEnergy(self): """Total kinetic energy as a rocketpy.Function of time.""" return self.rotationalEnergy + self.translationalEnergy @@ -2043,7 +2041,7 @@ def totalEnergy(self): return self.kineticEnergy + self.potentialEnergy # Thrust Power - @funcify_method("Time (s)", "Thrust Power (W)", "spline", "constant") + @funcify_method("Time (s)", "Thrust Power (W)", "spline", "zero") def thrustPower(self): """Thrust power as a rocketpy.Function of time.""" thrust = deepcopy(self.rocket.motor.thrust) @@ -2052,7 +2050,7 @@ def thrustPower(self): return thrustPower # Drag Power - @funcify_method("Time (s)", "Drag Power (W)", "spline", "constant") + @funcify_method("Time (s)", "Drag Power (W)", "spline", "zero") def dragPower(self): """Drag power as a rocketpy.Function of time.""" dragPower = self.R3 * self.speed @@ -2114,29 +2112,23 @@ def _frequency_analysis(self, func, samplingFrequency=1000): return FourierFrequencies, FourierAmplitude # Angular velocities frequency response - Fourier Analysis - @funcify_method( - "Frequency (Hz)", "ω1 Angle Fourier Amplitude", "spline", "constant" - ) + @funcify_method("Frequency (Hz)", "ω1 Angle Fourier Amplitude", "spline", "zero") def omega1FrequencyResponse(self): """Angular velocity 1 frequency response as a rocketpy.Function of frequency.""" return np.column_stack(self._frequency_analysis(self.w1)) - @funcify_method( - "Frequency (Hz)", "ω2 Angle Fourier Amplitude", "spline", "constant" - ) + @funcify_method("Frequency (Hz)", "ω2 Angle Fourier Amplitude", "spline", "zero") def omega2FrequencyResponse(self): """Angular velocity 2 frequency response as a rocketpy.Function of frequency.""" return np.column_stack(self._frequency_analysis(self.w2)) - @funcify_method( - "Frequency (Hz)", "ω3 Angle Fourier Amplitude", "spline", "constant" - ) + @funcify_method("Frequency (Hz)", "ω3 Angle Fourier Amplitude", "spline", "zero") def omega3FrequencyResponse(self): """Angular velocity 3 frequency response as a rocketpy.Function of frequency.""" return np.column_stack(self._frequency_analysis(self.w3)) @funcify_method( - "Frequency (Hz)", "Attitude Angle Fourier Amplitude", "spline", "constant" + "Frequency (Hz)", "Attitude Angle Fourier Amplitude", "spline", "zero" ) def attitudeFrequencyResponse(self): """Attitude frequency response as a rocketpy.Function of frequency.""" @@ -2285,7 +2277,7 @@ def bearing(self): bearing = np.column_stack((t, bearing)) return bearing - @funcify_method("Time (s)", "Latitude (°)", "linear") + @funcify_method("Time (s)", "Latitude (°)", "linear", "constant") def latitude(self): """Rocket latitude coordinate, in degrees, as a rocketpy.Function of time.""" lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians @@ -2301,7 +2293,7 @@ def latitude(self): ) return np.column_stack((self.time, latitude)) - @funcify_method("Time (s)", "Longitude (°)", "linear") + @funcify_method("Time (s)", "Longitude (°)", "linear", "constant") def longitude(self): """Rocket longitude coordinate, in degrees, as a rocketpy.Function of time.""" lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians From f112a3b77e651f33367ebc56e1753b8e6ec24f70 Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Wed, 16 Nov 2022 21:22:25 -0300 Subject: [PATCH 08/13] MAINT: minor improvements to frequency analysis --- rocketpy/Flight.py | 22 ++++++++++++++++------ 1 file changed, 16 insertions(+), 6 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index 649665a9a..c691f55dd 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -2086,7 +2086,7 @@ def angleOfAttack(self): return angleOfAttack # Stability and Control variables - def _frequency_analysis(self, func, samplingFrequency=1000): + def _frequency_analysis(self, func, samplingFrequency=100, sampleWindow=5): """Frequency analysis of any function of time between out of rail time and final simulation time. @@ -2094,6 +2094,11 @@ def _frequency_analysis(self, func, samplingFrequency=1000): ---------- func : rocketpy.Function, callable Function of time to be analyzed. + samplingFrequency : int, optional + Sampling frequency in Hz. The default is 100 Hz. + sampleWindow : float, optional + Sample window in seconds, starting from Flight.outOfRailTime. + The default is 5 seconds. Returns ------- @@ -2103,7 +2108,8 @@ def _frequency_analysis(self, func, samplingFrequency=1000): Array of amplitude values. """ samplingTimeStep = 1.0 / samplingFrequency - samplingRange = np.arange(self.outOfRailTime, self.tFinal, samplingTimeStep) + finalTime = self.outOfRailTime + sampleWindow + samplingRange = np.arange(self.outOfRailTime, finalTime, samplingTimeStep) sampledPoints = func(samplingRange) sampledPoints -= np.mean(sampledPoints) numberOfSamples = len(sampledPoints) @@ -2114,24 +2120,28 @@ def _frequency_analysis(self, func, samplingFrequency=1000): # Angular velocities frequency response - Fourier Analysis @funcify_method("Frequency (Hz)", "ω1 Angle Fourier Amplitude", "spline", "zero") def omega1FrequencyResponse(self): - """Angular velocity 1 frequency response as a rocketpy.Function of frequency.""" + """Angular velocity 1 frequency response as a rocketpy.Function of frequency, + as the rocket leaves the launch rail.""" return np.column_stack(self._frequency_analysis(self.w1)) @funcify_method("Frequency (Hz)", "ω2 Angle Fourier Amplitude", "spline", "zero") def omega2FrequencyResponse(self): - """Angular velocity 2 frequency response as a rocketpy.Function of frequency.""" + """Angular velocity 2 frequency response as a rocketpy.Function of frequency, + as the rocket leaves the launch rail.""" return np.column_stack(self._frequency_analysis(self.w2)) @funcify_method("Frequency (Hz)", "ω3 Angle Fourier Amplitude", "spline", "zero") def omega3FrequencyResponse(self): - """Angular velocity 3 frequency response as a rocketpy.Function of frequency.""" + """Angular velocity 3 frequency response as a rocketpy.Function of frequency, + as the rocket leaves the launch rail.""" return np.column_stack(self._frequency_analysis(self.w3)) @funcify_method( "Frequency (Hz)", "Attitude Angle Fourier Amplitude", "spline", "zero" ) def attitudeFrequencyResponse(self): - """Attitude frequency response as a rocketpy.Function of frequency.""" + """Attitude frequency response as a rocketpy.Function of frequency, as the + rocket leaves the launch rail.""" return np.column_stack(self._frequency_analysis(self.attitudeAngle)) # Static Margin From cf5028ae36c362bb4df711852072502c6b8deefd Mon Sep 17 00:00:00 2001 From: Lint Action Date: Thu, 17 Nov 2022 00:23:15 +0000 Subject: [PATCH 09/13] Fix code style issues with Black --- rocketpy/Flight.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index c691f55dd..6e85b84f7 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -2120,19 +2120,19 @@ def _frequency_analysis(self, func, samplingFrequency=100, sampleWindow=5): # Angular velocities frequency response - Fourier Analysis @funcify_method("Frequency (Hz)", "ω1 Angle Fourier Amplitude", "spline", "zero") def omega1FrequencyResponse(self): - """Angular velocity 1 frequency response as a rocketpy.Function of frequency, + """Angular velocity 1 frequency response as a rocketpy.Function of frequency, as the rocket leaves the launch rail.""" return np.column_stack(self._frequency_analysis(self.w1)) @funcify_method("Frequency (Hz)", "ω2 Angle Fourier Amplitude", "spline", "zero") def omega2FrequencyResponse(self): - """Angular velocity 2 frequency response as a rocketpy.Function of frequency, + """Angular velocity 2 frequency response as a rocketpy.Function of frequency, as the rocket leaves the launch rail.""" return np.column_stack(self._frequency_analysis(self.w2)) @funcify_method("Frequency (Hz)", "ω3 Angle Fourier Amplitude", "spline", "zero") def omega3FrequencyResponse(self): - """Angular velocity 3 frequency response as a rocketpy.Function of frequency, + """Angular velocity 3 frequency response as a rocketpy.Function of frequency, as the rocket leaves the launch rail.""" return np.column_stack(self._frequency_analysis(self.w3)) From 23f75609cd31cab1c0d54d7aae6f9743a9468ecc Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Thu, 17 Nov 2022 23:50:04 -0300 Subject: [PATCH 10/13] ENH: implement toFrequencyDomain --- rocketpy/Function.py | 54 ++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 54 insertions(+) diff --git a/rocketpy/Function.py b/rocketpy/Function.py index d369c86f6..c1c64d884 100644 --- a/rocketpy/Function.py +++ b/rocketpy/Function.py @@ -1081,6 +1081,60 @@ def __len__(self): """ return len(self.source) + # Define all conversion methods + def toFrequencyDomain(self, lower, upper, samplingFrequency, removeDC=True): + """Performs the conversion of the Function to the Frequency Domain and returns + the result. This is done by taking the Fourier transform of the Function. + The resulting frequency domain is symmetric, i.e., the negative frequencies are + included as well. + + Parameters + ---------- + lower : float + Lower bound of the time range. + upper : float + Upper bound of the time range. + samplingFrequency : float + Sampling frequency at which to perform the Fourier transform. + removeDC : bool, optional + If True, the DC component is removed from the Fourier transform. + + Returns + ------- + Function + The Function in the frequency domain. + + Examples + -------- + >>> from rocketpy import Function + >>> import numpy as np + >>> mainFrequency = 10 # Hz + >>> time = np.linspace(0, 10, 1000) + >>> signal = np.sin(2 * np.pi * mainFrequency * time) + >>> timeDomain = Function(np.array([time, signal]).T) + >>> frequencyDomain = timeDomain.toFrequencyDomain(lower=0, upper=10, samplingFrequency=100) + >>> peakFrequenciesIndex = np.where(frequencyDomain[:, 1] > 0.001) + >>> peakFrequencies = frequencyDomain[peakFrequenciesIndex, 0] + >>> print(peakFrequencies) + [[-10. 10.]] + """ + # Get the time domain data + samplingTimeStep = 1.0 / samplingFrequency + samplingRange = np.arange(lower, upper, samplingTimeStep) + numberOfSamples = len(samplingRange) + sampledPoints = self(samplingRange) + if removeDC: + sampledPoints -= np.mean(sampledPoints) + FourierAmplitude = np.abs(np.fft.fft(sampledPoints) / numberOfSamples) + FourierFrequencies = np.fft.fftfreq(numberOfSamples, samplingTimeStep) + return Function( + source=np.array([FourierFrequencies, FourierAmplitude]).T, + inputs="Frequency (Hz)", + outputs="Amplitude", + interpolation="linear", + extrapolation="zero", + ) + # Define all presentation methods def __call__(self, *args): """Plot the Function if no argument is given. If an From ba8b91d34830a647939a674ed411883d3d3f80b5 Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Thu, 17 Nov 2022 23:50:49 -0300 Subject: [PATCH 11/13] MAINT: replace _frequency_analysis with toFrequencyDomain --- rocketpy/Flight.py | 67 +++++++++++++++------------------------------- 1 file changed, 22 insertions(+), 45 deletions(-) diff --git a/rocketpy/Flight.py b/rocketpy/Flight.py index 6e85b84f7..f5ee97c9d 100644 --- a/rocketpy/Flight.py +++ b/rocketpy/Flight.py @@ -2085,66 +2085,43 @@ def angleOfAttack(self): return angleOfAttack - # Stability and Control variables - def _frequency_analysis(self, func, samplingFrequency=100, sampleWindow=5): - """Frequency analysis of any function of time between out of rail time - and final simulation time. - - Parameters - ---------- - func : rocketpy.Function, callable - Function of time to be analyzed. - samplingFrequency : int, optional - Sampling frequency in Hz. The default is 100 Hz. - sampleWindow : float, optional - Sample window in seconds, starting from Flight.outOfRailTime. - The default is 5 seconds. - - Returns - ------- - FourierFrequencies : numpy.ndarray - Array of frequencies. - FourierAmplitude : numpy.ndarray - Array of amplitude values. - """ - samplingTimeStep = 1.0 / samplingFrequency - finalTime = self.outOfRailTime + sampleWindow - samplingRange = np.arange(self.outOfRailTime, finalTime, samplingTimeStep) - sampledPoints = func(samplingRange) - sampledPoints -= np.mean(sampledPoints) - numberOfSamples = len(sampledPoints) - FourierAmplitude = np.abs(np.fft.fft(sampledPoints) / numberOfSamples) - FourierFrequencies = np.fft.fftfreq(numberOfSamples, samplingTimeStep) - return FourierFrequencies, FourierAmplitude - - # Angular velocities frequency response - Fourier Analysis - @funcify_method("Frequency (Hz)", "ω1 Angle Fourier Amplitude", "spline", "zero") + # Frequency response and stability variables + @funcify_method("Frequency (Hz)", "ω1 Fourier Amplitude", "spline", "zero") def omega1FrequencyResponse(self): """Angular velocity 1 frequency response as a rocketpy.Function of frequency, - as the rocket leaves the launch rail.""" - return np.column_stack(self._frequency_analysis(self.w1)) + as the rocket leaves the launch rail for 5 seconds of flight.""" + return self.w1.toFrequencyDomain( + self.outOfRailTime, self.outOfRailTime + 5, 100 + ) - @funcify_method("Frequency (Hz)", "ω2 Angle Fourier Amplitude", "spline", "zero") + @funcify_method("Frequency (Hz)", "ω2 Fourier Amplitude", "spline", "zero") def omega2FrequencyResponse(self): """Angular velocity 2 frequency response as a rocketpy.Function of frequency, - as the rocket leaves the launch rail.""" - return np.column_stack(self._frequency_analysis(self.w2)) + as the rocket leaves the launch rail for 5 seconds of flight.""" + return self.w2.toFrequencyDomain( + self.outOfRailTime, self.outOfRailTime + 5, 100 + ) - @funcify_method("Frequency (Hz)", "ω3 Angle Fourier Amplitude", "spline", "zero") + @funcify_method("Frequency (Hz)", "ω3 Fourier Amplitude", "spline", "zero") def omega3FrequencyResponse(self): """Angular velocity 3 frequency response as a rocketpy.Function of frequency, - as the rocket leaves the launch rail.""" - return np.column_stack(self._frequency_analysis(self.w3)) + as the rocket leaves the launch rail for 5 seconds of flight.""" + return self.w3.toFrequencyDomain( + self.outOfRailTime, self.outOfRailTime + 5, 100 + ) @funcify_method( "Frequency (Hz)", "Attitude Angle Fourier Amplitude", "spline", "zero" ) def attitudeFrequencyResponse(self): """Attitude frequency response as a rocketpy.Function of frequency, as the - rocket leaves the launch rail.""" - return np.column_stack(self._frequency_analysis(self.attitudeAngle)) + rocket leaves the launch rail for 5 seconds of flight.""" + return self.attitudeAngle.toFrequencyDomain( + lower=self.outOfRailTime, + upper=self.outOfRailTime + 5, + samplingFrequency=100, + ) - # Static Margin @cached_property def staticMargin(self): """Static margin of the rocket.""" From 351b9854e6a8ae61f28352667300108285c461de Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Thu, 17 Nov 2022 23:52:03 -0300 Subject: [PATCH 12/13] TST: delete temp test files after testing is done --- tests/test_flight.py | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/tests/test_flight.py b/tests/test_flight.py index 9c66a7b3a..640783d36 100644 --- a/tests/test_flight.py +++ b/tests/test_flight.py @@ -1,3 +1,4 @@ +import os import datetime from unittest.mock import patch @@ -503,6 +504,10 @@ def test_export_data(): test_1 = np.loadtxt("test_export_data_1.csv", delimiter=",") test_2 = np.loadtxt("test_export_data_2.csv", delimiter=",") + # Delete files + os.remove("test_export_data_1.csv") + os.remove("test_export_data_2.csv") + # Check if basic exported content matches data assert np.allclose(test_flight.x[:, 0], test_1[:, 0], atol=1e-5) == True assert np.allclose(test_flight.x[:, 1], test_1[:, 1], atol=1e-5) == True @@ -604,6 +609,10 @@ def test_export_KML(): lat.append(float(coords[i + 1])) z.append(float(coords[i + 2])) + # Delete temporary test file + test_1.close() + os.remove("test_export_data_1.kml") + assert np.allclose(test_flight.latitude[:, 1], lat, atol=1e-3) == True assert np.allclose(test_flight.longitude[:, 1], lon, atol=1e-3) == True assert np.allclose(test_flight.z[:, 1], z, atol=1e-3) == True From f9871e7ad446028d01945aef60da4c117c536cd3 Mon Sep 17 00:00:00 2001 From: Giovani Hidalgo Ceotto Date: Mon, 21 Nov 2022 08:37:10 -0300 Subject: [PATCH 13/13] FIX: update Fourier transform normalization factor Co-authored-by: MateusStano <69485049+MateusStano@users.noreply.github.com> --- rocketpy/Function.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/rocketpy/Function.py b/rocketpy/Function.py index c1c64d884..f0d64e60b 100644 --- a/rocketpy/Function.py +++ b/rocketpy/Function.py @@ -1125,7 +1125,7 @@ def toFrequencyDomain(self, lower, upper, samplingFrequency, removeDC=True): sampledPoints = self(samplingRange) if removeDC: sampledPoints -= np.mean(sampledPoints) - FourierAmplitude = np.abs(np.fft.fft(sampledPoints) / numberOfSamples) + FourierAmplitude = np.abs(np.fft.fft(sampledPoints) / (numberOfSamples / 2)) FourierFrequencies = np.fft.fftfreq(numberOfSamples, samplingTimeStep) return Function( source=np.array([FourierFrequencies, FourierAmplitude]).T,