diff --git a/docs/reference/classes/aero_surfaces/TubeFins.rst b/docs/reference/classes/aero_surfaces/TubeFins.rst new file mode 100644 index 000000000..27e9827a9 --- /dev/null +++ b/docs/reference/classes/aero_surfaces/TubeFins.rst @@ -0,0 +1,5 @@ +TubeFins Class +============== + +.. autoclass:: rocketpy.TubeFins + :members: diff --git a/docs/reference/classes/aero_surfaces/index.rst b/docs/reference/classes/aero_surfaces/index.rst index a3dad0417..1662c965e 100644 --- a/docs/reference/classes/aero_surfaces/index.rst +++ b/docs/reference/classes/aero_surfaces/index.rst @@ -12,6 +12,7 @@ AeroSurface Classes TrapezoidalFins EllipticalFins FreeFormFins + TubeFins Fin TrapezoidalFin EllipticalFin @@ -19,4 +20,4 @@ AeroSurface Classes RailButtons AirBrakes GenericSurface - LinearGenericSurface \ No newline at end of file + LinearGenericSurface diff --git a/docs/technical/aerodynamics/tube_fins.rst b/docs/technical/aerodynamics/tube_fins.rst new file mode 100644 index 000000000..46005d24b --- /dev/null +++ b/docs/technical/aerodynamics/tube_fins.rst @@ -0,0 +1,87 @@ +Tube Fin Aerodynamics +===================== + +RocketPy models a tube-fin set as a symmetric ring of uncanted ring airfoils. +The implementation follows the preliminary tube-fin model in OpenRocket's +`TubeFinSetCalc `_. +The normal-force derivative is based on Ribner's analysis of a ring airfoil in +nonaxial flow [1]_. + +Geometry +-------- + +Let :math:`n` be the number of tubes, :math:`R` the rocket-body radius, +:math:`r_i` the tube inner radius, :math:`r_o` the tube outer radius, and +:math:`L` the tube length. The tubes are distributed evenly around the rocket. +The current implementation requires each tube to touch its two neighbors: + +.. math:: + + r_o = R \frac{\sin(\pi / n)}{1 - \sin(\pi / n)}. + +This constraint also places each tube against the rocket body. Configurations +with fewer than three tubes, gaps between tubes, or overlapping tubes are +rejected. + +Normal Force +------------ + +The ring-airfoil aspect ratio and its modified form are + +.. math:: + + AR = \frac{2 r_i}{L}, \qquad AR' = \frac{2 AR}{\pi}. + +For a rocket reference area :math:`A_{ref} = \pi R^2`, the normal-force +coefficient derivative of the complete tube set is + +.. math:: + + C_{N_\alpha} = + \frac{n}{A_{ref}} + 2 \left(\frac{AR'}{1 + AR'}\right) \pi^2 r_i L. + +RocketPy applies this derivative symmetrically for positive and negative +angles of attack and caps the magnitude at 20 degrees: + +.. math:: + + C_N(\alpha) = C_{N_\alpha} + \operatorname{clip}\left(\alpha, -20^\circ, 20^\circ\right). + +Center of Pressure +------------------ + +For Mach numbers up to 0.5, the center of pressure is placed at the quarter +chord: + +.. math:: + + x_{CP} = \frac{L}{4}. + +The position is measured from the tube leading edge. OpenRocket moves this +position with Mach number above Mach 0.5; RocketPy does not yet implement that +correction. Simulations that exceed Mach 0.5 should use aerodynamic data from a +higher-fidelity source instead of this fixed-CP model. + +Model Limits +------------ + +The tube-fin surface contributes normal force and the corresponding pitch and +yaw moments about its center of pressure. It does not calculate: + +- friction or pressure drag from the tubes; +- roll forcing or damping; +- side-force or yaw behavior for asymmetric tube layouts; +- tube cant; or +- aerodynamic corrections for separated or overlapping tubes. + +Represent tube-fin drag in the rocket's power-on and power-off drag curves. +Use :class:`rocketpy.GenericSurface` when measured, wind-tunnel, or CFD +coefficients are available outside the limits above. + +References +---------- + +.. [1] Ribner, H. S. "The Ring Airfoil in Nonaxial Flow." *Journal of the + Aeronautical Sciences*, 14(9), 529--530, 1947. diff --git a/docs/technical/index.rst b/docs/technical/index.rst index 73583eba9..360d24c62 100644 --- a/docs/technical/index.rst +++ b/docs/technical/index.rst @@ -14,9 +14,10 @@ in their code. Equations of Motion v1 Elliptical Fins Individual Fin + Tube Fins Roll Moment Sensitivity Analysis References This section is still a work in progress, however, and not everything is documented yet. -If you have any questions, please contact the maintainers of RocketPy. \ No newline at end of file +If you have any questions, please contact the maintainers of RocketPy. diff --git a/docs/user/aerodynamics/surfaces.rst b/docs/user/aerodynamics/surfaces.rst index 38f266dca..11d7629c4 100644 --- a/docs/user/aerodynamics/surfaces.rst +++ b/docs/user/aerodynamics/surfaces.rst @@ -7,12 +7,11 @@ Aerodynamic Surfaces This page provides an overview of the aerodynamic surfaces available in RocketPy and explains how they connect to the rocket's simulation. -RocketPy models the aerodynamic forces and moments generated by three types -of surfaces: **nose cones**, **fins**, and **tails**. Each surface is -defined by its geometric parameters and, optionally, by an airfoil profile. -The aerodynamic coefficients are computed internally using the Barrowman -method and are used during the flight simulation to evaluate the rocket's -stability and control. +RocketPy models the aerodynamic forces and moments generated by **nose +cones**, **planar fins**, **tube fins**, and **tails**. Each surface is defined +by its geometric parameters. Planar fins can also use a measured airfoil lift +curve. The resulting aerodynamic coefficients are used during the flight +simulation to evaluate the rocket's stability and control. .. seealso:: @@ -112,6 +111,10 @@ RocketPy distinguishes between two levels of fin definition: separately. Individual fins are useful for canards, asymmetric configurations, or when you need fine-grained control. +- **Tube-fin sets** (:class:`rocketpy.TubeFins`): This class defines a + symmetric ring of cylindrical fins. Its normal-force derivative follows + Ribner's ring-airfoil model rather than the planar-fin Barrowman equations. + .. seealso:: For the mathematical model of individual fins, including the moment @@ -178,6 +181,31 @@ Parameters: - ``coordinates``: A list of ``(x, y)`` tuples defining the fin shape in the fin coordinate frame. +Tube Fins +~~~~~~~~~ + +Tube fins are defined by the number and length of the tubes, their inner and +outer radii, and the rocket-body radius at the mounting position. Add them to a +rocket with :meth:`rocketpy.Rocket.add_tube_fins` or create a +:class:`rocketpy.TubeFins` object and pass it to +:meth:`rocketpy.Rocket.add_surfaces`. + +The current implementation is a subsonic normal-force model with these +limits: + +- The center of pressure is fixed at one quarter of the tube length, measured + from the leading edge. Use the model only through Mach 0.5. +- Lift is capped at an absolute angle of attack of 20 degrees. +- At least three uncanted tubes must be distributed evenly around the body. + Every tube must touch the body and both adjacent tubes. Separated and + overlapping tube layouts are rejected. +- The model does not calculate tube-fin friction drag, pressure drag, roll, + side force, or yaw. Include tube-fin drag in the rocket's power-on and + power-off drag curves. + +For the equations and geometry constraint, see +:doc:`Tube Fin Aerodynamics `. + Common Fin Set Parameters ------------------------- @@ -260,6 +288,7 @@ Fins can be added to a rocket using the ``Rocket`` class methods: - :meth:`rocketpy.Rocket.add_trapezoidal_fins` - :meth:`rocketpy.Rocket.add_elliptical_fins` - :meth:`rocketpy.Rocket.add_free_form_fins` + - :meth:`rocketpy.Rocket.add_tube_fins` - :meth:`rocketpy.Rocket.add_surfaces` (for individual fins) Tail diff --git a/rocketpy/__init__.py b/rocketpy/__init__.py index 6008ff09b..27b7c034a 100644 --- a/rocketpy/__init__.py +++ b/rocketpy/__init__.py @@ -52,6 +52,7 @@ Tail, TrapezoidalFin, TrapezoidalFins, + TubeFins, ) from .sensitivity import SensitivityModel from .sensors import Accelerometer, Barometer, GnssReceiver, Gyroscope diff --git a/rocketpy/plots/aero_surface_plots.py b/rocketpy/plots/aero_surface_plots.py index eb97ce19b..a9f7b3a06 100644 --- a/rocketpy/plots/aero_surface_plots.py +++ b/rocketpy/plots/aero_surface_plots.py @@ -842,6 +842,41 @@ def draw(self, *, filename=None): show_or_save_plot(filename) +class _TubeFinsPlots(_AeroSurfacePlots): + """Class that contains all tube-fin plots.""" + + def draw(self, *, filename=None): + """Draw a side-view envelope of the tube-fin set.""" + axial, radial = self.aero_surface.shape_vec + _, ax = plt.subplots() + + ax.plot(axial, radial, color="#A60628", label="Tube-fin envelope") + ax.plot(axial, -radial, color="#A60628") + ax.plot( + [0, self.aero_surface.length], + [0, 0], + color="#7A68A6", + linestyle="--", + label="Rocket centerline", + ) + + cp_point = (self.aero_surface.cpz, 0) + ax.scatter(*cp_point, label="Center of Pressure", color="red", zorder=10) + ax.scatter(*cp_point, facecolors="none", edgecolors="red", s=300, zorder=10) + + limit = self.aero_surface.rocket_radius + 2 * self.aero_surface.outer_radius + ax.set_xlim(-0.02 * self.aero_surface.length, 1.02 * self.aero_surface.length) + ax.set_ylim(-1.05 * limit, 1.05 * limit) + ax.set_aspect("equal") + ax.set_xlabel("Length (m)") + ax.set_ylabel("Radius (m)") + ax.set_title("Tube Fin Set Side-View Envelope") + ax.grid(True, linestyle="--", linewidth=0.5) + ax.legend(bbox_to_anchor=(1.05, 1.0), loc="upper left") + plt.tight_layout() + show_or_save_plot(filename) + + class _TailPlots(_AeroSurfacePlots): """Class that contains all tail plots.""" diff --git a/rocketpy/prints/aero_surface_prints.py b/rocketpy/prints/aero_surface_prints.py index cc36f1b01..c3c02ea57 100644 --- a/rocketpy/prints/aero_surface_prints.py +++ b/rocketpy/prints/aero_surface_prints.py @@ -291,6 +291,21 @@ class _FreeFormFinPrints(_FinPrints): """Class that contains all free form fins prints.""" +class _TubeFinsPrints(_AeroSurfacePrints): + """Class that contains all tube-fin prints.""" + + def geometry(self): + """Print the geometric information of the tube-fin set.""" + print("Geometric information of the tube-fin set:") + print("------------------------------------------") + print(f"Number of tubes: {self.aero_surface.n}") + print(f"Tube length: {self.aero_surface.length:.3f} m") + print(f"Inner tube radius: {self.aero_surface.inner_radius:.3f} m") + print(f"Outer tube radius: {self.aero_surface.outer_radius:.3f} m") + print(f"Reference rocket radius: {self.aero_surface.rocket_radius:.3f} m") + print(f"Ring-airfoil aspect ratio: {self.aero_surface.aspect_ratio:.3f}\n") + + class _TailPrints(_AeroSurfacePrints): """Class that contains all tail prints.""" diff --git a/rocketpy/rocket/__init__.py b/rocketpy/rocket/__init__.py index afb7f0bb6..3fa907f2d 100644 --- a/rocketpy/rocket/__init__.py +++ b/rocketpy/rocket/__init__.py @@ -15,6 +15,7 @@ Tail, TrapezoidalFin, TrapezoidalFins, + TubeFins, ) from rocketpy.rocket.components import Components from rocketpy.rocket.parachute import Parachute diff --git a/rocketpy/rocket/aero_surface/__init__.py b/rocketpy/rocket/aero_surface/__init__.py index 7634d3500..ec8aed537 100644 --- a/rocketpy/rocket/aero_surface/__init__.py +++ b/rocketpy/rocket/aero_surface/__init__.py @@ -15,3 +15,4 @@ from rocketpy.rocket.aero_surface.nose_cone import NoseCone from rocketpy.rocket.aero_surface.rail_buttons import RailButtons from rocketpy.rocket.aero_surface.tail import Tail +from rocketpy.rocket.aero_surface.tube_fins import TubeFins diff --git a/rocketpy/rocket/aero_surface/tube_fins.py b/rocketpy/rocket/aero_surface/tube_fins.py new file mode 100644 index 000000000..ca585d08a --- /dev/null +++ b/rocketpy/rocket/aero_surface/tube_fins.py @@ -0,0 +1,291 @@ +import numbers + +import numpy as np + +from rocketpy.mathutils.function import Function +from rocketpy.plots.aero_surface_plots import _TubeFinsPlots +from rocketpy.prints.aero_surface_prints import _TubeFinsPrints + +from .aero_surface import AeroSurface + + +class TubeFins(AeroSurface): + """Defines a symmetric set of tube fins for subsonic flight. + + The aerodynamic model follows the Ribner ring-airfoil normal-force + derivative used by OpenRocket. It is limited to uncanted tube fins that + touch both the rocket body and their two neighboring tubes. The center of + pressure is fixed at the quarter chord, so the model is intended for + Mach numbers up to 0.5. + + Parameters + ---------- + n : int + Number of tubes. Must be at least 3. + length : int, float + Tube length along the rocket axis, in meters. + inner_radius : int, float + Inner radius of each tube, in meters. + outer_radius : int, float + Outer radius of each tube, in meters. For the supported touching + geometry, this must equal + ``rocket_radius * sin(pi / n) / (1 - sin(pi / n))``. + rocket_radius : int, float + Radius of the rocket body where the tube fins are mounted, in meters. + name : str, optional + Name of the tube-fin set. Default is ``"Tube Fins"``. + + Notes + ----- + This model calculates normal force only. Tube-fin friction and pressure + drag must be included in the rocket's power-on and power-off drag curves. + Cant, roll, side-force, yaw, separated tubes, and overlapping tubes are not + supported. + """ + + stall_angle = np.radians(20) + + def __init__( + self, + n, + length, + inner_radius, + outer_radius, + rocket_radius, + name="Tube Fins", + ): + self._n = n + self._length = length + self._inner_radius = inner_radius + self._outer_radius = outer_radius + self._rocket_radius = rocket_radius + + self._validate_geometry() + super().__init__( + name=name, + reference_area=np.pi * rocket_radius**2, + reference_length=2 * rocket_radius, + ) + + self._evaluate_all() + + self.prints = _TubeFinsPrints(self) + self.plots = _TubeFinsPlots(self) + + @staticmethod + def _touching_outer_radius(n, rocket_radius): + sin_half_angle = np.sin(np.pi / n) + return rocket_radius * sin_half_angle / (1 - sin_half_angle) + + def _validate_geometry(self): + if isinstance(self.n, bool) or not isinstance(self.n, numbers.Integral): + raise ValueError("'n' must be an integer greater than or equal to 3.") + if self.n < 3: + raise ValueError("'n' must be greater than or equal to 3.") + + dimensions = { + "length": self.length, + "inner_radius": self.inner_radius, + "outer_radius": self.outer_radius, + "rocket_radius": self.rocket_radius, + } + for parameter, value in dimensions.items(): + if isinstance(value, bool) or not isinstance(value, numbers.Real): + raise ValueError(f"'{parameter}' must be a positive real number.") + if not np.isfinite(value) or value <= 0: + raise ValueError(f"'{parameter}' must be finite and greater than zero.") + + if self.inner_radius >= self.outer_radius: + raise ValueError("'inner_radius' must be smaller than 'outer_radius'.") + + touching_radius = self._touching_outer_radius(self.n, self.rocket_radius) + tolerance = max(1e-12, touching_radius * 1e-6) + if not np.isclose( + self.outer_radius, touching_radius, rtol=1e-6, atol=tolerance + ): + geometry = ( + "separated" if self.outer_radius < touching_radius else "overlapping" + ) + raise ValueError( + f"The specified geometry produces {geometry} tube fins. " + "Only mutually tangent tubes are supported; for " + f"n={self.n} and rocket_radius={self.rocket_radius:g} m, " + f"outer_radius must be {touching_radius:g} m." + ) + + def _evaluate_all(self): + self.reference_area = np.pi * self.rocket_radius**2 + self.reference_length = 2 * self.rocket_radius + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_shape() + + def _set_geometry_attribute(self, attribute, value): + old_value = getattr(self, attribute) + setattr(self, attribute, value) + try: + self._validate_geometry() + except (TypeError, ValueError): + setattr(self, attribute, old_value) + raise + self._evaluate_all() + + @property + def n(self): + """Number of tubes in the set.""" + return self._n + + @n.setter + def n(self, value): + self._set_geometry_attribute("_n", value) + + @property + def length(self): + """Tube length along the rocket axis, in meters.""" + return self._length + + @length.setter + def length(self, value): + self._set_geometry_attribute("_length", value) + + @property + def inner_radius(self): + """Inner tube radius, in meters.""" + return self._inner_radius + + @inner_radius.setter + def inner_radius(self, value): + self._set_geometry_attribute("_inner_radius", value) + + @property + def outer_radius(self): + """Outer tube radius, in meters.""" + return self._outer_radius + + @outer_radius.setter + def outer_radius(self, value): + self._set_geometry_attribute("_outer_radius", value) + + @property + def rocket_radius(self): + """Reference rocket-body radius, in meters.""" + return self._rocket_radius + + @rocket_radius.setter + def rocket_radius(self, value): + self._set_geometry_attribute("_rocket_radius", value) + + @property + def rocket_diameter(self): + """Reference rocket-body diameter, in meters.""" + return 2 * self.rocket_radius + + def evaluate_geometrical_parameters(self): + """Evaluate the ring-airfoil aspect ratio and tube spacing.""" + self.aspect_ratio = 2 * self.inner_radius / self.length + self.touching_outer_radius = self._touching_outer_radius( + self.n, self.rocket_radius + ) + self.tube_separation = 2 * (self.touching_outer_radius - self.outer_radius) + + def evaluate_center_of_pressure(self): + """Set the subsonic center of pressure at the quarter chord.""" + self.cpx = 0 + self.cpy = 0 + self.cpz = self.length / 4 + self.cp = (self.cpx, self.cpy, self.cpz) + + def evaluate_lift_coefficient(self): + """Evaluate the Ribner normal-force derivative for the tube set.""" + modified_aspect_ratio = 2 * self.aspect_ratio / np.pi + single_tube_constant = ( + 2 + * (modified_aspect_ratio / (1 + modified_aspect_ratio)) + * np.pi**2 + * self.inner_radius + * self.length + ) + clalpha_value = self.n * single_tube_constant / self.reference_area + + self.clalpha = Function( + lambda mach: clalpha_value, + "Mach", + f"Lift coefficient derivative for {self.name}", + ) + self.cl = Function( + lambda alpha, mach: ( + self.clalpha(mach) * np.clip(alpha, -self.stall_angle, self.stall_angle) + ), + ["Alpha (rad)", "Mach"], + "Lift coefficient", + ) + return self.cl + + def evaluate_shape(self): + """Store a side-view outline for plotting the tube-fin envelope.""" + lower = self.rocket_radius + upper = self.rocket_radius + 2 * self.outer_radius + self.shape_vec = [ + np.array([0, self.length, self.length, 0, 0]), + np.array([lower, lower, upper, upper, lower]), + ] + + def info(self): + """Print tube-fin geometry and lift information.""" + self.prints.geometry() + self.prints.lift() + + def all_info(self): + """Print and plot all available tube-fin information.""" + self.prints.all() + self.plots.all() + + def draw(self, *, filename=None): + """Draw a side-view envelope of the tube-fin set.""" + return self.plots.draw(filename=filename) + + def to_dict(self, **kwargs): + data = { + "n": self.n, + "length": self.length, + "inner_radius": self.inner_radius, + "outer_radius": self.outer_radius, + "rocket_radius": self.rocket_radius, + "name": self.name, + } + + if kwargs.get("include_outputs", False): + clalpha = self.clalpha + cl = self.cl + if kwargs.get("discretize", False): + clalpha = clalpha.set_discrete(0, 0.5, 10, mutate_self=False) + cl = cl.set_discrete( + (-self.stall_angle, 0), + (self.stall_angle, 0.5), + (10, 10), + mutate_self=False, + ) + data.update( + { + "aspect_ratio": self.aspect_ratio, + "cp": self.cp, + "clalpha": clalpha, + "cl": cl, + "reference_area": self.reference_area, + "reference_length": self.reference_length, + } + ) + + return data + + @classmethod + def from_dict(cls, data): + return cls( + n=data["n"], + length=data["length"], + inner_radius=data["inner_radius"], + outer_radius=data["outer_radius"], + rocket_radius=data["rocket_radius"], + name=data["name"], + ) diff --git a/rocketpy/rocket/rocket.py b/rocketpy/rocket/rocket.py index 3f4748ac0..ba4bc9fcf 100644 --- a/rocketpy/rocket/rocket.py +++ b/rocketpy/rocket/rocket.py @@ -27,6 +27,7 @@ RailButtons, Tail, TrapezoidalFins, + TubeFins, ) from rocketpy.rocket.aero_surface.fins.elliptical_fin import EllipticalFin from rocketpy.rocket.aero_surface.fins.free_form_fin import FreeFormFin @@ -490,6 +491,11 @@ def fins(self): """A list containing all the fins currently added to the rocket.""" return self.aerodynamic_surfaces.get_by_type(Fins) + @property + def tube_fins(self): + """A list containing all tube-fin sets currently added to the rocket.""" + return self.aerodynamic_surfaces.get_by_type(TubeFins) + @property def tails(self): """A list with all the tails currently added to the rocket""" @@ -1190,6 +1196,7 @@ def add_surfaces(self, surfaces, positions): For Fins type, position refers to the z-coordinate of the root chord leading-edge point closest to the nose cone, before any cant-angle offset is considered. + For TubeFins type, position refers to the leading edge of the tubes. For Tail type, position is relative to the point belonging to the tail which is highest in the rocket coordinate system. For RailButtons type, position is relative to the lower rail button. @@ -1633,6 +1640,67 @@ def add_free_form_fins( self.add_surfaces(fin_set, position) return fin_set + def add_tube_fins( + self, + n, + length, + inner_radius, + outer_radius, + position, + radius=None, + name="Tube Fins", + ): + """Create and add a symmetric set of tube fins to the rocket. + + This first-order model uses the Ribner ring-airfoil normal-force slope + and a fixed quarter-chord center of pressure. It is intended for Mach + numbers up to 0.5 and angles of attack up to 20 degrees. + + Parameters + ---------- + n : int + Number of tubes. Must be at least 3. + length : int, float + Tube length along the rocket axis, in meters. + inner_radius : int, float + Inner radius of each tube, in meters. + outer_radius : int, float + Outer radius of each tube, in meters. The current model requires + neighboring tubes to touch, so this must equal + ``radius * sin(pi / n) / (1 - sin(pi / n))``. + position : int, float + Axial position of the tube leading edges in the user-defined rocket + coordinate system. + radius : int, float, optional + Rocket-body radius where the tubes are mounted. If ``None``, the + rocket radius is used. + name : str, optional + Name of the tube-fin set. Default is ``"Tube Fins"``. + + Returns + ------- + TubeFins + Tube-fin set created and added to the rocket. + + Notes + ----- + Only uncanted, mutually tangent tubes are supported. Component drag, + roll, side-force, yaw, separated tubes, and overlapping tubes are not + included in this model. Tube-fin drag must be represented in the + rocket's power-on and power-off drag curves. + """ + radius = self.radius if radius is None else radius + tube_fins = TubeFins( + n=n, + length=length, + inner_radius=inner_radius, + outer_radius=outer_radius, + rocket_radius=radius, + name=name, + ) + self.add_surfaces(tube_fins, position) + return tube_fins + def add_parachute( self, name, diff --git a/tests/unit/rocket/aero_surface/test_tube_fins.py b/tests/unit/rocket/aero_surface/test_tube_fins.py new file mode 100644 index 000000000..296a13ea1 --- /dev/null +++ b/tests/unit/rocket/aero_surface/test_tube_fins.py @@ -0,0 +1,141 @@ +import json + +import numpy as np +import pytest + +from rocketpy import TubeFins +from rocketpy._encoders import RocketPyDecoder, RocketPyEncoder + + +@pytest.fixture +def tube_fins(): + return TubeFins( + n=6, + length=0.1, + inner_radius=0.045, + outer_radius=0.05, + rocket_radius=0.05, + ) + + +def test_tube_fins_geometry_and_normal_force_slope(tube_fins): + aspect_ratio = 2 * tube_fins.inner_radius / tube_fins.length + modified_aspect_ratio = 2 * aspect_ratio / np.pi + expected_clalpha = ( + tube_fins.n + * 2 + * modified_aspect_ratio + / (1 + modified_aspect_ratio) + * np.pi**2 + * tube_fins.inner_radius + * tube_fins.length + / (np.pi * tube_fins.rocket_radius**2) + ) + + assert tube_fins.aspect_ratio == pytest.approx(aspect_ratio) + assert tube_fins.cp == pytest.approx((0, 0, tube_fins.length / 4)) + assert tube_fins.reference_area == pytest.approx(np.pi * tube_fins.rocket_radius**2) + assert tube_fins.reference_length == pytest.approx(2 * tube_fins.rocket_radius) + assert tube_fins.tube_separation == pytest.approx(0, abs=1e-12) + assert tube_fins.clalpha(0) == pytest.approx(expected_clalpha) + assert tube_fins.clalpha(0.5) == pytest.approx(expected_clalpha) + + +def test_tube_fins_lift_is_capped_at_twenty_degrees(tube_fins): + capped_lift = tube_fins.clalpha(0) * np.radians(20) + + assert tube_fins.cl(np.radians(10), 0) == pytest.approx(capped_lift / 2) + assert tube_fins.cl(np.radians(30), 0) == pytest.approx(capped_lift) + assert tube_fins.cl(np.radians(-30), 0) == pytest.approx(-capped_lift) + + +@pytest.mark.parametrize( + ("overrides", "message"), + [ + ({"n": 2}, "greater than or equal to 3"), + ({"n": 6.0}, "must be an integer"), + ({"length": 0}, "length.*greater than zero"), + ({"inner_radius": 0}, "inner_radius.*greater than zero"), + ({"inner_radius": 0.05}, "smaller than.*outer_radius"), + ({"outer_radius": 0}, "outer_radius.*greater than zero"), + ({"rocket_radius": 0}, "rocket_radius.*greater than zero"), + ({"length": np.nan}, "length.*finite"), + ({"outer_radius": 0.049}, "separated tube fins"), + ({"outer_radius": 0.06}, "overlapping tube fins"), + ], +) +def test_tube_fins_reject_unsupported_geometry(overrides, message): + parameters = { + "n": 6, + "length": 0.1, + "inner_radius": 0.045, + "outer_radius": 0.05, + "rocket_radius": 0.05, + } + parameters.update(overrides) + + with pytest.raises(ValueError, match=message): + TubeFins(**parameters) + + +def test_tube_fins_setters_update_dependent_values(tube_fins): + initial_clalpha = tube_fins.clalpha(0) + + tube_fins.length = 0.2 + assert tube_fins.cpz == pytest.approx(0.05) + assert tube_fins.aspect_ratio == pytest.approx(0.45) + assert tube_fins.clalpha(0) != pytest.approx(initial_clalpha) + + previous_outer_radius = tube_fins.outer_radius + with pytest.raises(ValueError, match="separated tube fins"): + tube_fins.outer_radius = 0.049 + assert tube_fins.outer_radius == previous_outer_radius + + +def test_tube_fins_add_to_rocket(calisto): + initial_clalpha = calisto.total_lift_coeff_der(0) + tube_fins = calisto.add_tube_fins( + n=6, + length=0.12, + inner_radius=0.055, + outer_radius=calisto.radius, + position=-1.1, + ) + + assert tube_fins in calisto.tube_fins + assert calisto.aerodynamic_surfaces[-1].component is tube_fins + assert calisto.aerodynamic_surfaces[-1].position.z == pytest.approx(-1.1) + assert calisto.total_lift_coeff_der(0) == pytest.approx( + initial_clalpha + tube_fins.clalpha(0) + ) + + +@pytest.mark.parametrize( + ("include_outputs", "discretize"), + [(False, False), (True, False), (True, True)], +) +def test_tube_fins_json_round_trip(tube_fins, include_outputs, discretize): + encoded = json.dumps( + tube_fins, + cls=RocketPyEncoder, + include_outputs=include_outputs, + discretize=discretize, + ) + decoded = json.loads(encoded, cls=RocketPyDecoder) + + assert isinstance(decoded, TubeFins) + assert decoded.n == tube_fins.n + assert decoded.length == pytest.approx(tube_fins.length) + assert decoded.inner_radius == pytest.approx(tube_fins.inner_radius) + assert decoded.outer_radius == pytest.approx(tube_fins.outer_radius) + assert decoded.rocket_radius == pytest.approx(tube_fins.rocket_radius) + assert decoded.cp == pytest.approx(tube_fins.cp) + assert decoded.clalpha(0) == pytest.approx(tube_fins.clalpha(0)) + + +def test_tube_fins_info_and_draw(tube_fins, capsys, monkeypatch): + monkeypatch.setattr("matplotlib.pyplot.show", lambda: None) + + assert tube_fins.info() is None + assert "Number of tubes: 6" in capsys.readouterr().out + assert tube_fins.draw() is None