diff --git a/.pylintrc b/.pylintrc index 4ae88a7c6..a0e9cdb7c 100644 --- a/.pylintrc +++ b/.pylintrc @@ -229,6 +229,11 @@ good-names=FlightPhases, center_of_mass_without_motor_to_CDM, motor_center_of_dry_mass_to_CDM, generic_motor_cesaroni_M1520, + R_phi, + R_delta, + R_pi, + R_uncanted, + R_body_to_fin, # Good variable names regexes, separated by a comma. If names match any regex, # they will always be accepted diff --git a/docs/reference/classes/aero_surfaces/EllipticalFin.rst b/docs/reference/classes/aero_surfaces/EllipticalFin.rst new file mode 100644 index 000000000..3f1403e8e --- /dev/null +++ b/docs/reference/classes/aero_surfaces/EllipticalFin.rst @@ -0,0 +1,5 @@ +EllipticalFin Class +------------------- + +.. autoclass:: rocketpy.EllipticalFin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/Fin.rst b/docs/reference/classes/aero_surfaces/Fin.rst new file mode 100644 index 000000000..b9c70c25e --- /dev/null +++ b/docs/reference/classes/aero_surfaces/Fin.rst @@ -0,0 +1,5 @@ +Fin Class +--------- + +.. autoclass:: rocketpy.Fin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/FreeFormFin.rst b/docs/reference/classes/aero_surfaces/FreeFormFin.rst new file mode 100644 index 000000000..636be2c9e --- /dev/null +++ b/docs/reference/classes/aero_surfaces/FreeFormFin.rst @@ -0,0 +1,5 @@ +FreeFormFin Class +----------------- + +.. autoclass:: rocketpy.FreeFormFin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/FreeFormFins.rst b/docs/reference/classes/aero_surfaces/FreeFormFins.rst new file mode 100644 index 000000000..c641ce156 --- /dev/null +++ b/docs/reference/classes/aero_surfaces/FreeFormFins.rst @@ -0,0 +1,5 @@ +FreeFormFins Class +------------------ + +.. autoclass:: rocketpy.FreeFormFins + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/TrapezoidalFin.rst b/docs/reference/classes/aero_surfaces/TrapezoidalFin.rst new file mode 100644 index 000000000..34fa631df --- /dev/null +++ b/docs/reference/classes/aero_surfaces/TrapezoidalFin.rst @@ -0,0 +1,5 @@ +TrapezoidalFin Class +-------------------- + +.. autoclass:: rocketpy.TrapezoidalFin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/_BaseFin.rst b/docs/reference/classes/aero_surfaces/_BaseFin.rst new file mode 100644 index 000000000..a04b02ee5 --- /dev/null +++ b/docs/reference/classes/aero_surfaces/_BaseFin.rst @@ -0,0 +1,5 @@ +_BaseFin Class +-------------- + +.. autoclass:: rocketpy.rocket.aero_surface.fins._base_fin._BaseFin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/_EllipticalMixin.rst b/docs/reference/classes/aero_surfaces/_EllipticalMixin.rst new file mode 100644 index 000000000..c552ae2bb --- /dev/null +++ b/docs/reference/classes/aero_surfaces/_EllipticalMixin.rst @@ -0,0 +1,5 @@ +_EllipticalMixin Class +---------------------- + +.. autoclass:: rocketpy.rocket.aero_surface.fins._elliptical_mixin._EllipticalMixin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/_FreeFormMixin.rst b/docs/reference/classes/aero_surfaces/_FreeFormMixin.rst new file mode 100644 index 000000000..bf7ff6a61 --- /dev/null +++ b/docs/reference/classes/aero_surfaces/_FreeFormMixin.rst @@ -0,0 +1,5 @@ +_FreeFormMixin Class +-------------------- + +.. autoclass:: rocketpy.rocket.aero_surface.fins._free_form_mixin._FreeFormMixin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/_TrapezoidalMixin.rst b/docs/reference/classes/aero_surfaces/_TrapezoidalMixin.rst new file mode 100644 index 000000000..ddb05ecf3 --- /dev/null +++ b/docs/reference/classes/aero_surfaces/_TrapezoidalMixin.rst @@ -0,0 +1,5 @@ +_TrapezoidalMixin Class +----------------------- + +.. autoclass:: rocketpy.rocket.aero_surface.fins._trapezoidal_mixin._TrapezoidalMixin + :members: \ No newline at end of file diff --git a/docs/reference/classes/aero_surfaces/index.rst b/docs/reference/classes/aero_surfaces/index.rst index c6e6a3efe..979fec555 100644 --- a/docs/reference/classes/aero_surfaces/index.rst +++ b/docs/reference/classes/aero_surfaces/index.rst @@ -11,7 +11,16 @@ AeroSurface Classes Fins TrapezoidalFins EllipticalFins + FreeFormFins + Fin + TrapezoidalFin + EllipticalFin + FreeFormFin RailButtons AirBrakes GenericSurface LinearGenericSurface + _EllipticalMixin + _TrapezoidalMixin + _FreeFormMixin + _BaseFin \ No newline at end of file diff --git a/docs/static/rocket/fin_forces.png b/docs/static/rocket/fin_forces.png new file mode 100644 index 000000000..41be9bfa2 Binary files /dev/null and b/docs/static/rocket/fin_forces.png differ diff --git a/docs/static/rocket/individual_fin_frame.png b/docs/static/rocket/individual_fin_frame.png new file mode 100644 index 000000000..67d417cdb Binary files /dev/null and b/docs/static/rocket/individual_fin_frame.png differ diff --git a/docs/technical/aerodynamics/individual_fins.rst b/docs/technical/aerodynamics/individual_fins.rst new file mode 100644 index 000000000..952772fad --- /dev/null +++ b/docs/technical/aerodynamics/individual_fins.rst @@ -0,0 +1,498 @@ +.. _individual_fins: + +==================== +Individual Fin Model +==================== + +:Author: Mateus Stano Junqueira +:Date: March 2025 + +Introduction +============ + +There are currently three ways to model the aerodynamic effects of fins in RocketPy: + +1. By use of the :class:`rocketpy.GenericSurface` or :class:`rocketpy.LinearGenericSurface` classes, + which are defined by receiving the aerodynamic coefficients directly; +2. By use of the :class:`rocketpy.Fins` classes (:class:`rocketpy.TrapezoidalFins`, :class:`rocketpy.EllipticalFins`, :class:`rocketpy.FreeFormFins`), + which are defined by receiving the geometric parameters of the **fin set** and + calculating the aerodynamic coefficients internally with the Barrowman method; +3. By use of the :class:`rocketpy.Fin` classes (:class:`rocketpy.TrapezoidalFin`, :class:`rocketpy.EllipticalFin`, :class:`rocketpy.FreeFormFin`) + (which are the base classes for the fin set classes), which are defined by receiving the + geometric parameters of **one single fin** and calculating the aerodynamic + coefficients internally with the Barrowman method, with exception of the + moment coefficients, which have been derived in this document. + +This document will focus on all the mathematical derivations and implementations +of the :class:`rocketpy.Fin` classes. + +Fin Coordinate Frame +==================== + +The fin coordinate frame is defined as follows: + +- The origin is at the leading edge of the fin, at the intersection of the + fin's root chord and the fin's leading edge; +- The z-axis is aligned with the fin's root chord, pointing towards the fin's + trailing edge; +- The y-axis is aligned with the fin's span, pointing towards the fin's tip chord; +- The x-axis completes the right-handed coordinate system. + +The fin is postioned at the rocket's body, given a ``rocket_radius`` :math:`r`, +an ``angular_position`` :math:`\Gamma` and a position :math:`p` along the +rocket's body. The angular position is given according to +:ref:`Angular Position Inputs `. + + +The fin can also be rotated by a cant angle :math:`\delta`. + +The image below shows the fin coordinate frame: + +.. image:: ../../static/rocket/individual_fin_frame.png + :width: 800 + :align: center + +.. note:: + A positive cant angle :math:`\delta` produces a negative body axis rolling + moment at zero angle of attack. + +The rotation matrix from the fin coordinate frame to the rocket's body frame is +define by, first a rotation around the y-axis by 180 degrees: + +.. math:: + \mathbf{R}_{y(\pi)} = \begin{bmatrix} + -1 & 0 & 0 \\ + 0 & 1 & 0 \\ + 0 & 0 & -1 + \end{bmatrix} + +Then a rotation around the z-axis by the angle :math:`\Gamma`: + +.. math:: + \mathbf{R}_{z(\Gamma)} = \begin{bmatrix} + \cos(\Gamma) & -\sin(\Gamma) & 0 \\ + \sin(\Gamma) & \cos(\Gamma) & 0 \\ + 0 & 0 & 1 + \end{bmatrix} + +Then a rotation around the y-axis by the cant angle :math:`\delta`: + +.. math:: + \mathbf{R}_{y(\delta)} = \begin{bmatrix} + \cos(\delta) & 0 & \sin(\delta) \\ + 0 & 1 & 0 \\ + -\sin(\delta) & 0 & \cos(\delta) + \end{bmatrix} + +The final rotation matrix is given by: + +.. math:: + \mathbf{R} = \mathbf{R}_{y(\delta)} \cdot \mathbf{R}_{z(\Gamma)} \cdot \mathbf{R}_{y(\pi)} + + +The position of the fin's coordinate frame origin in the rocket's body frame +is calculated by first assuming no a fin frame with no cant angle, then +calculating the position of the fin's leading edge (with cant angle) in this +frame, and finally translating this position to the fin's position in the +rocket's body frame. The position of the fin's real leading edge in this no cant +angle fin frame is given by the point :math:`\mathbf{P}^{\delta}_{le_f}`: + +.. math:: + \mathbf{P}^{\delta}_{le_f} = \begin{bmatrix} + -\frac{Cr}{2} \sin(\delta) \\ + 0 \\ + \frac{Cr}{2} (1 - \cos(\delta)) + \end{bmatrix} + +Then, describing this point to the rocket's body frame orientation (no +translation): + +.. math:: + \mathbf{P}^{\delta}_{le_b} = (\mathbf{R}_{z(\Gamma)} \cdot \mathbf{R}_{y(\pi)}) \cdot \mathbf{P}^{\delta}_{le_f} + +The position of the fin's leading edge with no cant angle in the rocket's body +frame is given by: + +.. math:: + \mathbf{P}^{\overline{\delta}}_{le_b} = \begin{bmatrix} + -r \sin(\Gamma) \\ + r \cos(\Gamma) \\ + p + \end{bmatrix} + +Finally, we add the position of the fin's leading edge with no cant angle to the +position of the fin's leading edge with cant angle in the rocket's body frame: + +.. math:: + \mathbf{P}_{le_b} = \mathbf{P}^{\overline{\delta}}_{le_b} + \mathbf{P}^{\delta}_{le_b} + + +Center of Pressure Position +=========================== + +In the Fin Coordinate Frame, the center of pressure is given by the Barrowman +method, and will here only be defined symbolically: + +.. math:: + \mathbf{cp}_f = \begin{bmatrix} + cp_x \\ + cp_y \\ + cp_z + \end{bmatrix} + +The center of pressure position in the rocket's body frame is given by: + +.. math:: + \mathbf{cp}_{rocket} = \mathbf{R} \cdot \mathbf{cp}_f + \mathbf{P}_{le_b} + +Aerodynamic Forces +================== + +.. note:: + The aerodynamic coefficients are defined according the Barrowman method. + +Given a stream velocity in the fin frame :math:`\mathbf{v}_{0f} = [v_{0x}, v_{0y}, v_{0z}]^{T}`, +the effective angle of attack of the fin is given by: + +.. math:: + \alpha_f = \arctan\left(\frac{v_{0x}}{v_{0z}}\right) + +This can also be seen as the angle between the fin's root chord and the stream +velocity vector in the fin frame. + +The aerodynamic force in the x-direction of the fin is given by: + +.. math:: + F_{x} = \frac{1}{2} \cdot \rho \cdot \|\mathbf{v}_{0f}\|^2 \cdot A_{r} \cdot C_{N}(\alpha_f, Ma) + +Where :math:`A_{r}` is the reference area of the fin, and :math:`C_{N}` is the +normal force coefficient, which is a function of the angle of attack and the +Mach number :math:`Ma`. +This force is then transformed to the rocket's body frame by the rotation matrix: + +.. math:: + \begin{bmatrix} + F_{x} \\ + F_{y} \\ + F_{z} + \end{bmatrix}_{rocket} = \mathbf{R} \cdot \begin{bmatrix} + F_{x} \\ + 0 \\ + 0 + \end{bmatrix}_{fin} + +Then, the moments are calculated by the cross product of the center of pressure +and the aerodynamic force: + +.. math:: + \begin{bmatrix} + M_{x} \\ + M_{y} \\ + M_{z} + \end{bmatrix}_{rocket} = \mathbf{cp}_{rocket} \times \begin{bmatrix} + F_{x} \\ + F_{y} \\ + F_{z} + \end{bmatrix}_{rocket} + +From the Barrowman method, the moment along the center axis of the rocket +(:math:`M_{z}`) is still missing the damping term, which is given by: + +.. math:: + M_{damp} = \frac{1}{2} \cdot \rho \cdot \|v_{0}\| \cdot A_{r} \cdot L_{r}^2 \cdot C_{ld\omega}(Ma) \cdot \frac{1}{2} \cdot \omega_z + +.. math:: + M_{z \, \text{final}} = M_{z} + M_{damp} + +Where :math:`C_{ld}` is the roll moment damping coefficient, :math:`L_{r}` +is the reference length, which is equal to the rocket diamete, and +:math:`\omega_z` is the angular velocity of the rocket around the z-axis. + +Adding Individual Fins to the Rocket +==================================== +.. jupyter-execute:: + :hide-code: + :hide-output: + + from rocketpy import * + env = Environment(latitude=32.990254, longitude=-106.974998, elevation=1400) + Pro75M1670 = SolidMotor( + thrust_source="../data/motors/cesaroni/Cesaroni_M1670.eng", + dry_mass=1.815, + dry_inertia=(0.125, 0.125, 0.002), + nozzle_radius=33 / 1000, + grain_number=5, + grain_density=1815, + grain_outer_radius=33 / 1000, + grain_initial_inner_radius=15 / 1000, + grain_initial_height=120 / 1000, + grain_separation=5 / 1000, + grains_center_of_mass_position=0.397, + center_of_dry_mass_position=0.317, + nozzle_position=0, + burn_time=3.9, + throat_radius=11 / 1000, + coordinate_system_orientation="nozzle_to_combustion_chamber", + ) + # IMPORTANT: modify the file paths below to match your own system + + example_rocket = Rocket( + radius=127 / 2000, + mass=14.426, + inertia=(6.321, 6.321, 0.034), + power_off_drag="../data/rockets/calisto/powerOffDragCurve.csv", + power_on_drag="../data/rockets/calisto/powerOnDragCurve.csv", + center_of_mass_without_motor=0, + coordinate_system_orientation="tail_to_nose", + ) + + rail_buttons = example_rocket.set_rail_buttons( + upper_button_position=0.0818, + lower_button_position=-0.618, + angular_position=45, + ) + example_rocket.add_motor(Pro75M1670, position=-1.255) + nose_cone = example_rocket.add_nose(length=0.55829, kind="vonKarman", position=1.278) + tail = example_rocket.add_tail( + top_radius=0.0635, bottom_radius=0.0435, length=0.060, position=-1.194656 + ) + example_rocket.add_trapezoidal_fins( + n=4, + root_chord=0.120, + tip_chord=0.060, + span=0.110, + cant_angle=0.0, + position=-1.04956, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + +Given a defined ``Rocket`` object, we can add individual fins to the rocket by +using the ``add_surfaces`` method. Here is an example of adding two canards +in the Calisto rocket from the :ref:`First Simulation ` example: + +.. jupyter-execute:: + + canard1 = TrapezoidalFin( + angular_position=0, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + canard2 = TrapezoidalFin( + angular_position=180, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + + # Position along the center axis of the rocket is specified here. + # If different positions are desired, the position can be specified as a list. + example_rocket.add_surfaces([canard1, canard2], positions = 0.35) + + example_rocket.draw(plane="yz") + +.. seealso:: + + There are three classes for defining fins in RocketPy given their geometry: + + - :class:`rocketpy.TrapezoidalFin` - For how to define a trapezoidal fin + - :class:`rocketpy.EllipticalFin` - For how to define an elliptical fin + - :class:`rocketpy.FreeFormFin` - For how to define a free form fin + + + +Fin Force Conventions +===================== + +.. - Explain positive cant angle resultant force +.. - if all fins have positive cant angle then they will all generate a force that will rotate the rocket in the same direction +.. which is negative roll +.. - show what a positive and negative cant angle on oposing fins looks like. (generate pitch moment -> pitch control) +.. - show what a positive and negative cant angle on oposing fins looks like. (generate yaw moment -> yaw control) +.. - example for 4 fins, show how to count the number of fins and how to access each of them + + + +Here we exemplify the fin force conventions relating the cant angle +(deflection angle) of the fins to the pitch, yaw and roll moments. We will +consider a rocket with four fins, to illustrate the concepts. The image below +show the sign convention for the forces acting on the fins, given positive cant +angles: + +.. image:: ../../static/rocket/fin_forces.png + :width: 800 + :align: center + + +Roll +^^^^ + +.. jupyter-execute:: + :hide-code: + :hide-output: + + example_rocket.aerodynamic_surfaces.pop() + example_rocket.aerodynamic_surfaces.pop() + + +A positive cant angle :math:`\delta` produces a negative roll moment at zero +angle of attack. Any fin with a positive cant angle will produce a negative roll +moment, and any fin with a negative cant angle will produce a positive roll +moment. + +Here is a flight of the calisto with canards defined with a positive cant angle: + +.. jupyter-execute:: + + + canard1 = TrapezoidalFin( + angular_position=0, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + canard2 = TrapezoidalFin( + angular_position=180, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + + example_rocket.add_surfaces([canard1, canard2], positions = 0.35) + + test_flight = Flight( + rocket=example_rocket, + environment=env, rail_length=5.2, inclination=85, heading=0, + terminate_on_apogee=True, + ) + + # Rolling Moment + test_flight.M3() + + # Rolling Speed + test_flight.w3() + + # Angle of attack + test_flight.partial_angle_of_attack.plot(test_flight.out_of_rail_time, 5) + + # Angle of sideslip + test_flight.angle_of_sideslip.plot(test_flight.out_of_rail_time, 5) + +Pitch +^^^^^ + +.. jupyter-execute:: + :hide-code: + :hide-output: + + example_rocket.aerodynamic_surfaces.pop() + example_rocket.aerodynamic_surfaces.pop() + +Given canards fins at 90 degrees and 270 degrees, having opposite cant angles, +a positive pitch moment will be generated. The following example shows the +effect of this configuration in the non-zero angle of attack flight of the +rocket: + +.. jupyter-execute:: + + + canard1 = TrapezoidalFin( + angular_position=90, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + canard2 = TrapezoidalFin( + angular_position=270, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=-0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + + example_rocket.add_surfaces([canard1, canard2], positions = 0.35) + + test_flight = Flight( + rocket=example_rocket, + environment=env, rail_length=5.2, inclination=85, heading=0, + terminate_on_apogee=True, + ) + + # Angle of attack + test_flight.partial_angle_of_attack.plot(test_flight.out_of_rail_time, 5) + + # Angle of sideslip + test_flight.angle_of_sideslip.plot(test_flight.out_of_rail_time, 5) + +Yaw +^^^ + +.. jupyter-execute:: + :hide-code: + :hide-output: + + example_rocket.aerodynamic_surfaces.pop() + example_rocket.aerodynamic_surfaces.pop() + +Given oposing canards at 0 degrees and 180 degrees, having opposite cant angles, +a positive yaw moment will be generated. The following example shows the +effect of this configuration in the non-zero angle of attack flight of the +rocket: + +.. jupyter-execute:: + + + canard1 = TrapezoidalFin( + angular_position=0, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + canard2 = TrapezoidalFin( + angular_position=180, + root_chord=0.060, + tip_chord=0.020, + span=0.03, + rocket_radius=example_rocket.radius, + cant_angle=-0.5, + airfoil=("../data/airfoils/NACA0012-radians.txt", "radians"), + ) + + example_rocket.add_surfaces([canard1, canard2], positions = 0.35) + + test_flight = Flight( + rocket=example_rocket, + environment=env, rail_length=5.2, inclination=85, heading=0, + terminate_on_apogee=True, + ) + + # Angle of attack + test_flight.partial_angle_of_attack.plot(test_flight.out_of_rail_time, 5) + + # Angle of sideslip + test_flight.angle_of_sideslip.plot(test_flight.out_of_rail_time, 5) + + + + + diff --git a/docs/technical/index.rst b/docs/technical/index.rst index b96e611ed..73583eba9 100644 --- a/docs/technical/index.rst +++ b/docs/technical/index.rst @@ -13,6 +13,7 @@ in their code. Equations of Motion v0 Equations of Motion v1 Elliptical Fins + Individual Fin Roll Moment Sensitivity Analysis References diff --git a/docs/user/rocket/rocket_axes.rst b/docs/user/rocket/rocket_axes.rst index 56fe17cad..aa99cb231 100644 --- a/docs/user/rocket/rocket_axes.rst +++ b/docs/user/rocket/rocket_axes.rst @@ -61,6 +61,8 @@ The following figure shows the two possibilities for the user input coordinate s the same as the **Body Axes Coordinate System**. The origin of the coordinate \ system may still be different. +.. _angular_position: + Angular Position Inputs ~~~~~~~~~~~~~~~~~~~~~~~ diff --git a/rocketpy/__init__.py b/rocketpy/__init__.py index f99a70f28..abd3d6bc1 100644 --- a/rocketpy/__init__.py +++ b/rocketpy/__init__.py @@ -28,8 +28,11 @@ AeroSurface, AirBrakes, Components, + EllipticalFin, EllipticalFins, + Fin, Fins, + FreeFormFin, FreeFormFins, GenericSurface, LinearGenericSurface, @@ -38,6 +41,7 @@ RailButtons, Rocket, Tail, + TrapezoidalFin, TrapezoidalFins, ) from .sensitivity import SensitivityModel diff --git a/rocketpy/plots/aero_surface_plots.py b/rocketpy/plots/aero_surface_plots.py index d501d1352..f599b8ffc 100644 --- a/rocketpy/plots/aero_surface_plots.py +++ b/rocketpy/plots/aero_surface_plots.py @@ -1,3 +1,5 @@ +# pylint: disable=too-many-statements + from abc import ABC, abstractmethod import matplotlib.pyplot as plt @@ -139,10 +141,6 @@ class _FinsPlots(_AeroSurfacePlots): """Abstract class that contains all fin plots. This class inherits from the _AeroSurfacePlots class.""" - @abstractmethod - def draw(self, *, filename=None): - pass - def airfoil(self): """Plots the airfoil information when the fin has an airfoil shape. If the fin does not have an airfoil shape, this method does nothing. @@ -164,7 +162,6 @@ def roll(self): None """ print("Roll parameters:") - # TODO: lacks a title in the plots self.aero_surface.roll_parameters[0]() self.aero_surface.roll_parameters[1]() @@ -197,10 +194,65 @@ def all(self): self.lift() +class _FinPlots(_AeroSurfacePlots): + """Abstract class that contains all fin plots. This class inherits from the + _AeroSurfacePlots class.""" + + def airfoil(self): + """Plots the airfoil information when the fin has an airfoil shape. If + the fin does not have an airfoil shape, this method does nothing. + + Returns + ------- + None + """ + + if self.aero_surface.airfoil: + print("Airfoil lift curve:") + self.aero_surface.airfoil_cl.plot_1d(force_data=True) + + def roll(self): + """Plots the roll parameters of the fin set. + + Returns + ------- + None + """ + print("Roll parameters:") + self.aero_surface.roll_parameters[0]() + self.aero_surface.roll_parameters[1]() + + def lift(self): + """Plots the lift coefficient of the aero surface as a function of Mach + and the angle of attack. A 3D plot is expected. See the rocketpy.Function + class for more information on how this plot is made. Also, this method + plots the lift coefficient considering a single fin and the lift + coefficient considering all fins. + + Returns + ------- + None + """ + print("Lift coefficient:") + self.aero_surface.cl() + self.aero_surface.clalpha_single_fin() + + def all(self): + """Plots all available fin plots. + + Returns + ------- + None + """ + self.draw() + self.airfoil() + self.roll() + self.lift() + + class _TrapezoidalFinsPlots(_FinsPlots): """Class that contains all trapezoidal fin plots.""" - # pylint: disable=too-many-statements def draw(self, *, filename=None): """Draw the fin shape along with some important information, including the center line, the quarter line and the center of pressure position. @@ -326,10 +378,129 @@ def draw(self, *, filename=None): show_or_save_plot(filename) +class _TrapezoidalFinPlots(_FinPlots): + """Class that contains all trapezoidal fin plots.""" + + def draw(self, *, filename=None): + """Draw the fin shape along with some important information, including + the center line, the quarter line and the center of pressure position. + + Returns + ------- + None + """ + # Color cycle [#348ABD, #A60628, #7A68A6, #467821, #D55E00, #CC79A7, + # #56B4E9, #009E73, #F0E442, #0072B2] + # Fin + leading_edge = plt.Line2D( + (0, self.aero_surface.sweep_length), + (0, self.aero_surface.span), + color="#A60628", + ) + tip = plt.Line2D( + ( + self.aero_surface.sweep_length, + self.aero_surface.sweep_length + self.aero_surface.tip_chord, + ), + (self.aero_surface.span, self.aero_surface.span), + color="#A60628", + ) + back_edge = plt.Line2D( + ( + self.aero_surface.sweep_length + self.aero_surface.tip_chord, + self.aero_surface.root_chord, + ), + (self.aero_surface.span, 0), + color="#A60628", + ) + root = plt.Line2D((self.aero_surface.root_chord, 0), (0, 0), color="#A60628") + + # Center and Quarter line + center_line = plt.Line2D( + ( + self.aero_surface.root_chord / 2, + self.aero_surface.sweep_length + self.aero_surface.tip_chord / 2, + ), + (0, self.aero_surface.span), + color="#7A68A6", + alpha=0.35, + linestyle="--", + label="Center Line", + ) + quarter_line = plt.Line2D( + ( + self.aero_surface.root_chord / 4, + self.aero_surface.sweep_length + self.aero_surface.tip_chord / 4, + ), + (0, self.aero_surface.span), + color="#7A68A6", + alpha=1, + linestyle="--", + label="Quarter Line", + ) + + # Center of pressure + cp_point = [self.aero_surface.cpz, self.aero_surface.Yma] + + # Mean Aerodynamic Chord + yma_start = ( + self.aero_surface.sweep_length + * (self.aero_surface.root_chord + 2 * self.aero_surface.tip_chord) + / (3 * (self.aero_surface.root_chord + self.aero_surface.tip_chord)) + ) + yma_end = ( + 2 * self.aero_surface.root_chord**2 + + self.aero_surface.root_chord * self.aero_surface.sweep_length + + 2 * self.aero_surface.root_chord * self.aero_surface.tip_chord + + 2 * self.aero_surface.sweep_length * self.aero_surface.tip_chord + + 2 * self.aero_surface.tip_chord**2 + ) / (3 * (self.aero_surface.root_chord + self.aero_surface.tip_chord)) + yma_line = plt.Line2D( + (yma_start, yma_end), + (self.aero_surface.Yma, self.aero_surface.Yma), + color="#467821", + linestyle="--", + label="Mean Aerodynamic Chord", + ) + + # Plotting + fig = plt.figure(figsize=(7, 4)) + with plt.style.context("bmh"): + ax = fig.add_subplot(111) + + # Fin + ax.add_line(leading_edge) + ax.add_line(tip) + ax.add_line(back_edge) + ax.add_line(root) + + ax.add_line(center_line) + ax.add_line(quarter_line) + ax.add_line(yma_line) + ax.scatter(*cp_point, label="Center of Pressure", color="red", s=100, zorder=10) + ax.scatter(*cp_point, facecolors="none", edgecolors="red", s=500, zorder=10) + + # Plot settings + xlim = ( + self.aero_surface.root_chord + if self.aero_surface.sweep_length + self.aero_surface.tip_chord + < self.aero_surface.root_chord + else self.aero_surface.sweep_length + self.aero_surface.tip_chord + ) + ax.set_xlim(0, xlim * 1.1) + ax.set_ylim(0, self.aero_surface.span * 1.1) + ax.set_xlabel("Root chord (m)") + ax.set_ylabel("Span (m)") + ax.set_title("Trapezoidal Fin Cross Section") + ax.legend(bbox_to_anchor=(1.05, 1.0), loc="upper left") + + plt.tight_layout() + plt.show() + + class _EllipticalFinsPlots(_FinsPlots): """Class that contains all elliptical fin plots.""" - # pylint: disable=too-many-statements def draw(self, *, filename=None): """Draw the fin shape along with some important information. These being: the center line and the center of pressure position. @@ -405,10 +576,153 @@ def draw(self, *, filename=None): show_or_save_plot(filename) +class _EllipticalFinPlots(_FinPlots): + """Class that contains all elliptical fin plots.""" + + def draw(self, *, filename=None): + """Draw the fin shape along with some important information. + These being: the center line and the center of pressure position. + + Returns + ------- + None + """ + # Ellipse + ellipse = Ellipse( + (self.aero_surface.root_chord / 2, 0), + self.aero_surface.root_chord, + self.aero_surface.span * 2, + fill=False, + edgecolor="#A60628", + linewidth=2, + ) + + # Mean Aerodynamic Chord # From Barrowman's theory + yma_length = 8 * self.aero_surface.root_chord / (3 * np.pi) + yma_start = (self.aero_surface.root_chord - yma_length) / 2 + yma_end = ( + self.aero_surface.root_chord + - (self.aero_surface.root_chord - yma_length) / 2 + ) + yma_line = plt.Line2D( + (yma_start, yma_end), + (self.aero_surface.Yma, self.aero_surface.Yma), + label="Mean Aerodynamic Chord", + color="#467821", + ) + + # Center Line + center_line = plt.Line2D( + (self.aero_surface.root_chord / 2, self.aero_surface.root_chord / 2), + (0, self.aero_surface.span), + color="#7A68A6", + alpha=0.35, + linestyle="--", + label="Center Line", + ) + + # Center of pressure + cp_point = [self.aero_surface.cpz, self.aero_surface.Yma] + + # Plotting + fig = plt.figure(figsize=(7, 4)) + with plt.style.context("bmh"): + ax = fig.add_subplot(111) + ax.add_patch(ellipse) + ax.add_line(yma_line) + ax.add_line(center_line) + ax.scatter(*cp_point, label="Center of Pressure", color="red", s=100, zorder=10) + ax.scatter(*cp_point, facecolors="none", edgecolors="red", s=500, zorder=10) + + # Plot settings + ax.set_xlim(0, self.aero_surface.root_chord) + ax.set_ylim(0, self.aero_surface.span * 1.1) + ax.set_xlabel("Root chord (m)") + ax.set_ylabel("Span (m)") + ax.set_title("Elliptical Fin Cross Section") + ax.legend(bbox_to_anchor=(1.05, 1.0), loc="upper left") + + plt.tight_layout() + plt.show() + + class _FreeFormFinsPlots(_FinsPlots): """Class that contains all free form fin plots.""" - # pylint: disable=too-many-statements + def draw(self, *, filename=None): + """Draw the fin shape along with some important information. + These being: the center line and the center of pressure position. + + Parameters + ---------- + filename : str | None, optional + The path the plot should be saved to. By default None, in which case + the plot will be shown instead of saved. Supported file endings are: + eps, jpg, jpeg, pdf, pgf, png, ps, raw, rgba, svg, svgz, tif, tiff + and webp (these are the formats supported by matplotlib). + + Returns + ------- + None + """ + # Color cycle [#348ABD, #A60628, #7A68A6, #467821, #D55E00, #CC79A7, + # #56B4E9, #009E73, #F0E442, #0072B2] + + # Center of pressure + cp_point = [self.aero_surface.cpz, self.aero_surface.Yma] + + # Mean Aerodynamic Chord + yma_line = plt.Line2D( + ( + self.aero_surface.mac_lead, + self.aero_surface.mac_lead + self.aero_surface.mac_length, + ), + (self.aero_surface.Yma, self.aero_surface.Yma), + color="#467821", + linestyle="--", + label="Mean Aerodynamic Chord", + ) + + # Plotting + fig = plt.figure(figsize=(7, 4)) + with plt.style.context("bmh"): + ax = fig.add_subplot(111) + + # Fin + ax.scatter( + self.aero_surface.shape_vec[0], + self.aero_surface.shape_vec[1], + color="#A60628", + ) + ax.plot( + self.aero_surface.shape_vec[0], + self.aero_surface.shape_vec[1], + color="#A60628", + ) + # line from the last point to the first point + ax.plot( + [self.aero_surface.shape_vec[0][-1], self.aero_surface.shape_vec[0][0]], + [self.aero_surface.shape_vec[1][-1], self.aero_surface.shape_vec[1][0]], + color="#A60628", + ) + + ax.add_line(yma_line) + ax.scatter(*cp_point, label="Center of Pressure", color="red", s=100, zorder=10) + ax.scatter(*cp_point, facecolors="none", edgecolors="red", s=500, zorder=10) + + # Plot settings + ax.set_xlabel("Root chord (m)") + ax.set_ylabel("Span (m)") + ax.set_title("Free Form Fin Cross Section") + ax.legend(bbox_to_anchor=(1.05, 1.0), loc="upper left") + + plt.tight_layout() + show_or_save_plot(filename) + + +class _FreeFormFinPlots(_FinPlots): + """Class that contains all free form fin plots.""" + def draw(self, *, filename=None): """Draw the fin shape along with some important information. These being: the center line and the center of pressure position. diff --git a/rocketpy/plots/rocket_plots.py b/rocketpy/plots/rocket_plots.py index 8eaaded16..46c4b2669 100644 --- a/rocketpy/plots/rocket_plots.py +++ b/rocketpy/plots/rocket_plots.py @@ -1,8 +1,9 @@ import matplotlib.pyplot as plt import numpy as np +from rocketpy.mathutils.vector_matrix import Vector from rocketpy.motors import EmptyMotor, HybridMotor, LiquidMotor, SolidMotor -from rocketpy.rocket.aero_surface import Fins, NoseCone, Tail +from rocketpy.rocket.aero_surface import Fin, Fins, NoseCone, Tail from rocketpy.rocket.aero_surface.generic_surface import GenericSurface from .plot_helpers import show_or_save_plot @@ -177,7 +178,7 @@ def draw(self, vis_args=None, plane="xz", *, filename=None): and webp (these are the formats supported by matplotlib). """ - self.__validate_aerodynamic_surfaces() + self.__validate_aerodynamic_surfaces(plane) if vis_args is None: vis_args = { @@ -197,9 +198,9 @@ def draw(self, vis_args=None, plane="xz", *, filename=None): csys = self.rocket._csys reverse = csys == 1 - self.rocket.aerodynamic_surfaces.sort_by_position(reverse=reverse) + surfaces = self.rocket.aerodynamic_surfaces.sort_by_position(reverse=reverse) - drawn_surfaces = self._draw_aerodynamic_surfaces(ax, vis_args, plane) + drawn_surfaces = self._draw_aerodynamic_surfaces(ax, vis_args, plane, surfaces) last_radius, last_x = self._draw_tubes(ax, drawn_surfaces, vis_args) self._draw_motor(last_radius, last_x, ax, vis_args) self._draw_rail_buttons(ax, vis_args) @@ -215,13 +216,15 @@ def draw(self, vis_args=None, plane="xz", *, filename=None): plt.tight_layout() show_or_save_plot(filename) - def __validate_aerodynamic_surfaces(self): + def __validate_aerodynamic_surfaces(self, plane): if not self.rocket.aerodynamic_surfaces: raise ValueError( "The rocket must have at least one aerodynamic surface to be drawn." ) + if plane not in ("xz", "yz"): + raise ValueError("The plane must be 'xz' or 'yz'. The default is 'xz'.") - def _draw_aerodynamic_surfaces(self, ax, vis_args, plane): + def _draw_aerodynamic_surfaces(self, ax, vis_args, plane, surfaces): """Draws the aerodynamic surfaces and saves the position of the points of interest for the tubes.""" # List of drawn surfaces with the position of points of interest @@ -234,13 +237,17 @@ def _draw_aerodynamic_surfaces(self, ax, vis_args, plane): # diameter changes. The final point of the last surface is the final # point of the last tube - for surface, position in self.rocket.aerodynamic_surfaces: + for surface, position in surfaces: if isinstance(surface, NoseCone): self._draw_nose_cone(ax, surface, position.z, drawn_surfaces, vis_args) elif isinstance(surface, Tail): self._draw_tail(ax, surface, position.z, drawn_surfaces, vis_args) elif isinstance(surface, Fins): - self._draw_fins(ax, surface, position.z, drawn_surfaces, vis_args) + self._draw_fins( + ax, surface, position.z, drawn_surfaces, vis_args, plane + ) + elif isinstance(surface, Fin): + self._draw_fin(ax, surface, position, drawn_surfaces, vis_args, plane) elif isinstance(surface, GenericSurface): self._draw_generic_surface( ax, surface, position, drawn_surfaces, vis_args, plane @@ -303,13 +310,15 @@ def _draw_tail(self, ax, surface, position, drawn_surfaces, vis_args): # Add the tail to the list of drawn surfaces drawn_surfaces.append((surface, position, surface.bottom_radius, x_tail[-1])) - def _draw_fins(self, ax, surface, position, drawn_surfaces, vis_args): + def _draw_fins(self, ax, surface, position, drawn_surfaces, vis_args, plane): """Draws the fins and saves the position of the points of interest for the tubes.""" num_fins = surface.n x_fin = -self.rocket._csys * surface.shape_vec[0] + position y_fin = surface.shape_vec[1] + surface.rocket_radius - rotation_angles = [2 * np.pi * i / num_fins for i in range(num_fins)] + rotation_angles = np.array([2 * np.pi * i / num_fins for i in range(num_fins)]) + if plane == "xz": + rotation_angles -= np.pi / 2 for angle in rotation_angles: # Create a rotation matrix for the current angle around the x-axis @@ -321,13 +330,6 @@ def _draw_fins(self, ax, surface, position, drawn_surfaces, vis_args): # Extract x and y coordinates of the rotated points x_rotated, y_rotated = rotated_points_2d - # Project points above the XY plane back into the XY plane (set z-coordinate to 0) - x_rotated = np.where( - rotated_points_2d[1] > 0, rotated_points_2d[0], x_rotated - ) - y_rotated = np.where( - rotated_points_2d[1] > 0, rotated_points_2d[1], y_rotated - ) ax.plot( x_rotated, y_rotated, @@ -337,6 +339,56 @@ def _draw_fins(self, ax, surface, position, drawn_surfaces, vis_args): drawn_surfaces.append((surface, position, surface.rocket_radius, x_rotated[-1])) + def _draw_fin(self, ax, surface, position, drawn_surfaces, vis_args, plane): + """Draws individual fins.""" + + # Get shape vec + xs = surface.shape_vec[0] + ys = surface.shape_vec[1] + zs = np.zeros_like(xs) + + # Define shape in fin coordinate system + x_fin = -zs + y_fin = ys + z_fin = xs + points = np.column_stack((x_fin, y_fin, z_fin)) + + # Move drawing coordinates to center of fin for cant angle rotation + xd = np.array([0, 0, max(xs) / 2]) + points -= xd + + # Rotate to body coordinate system + for i, p in enumerate(points): + points[i] = surface._rotation_fin_to_body @ Vector(p) + + rotated_xd = surface._rotation_fin_to_body @ Vector(xd) + points += np.array(rotated_xd) + + # Back to the drawing system + x_fin_rotated = points[:, 0] + y_fin_rotated = points[:, 1] + z_fin_rotated = points[:, 2] + + if plane == "xz": + x_rotated = self.rocket._csys * z_fin_rotated + position.z + y_rotated = x_fin_rotated + position.x + elif plane == "yz": + x_rotated = self.rocket._csys * z_fin_rotated + position.z + y_rotated = y_fin_rotated + position.y + else: # pragma: no cover + raise ValueError("Plane must be 'xz' or 'yz'.") + + ax.plot( + x_rotated, + y_rotated, + color=vis_args["fins"], + linewidth=vis_args["line_width"], + ) + + drawn_surfaces.append( + (surface, position.z, surface.rocket_radius, x_rotated[-1]) + ) + def _draw_generic_surface( self, ax, @@ -358,7 +410,7 @@ def _draw_generic_surface( x_pos = position[2] # y position of the surface is the y position in the plot y_pos = position[1] - else: # pragma: no cover + else: raise ValueError("Plane must be 'xz' or 'yz'.") ax.scatter( @@ -438,7 +490,7 @@ def _draw_motor(self, last_radius, last_x, ax, vis_args): self._draw_nozzle_tube(last_radius, last_x, nozzle_position, ax, vis_args) - def _generate_motor_patches(self, total_csys, ax): # pylint: disable=unused-argument + def _generate_motor_patches(self, total_csys, ax): """Generates motor patches for drawing""" motor_patches = [] @@ -646,7 +698,7 @@ def all(self): # Rocket draw if len(self.rocket.aerodynamic_surfaces) > 0: - print("\nRocket Draw") + print("\nRocket Drawing") print("-" * 40) self.draw() diff --git a/rocketpy/prints/aero_surface_prints.py b/rocketpy/prints/aero_surface_prints.py index 4eb42b08d..b2c81eaf9 100644 --- a/rocketpy/prints/aero_surface_prints.py +++ b/rocketpy/prints/aero_surface_prints.py @@ -156,9 +156,104 @@ def lift(self): "Lift Coefficient derivative (single fin) at Mach 0 and AoA 0: " f"{self.aero_surface.clalpha_single_fin(0):.3f}" ) + + def all(self): + """Prints all information of the fin set. + + Returns + ------- + None + """ + self.identity() + self.geometry() + self.airfoil() + self.roll() + self.lift() + + +class _FinPrints(_AeroSurfacePrints): + def geometry(self): + print("Geometric information of the fin set:") + print("-------------------------------------") + print(f"Reference rocket radius: {self.aero_surface.rocket_radius:.3f} m") + try: + print(f"Tip chord: {self.aero_surface.tip_chord:.3f} m") + except AttributeError: + pass # it isn't a trapezoidal fin, just don't worry about tip chord + print(f"Root chord: {self.aero_surface.root_chord:.3f} m") + print(f"Span: {self.aero_surface.span:.3f} m") print( - "Lift Coefficient derivative (fin set) at Mach 0 and AoA 0: " - f"{self.aero_surface.clalpha_multiple_fins(0):.3f}" + f"Cant angle: {self.aero_surface.cant_angle:.3f} ° or " + f"{self.aero_surface.cant_angle_rad:.3f} rad" + ) + print(f"Longitudinal section area: {self.aero_surface.Af:.3f} m²") + print(f"Aspect ratio: {self.aero_surface.AR:.3f} ") + print(f"Gamma_c: {self.aero_surface.gamma_c:.3f} m") + print(f"Mean aerodynamic chord: {self.aero_surface.Yma:.3f} m\n") + + def airfoil(self): + """Prints out airfoil related information of the fin set. + + Returns + ------- + None + """ + if self.aero_surface.airfoil: + print("Airfoil information:") + print("--------------------") + print( + "Number of points defining the lift curve: " + f"{len(self.aero_surface.airfoil_cl.x_array)}" + ) + print( + "Lift coefficient derivative at Mach 0 and AoA 0: " + f"{self.aero_surface.clalpha(0):.5f} 1/rad\n" + ) + + def roll(self): + """Prints out information about roll parameters + of the fin set. + + Returns + ------- + None + """ + print("Roll information of the fin set:") + print("--------------------------------") + print( + f"Geometric constant: {self.aero_surface.roll_geometrical_constant:.3f} m" + ) + print( + "Damping interference factor: " + f"{self.aero_surface.roll_damping_interference_factor:.3f} rad" + ) + print( + "Forcing interference factor: " + f"{self.aero_surface.roll_forcing_interference_factor:.3f} rad\n" + ) + + def lift(self): + """Prints out information about lift parameters + of the fin set. + + Returns + ------- + None + """ + print("Lift information of the fin set:") + print("--------------------------------") + print( + "Lift interference factor: " + f"{self.aero_surface.lift_interference_factor:.3f} m" + ) + print( + "Center of Pressure position in local coordinates: " + f"({self.aero_surface.cpx:.3f}, {self.aero_surface.cpy:.3f}, " + f"{self.aero_surface.cpz:.3f})" + ) + print( + "Lift Coefficient derivative (single fin) at Mach 0 and AoA 0: " + f"{self.aero_surface.clalpha_single_fin(0):.3f}" ) def all(self): @@ -179,14 +274,26 @@ class _TrapezoidalFinsPrints(_FinsPrints): """Class that contains all trapezoidal fins prints.""" +class _TrapezoidalFinPrints(_FinPrints): + """Class that contains all trapezoidal fin prints.""" + + class _EllipticalFinsPrints(_FinsPrints): """Class that contains all elliptical fins prints.""" +class _EllipticalFinPrints(_FinPrints): + """Class that contains all elliptical fin prints.""" + + class _FreeFormFinsPrints(_FinsPrints): """Class that contains all free form fins prints.""" +class _FreeFormFinPrints(_FinPrints): + """Class that contains all free form fins prints.""" + + class _TailPrints(_AeroSurfacePrints): """Class that contains all tail prints.""" diff --git a/rocketpy/rocket/__init__.py b/rocketpy/rocket/__init__.py index 0687b5ee5..87ea03c47 100644 --- a/rocketpy/rocket/__init__.py +++ b/rocketpy/rocket/__init__.py @@ -2,14 +2,18 @@ from rocketpy.rocket.aero_surface import ( AeroSurface, AirBrakes, + EllipticalFin, EllipticalFins, + Fin, Fins, + FreeFormFin, FreeFormFins, GenericSurface, LinearGenericSurface, NoseCone, RailButtons, Tail, + TrapezoidalFin, TrapezoidalFins, ) from rocketpy.rocket.components import Components diff --git a/rocketpy/rocket/aero_surface/__init__.py b/rocketpy/rocket/aero_surface/__init__.py index ad784f8d0..7634d3500 100644 --- a/rocketpy/rocket/aero_surface/__init__.py +++ b/rocketpy/rocket/aero_surface/__init__.py @@ -1,9 +1,13 @@ from rocketpy.rocket.aero_surface.aero_surface import AeroSurface from rocketpy.rocket.aero_surface.air_brakes import AirBrakes from rocketpy.rocket.aero_surface.fins import ( + EllipticalFin, EllipticalFins, + Fin, Fins, + FreeFormFin, FreeFormFins, + TrapezoidalFin, TrapezoidalFins, ) from rocketpy.rocket.aero_surface.generic_surface import GenericSurface diff --git a/rocketpy/rocket/aero_surface/aero_surface.py b/rocketpy/rocket/aero_surface/aero_surface.py index 15ca14f1d..6727476c6 100644 --- a/rocketpy/rocket/aero_surface/aero_surface.py +++ b/rocketpy/rocket/aero_surface/aero_surface.py @@ -2,6 +2,8 @@ import numpy as np +from rocketpy.mathutils.vector_matrix import Matrix + class AeroSurface(ABC): """Abstract class used to define aerodynamic surfaces.""" @@ -15,6 +17,14 @@ def __init__(self, name, reference_area, reference_length): self.cpy = 0 self.cpz = 0 + self._rotation_surface_to_body = Matrix( + [ + [-1, 0, 0], + [0, 1, 0], + [0, 0, -1], + ] + ) + @staticmethod def _beta(mach): """Defines a parameter that is often used in aerodynamic @@ -130,7 +140,7 @@ def compute_forces_and_moments( R1, R2, R3, M1, M2, M3 = 0, 0, 0, 0, 0, 0 cpz = cp[2] stream_vx, stream_vy, stream_vz = stream_velocity - if stream_vx**2 + stream_vy**2 != 0: # TODO: maybe try/except + if stream_vx**2 + stream_vy**2 != 0: # Normalize component stream velocity in body frame stream_vzn = stream_vz / stream_speed if -1 * stream_vzn < 1: diff --git a/rocketpy/rocket/aero_surface/fins/__init__.py b/rocketpy/rocket/aero_surface/fins/__init__.py index 941aa5465..dd678c625 100644 --- a/rocketpy/rocket/aero_surface/fins/__init__.py +++ b/rocketpy/rocket/aero_surface/fins/__init__.py @@ -1,4 +1,8 @@ +from rocketpy.rocket.aero_surface.fins.elliptical_fin import EllipticalFin from rocketpy.rocket.aero_surface.fins.elliptical_fins import EllipticalFins +from rocketpy.rocket.aero_surface.fins.fin import Fin from rocketpy.rocket.aero_surface.fins.fins import Fins +from rocketpy.rocket.aero_surface.fins.free_form_fin import FreeFormFin from rocketpy.rocket.aero_surface.fins.free_form_fins import FreeFormFins +from rocketpy.rocket.aero_surface.fins.trapezoidal_fin import TrapezoidalFin from rocketpy.rocket.aero_surface.fins.trapezoidal_fins import TrapezoidalFins diff --git a/rocketpy/rocket/aero_surface/fins/_base_fin.py b/rocketpy/rocket/aero_surface/fins/_base_fin.py new file mode 100644 index 000000000..a31dd230a --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/_base_fin.py @@ -0,0 +1,193 @@ +import math +from abc import abstractmethod + +import numpy as np + +from rocketpy.mathutils.function import Function + +from ..aero_surface import AeroSurface + + +class _BaseFin(AeroSurface): + """ + Base class for fins, shared by both Fin and Fins classes. + Inherits from AeroSurface. + + Handles shared initialization logic and common properties. + """ + + def __init__( + self, name, rocket_radius, root_chord, span, airfoil=None, cant_angle=0 + ): + """ + Initialize the base fin class. + + Parameters + ---------- + name : str + Name of the fin or fin set. + rocket_radius : float + Rocket radius in meters. + root_chord : float + Root chord of the fin in meters. + span : float + Span of the fin in meters. + airfoil : tuple, optional + Tuple containing airfoil data and unit ('degrees' or 'radians'). + cant_angle : float, optional + Cant angle in degrees. + """ + self.name = name + self._rocket_radius = rocket_radius + self._root_chord = root_chord + self._span = span + self._airfoil = airfoil + self._cant_angle = cant_angle + self._cant_angle_rad = math.radians(cant_angle) + + self.d = 2 * rocket_radius + self.ref_area = np.pi * rocket_radius**2 + + super().__init__(name, self.ref_area, self.d) + + @property + def rocket_radius(self): + return self._rocket_radius + + @rocket_radius.setter + def rocket_radius(self, value): + self._rocket_radius = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + @property + def root_chord(self): + return self._root_chord + + @root_chord.setter + def root_chord(self, value): + self._root_chord = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + @property + def span(self): + return self._span + + @span.setter + def span(self, value): + self._span = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + @property + def cant_angle(self): + return self._cant_angle + + @cant_angle.setter + def cant_angle(self, value): + self._cant_angle = value + self.cant_angle_rad = math.radians(value) + + @property + def cant_angle_rad(self): + return self._cant_angle_rad + + @cant_angle_rad.setter + def cant_angle_rad(self, value): + self._cant_angle_rad = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + @property + def airfoil(self): + return self._airfoil + + @airfoil.setter + def airfoil(self, value): + self._airfoil = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + def evaluate_single_fin_lift_coefficient(self): + if not self.airfoil: + # Defines clalpha2D as 2*pi for planar fins + clalpha2D_incompressible = 2 * np.pi + else: + # Defines clalpha2D as the derivative of the lift coefficient curve + # for the specific airfoil + self.airfoil_cl = Function( + self.airfoil[0], + title="Airfoil lift coefficient", + interpolation="linear", + ) + + # Differentiating at alpha = 0 to get cl_alpha + clalpha2D_incompressible = self.airfoil_cl.differentiate_complex_step( + x=1e-3, dx=1e-3 + ) + + # Convert to radians if needed + if self.airfoil[1] == "degrees": + clalpha2D_incompressible *= 180 / np.pi + + # Correcting for compressible flow (apply Prandtl-Glauert correction) + clalpha2D = Function(lambda mach: clalpha2D_incompressible / self._beta(mach)) + + # Diederich's Planform Correlation Parameter + planform_correlation_parameter = ( + 2 * np.pi * self.AR / (clalpha2D * np.cos(self.gamma_c)) + ) + + # Lift coefficient derivative for a single fin + def lift_source(mach): + return ( + clalpha2D(mach) + * planform_correlation_parameter(mach) + * (self.Af / self.ref_area) + * np.cos(self.gamma_c) + ) / ( + 2 + + planform_correlation_parameter(mach) + * np.sqrt(1 + (2 / planform_correlation_parameter(mach)) ** 2) + ) + + self.clalpha_single_fin = Function( + lift_source, + "Mach", + "Lift coefficient derivative for a single fin", + ) + + @abstractmethod + def evaluate_lift_coefficient(self): + pass + + @abstractmethod + def evaluate_roll_parameters(self): + pass + + @abstractmethod + def evaluate_center_of_pressure(self): + pass + + @abstractmethod + def evaluate_geometrical_parameters(self): + pass + + @abstractmethod + def evaluate_shape(self): + pass + + @abstractmethod + def draw(self): + pass diff --git a/rocketpy/rocket/aero_surface/fins/_elliptical_mixin.py b/rocketpy/rocket/aero_surface/fins/_elliptical_mixin.py new file mode 100644 index 000000000..5a5b58723 --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/_elliptical_mixin.py @@ -0,0 +1,151 @@ +import numpy as np + + +class _EllipticalMixin: + """Mixin class for elliptical fins. This class holds methods and properties + specific to elliptical fin shapes. It is designed to be used in conjunction + with other classes that define the overall behavior of the fins.""" + + def evaluate_geometrical_parameters(self): # pylint: disable=too-many-statements + """Calculates and saves fin set's geometrical parameters such as the + fins' area, aspect ratio and parameters for roll movement.""" + + # Compute auxiliary geometrical parameters + # pylint: disable=invalid-name + Af = (np.pi * self.root_chord / 2 * self.span) / 2 # Fin area + gamma_c = 0 # Zero for elliptical fins + AR = 2 * self.span**2 / Af # Fin aspect ratio + Yma = ( + self.span / (3 * np.pi) * np.sqrt(9 * np.pi**2 - 64) + ) # Span wise coord of mean aero chord + roll_geometrical_constant = ( + self.root_chord + * self.span + * ( + 3 * np.pi * self.span**2 + + 32 * self.rocket_radius * self.span + + 12 * np.pi * self.rocket_radius**2 + ) + / 48 + ) + + # Fin–body interference correction parameters + tau = (self.span + self.rocket_radius) / self.rocket_radius + lift_interference_factor = 1 + 1 / tau + if self.span > self.rocket_radius: + roll_damping_interference_factor = 1 + ( + (self.rocket_radius**2) + * ( + 2 + * (self.rocket_radius**2) + * np.sqrt(self.span**2 - self.rocket_radius**2) + * np.log( + ( + 2 + * self.span + * np.sqrt(self.span**2 - self.rocket_radius**2) + + 2 * self.span**2 + ) + / self.rocket_radius + ) + - 2 + * (self.rocket_radius**2) + * np.sqrt(self.span**2 - self.rocket_radius**2) + * np.log(2 * self.span) + + 2 * self.span**3 + - np.pi * self.rocket_radius * self.span**2 + - 2 * (self.rocket_radius**2) * self.span + + np.pi * self.rocket_radius**3 + ) + ) / ( + 2 + * (self.span**2) + * (self.span / 3 + np.pi * self.rocket_radius / 4) + * (self.span**2 - self.rocket_radius**2) + ) + elif self.span < self.rocket_radius: + roll_damping_interference_factor = 1 - ( + self.rocket_radius**2 + * ( + 2 * self.span**3 + - np.pi * self.span**2 * self.rocket_radius + - 2 * self.span * self.rocket_radius**2 + + np.pi * self.rocket_radius**3 + + 2 + * self.rocket_radius**2 + * np.sqrt(-(self.span**2) + self.rocket_radius**2) + * np.arctan( + (self.span) / (np.sqrt(-(self.span**2) + self.rocket_radius**2)) + ) + - np.pi + * self.rocket_radius**2 + * np.sqrt(-(self.span**2) + self.rocket_radius**2) + ) + ) / ( + 2 + * self.span + * (-(self.span**2) + self.rocket_radius**2) + * (self.span**2 / 3 + np.pi * self.span * self.rocket_radius / 4) + ) + else: + roll_damping_interference_factor = (28 - 3 * np.pi) / (4 + 3 * np.pi) + + roll_forcing_interference_factor = (1 / np.pi**2) * ( + (np.pi**2 / 4) * ((tau + 1) ** 2 / tau**2) + + ((np.pi * (tau**2 + 1) ** 2) / (tau**2 * (tau - 1) ** 2)) + * np.arcsin((tau**2 - 1) / (tau**2 + 1)) + - (2 * np.pi * (tau + 1)) / (tau * (tau - 1)) + + ((tau**2 + 1) ** 2) + / (tau**2 * (tau - 1) ** 2) + * (np.arcsin((tau**2 - 1) / (tau**2 + 1))) ** 2 + - (4 * (tau + 1)) + / (tau * (tau - 1)) + * np.arcsin((tau**2 - 1) / (tau**2 + 1)) + + (8 / (tau - 1) ** 2) * np.log((tau**2 + 1) / (2 * tau)) + ) + + # Store values + # pylint: disable=invalid-name + self.Af = Af # Fin area + self.AR = AR # Fin aspect ratio + self.gamma_c = gamma_c # Mid chord angle + self.Yma = Yma # Span wise coord of mean aero chord + self.roll_geometrical_constant = roll_geometrical_constant + self.tau = tau + self.lift_interference_factor = lift_interference_factor + self.roll_damping_interference_factor = roll_damping_interference_factor + self.roll_forcing_interference_factor = roll_forcing_interference_factor + + self.evaluate_shape() + + def evaluate_shape(self): + angles = np.arange(0, 180, 5) + x_array = self.root_chord / 2 + self.root_chord / 2 * np.cos(np.radians(angles)) + y_array = self.span * np.sin(np.radians(angles)) + self.shape_vec = [x_array, y_array] + + def info(self): + self.prints.geometry() + self.prints.lift() + + def all_info(self): + self.prints.all() + self.plots.all() + + def to_dict(self, include_outputs=False): + data = super().to_dict(include_outputs) + if include_outputs: + data.update( + { + "Af": self.Af, + "AR": self.AR, + "gamma_c": self.gamma_c, + "Yma": self.Yma, + "roll_geometrical_constant": self.roll_geometrical_constant, + "tau": self.tau, + "lift_interference_factor": self.lift_interference_factor, + "roll_damping_interference_factor": self.roll_damping_interference_factor, + "roll_forcing_interference_factor": self.roll_forcing_interference_factor, + } + ) + return data diff --git a/rocketpy/rocket/aero_surface/fins/_free_form_mixin.py b/rocketpy/rocket/aero_surface/fins/_free_form_mixin.py new file mode 100644 index 000000000..71df52c52 --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/_free_form_mixin.py @@ -0,0 +1,228 @@ +import warnings + +import numpy as np + + +class _FreeFormMixin: + """Mixin class for free form fins. This class holds methods and properties + specific to free form fin shapes. It is designed to be used in conjunction + with other classes that define the overall behavior of the fins.""" + + def _initialize(self, shape_points): + self.shape_points = shape_points + + down = False + for i in range(1, len(shape_points)): + if shape_points[i][1] > shape_points[i - 1][1] and down: + warnings.warn( + "Jagged fin shape detected. This may cause small inaccuracies " + "center of pressure and pitch moment calculations." + ) + break + if shape_points[i][1] < shape_points[i - 1][1]: + down = True + i += 1 + + root_chord = abs(shape_points[0][0] - shape_points[-1][0]) + ys = [point[1] for point in shape_points] + span = max(ys) - min(ys) + return ( + root_chord, + span, + ) + + def evaluate_geometrical_parameters(self): # pylint: disable=too-many-statements + """Calculates and saves the fin set's geometrical parameters such as the + fin area, aspect ratio, and parameters related to roll movement. This + method uses the same calculations to those in OpenRocket for free-form + fin shapes.""" + # pylint: disable=invalid-name + # pylint: disable=too-many-locals + # Calculate the fin area (Af) using the Shoelace theorem (polygon area formula) + Af = 0 + for i in range(len(self.shape_points) - 1): + x1, y1 = self.shape_points[i] + x2, y2 = self.shape_points[i + 1] + Af += (y1 + y2) * (x1 - x2) + Af = abs(Af) / 2 + if Af < 1e-6: + raise ValueError("Fin area is too small. Check the shape_points.") + + # Calculate aspect ratio (AR) and lift interference factors + AR = 2 * self.span**2 / Af # Aspect ratio + tau = (self.span + self.rocket_radius) / self.rocket_radius + lift_interference_factor = 1 + 1 / tau + + # Calculate roll forcing interference factor using OpenRocket's approach + roll_forcing_interference_factor = (1 / np.pi**2) * ( + (np.pi**2 / 4) * ((tau + 1) ** 2 / tau**2) + + ((np.pi * (tau**2 + 1) ** 2) / (tau**2 * (tau - 1) ** 2)) + * np.arcsin((tau**2 - 1) / (tau**2 + 1)) + - (2 * np.pi * (tau + 1)) / (tau * (tau - 1)) + + ((tau**2 + 1) ** 2 / (tau**2 * (tau - 1) ** 2)) + * (np.arcsin((tau**2 - 1) / (tau**2 + 1))) ** 2 + - (4 * (tau + 1)) + / (tau * (tau - 1)) + * np.arcsin((tau**2 - 1) / (tau**2 + 1)) + + (8 / (tau - 1) ** 2) * np.log((tau**2 + 1) / (2 * tau)) + ) + + # Define number of interpolation points along the span of the fin + points_per_line = 40 # Same as OpenRocket + + # Initialize arrays for leading/trailing edge and chord lengths + chord_lead = np.ones(points_per_line) * np.inf # Leading edge x coordinates + chord_trail = np.ones(points_per_line) * -np.inf # Trailing edge x coordinates + chord_length = np.zeros( + points_per_line + ) # Chord length for each spanwise section + + # Discretize fin shape and calculate chord length, leading, and trailing edges + for p in range(1, len(self.shape_points)): + x1, y1 = self.shape_points[p - 1] + x2, y2 = self.shape_points[p] + + # Compute corresponding points along the fin span (clamp to valid range) + prev_idx = int(y1 / self.span * (points_per_line - 1)) + curr_idx = int(y2 / self.span * (points_per_line - 1)) + prev_idx = np.clip(prev_idx, 0, points_per_line - 1) + curr_idx = np.clip(curr_idx, 0, points_per_line - 1) + + if prev_idx > curr_idx: + prev_idx, curr_idx = curr_idx, prev_idx + + # Compute intersection of fin edge with each spanwise section + for i in range(prev_idx, curr_idx + 1): + y = i * self.span / (points_per_line - 1) + if y1 != y2: + x = np.clip( + (y - y2) / (y1 - y2) * x1 + (y1 - y) / (y1 - y2) * x2, + min(x1, x2), + max(x1, x2), + ) + else: + x = x1 # Handle horizontal segments + + # Update leading and trailing edge positions + chord_lead[i] = min(chord_lead[i], x) + chord_trail[i] = max(chord_trail[i], x) + + # Update chord length + if y1 < y2: + chord_length[i] -= x + else: + chord_length[i] += x + + # Replace infinities and handle invalid values in chord_lead and chord_trail + for i in range(points_per_line): + if ( + np.isinf(chord_lead[i]) + or np.isinf(chord_trail[i]) + or np.isnan(chord_lead[i]) + or np.isnan(chord_trail[i]) + ): + chord_lead[i] = 0 + chord_trail[i] = 0 + if chord_length[i] < 0 or np.isnan(chord_length[i]): + chord_length[i] = 0 + if chord_length[i] > chord_trail[i] - chord_lead[i]: + chord_length[i] = chord_trail[i] - chord_lead[i] + + # Initialize integration variables for various aerodynamic and roll properties + radius = self.rocket_radius + total_area = 0 + mac_length = 0 # Mean aerodynamic chord length + mac_lead = 0 # Mean aerodynamic chord leading edge + mac_span = 0 # Mean aerodynamic chord spanwise position (Yma) + cos_gamma_sum = 0 # Sum of cosine of the sweep angle + roll_geometrical_constant = 0 + roll_damping_numerator = 0 + roll_damping_denominator = 0 + + # Perform integration over spanwise sections + dy = self.span / (points_per_line - 1) + for i in range(points_per_line): + chord = chord_trail[i] - chord_lead[i] + y = i * dy + + # Update integration variables + mac_length += chord * chord + mac_span += y * chord + mac_lead += chord_lead[i] * chord + total_area += chord + roll_geometrical_constant += chord_length[i] * (radius + y) ** 2 + roll_damping_numerator += radius**3 * chord / (radius + y) ** 2 + roll_damping_denominator += (radius + y) * chord + + # Update cosine of sweep angle (cos_gamma) + if i > 0: + dx = (chord_trail[i] + chord_lead[i]) / 2 - ( + chord_trail[i - 1] + chord_lead[i - 1] + ) / 2 + cos_gamma_sum += dy / np.hypot(dx, dy) + + # Finalize mean aerodynamic chord properties + mac_length *= dy + mac_span *= dy + mac_lead *= dy + total_area *= dy + roll_geometrical_constant *= dy + roll_damping_numerator *= dy + roll_damping_denominator *= dy + + mac_length /= total_area + mac_span /= total_area + mac_lead /= total_area + cos_gamma = cos_gamma_sum / (points_per_line - 1) + + # Store computed values + self.Af = Af # Fin area + self.AR = AR # Aspect ratio + self.gamma_c = np.arccos(cos_gamma) # Sweep angle + self.Yma = mac_span # Mean aerodynamic chord spanwise position + self.mac_length = mac_length + self.mac_lead = mac_lead + self.tau = tau + self.roll_geometrical_constant = roll_geometrical_constant + self.lift_interference_factor = lift_interference_factor + self.roll_forcing_interference_factor = roll_forcing_interference_factor + self.roll_damping_interference_factor = 1 + ( + roll_damping_numerator / roll_damping_denominator + ) + + # Evaluate the shape and finalize geometry + self.evaluate_shape() + + def evaluate_shape(self): + x_array, y_array = zip(*self.shape_points) + self.shape_vec = [np.array(x_array), np.array(y_array)] + + def to_dict(self, include_outputs=False): + data = super().to_dict(include_outputs) + data["shape_points"] = self.shape_points + + if include_outputs: + data.update( + { + "Af": self.Af, + "AR": self.AR, + "gamma_c": self.gamma_c, + "Yma": self.Yma, + "mac_length": self.mac_length, + "mac_lead": self.mac_lead, + "roll_geometrical_constant": self.roll_geometrical_constant, + "tau": self.tau, + "lift_interference_factor": self.lift_interference_factor, + "roll_forcing_interference_factor": self.roll_forcing_interference_factor, + "roll_damping_interference_factor": self.roll_damping_interference_factor, + } + ) + return data + + def info(self): + self.prints.geometry() + self.prints.lift() + + def all_info(self): + self.prints.all() + self.plots.all() diff --git a/rocketpy/rocket/aero_surface/fins/_trapezoidal_mixin.py b/rocketpy/rocket/aero_surface/fins/_trapezoidal_mixin.py new file mode 100644 index 000000000..daa9e9144 --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/_trapezoidal_mixin.py @@ -0,0 +1,179 @@ +import numpy as np + + +class _TrapezoidalMixin: + """Mixin class for trapezoidal fins. This class holds methods and properties + specific to trapezoidal fin shapes. It is designed to be used in conjunction + with other classes that define the overall behavior of the fins.""" + + @property + def tip_chord(self): + return self._tip_chord + + @tip_chord.setter + def tip_chord(self, value): + self._tip_chord = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + @property + def sweep_angle(self): + return self._sweep_angle + + @sweep_angle.setter + def sweep_angle(self, value): + self._sweep_angle = value + self._sweep_length = np.tan(value * np.pi / 180) * self.span + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + @property + def sweep_length(self): + return self._sweep_length + + @sweep_length.setter + def sweep_length(self, value): + self._sweep_length = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + def _initialize(self, sweep_length, sweep_angle, root_chord, tip_chord, span): + # Check if sweep angle or sweep length is given + if sweep_length is not None and sweep_angle is not None: + raise ValueError("Cannot use sweep_length and sweep_angle together") + elif sweep_angle is not None: + sweep_length = np.tan(sweep_angle * np.pi / 180) * span + elif sweep_length is None: + sweep_length = root_chord - tip_chord + else: # Sweep length is given + pass + + self._tip_chord = tip_chord + self._sweep_length = sweep_length + self._sweep_angle = sweep_angle + + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + def evaluate_geometrical_parameters(self): + """Calculates and saves fin set's geometrical parameters such as the + fins' area, aspect ratio and parameters for roll movement.""" + # pylint: disable=invalid-name + Yr = self.root_chord + self.tip_chord + Af = Yr * self.span / 2 # Fin area + AR = 2 * self.span**2 / Af # Fin aspect ratio + gamma_c = np.arctan( + (self.sweep_length + 0.5 * self.tip_chord - 0.5 * self.root_chord) + / (self.span) + ) + Yma = ( + (self.span / 3) * (self.root_chord + 2 * self.tip_chord) / Yr + ) # Span wise coord of mean aero chord + + # Fin–body interference correction parameters + tau = (self.span + self.rocket_radius) / self.rocket_radius + lift_interference_factor = 1 + 1 / tau + lambda_ = self.tip_chord / self.root_chord + + # Parameters for Roll Moment. + # Documented at: https://docs.rocketpy.org/en/latest/technical/ + roll_geometrical_constant = ( + (self.root_chord + 3 * self.tip_chord) * self.span**3 + + 4 + * (self.root_chord + 2 * self.tip_chord) + * self.rocket_radius + * self.span**2 + + 6 * (self.root_chord + self.tip_chord) * self.span * self.rocket_radius**2 + ) / 12 + roll_damping_interference_factor = 1 + ( + ((tau - lambda_) / (tau)) - ((1 - lambda_) / (tau - 1)) * np.log(tau) + ) / ( + ((tau + 1) * (tau - lambda_)) / (2) + - ((1 - lambda_) * (tau**3 - 1)) / (3 * (tau - 1)) + ) + roll_forcing_interference_factor = (1 / np.pi**2) * ( + (np.pi**2 / 4) * ((tau + 1) ** 2 / tau**2) + + ((np.pi * (tau**2 + 1) ** 2) / (tau**2 * (tau - 1) ** 2)) + * np.arcsin((tau**2 - 1) / (tau**2 + 1)) + - (2 * np.pi * (tau + 1)) / (tau * (tau - 1)) + + ((tau**2 + 1) ** 2) + / (tau**2 * (tau - 1) ** 2) + * (np.arcsin((tau**2 - 1) / (tau**2 + 1))) ** 2 + - (4 * (tau + 1)) + / (tau * (tau - 1)) + * np.arcsin((tau**2 - 1) / (tau**2 + 1)) + + (8 / (tau - 1) ** 2) * np.log((tau**2 + 1) / (2 * tau)) + ) + + # Store values + self.Yr = Yr + self.Af = Af # Fin area + self.AR = AR # Aspect Ratio + self.gamma_c = gamma_c # Mid chord angle + self.Yma = Yma # Span wise coord of mean aero chord + self.roll_geometrical_constant = roll_geometrical_constant + self.tau = tau + self.lift_interference_factor = lift_interference_factor + self.λ = lambda_ # pylint: disable=non-ascii-name + self.roll_damping_interference_factor = roll_damping_interference_factor + self.roll_forcing_interference_factor = roll_forcing_interference_factor + + self.evaluate_shape() + + def evaluate_shape(self): + if self.sweep_length: + points = [ + (0, 0), + (self.sweep_length, self.span), + (self.sweep_length + self.tip_chord, self.span), + (self.root_chord, 0), + ] + else: + points = [ + (0, 0), + (self.root_chord - self.tip_chord, self.span), + (self.root_chord, self.span), + (self.root_chord, 0), + ] + + x_array, y_array = zip(*points) + self.shape_vec = [np.array(x_array), np.array(y_array)] + + def info(self): + self.prints.geometry() + self.prints.lift() + + def all_info(self): + self.prints.all() + self.plots.all() + + def to_dict(self, include_outputs=False): + data = super().to_dict(include_outputs) + data["tip_chord"] = self.tip_chord + + if include_outputs: + data.update( + { + "sweep_length": self.sweep_length, + "sweep_angle": self.sweep_angle, + "shape_vec": self.shape_vec, + "Af": self.Af, + "AR": self.AR, + "gamma_c": self.gamma_c, + "Yma": self.Yma, + "roll_geometrical_constant": self.roll_geometrical_constant, + "tau": self.tau, + "lift_interference_factor": self.lift_interference_factor, + "roll_damping_interference_factor": self.roll_damping_interference_factor, + "roll_forcing_interference_factor": self.roll_forcing_interference_factor, + } + ) + return data diff --git a/rocketpy/rocket/aero_surface/fins/elliptical_fin.py b/rocketpy/rocket/aero_surface/fins/elliptical_fin.py new file mode 100644 index 000000000..9b7293a2f --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/elliptical_fin.py @@ -0,0 +1,212 @@ +from rocketpy.plots.aero_surface_plots import _EllipticalFinPlots +from rocketpy.prints.aero_surface_prints import _EllipticalFinPrints +from rocketpy.rocket.aero_surface.fins._elliptical_mixin import _EllipticalMixin +from rocketpy.rocket.aero_surface.fins.fin import Fin + + +class EllipticalFin(_EllipticalMixin, Fin): + """Class that defines and holds information for an elliptical fin set. + + This class inherits from the Fin class. + + Note + ---- + Local coordinate system: + - Origin located at the top of the root chord. + - Z axis along the longitudinal axis of symmetry, positive downwards (top -> bottom). + - Y axis perpendicular to the Z axis, in the span direction, positive upwards. + - X axis completes the right-handed coordinate system. + + See Also + -------- + Fin + + Attributes + ---------- + EllipticalFin.rocket_radius : float + The reference rocket radius used for lift coefficient normalization, in + meters. + EllipticalFin.airfoil : tuple + Tuple of two items. First is the airfoil lift curve. + Second is the unit of the curve (radians or degrees) + EllipticalFin.cant_angle : float + Fins cant angle with respect to the rocket centerline, in degrees. + EllipticalFin.cant_angle_rad : float + Fins cant angle with respect to the rocket centerline, in radians. + EllipticalFin.root_chord : float + Fin root chord in meters. + EllipticalFin.span : float + Fin span in meters. + EllipticalFin.name : string + Name of fin set. + EllipticalFin.sweep_length : float + Fins sweep length in meters. By sweep length, understand the axial + distance between the fin root leading edge and the fin tip leading edge + measured parallel to the rocket centerline. + EllipticalFin.sweep_angle : float + Fins sweep angle with respect to the rocket centerline. Must + be given in degrees. + EllipticalFin.d : float + Reference diameter of the rocket, in meters. + EllipticalFin.ref_area : float + Reference area of the rocket. + EllipticalFin.Af : float + Area of the longitudinal section of each fin in the set. + EllipticalFin.AR : float + Aspect ratio of each fin in the set. + EllipticalFin.gamma_c : float + Fin mid-chord sweep angle. + EllipticalFin.Yma : float + Span wise position of the mean aerodynamic chord. + EllipticalFin.roll_geometrical_constant : float + Geometrical constant used in roll calculations. + EllipticalFin.tau : float + Geometrical relation used to simplify lift and roll calculations. + EllipticalFin.lift_interference_factor : float + Factor of Fin-Body interference in the lift coefficient. + EllipticalFin.cp : tuple + Tuple with the x, y and z local coordinates of the fin set center of + pressure. Has units of length and is given in meters. + EllipticalFin.cpx : float + Fin set local center of pressure x coordinate. Has units of length and + is given in meters. + EllipticalFin.cpy : float + Fin set local center of pressure y coordinate. Has units of length and + is given in meters. + EllipticalFin.cpz : float + Fin set local center of pressure z coordinate. Has units of length and + is given in meters. + EllipticalFin.cl : Function + Function which defines the lift coefficient as a function of the angle + of attack and the Mach number. Takes as input the angle of attack in + radians and the Mach number. Returns the lift coefficient. + EllipticalFin.clalpha : float + Lift coefficient slope. Has units of 1/rad. + """ + + def __init__( + self, + angular_position, + root_chord, + span, + rocket_radius, + cant_angle=0, + airfoil=None, + name="Fins", + ): + """Initialize EllipticalFin class. + + Parameters + ---------- + angular_position : float + Angular position of the fin in degrees measured as the rotation + around the symmetry axis of the rocket relative to one of the other + principal axis. See :ref:`Angular Position Inputs ` + root_chord : int, float + Fin root chord in meters. + span : int, float + Fin span in meters. + rocket_radius : int, float + Reference radius to calculate lift coefficient, in meters. + cant_angle : int, float, optional + Fins cant angle with respect to the rocket centerline. Must + be given in degrees. + sweep_length : int, float, optional + Fins sweep length in meters. By sweep length, understand the axial + distance between the fin root leading edge and the fin tip leading + edge measured parallel to the rocket centerline. If not given, the + sweep length is assumed to be equal the root chord minus the tip + chord, in which case the fin is a right trapezoid with its base + perpendicular to the rocket's axis. Cannot be used in conjunction + with sweep_angle. + sweep_angle : int, float, optional + Fins sweep angle with respect to the rocket centerline. Must + be given in degrees. If not given, the sweep angle is automatically + calculated, in which case the fin is assumed to be a right trapezoid + with its base perpendicular to the rocket's axis. + Cannot be used in conjunction with sweep_length. + airfoil : tuple, optional + Default is null, in which case fins will be treated as flat plates. + Otherwise, if tuple, fins will be considered as airfoils. The + tuple's first item specifies the airfoil's lift coefficient + by angle of attack and must be either a .csv, .txt, ndarray + or callable. The .csv and .txt files can contain a single line + header and the first column must specify the angle of attack, while + the second column must specify the lift coefficient. The + ndarray should be as [(x0, y0), (x1, y1), (x2, y2), ...] + where x0 is the angle of attack and y0 is the lift coefficient. + If callable, it should take an angle of attack as input and + return the lift coefficient at that angle of attack. + The tuple's second item is the unit of the angle of attack, + accepting either "radians" or "degrees". + name : str + Name of fin set. + + Returns + ------- + None + """ + + super().__init__( + angular_position, + root_chord, + span, + rocket_radius, + cant_angle, + airfoil, + name, + ) + + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + self.prints = _EllipticalFinPrints(self) + self.plots = _EllipticalFinPlots(self) + + def evaluate_center_of_pressure(self): + """Calculates and returns the center of pressure of the fin set in local + coordinates. The center of pressure position is saved and stored as a + tuple. + + Returns + ------- + None + """ + # Center of pressure position in local coordinates + cpz = 0.288 * self.root_chord + self.cpx = 0 + self.cpy = self.Yma + self.cpz = cpz + self.cp = (self.cpx, self.cpy, self.cpz) + + def to_dict(self, include_outputs=False): + data = super().to_dict(include_outputs) + if include_outputs: + data.update( + { + "Af": self.Af, + "AR": self.AR, + "gamma_c": self.gamma_c, + "Yma": self.Yma, + "roll_geometrical_constant": self.roll_geometrical_constant, + "tau": self.tau, + "lift_interference_factor": self.lift_interference_factor, + "roll_damping_interference_factor": self.roll_damping_interference_factor, + "roll_forcing_interference_factor": self.roll_forcing_interference_factor, + } + ) + return data + + @classmethod + def from_dict(cls, data): + return cls( + angular_position=data["angular_position"], + root_chord=data["root_chord"], + span=data["span"], + rocket_radius=data["rocket_radius"], + cant_angle=data["cant_angle"], + airfoil=data["airfoil"], + name=data["name"], + ) diff --git a/rocketpy/rocket/aero_surface/fins/elliptical_fins.py b/rocketpy/rocket/aero_surface/fins/elliptical_fins.py index 08c51ab64..05989f71c 100644 --- a/rocketpy/rocket/aero_surface/fins/elliptical_fins.py +++ b/rocketpy/rocket/aero_surface/fins/elliptical_fins.py @@ -1,12 +1,11 @@ -import numpy as np - from rocketpy.plots.aero_surface_plots import _EllipticalFinsPlots from rocketpy.prints.aero_surface_prints import _EllipticalFinsPrints +from rocketpy.rocket.aero_surface.fins._elliptical_mixin import _EllipticalMixin from .fins import Fins -class EllipticalFins(Fins): +class EllipticalFins(_EllipticalMixin, Fins): """Class that defines and holds information for an elliptical fin set. This class inherits from the Fins class. @@ -35,9 +34,6 @@ class EllipticalFins(Fins): Second is the unit of the curve (radians or degrees) EllipticalFins.cant_angle : float Fins cant angle with respect to the rocket centerline, in degrees. - EllipticalFins.changing_attribute_dict : dict - Dictionary that stores the name and the values of the attributes that - may be changed during a simulation. Useful for control systems. EllipticalFins.cant_angle_rad : float Fins cant angle with respect to the rocket centerline, in radians. EllipticalFins.root_chord : float @@ -186,155 +182,6 @@ def evaluate_center_of_pressure(self): self.cpz = cpz self.cp = (self.cpx, self.cpy, self.cpz) - def evaluate_geometrical_parameters(self): # pylint: disable=too-many-statements - """Calculates and saves fin set's geometrical parameters such as the - fins' area, aspect ratio and parameters for roll movement. - - Returns - ------- - None - """ - - # Compute auxiliary geometrical parameters - # pylint: disable=invalid-name - Af = (np.pi * self.root_chord / 2 * self.span) / 2 # Fin area - gamma_c = 0 # Zero for elliptical fins - AR = 2 * self.span**2 / Af # Fin aspect ratio - Yma = ( - self.span / (3 * np.pi) * np.sqrt(9 * np.pi**2 - 64) - ) # Span wise coord of mean aero chord - roll_geometrical_constant = ( - self.root_chord - * self.span - * ( - 3 * np.pi * self.span**2 - + 32 * self.rocket_radius * self.span - + 12 * np.pi * self.rocket_radius**2 - ) - / 48 - ) - - # Fin–body interference correction parameters - tau = (self.span + self.rocket_radius) / self.rocket_radius - lift_interference_factor = 1 + 1 / tau - if self.span > self.rocket_radius: - roll_damping_interference_factor = 1 + ( - (self.rocket_radius**2) - * ( - 2 - * (self.rocket_radius**2) - * np.sqrt(self.span**2 - self.rocket_radius**2) - * np.log( - ( - 2 - * self.span - * np.sqrt(self.span**2 - self.rocket_radius**2) - + 2 * self.span**2 - ) - / self.rocket_radius - ) - - 2 - * (self.rocket_radius**2) - * np.sqrt(self.span**2 - self.rocket_radius**2) - * np.log(2 * self.span) - + 2 * self.span**3 - - np.pi * self.rocket_radius * self.span**2 - - 2 * (self.rocket_radius**2) * self.span - + np.pi * self.rocket_radius**3 - ) - ) / ( - 2 - * (self.span**2) - * (self.span / 3 + np.pi * self.rocket_radius / 4) - * (self.span**2 - self.rocket_radius**2) - ) - elif self.span < self.rocket_radius: - roll_damping_interference_factor = 1 - ( - self.rocket_radius**2 - * ( - 2 * self.span**3 - - np.pi * self.span**2 * self.rocket_radius - - 2 * self.span * self.rocket_radius**2 - + np.pi * self.rocket_radius**3 - + 2 - * self.rocket_radius**2 - * np.sqrt(-(self.span**2) + self.rocket_radius**2) - * np.arctan( - (self.span) / (np.sqrt(-(self.span**2) + self.rocket_radius**2)) - ) - - np.pi - * self.rocket_radius**2 - * np.sqrt(-(self.span**2) + self.rocket_radius**2) - ) - ) / ( - 2 - * self.span - * (-(self.span**2) + self.rocket_radius**2) - * (self.span**2 / 3 + np.pi * self.span * self.rocket_radius / 4) - ) - else: - roll_damping_interference_factor = (28 - 3 * np.pi) / (4 + 3 * np.pi) - - roll_forcing_interference_factor = (1 / np.pi**2) * ( - (np.pi**2 / 4) * ((tau + 1) ** 2 / tau**2) - + ((np.pi * (tau**2 + 1) ** 2) / (tau**2 * (tau - 1) ** 2)) - * np.arcsin((tau**2 - 1) / (tau**2 + 1)) - - (2 * np.pi * (tau + 1)) / (tau * (tau - 1)) - + ((tau**2 + 1) ** 2) - / (tau**2 * (tau - 1) ** 2) - * (np.arcsin((tau**2 - 1) / (tau**2 + 1))) ** 2 - - (4 * (tau + 1)) - / (tau * (tau - 1)) - * np.arcsin((tau**2 - 1) / (tau**2 + 1)) - + (8 / (tau - 1) ** 2) * np.log((tau**2 + 1) / (2 * tau)) - ) - - # Store values - # pylint: disable=invalid-name - self.Af = Af # Fin area - self.AR = AR # Fin aspect ratio - self.gamma_c = gamma_c # Mid chord angle - self.Yma = Yma # Span wise coord of mean aero chord - self.roll_geometrical_constant = roll_geometrical_constant - self.tau = tau - self.lift_interference_factor = lift_interference_factor - self.roll_damping_interference_factor = roll_damping_interference_factor - self.roll_forcing_interference_factor = roll_forcing_interference_factor - - self.evaluate_shape() - - def evaluate_shape(self): - angles = np.arange(0, 180, 5) - x_array = self.root_chord / 2 + self.root_chord / 2 * np.cos(np.radians(angles)) - y_array = self.span * np.sin(np.radians(angles)) - self.shape_vec = [x_array, y_array] - - def info(self): - self.prints.geometry() - self.prints.lift() - - def all_info(self): - self.prints.all() - self.plots.all() - - def to_dict(self, include_outputs=False): - data = super().to_dict(include_outputs) - if include_outputs: - data.update( - { - "Af": self.Af, - "AR": self.AR, - "gamma_c": self.gamma_c, - "Yma": self.Yma, - "roll_geometrical_constant": self.roll_geometrical_constant, - "tau": self.tau, - "lift_interference_factor": self.lift_interference_factor, - "roll_damping_interference_factor": self.roll_damping_interference_factor, - "roll_forcing_interference_factor": self.roll_forcing_interference_factor, - } - ) - return data - @classmethod def from_dict(cls, data): return cls( diff --git a/rocketpy/rocket/aero_surface/fins/fin.py b/rocketpy/rocket/aero_surface/fins/fin.py new file mode 100644 index 000000000..5bdfa996d --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/fin.py @@ -0,0 +1,452 @@ +import math + +import numpy as np + +from rocketpy.mathutils.function import Function +from rocketpy.mathutils.vector_matrix import Matrix, Vector +from rocketpy.rocket.aero_surface.fins._base_fin import _BaseFin + + +class Fin(_BaseFin): + """Abstract class that holds common methods for the individual fin classes. + Cannot be instantiated. + + Note + ---- + Local coordinate system: + - Origin located at the top of the root chord. + - Z axis along the longitudinal axis of symmetry, positive downwards (top -> bottom). + - Y axis perpendicular to the Z axis, in the span direction, positive upwards. + - X axis completes the right-handed coordinate system. + + Attributes + ---------- + Fin.rocket_radius : float + The reference rocket radius used for lift coefficient normalization, + in meters. + Fin.airfoil : tuple + Tuple of two items. First is the airfoil lift curve. + Second is the unit of the curve (radians or degrees). + Fin.cant_angle : float + Fin cant angle with respect to the rocket centerline, in degrees. + Fin.changing_attribute_dict : dict + Dictionary that stores the name and the values of the attributes that + may be changed during a simulation. Useful for control systems. + Fin.cant_angle_rad : float + Fin cant angle with respect to the rocket centerline, in radians. + Fin.root_chord : float + Fin root chord in meters. + Fin.tip_chord : float + Fin tip chord in meters. + Fin.span : float + Fin span in meters. + Fin.name : string + Name of fin set. + Fin.sweep_length : float + Fin sweep length in meters. By sweep length, understand the axial + distance between the fin root leading edge and the fin tip leading edge + measured parallel to the rocket centerline. + Fin.sweep_angle : float + Fin sweep angle with respect to the rocket centerline. Must + be given in degrees. + Fin.d : float + Reference diameter of the rocket. Has units of length and is given + in meters. + Fin.ref_area : float + Reference area of the rocket. + Fin.Af : float + Area of the longitudinal section of each fin in the set. + Fin.AR : float + Aspect ratio of each fin in the set. + Fin.gamma_c : float + Fin mid-chord sweep angle. + Fin.Yma : float + Span wise position of the mean aerodynamic chord. + Fin.roll_geometrical_constant : float + Geometrical constant used in roll calculations. + Fin.tau : float + Geometrical relation used to simplify lift and roll calculations. + Fin.lift_interference_factor : float + Factor of Fin-Body interference in the lift coefficient. + Fin.cp : tuple + Tuple with the x, y and z local coordinates of the fin set center of + pressure. Has units of length and is given in meters. + Fin.cpx : float + Fin set local center of pressure x coordinate. Has units of length and + is given in meters. + Fin.cpy : float + Fin set local center of pressure y coordinate. Has units of length and + is given in meters. + Fin.cpz : float + Fin set local center of pressure z coordinate. Has units of length and + is given in meters. + Fin.cl : Function + Function which defines the lift coefficient as a function of the angle + of attack and the Mach number. Takes as input the angle of attack in + radians and the Mach number. Returns the lift coefficient. + Fin.clalpha : float + Lift coefficient slope. Has units of 1/rad. + Fin.roll_parameters : list + List containing the roll moment lift coefficient, the roll moment + damping coefficient and the cant angle in radians. + """ + + def __init__( + self, + angular_position, + root_chord, + span, + rocket_radius, + cant_angle=0, + airfoil=None, + name="Fin", + ): + """Initialize Fin class. + + Parameters + ---------- + angular_position : float + Angular position of the fin in degrees measured as the rotation + around the symmetry axis of the rocket relative to one of the other + principal axis. See :ref:`Angular Position Inputs ` + root_chord : int, float + Fin root chord in meters. + span : int, float + Fin span in meters. + rocket_radius : int, float + Reference rocket radius used for lift coefficient normalization. + cant_angle : int, float, optional + Fin cant angle with respect to the rocket centerline. Must + be given in degrees. + airfoil : tuple, optional + Default is null, in which case fins will be treated as flat plates. + Otherwise, if tuple, fins will be considered as airfoils. The + tuple's first item specifies the airfoil's lift coefficient + by angle of attack and must be either a .csv, .txt, ndarray + or callable. The .csv and .txt files can contain a single line + header and the first column must specify the angle of attack, while + the second column must specify the lift coefficient. The + ndarray should be as [(x0, y0), (x1, y1), (x2, y2), ...] + where x0 is the angle of attack and y0 is the lift coefficient. + If callable, it should take an angle of attack as input and + return the lift coefficient at that angle of attack. + The tuple's second item is the unit of the angle of attack, + accepting either "radians" or "degrees". + name : str + Name of fin. + """ + super().__init__( + name=name, + rocket_radius=rocket_radius, + root_chord=root_chord, + span=span, + airfoil=airfoil, + cant_angle=cant_angle, + ) + + # Store values + self._angular_position = angular_position + self._angular_position_rad = math.radians(angular_position) + + @property + def cant_angle(self): + return self._cant_angle + + @cant_angle.setter + def cant_angle(self, value): + self._cant_angle = value + self.cant_angle_rad = math.radians(value) + + @property + def cant_angle_rad(self): + return self._cant_angle_rad + + @cant_angle_rad.setter + def cant_angle_rad(self, value): + self._cant_angle_rad = value + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + self.evaluate_rotation_matrix() + + @property + def angular_position(self): + return self._angular_position + + @angular_position.setter + def angular_position(self, value): + self._angular_position = value + self.angular_position_rad = math.radians(value) + + @property + def angular_position_rad(self): + return self._angular_position_rad + + @angular_position_rad.setter + def angular_position_rad(self, value): + self._angular_position_rad = value + self.evaluate_rotation_matrix() + + def evaluate_lift_coefficient(self): + """Calculates and returns the fin set's lift coefficient. + The lift coefficient is saved and returned. This function + also calculates and saves the lift coefficient derivative + for a single fin and the lift coefficient derivative for + a number of n fins corrected for Fin-Body interference. + + Returns + ------- + None + """ + self.evaluate_single_fin_lift_coefficient() + + self.clalpha = self.clalpha_single_fin * self.lift_interference_factor + + # Cl = clalpha * alpha + self.cl = Function( + lambda alpha, mach: alpha * self.clalpha(mach), + ["Alpha (rad)", "Mach"], + "Lift coefficient", + ) + + return self.cl + + def evaluate_roll_parameters(self): + """Calculates and returns the fin set's roll coefficients. + The roll coefficients are saved in a list. + + Returns + ------- + self.roll_parameters : list + List containing the roll moment lift coefficient, the + roll moment damping coefficient and the cant angle in + radians + """ + clf_delta = ( + self.roll_forcing_interference_factor + * (self.Yma + self.rocket_radius) + * self.clalpha_single_fin + / self.reference_length + ) # Function of mach number + clf_delta.set_inputs("Mach") + clf_delta.set_outputs("Roll moment forcing coefficient derivative") + clf_delta.set_title( + "Roll moment forcing coefficient derivative vs. Mach number" + ) + cld_omega = -( + 2 + * self.roll_damping_interference_factor + * self.clalpha_single_fin + * np.cos(self.cant_angle_rad) + * self.roll_geometrical_constant + / (self.ref_area * self.d**2) + ) # Function of mach number + cld_omega.set_inputs("Mach") + cld_omega.set_outputs("Roll moment damping coefficient derivative") + cld_omega.set_title( + "Roll moment damping coefficient derivative vs. Mach number" + ) + self.roll_parameters = [clf_delta, cld_omega, self.cant_angle_rad] + return self.roll_parameters + + def evaluate_rotation_matrix(self): + """Calculates and returns the rotation matrix from the rocket body frame + to the fin frame. + + Note + ---- + Local coordinate system: + + - Origin located at the leading edge of the root chord. + - Z axis along the longitudinal axis of the fin, positive + downwards (leading edge -> trailing edge). + - Y axis perpendicular to the Z axis, in the span direction, + positive upwards (root chord -> tip chord). + - X axis completes the right-handed coordinate system. + + + Returns + ------- + None + + References + ---------- + :ref:`Individual Fin Model ` + """ + phi = self.angular_position_rad + delta = self.cant_angle_rad + + # Rotation about body Z by angular position + R_phi = Matrix( + [ + [np.cos(phi), -np.sin(phi), 0], + [np.sin(phi), np.cos(phi), 0], + [0, 0, 1], + ] + ) + + # Cant rotation about body Y + R_delta = Matrix( + [ + [np.cos(delta), 0, -np.sin(delta)], + [0, 1, 0], + [np.sin(delta), 0, np.cos(delta)], + ] + ) + + # 180 flip about Y to align fin leading/trailing edge + R_pi = Matrix( + [ + [-1, 0, 0], + [0, 1, 0], + [0, 0, -1], + ] + ) + + # Uncanted body to fin, then apply cant + R_uncanted = R_phi @ R_pi + R_body_to_fin = R_delta @ R_uncanted + + # Store for downstream transforms + self._rotation_fin_to_body_uncanted = R_uncanted.transpose + self._rotation_body_to_fin = R_body_to_fin + self._rotation_fin_to_body = R_body_to_fin.transpose + self._rotation_surface_to_body = self._rotation_fin_to_body + + def compute_forces_and_moments( + self, + stream_velocity, + stream_speed, + stream_mach, + rho, + cp, + omega, + *args, + ): # pylint: disable=arguments-differ + """Computes the forces and moments acting on the aerodynamic surface. + + Parameters + ---------- + stream_velocity : tuple of float + The velocity of the airflow relative to the surface. + stream_speed : float + The magnitude of the airflow speed. + stream_mach : float + The Mach number of the airflow. + rho : float + Air density. + cp : Vector + Center of pressure coordinates in the body frame. + omega: tuple[float, float, float] + Tuple containing angular velocities around the x, y, z axes. + + Returns + ------- + tuple of float + The aerodynamic forces (lift, side_force, drag) and moments + (pitch, yaw, roll) in the body frame. + """ + R1, R2, R3, M1, M2, M3 = 0, 0, 0, 0, 0, 0 + + # stream velocity in fin frame + stream_velocity_f = self._rotation_body_to_fin @ stream_velocity + + attack_angle = np.arctan2(stream_velocity_f[0], stream_velocity_f[2]) + # Force in the X direction of the fin + X = ( + 0.5 + * rho + * stream_speed**2 + * self.reference_area + * self.cl.get_value_opt(attack_angle, stream_mach) + ) + # Force in body frame + R1, R2, R3 = self._rotation_fin_to_body @ Vector([X, 0, 0]) + # Moments + M1, M2, M3 = cp ^ Vector([R1, R2, R3]) + # Apply roll interference factor, disregarding lift interference factor + M3 *= self.roll_forcing_interference_factor / self.lift_interference_factor + + # Roll damping + _, cld_omega, _ = self.roll_parameters + M3_damping = ( + (1 / 2 * rho * stream_speed) + * self.reference_area + * (self.reference_length) ** 2 + * cld_omega.get_value_opt(stream_mach) + * omega[2] # omega3 + / 2 + ) + M3 += M3_damping + return R1, R2, R3, M1, M2, M3 + + def _compute_leading_edge_position(self, position, _csys): + """Computes the position of the fin leading edge in a rocket's user, + given its position in a rocket.""" + # Point from deflection from cant angle in the plane perpendicular to + # the fuselage where the fin is located in the fin frame + p = Vector( + [ + -self.root_chord / 2 * np.sin(self.cant_angle_rad), + 0, + self.root_chord / 2 * (1 - np.cos(self.cant_angle_rad)), + ] + ) + # Rotate the point to the body frame orientation + p = self._rotation_fin_to_body_uncanted @ p + + # Rotate the point to the user-defined coordinate system + p = Vector([p.x * _csys, p.y, p.z * _csys]) + + # Calculate the position of the fin leading edge in the user frame + # as if no cant angle was applied + position = Vector( + [ + -self.rocket_radius * math.sin(self.angular_position_rad) * _csys, + self.rocket_radius * math.cos(self.angular_position_rad), + position, + ] + ) + + # Translate the position of the fin leading edge to the position of the + # fin leading edge with cant angle + position += p + return position + + def to_dict(self, include_outputs=False): + data = { + "angular_position": self.angular_position, + "root_chord": self.root_chord, + "span": self.span, + "rocket_radius": self.rocket_radius, + "cant_angle": self.cant_angle, + "airfoil": self.airfoil, + "name": self.name, + } + + if include_outputs: + data.update( + { + "cp": self.cp, + "cl": self.cl, + "roll_parameters": self.roll_parameters, + "d": self.d, + "ref_area": self.ref_area, + } + ) + return data + + def draw(self, *, filename=None): + """Draw the fin shape along with some important information, including + the center line, the quarter line and the center of pressure position. + + Parameters + ---------- + filename : str | None, optional + The path the plot should be saved to. By default None, in which case + the plot will be shown instead of saved. Supported file endings are: + eps, jpg, jpeg, pdf, pgf, png, ps, raw, rgba, svg, svgz, tif, tiff + and webp (these are the formats supported by matplotlib). + """ + self.plots.draw(filename=filename) diff --git a/rocketpy/rocket/aero_surface/fins/fins.py b/rocketpy/rocket/aero_surface/fins/fins.py index 26a4a7111..5efec62e5 100644 --- a/rocketpy/rocket/aero_surface/fins/fins.py +++ b/rocketpy/rocket/aero_surface/fins/fins.py @@ -1,11 +1,10 @@ import numpy as np from rocketpy.mathutils.function import Function +from rocketpy.rocket.aero_surface.fins._base_fin import _BaseFin -from ..aero_surface import AeroSurface - -class Fins(AeroSurface): +class Fins(_BaseFin): """Abstract class that holds common methods for the fin classes. Cannot be instantiated. @@ -132,27 +131,18 @@ def __init__( accepting either "radians" or "degrees". name : str Name of fin set. - - Returns - ------- - None """ - # Compute auxiliary geometrical parameters - d = 2 * rocket_radius - ref_area = np.pi * rocket_radius**2 # Reference area - - super().__init__(name, ref_area, d) + super().__init__( + name=name, + rocket_radius=rocket_radius, + root_chord=root_chord, + span=span, + airfoil=airfoil, + cant_angle=-cant_angle, + ) # Store values self._n = n - self._rocket_radius = rocket_radius - self._airfoil = airfoil - self._cant_angle = cant_angle - self._root_chord = root_chord - self._span = span - self.name = name - self.d = d - self.ref_area = ref_area # Reference area @property def n(self): @@ -166,134 +156,26 @@ def n(self, value): self.evaluate_lift_coefficient() self.evaluate_roll_parameters() - @property - def root_chord(self): - return self._root_chord - - @root_chord.setter - def root_chord(self, value): - self._root_chord = value - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - - @property - def span(self): - return self._span - - @span.setter - def span(self, value): - self._span = value - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - - @property - def rocket_radius(self): - return self._rocket_radius - - @rocket_radius.setter - def rocket_radius(self, value): - self._rocket_radius = value - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - - @property - def cant_angle(self): - return self._cant_angle - - @cant_angle.setter - def cant_angle(self, value): - self._cant_angle = value - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - - @property - def airfoil(self): - return self._airfoil - - @airfoil.setter - def airfoil(self, value): - self._airfoil = value - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - def evaluate_lift_coefficient(self): """Calculates and returns the fin set's lift coefficient. The lift coefficient is saved and returned. This function also calculates and saves the lift coefficient derivative for a single fin and the lift coefficient derivative for a number of n fins corrected for Fin-Body interference. - - Returns - ------- - None """ - if not self.airfoil: - # Defines clalpha2D as 2*pi for planar fins - clalpha2D_incompressible = 2 * np.pi - else: - # Defines clalpha2D as the derivative of the lift coefficient curve - # for the specific airfoil - self.airfoil_cl = Function( - self.airfoil[0], - interpolation="linear", - ) - - # Differentiating at alpha = 0 to get cl_alpha - clalpha2D_incompressible = self.airfoil_cl.differentiate_complex_step( - x=1e-3, dx=1e-3 - ) - - # Convert to radians if needed - if self.airfoil[1] == "degrees": - clalpha2D_incompressible *= 180 / np.pi - - # Correcting for compressible flow (apply Prandtl-Glauert correction) - clalpha2D = Function(lambda mach: clalpha2D_incompressible / self._beta(mach)) - - # Diederich's Planform Correlation Parameter - planform_correlation_parameter = ( - 2 * np.pi * self.AR / (clalpha2D * np.cos(self.gamma_c)) - ) - - # Lift coefficient derivative for a single fin - def lift_source(mach): - return ( - clalpha2D(mach) - * planform_correlation_parameter(mach) - * (self.Af / self.ref_area) - * np.cos(self.gamma_c) - ) / ( - 2 - + planform_correlation_parameter(mach) - * np.sqrt(1 + (2 / planform_correlation_parameter(mach)) ** 2) - ) - - self.clalpha_single_fin = Function( - lift_source, - "Mach", - "Lift coefficient derivative for a single fin", - ) + self.evaluate_single_fin_lift_coefficient() # Lift coefficient derivative for n fins corrected with Fin-Body interference self.clalpha_multiple_fins = ( - self.lift_interference_factor - * self.fin_num_correction(self.n) + self.fin_num_correction(self.n) + * self.lift_interference_factor * self.clalpha_single_fin ) # Function of mach number self.clalpha_multiple_fins.set_inputs("Mach") self.clalpha_multiple_fins.set_outputs( f"Lift coefficient derivative for {self.n:.0f} fins" ) + self.clalpha = self.clalpha_multiple_fins # Cl = clalpha * alpha @@ -316,19 +198,19 @@ def evaluate_roll_parameters(self): roll moment damping coefficient and the cant angle in radians """ - - self.cant_angle_rad = np.radians(self.cant_angle) - clf_delta = ( self.roll_forcing_interference_factor - * self.n + * self.fin_num_correction(self.n) * (self.Yma + self.rocket_radius) * self.clalpha_single_fin - / self.d + / self.reference_length ) # Function of mach number clf_delta.set_inputs("Mach") clf_delta.set_outputs("Roll moment forcing coefficient derivative") - cld_omega = ( + clf_delta.set_title( + "Roll moment forcing coefficient derivative vs. Mach number" + ) + cld_omega = -( 2 * self.roll_damping_interference_factor * self.n @@ -339,6 +221,9 @@ def evaluate_roll_parameters(self): ) # Function of mach number cld_omega.set_inputs("Mach") cld_omega.set_outputs("Roll moment damping coefficient derivative") + cld_omega.set_title( + "Roll moment damping coefficient derivative vs. Mach number" + ) self.roll_parameters = [clf_delta, cld_omega, self.cant_angle_rad] return self.roll_parameters @@ -423,7 +308,7 @@ def compute_forces_and_moments( * omega[2] / 2 ) - M3 = M3_forcing - M3_damping + M3 = M3_forcing + M3_damping return R1, R2, R3, M1, M2, M3 def to_dict(self, include_outputs=False): @@ -461,9 +346,5 @@ def draw(self, *, filename=None): the plot will be shown instead of saved. Supported file endings are: eps, jpg, jpeg, pdf, pgf, png, ps, raw, rgba, svg, svgz, tif, tiff and webp (these are the formats supported by matplotlib). - - Returns - ------- - None """ self.plots.draw(filename=filename) diff --git a/rocketpy/rocket/aero_surface/fins/free_form_fin.py b/rocketpy/rocket/aero_surface/fins/free_form_fin.py new file mode 100644 index 000000000..d36962835 --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/free_form_fin.py @@ -0,0 +1,181 @@ +from rocketpy.plots.aero_surface_plots import _FreeFormFinPlots +from rocketpy.prints.aero_surface_prints import _FreeFormFinPrints +from rocketpy.rocket.aero_surface.fins._free_form_mixin import _FreeFormMixin +from rocketpy.rocket.aero_surface.fins.fin import Fin + + +class FreeFormFin(_FreeFormMixin, Fin): + """Class that defines and holds information for a free form fin set. + + This class inherits from the Fin class. + + Note + ---- + Local coordinate system: + - Origin located at the top of the root chord. + - Z axis along the longitudinal axis of symmetry, positive downwards (top -> bottom). + - Y axis perpendicular to the Z axis, in the span direction, positive upwards. + - X axis completes the right-handed coordinate system. + + See Also + -------- + Fin + + Attributes + ---------- + FreeFormFin.n : int + Number of fins in fin set. + FreeFormFin.rocket_radius : float + The reference rocket radius used for lift coefficient normalization, in + meters. + FreeFormFin.airfoil : tuple + Tuple of two items. First is the airfoil lift curve. + Second is the unit of the curve (radians or degrees). + FreeFormFin.cant_angle : float + Fins cant angle with respect to the rocket centerline, in degrees. + FreeFormFin.cant_angle_rad : float + Fins cant angle with respect to the rocket centerline, in radians. + FreeFormFin.root_chord : float + Fin root chord in meters. + FreeFormFin.span : float + Fin span in meters. + FreeFormFin.name : string + Name of fin set. + FreeFormFin.d : float + Reference diameter of the rocket, in meters. + FreeFormFin.ref_area : float + Reference area of the rocket, in m². + FreeFormFin.Af : float + Area of the longitudinal section of each fin in the set. + FreeFormFin.AR : float + Aspect ratio of each fin in the set + FreeFormFin.gamma_c : float + Fin mid-chord sweep angle. + FreeFormFin.Yma : float + Span wise position of the mean aerodynamic chord. + FreeFormFin.roll_geometrical_constant : float + Geometrical constant used in roll calculations. + FreeFormFin.tau : float + Geometrical relation used to simplify lift and roll calculations. + FreeFormFin.lift_interference_factor : float + Factor of Fin-Body interference in the lift coefficient. + FreeFormFin.cp : tuple + Tuple with the x, y and z local coordinates of the fin set center of + pressure. Has units of length and is given in meters. + FreeFormFin.cpx : float + Fin set local center of pressure x coordinate. Has units of length and + is given in meters. + FreeFormFin.cpy : float + Fin set local center of pressure y coordinate. Has units of length and + is given in meters. + FreeFormFin.cpz : float + Fin set local center of pressure z coordinate. Has units of length and + is given in meters. + FreeFormFin.cl : Function + Function which defines the lift coefficient as a function of the angle + of attack and the Mach number. Takes as input the angle of attack in + radians and the Mach number. Returns the lift coefficient. + FreeFormFin.clalpha : float + Lift coefficient slope. Has units of 1/rad. + FreeFormFin.mac_length : float + Mean aerodynamic chord length of the fin set. + FreeFormFin.mac_lead : float + Mean aerodynamic chord leading edge x coordinate. + """ + + def __init__( + self, + angular_position, + shape_points, + rocket_radius, + cant_angle=0, + airfoil=None, + name="Fins", + ): + """Initialize FreeFormFin class. + + Parameters + ---------- + angular_position : float + Angular position of the fin in degrees measured as the rotation + around the symmetry axis of the rocket relative to one of the other + principal axis. See :ref:`Angular Position Inputs ` + shape_points : list + List of tuples (x, y) containing the coordinates of the fin's + geometry defining points. The point (0, 0) is the root leading edge. + Positive x is rearwards, positive y is upwards (span direction). + The shape will be interpolated between the points, in the order + they are given. The last point connects to the first point, and + represents the trailing edge. + rocket_radius : int, float + Reference radius to calculate lift coefficient, in meters. + cant_angle : int, float, optional + Fins cant angle with respect to the rocket centerline. Must + be given in degrees. + airfoil : tuple, optional + Default is null, in which case fins will be treated as flat plates. + Otherwise, if tuple, fins will be considered as airfoils. The + tuple's first item specifies the airfoil's lift coefficient + by angle of attack and must be either a .csv, .txt, ndarray + or callable. The .csv and .txt files can contain a single line + header and the first column must specify the angle of attack, while + the second column must specify the lift coefficient. The + ndarray should be as [(x0, y0), (x1, y1), (x2, y2), ...] + where x0 is the angle of attack and y0 is the lift coefficient. + If callable, it should take an angle of attack as input and + return the lift coefficient at that angle of attack. + The tuple's second item is the unit of the angle of attack, + accepting either "radians" or "degrees". + name : str + Name of fin set. + + Returns + ------- + None + """ + root_chord, span = self._initialize(shape_points) + + super().__init__( + angular_position, + root_chord, + span, + rocket_radius, + cant_angle, + airfoil, + name, + ) + + self.evaluate_geometrical_parameters() + self.evaluate_center_of_pressure() + self.evaluate_lift_coefficient() + self.evaluate_roll_parameters() + + self.prints = _FreeFormFinPrints(self) + self.plots = _FreeFormFinPlots(self) + + def evaluate_center_of_pressure(self): + """Calculates and returns the center of pressure of the fin set in local + coordinates. The center of pressure position is saved and stored as a + tuple. + + Returns + ------- + None + """ + # Center of pressure position in local coordinates + cpz = self.mac_lead + 0.25 * self.mac_length + self.cpx = 0 + self.cpy = self.Yma + self.cpz = cpz + self.cp = (self.cpx, self.cpy, self.cpz) + + @classmethod + def from_dict(cls, data): + return cls( + angular_position=data["angular_position"], + shape_points=data["shape_points"], + rocket_radius=data["rocket_radius"], + cant_angle=data["cant_angle"], + airfoil=data["airfoil"], + name=data["name"], + ) diff --git a/rocketpy/rocket/aero_surface/fins/free_form_fins.py b/rocketpy/rocket/aero_surface/fins/free_form_fins.py index 7cae4e556..e16021d06 100644 --- a/rocketpy/rocket/aero_surface/fins/free_form_fins.py +++ b/rocketpy/rocket/aero_surface/fins/free_form_fins.py @@ -1,14 +1,11 @@ -import warnings - -import numpy as np - from rocketpy.plots.aero_surface_plots import _FreeFormFinsPlots from rocketpy.prints.aero_surface_prints import _FreeFormFinsPrints +from rocketpy.rocket.aero_surface.fins._free_form_mixin import _FreeFormMixin from .fins import Fins -class FreeFormFins(Fins): +class FreeFormFins(_FreeFormMixin, Fins): """Class that defines and holds information for a free form fin set. This class inherits from the Fins class. @@ -135,23 +132,7 @@ def __init__( ------- None """ - self.shape_points = shape_points - - down = False - for i in range(1, len(shape_points)): - if shape_points[i][1] > shape_points[i - 1][1] and down: - warnings.warn( - "Jagged fin shape detected. This may cause small inaccuracies " - "center of pressure and pitch moment calculations." - ) - break - if shape_points[i][1] < shape_points[i - 1][1]: - down = True - i += 1 - - root_chord = abs(shape_points[0][0] - shape_points[-1][0]) - ys = [point[1] for point in shape_points] - span = max(ys) - min(ys) + root_chord, span = self._initialize(shape_points) super().__init__( n, @@ -187,215 +168,13 @@ def evaluate_center_of_pressure(self): self.cpz = cpz self.cp = (self.cpx, self.cpy, self.cpz) - def evaluate_geometrical_parameters(self): # pylint: disable=too-many-statements - """ - Calculates and saves the fin set's geometrical parameters such as the - fin area, aspect ratio, and parameters related to roll movement. This - method uses the same calculations to those in OpenRocket for free-form - fin shapes. - - Returns - ------- - None - """ - # pylint: disable=invalid-name - # pylint: disable=too-many-locals - # Calculate the fin area (Af) using the Shoelace theorem (polygon area formula) - Af = 0 - for i in range(len(self.shape_points) - 1): - x1, y1 = self.shape_points[i] - x2, y2 = self.shape_points[i + 1] - Af += (y1 + y2) * (x1 - x2) - Af = abs(Af) / 2 - if Af < 1e-6: - raise ValueError("Fin area is too small. Check the shape_points.") - - # Calculate aspect ratio (AR) and lift interference factors - AR = 2 * self.span**2 / Af # Aspect ratio - tau = (self.span + self.rocket_radius) / self.rocket_radius - lift_interference_factor = 1 + 1 / tau - - # Calculate roll forcing interference factor using OpenRocket's approach - roll_forcing_interference_factor = (1 / np.pi**2) * ( - (np.pi**2 / 4) * ((tau + 1) ** 2 / tau**2) - + ((np.pi * (tau**2 + 1) ** 2) / (tau**2 * (tau - 1) ** 2)) - * np.arcsin((tau**2 - 1) / (tau**2 + 1)) - - (2 * np.pi * (tau + 1)) / (tau * (tau - 1)) - + ((tau**2 + 1) ** 2 / (tau**2 * (tau - 1) ** 2)) - * (np.arcsin((tau**2 - 1) / (tau**2 + 1))) ** 2 - - (4 * (tau + 1)) - / (tau * (tau - 1)) - * np.arcsin((tau**2 - 1) / (tau**2 + 1)) - + (8 / (tau - 1) ** 2) * np.log((tau**2 + 1) / (2 * tau)) - ) - - # Define number of interpolation points along the span of the fin - points_per_line = 40 # Same as OpenRocket - - # Initialize arrays for leading/trailing edge and chord lengths - chord_lead = np.ones(points_per_line) * np.inf # Leading edge x coordinates - chord_trail = np.ones(points_per_line) * -np.inf # Trailing edge x coordinates - chord_length = np.zeros( - points_per_line - ) # Chord length for each spanwise section - - # Discretize fin shape and calculate chord length, leading, and trailing edges - for p in range(1, len(self.shape_points)): - x1, y1 = self.shape_points[p - 1] - x2, y2 = self.shape_points[p] - - # Compute corresponding points along the fin span (clamp to valid range) - prev_idx = int(y1 / self.span * (points_per_line - 1)) - curr_idx = int(y2 / self.span * (points_per_line - 1)) - prev_idx = np.clip(prev_idx, 0, points_per_line - 1) - curr_idx = np.clip(curr_idx, 0, points_per_line - 1) - - if prev_idx > curr_idx: - prev_idx, curr_idx = curr_idx, prev_idx - - # Compute intersection of fin edge with each spanwise section - for i in range(prev_idx, curr_idx + 1): - y = i * self.span / (points_per_line - 1) - if y1 != y2: - x = np.clip( - (y - y2) / (y1 - y2) * x1 + (y1 - y) / (y1 - y2) * x2, - min(x1, x2), - max(x1, x2), - ) - else: - x = x1 # Handle horizontal segments - - # Update leading and trailing edge positions - chord_lead[i] = min(chord_lead[i], x) - chord_trail[i] = max(chord_trail[i], x) - - # Update chord length - if y1 < y2: - chord_length[i] -= x - else: - chord_length[i] += x - - # Replace infinities and handle invalid values in chord_lead and chord_trail - for i in range(points_per_line): - if ( - np.isinf(chord_lead[i]) - or np.isinf(chord_trail[i]) - or np.isnan(chord_lead[i]) - or np.isnan(chord_trail[i]) - ): - chord_lead[i] = 0 - chord_trail[i] = 0 - if chord_length[i] < 0 or np.isnan(chord_length[i]): - chord_length[i] = 0 - if chord_length[i] > chord_trail[i] - chord_lead[i]: - chord_length[i] = chord_trail[i] - chord_lead[i] - - # Initialize integration variables for various aerodynamic and roll properties - radius = self.rocket_radius - total_area = 0 - mac_length = 0 # Mean aerodynamic chord length - mac_lead = 0 # Mean aerodynamic chord leading edge - mac_span = 0 # Mean aerodynamic chord spanwise position (Yma) - cos_gamma_sum = 0 # Sum of cosine of the sweep angle - roll_geometrical_constant = 0 - roll_damping_numerator = 0 - roll_damping_denominator = 0 - - # Perform integration over spanwise sections - dy = self.span / (points_per_line - 1) - for i in range(points_per_line): - chord = chord_trail[i] - chord_lead[i] - y = i * dy - - # Update integration variables - mac_length += chord * chord - mac_span += y * chord - mac_lead += chord_lead[i] * chord - total_area += chord - roll_geometrical_constant += chord_length[i] * (radius + y) ** 2 - roll_damping_numerator += radius**3 * chord / (radius + y) ** 2 - roll_damping_denominator += (radius + y) * chord - - # Update cosine of sweep angle (cos_gamma) - if i > 0: - dx = (chord_trail[i] + chord_lead[i]) / 2 - ( - chord_trail[i - 1] + chord_lead[i - 1] - ) / 2 - cos_gamma_sum += dy / np.hypot(dx, dy) - - # Finalize mean aerodynamic chord properties - mac_length *= dy - mac_span *= dy - mac_lead *= dy - total_area *= dy - roll_geometrical_constant *= dy - roll_damping_numerator *= dy - roll_damping_denominator *= dy - - mac_length /= total_area - mac_span /= total_area - mac_lead /= total_area - cos_gamma = cos_gamma_sum / (points_per_line - 1) - - # Store computed values - self.Af = Af # Fin area - self.AR = AR # Aspect ratio - self.gamma_c = np.arccos(cos_gamma) # Sweep angle - self.Yma = mac_span # Mean aerodynamic chord spanwise position - self.mac_length = mac_length - self.mac_lead = mac_lead - self.tau = tau - self.roll_geometrical_constant = roll_geometrical_constant - self.lift_interference_factor = lift_interference_factor - self.roll_forcing_interference_factor = roll_forcing_interference_factor - self.roll_damping_interference_factor = 1 + ( - roll_damping_numerator / roll_damping_denominator - ) - - # Evaluate the shape and finalize geometry - self.evaluate_shape() - - def evaluate_shape(self): - x_array, y_array = zip(*self.shape_points) - self.shape_vec = [np.array(x_array), np.array(y_array)] - - def to_dict(self, include_outputs=False): - data = super().to_dict(include_outputs) - data["shape_points"] = self.shape_points - - if include_outputs: - data.update( - { - "Af": self.Af, - "AR": self.AR, - "gamma_c": self.gamma_c, - "Yma": self.Yma, - "mac_length": self.mac_length, - "mac_lead": self.mac_lead, - "roll_geometrical_constant": self.roll_geometrical_constant, - "tau": self.tau, - "lift_interference_factor": self.lift_interference_factor, - "roll_forcing_interference_factor": self.roll_forcing_interference_factor, - "roll_damping_interference_factor": self.roll_damping_interference_factor, - } - ) - return data - @classmethod def from_dict(cls, data): return cls( - data["n"], - data["shape_points"], - data["rocket_radius"], - data["cant_angle"], - data["airfoil"], - data["name"], + n=data["n"], + shape_points=data["shape_points"], + rocket_radius=data["rocket_radius"], + cant_angle=data["cant_angle"], + airfoil=data["airfoil"], + name=data["name"], ) - - def info(self): - self.prints.geometry() - self.prints.lift() - - def all_info(self): - self.prints.all() - self.plots.all() diff --git a/rocketpy/rocket/aero_surface/fins/trapezoidal_fin.py b/rocketpy/rocket/aero_surface/fins/trapezoidal_fin.py new file mode 100644 index 000000000..f5aed6688 --- /dev/null +++ b/rocketpy/rocket/aero_surface/fins/trapezoidal_fin.py @@ -0,0 +1,197 @@ +from rocketpy.plots.aero_surface_plots import _TrapezoidalFinPlots +from rocketpy.prints.aero_surface_prints import _TrapezoidalFinPrints +from rocketpy.rocket.aero_surface.fins._trapezoidal_mixin import _TrapezoidalMixin + +from .fin import Fin + + +class TrapezoidalFin(_TrapezoidalMixin, Fin): + """A class used to represent a single trapezoidal fin. + + This class inherits from the Fin class. + + Note + ---- + Local coordinate system: + - Origin located at the top of the root chord. + - Z axis along the longitudinal axis of symmetry, positive downwards (top -> bottom). + - Y axis perpendicular to the Z axis, in the span direction, positive upwards. + - X axis completes the right-handed coordinate system. + + See Also + -------- + Fin : Parent class + + Attributes + ---------- + TrapezoidalFin.angular_position : float + Angular position of the fin set with respect to the rocket centerline, + in degrees. + TrapezoidalFin.rocket_radius : float + The reference rocket radius used for lift coefficient normalization, in + meters. + TrapezoidalFin.airfoil : tuple + Tuple of two items. First is the airfoil lift curve. + Second is the unit of the curve (radians or degrees). + TrapezoidalFin.cant_angle : float + Fins cant angle with respect to the rocket centerline, in degrees. + TrapezoidalFin.cant_angle_rad : float + Fins cant angle with respect to the rocket centerline, in radians. + TrapezoidalFin.root_chord : float + Fin root chord in meters. + TrapezoidalFin.tip_chord : float + Fin tip chord in meters. + TrapezoidalFin.span : float + Fin span in meters. + TrapezoidalFin.name : string + Name of fin set. + TrapezoidalFin.sweep_length : float + Fins sweep length in meters. By sweep length, understand the axial + distance between the fin root leading edge and the fin tip leading edge + measured parallel to the rocket centerline. + TrapezoidalFin.sweep_angle : float + Fins sweep angle with respect to the rocket centerline. Must + be given in degrees. + TrapezoidalFin.d : float + Reference diameter of the rocket, in meters. + TrapezoidalFins.fin_area : float + Area of the longitudinal section of each fin in the set. + TrapezoidalFins.f_ar : float + Aspect ratio of each fin in the set + TrapezoidalFin.gamma_c : float + Fin mid-chord sweep angle. + TrapezoidalFin.yma : float + Span wise position of the mean aerodynamic chord. + TrapezoidalFin.roll_geometrical_constant : float + Geometrical constant used in roll calculations. + TrapezoidalFin.tau : float + Geometrical relation used to simplify lift and roll calculations. + TrapezoidalFin.lift_interference_factor : float + Factor of Fin-Body interference in the lift coefficient. + TrapezoidalFin.cp : tuple + Tuple with the x, y and z local coordinates of the fin set center of + pressure. Has units of length and is given in meters. + TrapezoidalFin.cpx : float + Fin set local center of pressure x coordinate. Has units of length and + is given in meters. + TrapezoidalFin.cpy : float + Fin set local center of pressure y coordinate. Has units of length and + is given in meters. + TrapezoidalFin.cpz : float + Fin set local center of pressure z coordinate. Has units of length and + is given in meters. + """ + + def __init__( + self, + angular_position, + root_chord, + tip_chord, + span, + rocket_radius, + cant_angle=0, + sweep_length=None, + sweep_angle=None, + airfoil=None, + name="Fins", + ): + """Initializes the TrapezoidalFin class. + + Parameters + ---------- + angular_position : float + Angular position of the fin in degrees measured as the rotation + around the symmetry axis of the rocket relative to one of the other + principal axis. See :ref:`Angular Position Inputs ` + root_chord : int, float + Fin root chord in meters. + tip_chord : int, float + Fin tip chord in meters. + span : int, float + Fin span in meters. + rocket_radius : int, float + Reference radius to calculate lift coefficient, in meters. + cant_angle : int, float, optional + Fins cant angle with respect to the rocket centerline. Must + be given in degrees. + sweep_length : int, float, optional + Fins sweep length in meters. By sweep length, understand the axial + distance between the fin root leading edge and the fin tip leading + edge measured parallel to the rocket centerline. If not given, the + sweep length is assumed to be equal the root chord minus the tip + chord, in which case the fin is a right trapezoid with its base + perpendicular to the rocket's axis. Cannot be used in conjunction + with sweep_angle. + sweep_angle : int, float, optional + Fins sweep angle with respect to the rocket centerline. Must + be given in degrees. If not given, the sweep angle is automatically + calculated, in which case the fin is assumed to be a right trapezoid + with its base perpendicular to the rocket's axis. + Cannot be used in conjunction with sweep_length. + airfoil : tuple, optional + Default is null, in which case fins will be treated as flat plates. + Otherwise, if tuple, fins will be considered as airfoils. The + tuple's first item specifies the airfoil's lift coefficient + by angle of attack and must be either a .csv, .txt, ndarray + or callable. The .csv and .txt files can contain a single line + header and the first column must specify the angle of attack, while + the second column must specify the lift coefficient. The + ndarray should be as [(x0, y0), (x1, y1), (x2, y2), ...] + where x0 is the angle of attack and y0 is the lift coefficient. + If callable, it should take an angle of attack as input and + return the lift coefficient at that angle of attack. + The tuple's second item is the unit of the angle of attack, + accepting either "radians" or "degrees". + name : str + Name of fin set. + """ + super().__init__( + angular_position, + root_chord, + span, + rocket_radius, + cant_angle, + airfoil, + name, + ) + + self._initialize(sweep_length, sweep_angle, root_chord, tip_chord, span) + self.evaluate_rotation_matrix() + + self.prints = _TrapezoidalFinPrints(self) + self.plots = _TrapezoidalFinPlots(self) + + def evaluate_center_of_pressure(self): + """Calculates and returns the center of pressure of the fin set in local + coordinates. The center of pressure position is saved and stored as a + tuple. + + Returns + ------- + None + """ + # Center of pressure position in local coordinates + cpz = (self.sweep_length / 3) * ( + (self.root_chord + 2 * self.tip_chord) / (self.root_chord + self.tip_chord) + ) + (1 / 6) * ( + self.root_chord + + self.tip_chord + - self.root_chord * self.tip_chord / (self.root_chord + self.tip_chord) + ) + self.cpx = 0 + self.cpy = self.Yma + self.cpz = cpz + self.cp = (self.cpx, self.cpy, self.cpz) + + @classmethod + def from_dict(cls, data): + return cls( + angular_position=data["angular_position"], + root_chord=data["root_chord"], + tip_chord=data["tip_chord"], + span=data["span"], + rocket_radius=data["rocket_radius"], + cant_angle=data["cant_angle"], + airfoil=data["airfoil"], + name=data["name"], + ) diff --git a/rocketpy/rocket/aero_surface/fins/trapezoidal_fins.py b/rocketpy/rocket/aero_surface/fins/trapezoidal_fins.py index 61e2b78fc..2cf9686de 100644 --- a/rocketpy/rocket/aero_surface/fins/trapezoidal_fins.py +++ b/rocketpy/rocket/aero_surface/fins/trapezoidal_fins.py @@ -1,12 +1,11 @@ -import numpy as np - from rocketpy.plots.aero_surface_plots import _TrapezoidalFinsPlots from rocketpy.prints.aero_surface_prints import _TrapezoidalFinsPrints +from rocketpy.rocket.aero_surface.fins._trapezoidal_mixin import _TrapezoidalMixin from .fins import Fins -class TrapezoidalFins(Fins): +class TrapezoidalFins(_TrapezoidalMixin, Fins): """Class that defines and holds information for a trapezoidal fin set. This class inherits from the Fins class. @@ -35,9 +34,6 @@ class TrapezoidalFins(Fins): Second is the unit of the curve (radians or degrees). TrapezoidalFins.cant_angle : float Fins cant angle with respect to the rocket centerline, in degrees. - TrapezoidalFins.changing_attribute_dict : dict - Dictionary that stores the name and the values of the attributes that - may be changed during a simulation. Useful for control systems. TrapezoidalFins.cant_angle_rad : float Fins cant angle with respect to the rocket centerline, in radians. TrapezoidalFins.root_chord : float @@ -169,66 +165,11 @@ def __init__( name, ) - # Check if sweep angle or sweep length is given - if sweep_length is not None and sweep_angle is not None: - raise ValueError("Cannot use sweep_length and sweep_angle together") - elif sweep_angle is not None: - sweep_length = np.tan(sweep_angle * np.pi / 180) * span - elif sweep_length is None: - sweep_length = root_chord - tip_chord - else: - # Sweep length is given - pass - - self._tip_chord = tip_chord - self._sweep_length = sweep_length - self._sweep_angle = sweep_angle - - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() + self._initialize(sweep_length, sweep_angle, root_chord, tip_chord, span) self.prints = _TrapezoidalFinsPrints(self) self.plots = _TrapezoidalFinsPlots(self) - @property - def tip_chord(self): - return self._tip_chord - - @tip_chord.setter - def tip_chord(self, value): - self._tip_chord = value - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - - @property - def sweep_angle(self): - return self._sweep_angle - - @sweep_angle.setter - def sweep_angle(self, value): - self._sweep_angle = value - self._sweep_length = np.tan(value * np.pi / 180) * self.span - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - - @property - def sweep_length(self): - return self._sweep_length - - @sweep_length.setter - def sweep_length(self, value): - self._sweep_length = value - self.evaluate_geometrical_parameters() - self.evaluate_center_of_pressure() - self.evaluate_lift_coefficient() - self.evaluate_roll_parameters() - def evaluate_center_of_pressure(self): """Calculates and returns the center of pressure of the fin set in local coordinates. The center of pressure position is saved and stored as a @@ -251,126 +192,6 @@ def evaluate_center_of_pressure(self): self.cpz = cpz self.cp = (self.cpx, self.cpy, self.cpz) - def evaluate_geometrical_parameters(self): # pylint: disable=too-many-statements - """Calculates and saves fin set's geometrical parameters such as the - fins' area, aspect ratio and parameters for roll movement. - - Returns - ------- - None - """ - # pylint: disable=invalid-name - Yr = self.root_chord + self.tip_chord - Af = Yr * self.span / 2 # Fin area - AR = 2 * self.span**2 / Af # Fin aspect ratio - gamma_c = np.arctan( - (self.sweep_length + 0.5 * self.tip_chord - 0.5 * self.root_chord) - / (self.span) - ) - Yma = ( - (self.span / 3) * (self.root_chord + 2 * self.tip_chord) / Yr - ) # Span wise coord of mean aero chord - - # Fin–body interference correction parameters - tau = (self.span + self.rocket_radius) / self.rocket_radius - lift_interference_factor = 1 + 1 / tau - lambda_ = self.tip_chord / self.root_chord - - # Parameters for Roll Moment. - # Documented at: https://docs.rocketpy.org/en/latest/technical/ - roll_geometrical_constant = ( - (self.root_chord + 3 * self.tip_chord) * self.span**3 - + 4 - * (self.root_chord + 2 * self.tip_chord) - * self.rocket_radius - * self.span**2 - + 6 * (self.root_chord + self.tip_chord) * self.span * self.rocket_radius**2 - ) / 12 - roll_damping_interference_factor = 1 + ( - ((tau - lambda_) / (tau)) - ((1 - lambda_) / (tau - 1)) * np.log(tau) - ) / ( - ((tau + 1) * (tau - lambda_)) / (2) - - ((1 - lambda_) * (tau**3 - 1)) / (3 * (tau - 1)) - ) - roll_forcing_interference_factor = (1 / np.pi**2) * ( - (np.pi**2 / 4) * ((tau + 1) ** 2 / tau**2) - + ((np.pi * (tau**2 + 1) ** 2) / (tau**2 * (tau - 1) ** 2)) - * np.arcsin((tau**2 - 1) / (tau**2 + 1)) - - (2 * np.pi * (tau + 1)) / (tau * (tau - 1)) - + ((tau**2 + 1) ** 2) - / (tau**2 * (tau - 1) ** 2) - * (np.arcsin((tau**2 - 1) / (tau**2 + 1))) ** 2 - - (4 * (tau + 1)) - / (tau * (tau - 1)) - * np.arcsin((tau**2 - 1) / (tau**2 + 1)) - + (8 / (tau - 1) ** 2) * np.log((tau**2 + 1) / (2 * tau)) - ) - - # Store values - self.Yr = Yr - self.Af = Af # Fin area - self.AR = AR # Aspect Ratio - self.gamma_c = gamma_c # Mid chord angle - self.Yma = Yma # Span wise coord of mean aero chord - self.roll_geometrical_constant = roll_geometrical_constant - self.tau = tau - self.lift_interference_factor = lift_interference_factor - self.λ = lambda_ # pylint: disable=non-ascii-name - self.roll_damping_interference_factor = roll_damping_interference_factor - self.roll_forcing_interference_factor = roll_forcing_interference_factor - - self.evaluate_shape() - - def evaluate_shape(self): - if self.sweep_length: - points = [ - (0, 0), - (self.sweep_length, self.span), - (self.sweep_length + self.tip_chord, self.span), - (self.root_chord, 0), - ] - else: - points = [ - (0, 0), - (self.root_chord - self.tip_chord, self.span), - (self.root_chord, self.span), - (self.root_chord, 0), - ] - - x_array, y_array = zip(*points) - self.shape_vec = [np.array(x_array), np.array(y_array)] - - def info(self): - self.prints.geometry() - self.prints.lift() - - def all_info(self): - self.prints.all() - self.plots.all() - - def to_dict(self, include_outputs=False): - data = super().to_dict(include_outputs) - data["tip_chord"] = self.tip_chord - - if include_outputs: - data.update( - { - "sweep_length": self.sweep_length, - "sweep_angle": self.sweep_angle, - "shape_vec": self.shape_vec, - "Af": self.Af, - "AR": self.AR, - "gamma_c": self.gamma_c, - "Yma": self.Yma, - "roll_geometrical_constant": self.roll_geometrical_constant, - "tau": self.tau, - "lift_interference_factor": self.lift_interference_factor, - "roll_damping_interference_factor": self.roll_damping_interference_factor, - "roll_forcing_interference_factor": self.roll_forcing_interference_factor, - } - ) - return data - @classmethod def from_dict(cls, data): return cls( diff --git a/rocketpy/rocket/aero_surface/generic_surface.py b/rocketpy/rocket/aero_surface/generic_surface.py index d1982ae04..8ca6afd22 100644 --- a/rocketpy/rocket/aero_surface/generic_surface.py +++ b/rocketpy/rocket/aero_surface/generic_surface.py @@ -84,6 +84,8 @@ def __init__( self.cpz = center_of_pressure[2] self.name = name + self._rotation_surface_to_body = Matrix([[1, 0, 0], [0, 1, 0], [0, 0, 1]]) + default_coefficients = self._get_default_coefficients() self._check_coefficients(coefficients, default_coefficients) coefficients = self._complete_coefficients(coefficients, default_coefficients) diff --git a/rocketpy/rocket/aero_surface/rail_buttons.py b/rocketpy/rocket/aero_surface/rail_buttons.py index 879f013d7..470e55306 100644 --- a/rocketpy/rocket/aero_surface/rail_buttons.py +++ b/rocketpy/rocket/aero_surface/rail_buttons.py @@ -17,6 +17,7 @@ class RailButtons(AeroSurface): Angular position of the rail buttons in degrees measured as the rotation around the symmetry axis of the rocket relative to one of the other principal axis. + See :ref:`Angular Position Inputs ` RailButtons.angular_position_rad : float Angular position of the rail buttons in radians. """ diff --git a/rocketpy/rocket/components.py b/rocketpy/rocket/components.py index 66448eb69..a604a997f 100644 --- a/rocketpy/rocket/components.py +++ b/rocketpy/rocket/components.py @@ -1,4 +1,5 @@ from collections import namedtuple +from copy import deepcopy class Components: @@ -186,7 +187,7 @@ def clear(self): self._components.clear() def sort_by_position(self, reverse=False): - """Sort the list of components by z axis position. + """Returns a new Components object sorted components by z axis position. Parameters ---------- @@ -196,9 +197,12 @@ def sort_by_position(self, reverse=False): Returns ------- - None + Components + A new Components object sorted by component position. """ - self._components.sort(key=lambda x: x.position.z, reverse=reverse) + components = deepcopy(self) + components._components.sort(key=lambda x: x.position.z, reverse=reverse) + return components def to_dict(self, include_outputs=False): # pylint: disable=unused-argument return { diff --git a/rocketpy/rocket/rocket.py b/rocketpy/rocket/rocket.py index bf938d4be..93a1aa9df 100644 --- a/rocketpy/rocket/rocket.py +++ b/rocketpy/rocket/rocket.py @@ -18,7 +18,10 @@ Tail, TrapezoidalFins, ) +from rocketpy.rocket.aero_surface.fins.elliptical_fin import EllipticalFin +from rocketpy.rocket.aero_surface.fins.free_form_fin import FreeFormFin from rocketpy.rocket.aero_surface.fins.free_form_fins import FreeFormFins +from rocketpy.rocket.aero_surface.fins.trapezoidal_fin import TrapezoidalFin from rocketpy.rocket.aero_surface.generic_surface import GenericSurface from rocketpy.rocket.components import Components from rocketpy.rocket.parachute import Parachute @@ -595,14 +598,20 @@ def __evaluate_single_surface_cp_to_cdm(self, surface, position): """Calculates the relative position of each aerodynamic surface center of pressure to the rocket's center of dry mass in Body Axes Coordinate System.""" - pos = Vector( + # position of the surfaces coordinate system origin in body frame + pos_origin = Vector( [ - (position.x - self.cm_eccentricity_x) * self._csys - surface.cpx, - (position.y - self.cm_eccentricity_y) - surface.cpy, - (position.z - self.center_of_dry_mass_position) * self._csys - - surface.cpz, + (position.x - self.cm_eccentricity_x) * self._csys, + (position.y - self.cm_eccentricity_y), + (position.z - self.center_of_dry_mass_position) * self._csys, ] ) + # position of the center of pressure in body frame + pos = ( + surface._rotation_surface_to_body + @ Vector([surface.cpx, surface.cpy, surface.cpz]) + + pos_origin + ) # TODO: this should be recomputed whenever cant angle changes for fin self.surfaces_cp_to_cdm[surface] = pos def evaluate_stability_margin(self): @@ -996,11 +1005,17 @@ def __add_single_surface(self, surface, position): """Adds a single aerodynamic surface to the rocket. Makes checks for rail buttons case, and position type. """ - position = ( - Vector([0, 0, position]) - if not isinstance(position, (Vector, tuple, list)) - else Vector(position) - ) + if isinstance(surface, (TrapezoidalFin, EllipticalFin, FreeFormFin)): + # TODO: this depends on cant angle, so it should somehow be + # recalculated whenever the cant angle of the fin changes + position = surface._compute_leading_edge_position(position, self._csys) + else: + position = ( + Vector([0, 0, position]) + if not isinstance(position, (Vector, tuple, list)) + else Vector(position) + ) + if isinstance(surface, RailButtons): self.rail_buttons = Components() self.rail_buttons.add(surface, position) @@ -1018,14 +1033,20 @@ def add_surfaces(self, surfaces, positions): surfaces : list, AeroSurface, NoseCone, TrapezoidalFins, EllipticalFins, Tail, RailButtons Aerodynamic surface to be added to the rocket. Can be a list of AeroSurface if more than one surface is to be added. - positions : int, float, list, tuple, Vector - Position, in m, of the aerodynamic surface's center of pressure - relative to the user defined rocket coordinate system. - If a list is passed, it will correspond to the position of each item - in the surfaces list. - For NoseCone type, position is relative to the nose cone tip. - For Fins type, position is relative to the point belonging to - the root chord which is highest in the rocket coordinate system. + positions : int, float, tuple, list, Vector + Position(s) of the aerodynamic surface's reference point. Can be: + + - a single number (int or float) giving the z-coordinate along + the rocket axis. + - a sequence of three numbers (x, y, z) representing the full + position in the user-defined coordinate system. + + If passing multiple surfaces, provide a list of positions matching + each surface in order. + For NoseCone type, position is the tip coordinate along the axis. + 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 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. @@ -1038,10 +1059,18 @@ def add_surfaces(self, surfaces, positions): ------- None """ - try: + if isinstance(surfaces, list): + if isinstance(positions, list): + if len(surfaces) != len(positions): + raise ValueError( + "The number of surfaces and positions must be the same." + ) + else: + positions = [positions] * len(surfaces) + for surface, position in zip(surfaces, positions): self.__add_single_surface(surface, position) - except TypeError: + else: self.__add_single_surface(surfaces, positions) self.evaluate_center_of_pressure() @@ -1215,10 +1244,10 @@ def add_trapezoidal_fins( tip_chord : int, float Fin tip chord in meters. position : int, float - Fin set position relative to the rocket's coordinate system. - By fin set position, understand the point belonging to the root - chord which is highest in the rocket coordinate system (i.e. - the point closest to the nose cone tip). + Fin set position in the z coordinate of the user defined rocket + coordinate system. By fin set position, understand the point + belonging to the root chord which is highest in the rocket + coordinate system (i.e. the point closest to the nose cone tip). See Also -------- @@ -1264,6 +1293,12 @@ def add_trapezoidal_fins( fin_set : TrapezoidalFins Fin set object created. """ + if n <= 2: + raise ValueError( + "Number of fins must be greater than 2." + "If you want to add 2 or 1 fins, create a TrapezoidalFin object " + " and add it to the rocket using the add_surfaces method." + ) # Modify radius if not given, use rocket radius, otherwise use given. radius = radius if radius is not None else self.radius @@ -1311,10 +1346,10 @@ def add_elliptical_fins( span : int, float Fin span in meters. position : int, float - Fin set position relative to the rocket's coordinate system. By fin - set position, understand the point belonging to the root chord which - is highest in the rocket coordinate system (i.e. the point - closest to the nose cone tip). + Fin set position in the z coordinate of the user defined rocket + coordinate system. By fin set position, understand the point + belonging to the root chord which is highest in the rocket + coordinate system (i.e. the point closest to the nose cone tip). See Also -------- @@ -1350,6 +1385,13 @@ def add_elliptical_fins( fin_set : EllipticalFins Fin set object created. """ + if n <= 2: + raise ValueError( + "Number of fins must be greater than 2." + "If you want to add 2 or 1 fins, create a EllipticalFin object " + " and add it to the rocket using the add_surfaces method." + ) + radius = radius if radius is not None else self.radius fin_set = EllipticalFins(n, root_chord, span, radius, cant_angle, airfoil, name) self.add_surfaces(fin_set, position) @@ -1381,10 +1423,10 @@ def add_free_form_fins( The shape will be interpolated between the points, in the order they are given. The last point connects to the first point. position : int, float - Fin set position relative to the rocket's coordinate system. - By fin set position, understand the point belonging to the root - chord which is highest in the rocket coordinate system (i.e. - the point closest to the nose cone tip). + Fin set position in the z coordinate of the user defined rocket + coordinate system. By fin set position, understand the point + belonging to the root chord which is highest in the rocket + coordinate system (i.e. the point closest to the nose cone tip). See Also -------- @@ -1416,6 +1458,12 @@ def add_free_form_fins( fin_set : FreeFormFins Fin set object created. """ + if n <= 2: + raise ValueError( + "Number of fins must be greater than 2." + "If you want to add 2 or 1 fins, create a FreeFormFin object " + " and add it to the rocket using the add_surfaces method." + ) # Modify radius if not given, use rocket radius, otherwise use given. radius = radius if radius is not None else self.radius @@ -1518,8 +1566,7 @@ def add_sensor(self, sensor, position): must be in the format (x, y, z) where x, y, and z are defined in the rocket's user defined coordinate system. If a single value is passed, it is assumed to be along the z-axis (centerline) of the - rocket's user defined coordinate system and angular_position and - radius must be given. + rocket's user defined coordinate system. Returns ------- @@ -1603,9 +1650,9 @@ def add_air_brakes( listed in the same order as they are provided in the `interactive_objects` 7. `sensors` (list): A list of sensors that are attached to the - rocket. The most recent measurements of the sensors are provided - with the ``sensor.measurement`` attribute. The sensors are - listed in the same order as they are added to the rocket + rocket. The most recent measurements of the sensors are provided + with the ``sensor.measurement`` attribute. The sensors are + listed in the same order as they are added to the rocket ``interactive_objects`` This function will be called during the simulation at the specified @@ -1710,7 +1757,7 @@ def set_rail_buttons( as the rotation around the symmetry axis of the rocket relative to one of the other principal axis. Default value is 45 degrees, generally used in rockets with - 4 fins. + 4 fins. See :ref:`Angular Position Inputs ` radius : int, float, optional Fuselage radius where the rail buttons are located. diff --git a/tests/integration/test_rocket.py b/tests/integration/test_rocket.py index c47096617..e36dd933e 100644 --- a/tests/integration/test_rocket.py +++ b/tests/integration/test_rocket.py @@ -19,7 +19,7 @@ def test_airfoil( calisto.add_surfaces(calisto_tail, -1.313) test_rocket.add_trapezoidal_fins( - 2, + 3, span=0.100, root_chord=0.120, tip_chord=0.040, @@ -28,7 +28,7 @@ def test_airfoil( name="NACA0012", ) test_rocket.add_trapezoidal_fins( - 2, + 3, span=0.100, root_chord=0.120, tip_chord=0.040,