ENH: added radius and porosity to parachute

Author: ArthurJWH

rocketpy/rocket/parachute.py

@@ -102,6 +102,8 @@ def __init__(
         sampling_rate,
         lag=0,
         noise=(0, 0, 0),
+        parachute_radius=1.5,
+        porosity=0.0432,
     ):
         """Initializes Parachute class.
 
@@ -161,6 +163,8 @@ def __init__(
         self.sampling_rate = sampling_rate
         self.lag = lag
         self.noise = noise
+        self.parachute_radius = parachute_radius
+        self.porosity = porosity
         self.noise_signal = [[-1e-6, np.random.normal(noise[0], noise[1])]]
         self.noisy_pressure_signal = []
         self.clean_pressure_signal = []
@@ -264,6 +268,8 @@ def to_dict(self, include_outputs=False):
             "sampling_rate": self.sampling_rate,
             "lag": self.lag,
             "noise": self.noise,
+            "parachute_radius": self.parachute_radius,
+            "porosity": self.porosity,
         }
 
         if include_outputs:
@@ -290,6 +296,8 @@ def from_dict(cls, data):
             sampling_rate=data["sampling_rate"],
             lag=data["lag"],
             noise=data["noise"],
+            parachute_radius=data["parachute_radius"],
+            porosity=data["porosity"],
         )
 
         return parachute

rocketpy/simulation/flight.py

@@ -767,7 +767,13 @@ def __simulate(self, verbose):
                         callbacks = [
                             lambda self, parachute_cd_s=parachute.cd_s: setattr(
                                 self, "parachute_cd_s", parachute_cd_s
-                            )
+                            ),
+                            lambda self, parachute_radius=parachute.radius: setattr(
+                                self, "parachute_radius", parachute_radius
+                            ),
+                            lambda self, parachute_porosity=parachute.porosity: setattr(
+                                self, "parachute_porosity", parachute_porosity
+                            ),
                         ]
                         self.flight_phases.add_phase(
                             node.t + parachute.lag,
@@ -1013,7 +1019,13 @@ def __simulate(self, verbose):
                                             lambda self,
                                             parachute_cd_s=parachute.cd_s: setattr(
                                                 self, "parachute_cd_s", parachute_cd_s
-                                            )
+                                            ),
+                                            lambda self, parachute_radius=parachute.radius: setattr(
+                                                self, "parachute_radius", parachute_radius
+                                            ),
+                                            lambda self, parachute_porosity=parachute.porosity: setattr(
+                                                self, "parachute_porosity", parachute_porosity
+                                            ),
                                         ]
                                         self.flight_phases.add_phase(
                                             overshootable_node.t + parachute.lag,
@@ -1961,22 +1973,26 @@ def u_dot_parachute(self, t, u, post_processing=False):
 
         # Get Parachute data
         cd_s = self.parachute_cd_s
+        R = self.parachute_radius
+        porosity = self.parachute_porosity
+
 
         # Get the mass of the rocket
         mp = self.rocket.dry_mass
 
         # Define constants
-        ka = 1  # Added mass coefficient (depends on parachute's porosity)
-        R = 1.5  # Parachute radius
+        ka = 1.068 * (1 - 1.465 * porosity - 0.25975 * porosity**2 + 1.2626 * porosity**3)
         # to = 1.2
         # eta = 1
         # Rdot = (6 * R * (1 - eta) / (1.2**6)) * (
         #     (1 - eta) * t**5 + eta * (to**3) * (t**2)
         # )
         # Rdot = 0
 
+        # tf = 8 * nominal diameter / velocity at line stretch
+
         # Calculate added mass
-        ma = ka * rho * (4 / 3) * np.pi * R**3
+        ma = ka * rho * (4 / 3) * np.pi * R**3 # ma = ka * rho * (4 / 3) * np.pi * R**2 * height
 
         # Calculate freestream speed
         freestream_x = vx - wind_velocity_x
@@ -1985,9 +2001,9 @@ def u_dot_parachute(self, t, u, post_processing=False):
         free_stream_speed = (freestream_x**2 + freestream_y**2 + freestream_z**2) ** 0.5
 
         # Determine drag force
-        pseudo_drag = -0.5 * rho * cd_s * free_stream_speed
+        pseudo_drag = -0.5 * rho * cd_s * free_stream_speed # * Area
         # pseudo_drag = pseudo_drag - ka * rho * 4 * np.pi * (R**2) * Rdot
-        Dx = pseudo_drag * freestream_x
+        Dx = pseudo_drag * freestream_x # add eta efficiency for wake
         Dy = pseudo_drag * freestream_y
         Dz = pseudo_drag * freestream_z
         ax = Dx / (mp + ma)