Source code for pterasoftware.problems

"""Contains the SteadyProblem and UnsteadyProblem classes.

**Contains the following classes:**

SteadyProblem: A class used to contain steady aerodynamics problems.

UnsteadyProblem: A class used to contain unsteady aerodynamics problems.

**Contains the following functions:**

None
"""

from __future__ import annotations

import math

import numpy as np

from . import _parameter_validation, _transformations, geometry, movements
from . import operating_point as operating_point_mod


[docs] class SteadyProblem: """A class used to contain steady aerodynamics problems. **Contains the following methods:** reynolds_numbers: A tuple of Reynolds numbers, one for each Airplane in the SteadyProblem. """ __slots__ = ( "_airplanes", "_operating_point", "_reynolds_numbers", ) def __init__( self, airplanes: list[geometry.airplane.Airplane], operating_point: operating_point_mod.OperatingPoint, ) -> None: """The initialization method. :param airplanes: The list of the Airplanes for this SteadyProblem. :param operating_point: The OperatingPoint for this SteadyProblem. :return: None """ # Validate and store immutable attributes. if not isinstance(airplanes, list): raise TypeError("airplanes must be a list.") if len(airplanes) < 1: raise ValueError("airplanes must have at least one element.") for airplane in airplanes: if not isinstance(airplane, geometry.airplane.Airplane): raise TypeError("Every element in airplanes must be an Airplane.") # Store as tuple to prevent external mutation via .append(), .pop(), etc. self._airplanes: tuple[geometry.airplane.Airplane, ...] = tuple(airplanes) if not isinstance(operating_point, operating_point_mod.OperatingPoint): raise TypeError("operating_point must be an OperatingPoint.") self._operating_point = operating_point # Initialize the caches for the properties derived from the immutable # attributes. self._reynolds_numbers: tuple[float, ...] | None = None # Validate that the first Airplane has Cg_GP1_CgP1 set to zeros. self._airplanes[0].validate_first_airplane_constraints() # Populate GP1_CgP1 coordinates for all Airplanes' Panels. This finds the # Panels' positions in the first Airplane's geometry axes, relative to the # first Airplane's CG based on their locally defined positions. for airplane in self._airplanes: # Compute the passive transformation matrix from this Airplane's local # geometry axes, relative to its CG, to the first Airplane's geometry axes, # relative to the first Airplane's CG. T_pas_G_Cg_to_GP1_CgP1 = airplane.T_pas_G_Cg_to_GP1_CgP1 for wing in airplane.wings: assert wing.panels is not None for panel in np.ravel(wing.panels): panel.Frpp_GP1_CgP1 = _transformations.apply_T_to_vectors( T_pas_G_Cg_to_GP1_CgP1, panel.Frpp_G_Cg, has_point=True ) panel.Flpp_GP1_CgP1 = _transformations.apply_T_to_vectors( T_pas_G_Cg_to_GP1_CgP1, panel.Flpp_G_Cg, has_point=True ) panel.Blpp_GP1_CgP1 = _transformations.apply_T_to_vectors( T_pas_G_Cg_to_GP1_CgP1, panel.Blpp_G_Cg, has_point=True ) panel.Brpp_GP1_CgP1 = _transformations.apply_T_to_vectors( T_pas_G_Cg_to_GP1_CgP1, panel.Brpp_G_Cg, has_point=True ) # --- Immutable: read only properties --- @property def airplanes(self) -> tuple[geometry.airplane.Airplane, ...]: return self._airplanes @property def operating_point(self) -> operating_point_mod.OperatingPoint: return self._operating_point # --- Immutable derived: manual lazy caching --- @property def reynolds_numbers(self) -> tuple[float, ...]: """A tuple of Reynolds numbers, one for each Airplane in the SteadyProblem. **Notes:** The Reynolds number is calculated as: Re = (V x L) / nu, where V is the freestream speed, observed from the Earth frame (vCg__E from OperatingPoint, m/s), L is the characteristic length (c_ref from Airplane, m), and nu is the kinematic viscosity (nu from OperatingPoint, m^2/s). These Reynolds numbers only consider the freestream speed, not any apparent velocity due to prescribed motion, so be careful interpreting it for cases where this SteadyProblem corresponds to one time step in an UnsteadyProblem. :return: A tuple of Reynolds numbers, one for each Airplane. """ if self._reynolds_numbers is None: v = self._operating_point.vCg__E nu = self._operating_point.nu reynolds_list = [] for airplane in self._airplanes: c_ref = airplane.c_ref assert c_ref is not None, "Airplane c_ref must be set to calculate Re" re = (v * c_ref) / nu reynolds_list.append(re) # Store as tuple to prevent external mutation. self._reynolds_numbers = tuple(reynolds_list) return self._reynolds_numbers
[docs] class UnsteadyProblem: """A class used to contain unsteady aerodynamics problems. **Contains the following methods:** None """ __slots__ = ( "_movement", "_only_final_results", "_num_steps", "_delta_time", "_max_wake_rows", "_first_averaging_step", "_first_results_step", "finalForces_W", "finalForceCoefficients_W", "finalMoments_W_CgP1", "finalMomentCoefficients_W_CgP1", "finalMeanForces_W", "finalMeanForceCoefficients_W", "finalMeanMoments_W_CgP1", "finalMeanMomentCoefficients_W_CgP1", "finalRmsForces_W", "finalRmsForceCoefficients_W", "finalRmsMoments_W_CgP1", "finalRmsMomentCoefficients_W_CgP1", "_steady_problems", ) def __init__( self, movement: movements.movement.Movement, only_final_results: bool | np.bool_ = False, ) -> None: """The initialization method. :param movement: The Movement that contains this UnsteadyProblem's OperatingPointMovement and AirplaneMovements. :param only_final_results: Determines whether the Solver will only calculate loads for the final time step (for static Movements) or (for non static Movements) for will only calculate loads for the time steps in the final complete motion cycle (of the Movement's sub Movement with the longest period), which increases simulation speed. Can be a bool or a numpy bool and will be converted internally to a bool. The default is False. :return: None """ # Validate and store immutable attributes. if not isinstance(movement, movements.movement.Movement): raise TypeError("movement must be a Movement.") self._movement = movement self._only_final_results = _parameter_validation.boolLike_return_bool( only_final_results, "only_final_results" ) self._num_steps: int = self._movement.num_steps self._delta_time: float = self._movement.delta_time self._max_wake_rows: int | None = self._movement.max_wake_rows # For UnsteadyProblems with a static Movement, we are typically interested in # the final time step's forces and moments, which, assuming convergence, will be # the most accurate. For UnsteadyProblems with cyclic movement, (e.g. flapping # wings) we are typically interested in the forces and moments averaged over the # last cycle simulated. Use the LCM of all motion periods to ensure we average # over a complete cycle of all motions. _movement_lcm_period = self._movement.lcm_period self._first_averaging_step: int if _movement_lcm_period == 0: self._first_averaging_step = self._num_steps - 1 else: self._first_averaging_step = max( 0, math.floor(self._num_steps - (_movement_lcm_period / self._delta_time)), ) # If we only wants to calculate forces and moments for the final cycle (for a # cyclic Movement) or for the final time step (for a static Movement) set the # first step to calculate results to the first averaging step. Otherwise, set it # to the zero, which is the first time step. self._first_results_step: int if self._only_final_results: self._first_results_step = self._first_averaging_step else: self._first_results_step = 0 # Initialize empty lists to hold the final loads and load coefficients each # Airplane experiences. These will only be populated if this UnsteadyProblem's # Movement is static. These are mutable and populated by the solver. self.finalForces_W: list[np.ndarray] = [] self.finalForceCoefficients_W: list[np.ndarray] = [] self.finalMoments_W_CgP1: list[np.ndarray] = [] self.finalMomentCoefficients_W_CgP1: list[np.ndarray] = [] # Initialize empty lists to hold the final cycle-averaged loads and load # coefficients each Airplane experiences. These will only be populated if this # UnsteadyProblem's Movement is cyclic. These are mutable and populated by the # solver. self.finalMeanForces_W: list[np.ndarray] = [] self.finalMeanForceCoefficients_W: list[np.ndarray] = [] self.finalMeanMoments_W_CgP1: list[np.ndarray] = [] self.finalMeanMomentCoefficients_W_CgP1: list[np.ndarray] = [] # Initialize empty lists to hold the final cycle-root-mean-squared loads and # load coefficients each airplane object experiences. These will only be # populated for variable geometry problems. These are mutable and populated by # the solver. self.finalRmsForces_W: list[np.ndarray] = [] self.finalRmsForceCoefficients_W: list[np.ndarray] = [] self.finalRmsMoments_W_CgP1: list[np.ndarray] = [] self.finalRmsMomentCoefficients_W_CgP1: list[np.ndarray] = [] # Initialize an empty list to hold the SteadyProblems as they are generated. steady_problems_temp: list[SteadyProblem] = [] # Iterate through the UnsteadyProblem's time steps. for step_id in range(self._num_steps): # Get the Airplanes and the OperatingPoint associated with this time step. these_airplanes = [] for this_base_airplane in movement.airplanes: these_airplanes.append(this_base_airplane[step_id]) this_operating_point = movement.operating_points[step_id] # Initialize the SteadyProblem at this time step. this_steady_problem = SteadyProblem( airplanes=these_airplanes, operating_point=this_operating_point ) # Append this SteadyProblem to the temporary list. steady_problems_temp.append(this_steady_problem) # Store as tuple to prevent external mutation via .append(), .pop(), etc. self._steady_problems: tuple[SteadyProblem, ...] = tuple(steady_problems_temp) # --- Immutable: read only properties --- @property def movement(self) -> movements.movement.Movement: return self._movement @property def only_final_results(self) -> bool: return self._only_final_results @property def num_steps(self) -> int: return self._num_steps @property def delta_time(self) -> float: return self._delta_time @property def first_averaging_step(self) -> int: return self._first_averaging_step @property def first_results_step(self) -> int: return self._first_results_step @property def max_wake_rows(self) -> int | None: return self._max_wake_rows @property def steady_problems(self) -> tuple[SteadyProblem, ...]: return self._steady_problems