diff --git a/CHANGELOG.md b/CHANGELOG.md index e37e9f030..e0695252f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -32,6 +32,7 @@ Attention: The newest changes should be on top --> ### Added +- ENH: `StochasticFreeFormFins`, so free-form fin sets can be used in Monte Carlo simulations. The outline is randomized as a block, since a shape is only meaningful as a complete set of points: every coordinate is perturbed by its own draw, the fin root is held on the body line, and a list of candidate outlines can have a different number of points in each. [#953](https://github.com/RocketPy-Team/RocketPy/issues/953) - ENH: Support for Open Meteo API in the `Environment` class, adding the `open_meteo` and `open_meteo_ensemble` atmospheric models. Pressure-level forecasts, past forecasts (from 2021 onwards) and ensembles are read straight from a keyless JSON API, with no external files and no netCDF/OPeNDAP dependency. [#520](https://github.com/RocketPy-Team/RocketPy/issues/520) [#1119](https://github.com/RocketPy-Team/RocketPy/pull/1119) - ENH: Add simplified opening shock force estimation [#1092](https://github.com/RocketPy-Team/RocketPy/pull/1092) - ENH: Add Qodo PR-Agent workflow using Google Gemini [#1089](https://github.com/RocketPy-Team/RocketPy/pull/1089) @@ -47,9 +48,11 @@ Attention: The newest changes should be on top --> ### Fixed +- BUG: Pick between the candidate values of a list input with the stochastic model's own seeded generator. `random.choice` was used, which draws from the interpreter-wide stream that `_set_stochastic` does not reseed, so a fixed seed did not reproduce the values chosen from a list, and Monte Carlo workers forked from one process walked a single shared stream instead of sampling independently. Fixed-seed baselines that vary a list input change. [#953](https://github.com/RocketPy-Team/RocketPy/issues/953) - BUG: Report the atmospheric model time period and ensemble member count for lower-case model types. `set_atmospheric_model` documents `type` as case-insensitive, but `Environment.info()` and `all_info()` compared against capitalised literals, so `type="ensemble"` printed no time period and no member count, and skipped the ensemble comparison plot. [#520](https://github.com/RocketPy-Team/RocketPy/issues/520) - BUG: Give each `CustomSampler` input its own deterministic stream, and seed samplers sharing one generator once as a group. Existing fixed-seed `CustomSampler` baselines change, and samplers built on the legacy `RandomState` must move to `default_rng` because seeds now carry the full 128 bits. [#1102](https://github.com/RocketPy-Team/RocketPy/pull/1102) - BUG: Accept the callable parachute triggers `StochasticParachute` documents, and reject the ones it cannot mean. An invalid string, an empty list or a boolean now fails during validation instead of reaching `Parachute` or becoming a one-metre height trigger. [#1103](https://github.com/RocketPy-Team/RocketPy/pull/1103) +- BUG: Accept a deterministic aerodynamic surface in `StochasticRocket.add_nose`, `add_trapezoidal_fins`, `add_elliptical_fins` and `add_tail`. Each wrapped the surface with a `component=` keyword none of the stochastic classes accept, so passing anything other than an already-stochastic surface raised a `TypeError`. [#953](https://github.com/RocketPy-Team/RocketPy/issues/953) - BUG: rocket with a late-starting thrust curve never leaves the rail [#1085](https://github.com/RocketPy-Team/RocketPy/pull/1085) ## [v1.13.0] - 2026-07-21 diff --git a/docs/reference/classes/monte_carlo/stochastic_models/index.rst b/docs/reference/classes/monte_carlo/stochastic_models/index.rst index ca8b2b1e2..d3c9bb8a1 100644 --- a/docs/reference/classes/monte_carlo/stochastic_models/index.rst +++ b/docs/reference/classes/monte_carlo/stochastic_models/index.rst @@ -19,6 +19,7 @@ input parameters, enabling robust Monte Carlo simulations. stochastic_nose_cone stochastic_trapezoidal_fins stochastic_elliptical_fins + stochastic_free_form_fins stochastic_tail stochastic_rail_buttons stochastic_rocket diff --git a/docs/reference/classes/monte_carlo/stochastic_models/stochastic_free_form_fins.rst b/docs/reference/classes/monte_carlo/stochastic_models/stochastic_free_form_fins.rst new file mode 100644 index 000000000..a1c8391c4 --- /dev/null +++ b/docs/reference/classes/monte_carlo/stochastic_models/stochastic_free_form_fins.rst @@ -0,0 +1,5 @@ +Stochastic Free Form Fins +------------------------- + +.. autoclass:: rocketpy.stochastic.StochasticFreeFormFins + :members: diff --git a/docs/user/stochastic.rst b/docs/user/stochastic.rst index 062b034e2..6e3376236 100644 --- a/docs/user/stochastic.rst +++ b/docs/user/stochastic.rst @@ -93,6 +93,43 @@ passed in a few different ways: gives you the full control of how the samples are generated. See :ref:`custom_sampler` for more details. +.. note:: + The formats above assume each argument holds a single number. The + ``shape_points`` of :class:`rocketpy.stochastic.StochasticFreeFormFins` is + the exception: a fin outline is only meaningful as a complete set of points, + so the deviation given applies to the outline as a block, with every + coordinate of every point perturbed by its own draw. The fin root is held on + the body line, so a point nominally at ``y = 0`` keeps that value and no + point ends up inside the airframe. + + Because the deviation has to centre on the nominal coordinate, only the + distributions that read their first argument as that centre can be used here: + *"normal"*, *"gumbel"*, *"laplace"* and *"logistic"*. The others take bounds + (*"uniform"*) or shape parameters (*"wald"*, *"gamma"*, ...), which a set of + coordinates cannot be, and are rejected when the object is created. + + A list means either one fixed outline, used as given, or a list of candidate + outlines to choose between, which do not have to have the same number of + points:: + + # One millimetre of deviation on every coordinate + StochasticFreeFormFins(free_form_fins=fins, shape_points=0.001) + + # One fixed outline, not randomized + StochasticFreeFormFins( + free_form_fins=fins, + shape_points=[(0, 0), (0.08, 0.1), (0.12, 0)], + ) + + # Choose between two outlines + StochasticFreeFormFins( + free_form_fins=fins, + shape_points=[[(0, 0), (0.08, 0.1), (0.12, 0)], [(0, 0), (0.06, 0.12), (0.12, 0)]], + ) + + A ``CustomSampler`` given for this argument has to yield a whole outline per + sample, since what it returns replaces the outline instead of perturbing it. + .. note:: In statistics, the terms "Normal" and "Gaussian" refer to the same type of \ distribution. This distribution is commonly used and is the default for the \ diff --git a/rocketpy/__init__.py b/rocketpy/__init__.py index d8720db4c..6008ff09b 100644 --- a/rocketpy/__init__.py +++ b/rocketpy/__init__.py @@ -62,6 +62,7 @@ StochasticEllipticalFins, StochasticEnvironment, StochasticFlight, + StochasticFreeFormFins, StochasticNoseCone, StochasticParachute, StochasticRocket, diff --git a/rocketpy/stochastic/__init__.py b/rocketpy/stochastic/__init__.py index ffadfaaaf..0045baefd 100644 --- a/rocketpy/stochastic/__init__.py +++ b/rocketpy/stochastic/__init__.py @@ -9,6 +9,7 @@ from .stochastic_aero_surfaces import ( StochasticAirBrakes, StochasticEllipticalFins, + StochasticFreeFormFins, StochasticNoseCone, StochasticRailButtons, StochasticTail, diff --git a/rocketpy/stochastic/stochastic_aero_surfaces.py b/rocketpy/stochastic/stochastic_aero_surfaces.py index 27d3d89a9..137a00385 100644 --- a/rocketpy/stochastic/stochastic_aero_surfaces.py +++ b/rocketpy/stochastic/stochastic_aero_surfaces.py @@ -1,17 +1,22 @@ """ Defines the StochasticNoseCone, StochasticTrapezoidalFins, -StochasticEllipticalFins, StochasticTail and StochasticRailButtons classes. +StochasticEllipticalFins, StochasticFreeFormFins, StochasticTail and +StochasticRailButtons classes. """ +import numpy as np + from rocketpy.rocket.aero_surface import ( AirBrakes, EllipticalFins, + FreeFormFins, NoseCone, RailButtons, Tail, TrapezoidalFins, ) +from .custom_sampler import CustomSampler from .stochastic_model import StochasticModel @@ -305,6 +310,362 @@ def create_object(self): return EllipticalFins(**generated_dict) +class StochasticFreeFormFins(StochasticModel): + """A Stochastic Free Form Fins class that inherits from StochasticModel. + + See Also + -------- + :ref:`stochastic_model` and + :class:`FreeFormFins ` + + Attributes + ---------- + object : FreeFormFins + FreeFormFins object to be used for validation. + n : list[int] + List with an integer representing the number of fins. This attribute + can be randomized. + shape_points : tuple, list, numpy.ndarray, int, float + The (x, y) points defining the fin outline, in meters. Unlike the other + fin sets, this geometry is a whole list of points rather than a single + scalar, so the deviation given applies to the outline as a block: every + coordinate of every point is perturbed, each by its own draw. See the + ``shape_points`` parameter of :meth:`__init__` for the accepted formats. + rocket_radius : tuple, list, int, float + Rocket radius of the fins in meters. + cant_angle : tuple, list, int, float + Cant angle of the fins in degrees. + airfoil : list + List of tuples in the form of (airfoil file path, airfoil name). + name : list[str] + List with the fins object name. This attribute can not be randomized. + """ + + # The outline is the whole nominal value of this input, not a single number. + array_valued_inputs = ("shape_points",) + + # The distributions that can mean a deviation around a nominal coordinate, + # which is what perturbing an outline asks of them. The rest of what + # ``get_distribution`` offers reads its arguments as bounds (``uniform``) or + # as shape parameters (``wald``, ``gamma``, ``poisson``, ...), neither of + # which an outline of coordinates can be. + _outline_distributions = ("normal", "gumbel", "laplace", "logistic") + + def __init__( + self, + free_form_fins=None, + n=None, + shape_points=None, + rocket_radius=None, + cant_angle=None, + airfoil=None, + ): + """Initializes the Stochastic Free Form Fins class. + + See Also + -------- + :ref:`stochastic_model` + + Parameters + ---------- + free_form_fins : FreeFormFins + FreeFormFins object to be used for validation. + shape_points : tuple, list, numpy.ndarray, int, float, optional + The (x, y) points defining the fin outline, in meters. The whole + outline is perturbed as a block, since a fin shape is only + meaningful as a complete set of points: the deviation given applies + to every coordinate of every point, each drawn independently of the + others. The fin root is held on the body line, so a point nominally + at ``y = 0`` keeps that value and no point is moved inside the + airframe. The accepted formats are: + + - ``int`` or ``float``: standard deviation applied to every + coordinate of the nominal outline, drawn from a normal + distribution. + - ``tuple``: ``(standard deviation, distribution name)``, or + ``(nominal outline, standard deviation[, distribution name])``. + The distribution must be one of ``"normal"``, ``"gumbel"``, + ``"laplace"`` or ``"logistic"``, the ones that take the nominal + coordinate as their centre. + - ``list`` or ``numpy.ndarray``: either one fixed outline, e.g. + ``[(0, 0), (0.1, 0.1), (0.1, 0)]``, which is used as given and + not randomized; or a list of candidate outlines, e.g. + ``[[(0, 0), (0.1, 0.1), (0.1, 0)], [(0, 0), (0.1, 0.12), (0.1, 0)]]``, + one of which is chosen per simulation. The candidates need not + all have the same number of points. An empty list means the + nominal outline of the object passed, unrandomized, as it does + for every other argument. + - ``CustomSampler``: has to yield a whole outline per sample, since + the value it returns replaces the outline instead of perturbing + it. + rocket_radius : tuple, list, int, float, optional + Rocket radius of the fins in meters. + cant_angle : tuple, list, int, float, optional + Cant angle of the fins in degrees. + airfoil : list[tuple], optional + List of tuples in the form of (airfoil file path, airfoil name). + """ + # TODO: never vary the number of fins. It is a fixed parameter. + self._validate_positive_int_list("n", n) + self._validate_airfoil(airfoil) + shape_points = self._validate_shape_points(shape_points) + super().__init__( + free_form_fins, + n=n, + shape_points=shape_points, + rocket_radius=rocket_radius, + cant_angle=cant_angle, + airfoil=airfoil, + name=None, + ) + + def _validate_shape_points(self, shape_points): + """Validate the ``shape_points`` input and normalize it to a form the + base class can randomize. + + A fin outline is a sequence of points, so it does not fit the + scalar-per-input assumption the base class makes. Two formats would be + silently misread if passed straight through, and both are the natural + thing for a user to write: + + - a bare outline ``[(0, 0), (0.1, 0.1), (0.1, 0)]`` is a ``list``, which + the base class reads as a list of candidate values and would sample a + single ``(x, y)`` point from. It is wrapped here so it is treated as + the one candidate outline it is. + - a ``(nominal outline, standard deviation)`` tuple carries an outline + where the base class reads a distribution argument, so the outline is + checked here and the rest is left to the base class. + + Parameters + ---------- + shape_points : tuple, list, numpy.ndarray, int, float, optional + Value of the ``shape_points`` input argument. + + Returns + ------- + tuple, list, int, float or None + The input, normalized so the base class randomizes the outline as a + block. Outlines come back as ``(n, 2)`` arrays of floats. + + Raises + ------ + AssertionError + If the input is not in a valid format. + """ + if shape_points is None or isinstance( + shape_points, (int, float, CustomSampler) + ): + # A number is a standard deviation around the nominal outline, which + # the base class looks up and hands to the distribution as an array. + # A sampler yields whole outlines, so it replaces that machinery. + return shape_points + + if isinstance(shape_points, tuple): + return self._validate_shape_points_tuple(shape_points) + + if isinstance(shape_points, (list, np.ndarray)): + if len(shape_points) == 0: + # An empty list means the nominal value everywhere else, and + # nothing about this argument makes it mean something else. + return [] + if self._is_outline(shape_points): + # A bare outline is the one candidate it describes. Left as a + # list of points it would be read as a list of candidates and + # sampled down to a single (x, y) point. + return [self._validate_outline(shape_points)] + return [self._validate_outline(outline) for outline in shape_points] + + raise AssertionError( + "`shape_points` must be a tuple, list, numpy array, int, or float " + "or a custom sampler" + ) + + def _validate_shape_points_tuple(self, shape_points): + """Validate a ``shape_points`` tuple. + + Accepts ``(standard deviation, distribution name)``, in which case the + nominal outline comes from the object passed, and + ``(nominal outline, standard deviation[, distribution name])``. + + Parameters + ---------- + shape_points : tuple + Value of the ``shape_points`` input argument. + + Returns + ------- + tuple + The input tuple, with any nominal outline converted to an ``(n, 2)`` + array of floats so the standard deviation broadcasts over every + coordinate. + + Raises + ------ + AssertionError + If the input is not in a valid format. + """ + if len(shape_points) not in (2, 3): + raise AssertionError("'shape_points': tuple must have length 2 or 3") + + if isinstance(shape_points[0], (int, float)): + # (standard deviation, distribution name), the nominal outline being + # taken from the object passed. A number in the second item would + # make the first one the nominal value, which for this argument is + # an outline rather than a number. + if not isinstance(shape_points[1], str): + raise AssertionError( + "'shape_points': when the first item of a tuple is a " + "standard deviation, the second must be a string naming a " + "valid numpy.random distribution function." + ) + self._validate_outline_distribution(shape_points[1]) + return shape_points + + # (nominal outline, standard deviation[, distribution name]). The second + # item is checked here rather than left to the base class, which also + # accepts a string there and would read the outline as the deviation. + outline = self._validate_outline(shape_points[0]) + if not isinstance(shape_points[1], (int, float)): + raise AssertionError( + "'shape_points': second item of tuple must be an int or float " + "standard deviation." + ) + if len(shape_points) == 3: + if not isinstance(shape_points[2], str): + raise AssertionError( + "'shape_points': Third item of tuple must be a string containing " + "the name of a valid numpy.random distribution function." + ) + self._validate_outline_distribution(shape_points[2]) + return (outline,) + tuple(shape_points[1:]) + + @classmethod + def _validate_outline_distribution(cls, distribution_name): + """Reject distributions that cannot mean a deviation around a coordinate. + + The distribution is called as ``dist_func(nominal_outline, std_dev)``, so + only the ones that read the first argument as the centre of the draw can + perturb an outline. ``uniform`` would read the outline as its lower bound + and the deviation as a single upper bound, leaving an empty range for + every coordinate above it, and ``wald`` and the shape-parameter + distributions reject the zeros that a root point has. + + Raises + ------ + AssertionError + If the distribution cannot be applied to an outline. + """ + if distribution_name not in cls._outline_distributions: + accepted = ", ".join(repr(name) for name in cls._outline_distributions) + raise AssertionError( + f"'shape_points': the '{distribution_name}' distribution cannot " + f"be applied to an outline. Use one of {accepted}, which take " + "the nominal coordinate as the centre of the deviation." + ) + + @staticmethod + def _is_outline(value): + """Return True if ``value`` is a single (x, y) outline. + + Used to tell a bare outline apart from a list of candidate outlines, + which are the two things a list input can mean. The conversion is what + decides it: numpy refuses a ragged or non-numeric sequence, and a list + of candidates whose outlines have different numbers of points is exactly + that, so those are left for the caller to check one at a time. + """ + try: + array = np.asarray(value, dtype=float) + except (ValueError, TypeError): + return False + return array.ndim == 2 and array.shape[1] == 2 + + @staticmethod + def _validate_outline(outline): + """Validate a single (x, y) fin outline. + + Returns + ------- + numpy.ndarray + The outline as an ``(n, 2)`` array of floats. + + Raises + ------ + AssertionError + If the outline is not a sequence of at least three (x, y) numbers. + """ + if not StochasticFreeFormFins._is_outline(outline): + raise AssertionError( + "`shape_points` outlines must be sequences of (x, y) numbers, " + "i.e. have shape (n, 2)." + ) + array = np.asarray(outline, dtype=float) + if array.shape[0] < 3: + raise AssertionError( + "`shape_points` outlines must have at least 3 points to " + "enclose an area." + ) + return array + + # pylint: disable=stop-iteration-return + def dict_generator(self): + """Generate the input arguments, with the fin root kept on the body line. + + Yields + ------ + dict + Dictionary with the randomly generated input arguments. + """ + generated_dict = next(super().dict_generator()) + if isinstance(self.shape_points, tuple): + # Only a perturbed outline can have drifted off the body line. One + # chosen from a list of candidates, or one a sampler produced, is + # used exactly as it was given. + generated_dict["shape_points"] = self._keep_root_on_body_line( + self.shape_points[0], generated_dict["shape_points"] + ) + yield generated_dict + + @staticmethod + def _keep_root_on_body_line(nominal_outline, sampled_outline): + """Hold the root of a perturbed outline on the body line. + + :class:`FreeFormFins ` measures the span from + ``y = 0`` and slices the chords over that interval, so a root point that + drifts off the line puts part of the fin inside the airframe and inflates + the span those chords are measured against. Points nominally on the line + are kept there, and no other point is allowed to cross it. + + Parameters + ---------- + nominal_outline : numpy.ndarray + The unperturbed outline, which says which points are on the line. + sampled_outline : numpy.ndarray + The perturbed outline. + + Returns + ------- + numpy.ndarray + The perturbed outline, with its root back on the body line. + """ + nominal_outline = np.asarray(nominal_outline, dtype=float) + sampled_outline = np.array(sampled_outline, dtype=float) + sampled_outline[nominal_outline[:, 1] == 0, 1] = 0.0 + sampled_outline[:, 1] = np.maximum(sampled_outline[:, 1], 0.0) + return sampled_outline + + def create_object(self): + """Creates and returns a FreeFormFins object from the randomly + generated input arguments. + + Returns + ------- + fins : FreeFormFins + FreeFormFins object with the randomly generated input arguments. + """ + generated_dict = next(self.dict_generator()) + return FreeFormFins(**generated_dict) + + class StochasticTail(StochasticModel): """A Stochastic Tail class that inherits from StochasticModel. diff --git a/rocketpy/stochastic/stochastic_model.py b/rocketpy/stochastic/stochastic_model.py index be2438a0c..79829beff 100644 --- a/rocketpy/stochastic/stochastic_model.py +++ b/rocketpy/stochastic/stochastic_model.py @@ -3,8 +3,6 @@ Stochastic classes. """ -from random import choice - import numpy as np from rocketpy.mathutils.function import Function @@ -31,6 +29,18 @@ def _names_as_spawn_key(input_names): ) +def _format_number(value): + """Format a nominal value or a standard deviation for the attribute report. + + An array-valued input, such as a fin outline, has no single number to show, + and a fixed-width format raises a ``TypeError`` on it, so its shape stands + in for the coordinates. + """ + if np.ndim(value) == 0: + return f"{value:.5f}" + return f"array of shape {np.shape(value)}" + + def _sampler_seed(seed, input_names): """Derive a seed for one sampler, or for one group that shares a generator. @@ -79,6 +89,12 @@ class StochasticModel: "ensemble_member", ] + # Arguments whose nominal value is an array of numbers rather than a single + # number, such as the outline of a free-form fin. Declared by name so that + # validation, sampling and the attribute report all agree on which ones they + # are, instead of each deciding for itself. + array_valued_inputs = () + def __init__(self, obj, seed=None, **kwargs): """ Initialize the StochasticModel class with validated input arguments. @@ -118,6 +134,13 @@ def _set_stochastic(self, seed=None): Seed for the random number generator. """ self.__random_number_generator = np.random.default_rng(seed) + # A stream of its own, derived from the same seed, for picking between + # the candidate values of a list input. Kept apart from the one above so + # that declaring a list input does not shift the numbers every other + # input draws, which would move each existing fixed-seed baseline. + self.__choice_generator = np.random.default_rng( + _sampler_seed(seed, ("__list_choice__",)) + ) self.last_rnd_dict = {} self._reset_custom_samplers(seed) @@ -153,6 +176,53 @@ def _set_stochastic(self, seed=None): def __repr__(self): return f"'{self.__class__.__name__}() object'" + def _choose(self, values): + """Pick one of the candidate values of a list input. + + ``random.choice`` was used here, which draws from the interpreter-wide + stream that ``_set_stochastic`` does not reseed: the same seed did not + reproduce the same choices, and Monte Carlo workers forked from one + process inherited a single stream and walked it together instead of + sampling independently. + + Parameters + ---------- + values : list + Candidate values of the input. + + Returns + ------- + object + One of the candidates, or ``values`` itself when there are none. + """ + if len(values) == 0: + return values + return values[self.__choice_generator.integers(len(values))] + + def _nominal_value(self, input_name, value): + """Return the nominal value of an input as the distribution needs it. + + The distributions are called as ``dist_func(nominal, std_dev)``, so an + array-valued input has to arrive as an array for the deviation to + broadcast over its entries. A list of ``(x, y)`` tuples, which is how a + fin outline is written, would not. + + Parameters + ---------- + input_name : str + Name of the input argument. + value : object + Nominal value of the input argument. + + Returns + ------- + object + The value, as an array of floats for the array-valued inputs. + """ + if input_name in self.array_valued_inputs: + return np.asarray(value, dtype=float) + return value + def _validate_tuple(self, input_name, input_value, getattr=getattr): # pylint: disable=redefined-builtin """ Validate tuple arguments. @@ -183,8 +253,15 @@ def _validate_tuple(self, input_name, input_value, getattr=getattr): # pylint: ]: raise AssertionError(f"'{input_name}': tuple must have length 2 or 3") if not isinstance(input_value[0], (int, float)): - raise AssertionError( - f"'{input_name}': First item of tuple must be an int or float" + if input_name not in self.array_valued_inputs: + raise AssertionError( + f"'{input_name}': First item of tuple must be an int or float" + ) + # An array-valued input carries its whole nominal value here, so the + # single number the others require is not what to expect. The child + # class that declared it has already checked the value itself. + input_value = (self._nominal_value(input_name, input_value[0]),) + tuple( + input_value[1:] ) if len(input_value) == 2: @@ -227,7 +304,11 @@ def _validate_tuple_length_two(self, input_name, input_value, getattr=getattr): # function. In this case, the nominal value will be taken from the # object passed. dist_func = get_distribution(input_value[1], self.__random_number_generator) - return (getattr(self.obj, input_name), input_value[0], dist_func) + return ( + self._nominal_value(input_name, getattr(self.obj, input_name)), + input_value[0], + dist_func, + ) else: # if second item is an int or float, then it is assumed that the # first item is the nominal value and the second item is the @@ -326,7 +407,7 @@ def _validate_scalar(self, input_name, input_value, getattr=getattr): # pylint: distribution function). """ return ( - getattr(self.obj, input_name), + self._nominal_value(input_name, getattr(self.obj, input_name)), input_value, get_distribution("normal", self.__random_number_generator), ) @@ -632,7 +713,7 @@ def dict_generator(self): dist_sampler = value[-1] generated_dict[arg] = dist_sampler(value[0], value[1]) elif isinstance(value, list): - generated_dict[arg] = choice(value) if value else value + generated_dict[arg] = self._choose(value) elif isinstance(value, CustomSampler): try: generated_dict[arg] = value.sample(n_samples=1)[0] @@ -671,13 +752,14 @@ def format_attribute(attr, value): upper_bound = std_dev return ( f"\t{attr.ljust(max_str_length)} " - f"{lower_bound:.5f}, {upper_bound:.5f} ({dist_func.__name__})" + f"{_format_number(lower_bound)}, " + f"{_format_number(upper_bound)} ({dist_func.__name__})" ) else: return ( f"\t{attr.ljust(max_str_length)} " - f"{nominal_value:.5f} ± " - f"{std_dev:.5f} ({dist_func.__name__})" + f"{_format_number(nominal_value)} ± " + f"{_format_number(std_dev)} ({dist_func.__name__})" ) elif isinstance(value, CustomSampler): sampler_name = type(value).__name__ diff --git a/rocketpy/stochastic/stochastic_rocket.py b/rocketpy/stochastic/stochastic_rocket.py index 33a364f18..895e9a2a4 100644 --- a/rocketpy/stochastic/stochastic_rocket.py +++ b/rocketpy/stochastic/stochastic_rocket.py @@ -1,7 +1,6 @@ """Defines the StochasticRocket class.""" import warnings -from random import choice from rocketpy.control import _Controller from rocketpy.mathutils.vector_matrix import Vector @@ -11,6 +10,7 @@ from rocketpy.rocket.aero_surface import ( AirBrakes, EllipticalFins, + FreeFormFins, NoseCone, RailButtons, Tail, @@ -25,6 +25,7 @@ from .stochastic_aero_surfaces import ( StochasticAirBrakes, StochasticEllipticalFins, + StochasticFreeFormFins, StochasticNoseCone, StochasticRailButtons, StochasticTail, @@ -273,7 +274,10 @@ def _add_surfaces(self, surfaces, positions, type_, stochastic_type, error_messa if not isinstance(surfaces, (type_, stochastic_type)): raise AssertionError(error_message) if isinstance(surfaces, type_): - surfaces = stochastic_type(component=surfaces) + # Positionally: the stochastic classes each name this first + # parameter after their own surface (`nosecone`, `tail`, ...), so + # there is no one keyword that reaches all of them. + surfaces = stochastic_type(surfaces) self.__components_map[surfaces] = positions self.aerodynamic_surfaces.add( surfaces, self._validate_position(surfaces, positions) @@ -333,6 +337,24 @@ def add_elliptical_fins(self, fins, position=None): "`fins` must be of EllipticalFins or StochasticEllipticalFins type", ) + def add_free_form_fins(self, fins, position=None): + """Adds a stochastic free form fins to the stochastic rocket. + + Parameters + ---------- + fins : StochasticFreeFormFins or FreeFormFins + The free form fins to be added to the stochastic rocket. + position : tuple, list, int, float, optional + The position of the free form fins. + """ + self._add_surfaces( + fins, + position, + FreeFormFins, + StochasticFreeFormFins, + "`fins` must be of FreeFormFins or StochasticFreeFormFins type", + ) + def add_tail(self, tail, position=None): """Adds a stochastic tail to the stochastic rocket. @@ -628,7 +650,7 @@ def _randomize_position(self, position): return position[-1](position[0].z, position[1]) return position[-1](position[0], position[1]) elif isinstance(position, list): - return choice(position) if position else position + return self._choose(position) # pylint: disable=stop-iteration-return def dict_generator(self): @@ -638,8 +660,8 @@ def dict_generator(self): all attributes of the class and generating a random value for each attribute. The random values are generated according to the format of each attribute. Tuples are generated using the distribution function - specified in the tuple. Lists are generated using the random.choice - function. + specified in the tuple. Lists are generated by picking one of their + values with this model's own seeded generator. Parameters ---------- diff --git a/tests/fixtures/monte_carlo/stochastic_fixtures.py b/tests/fixtures/monte_carlo/stochastic_fixtures.py index 6610666cf..45c4538e1 100644 --- a/tests/fixtures/monte_carlo/stochastic_fixtures.py +++ b/tests/fixtures/monte_carlo/stochastic_fixtures.py @@ -7,6 +7,7 @@ from rocketpy.stochastic import ( StochasticEnvironment, StochasticFlight, + StochasticFreeFormFins, StochasticNoseCone, StochasticParachute, StochasticRailButtons, @@ -117,6 +118,28 @@ def stochastic_trapezoidal_fins(calisto_trapezoidal_fins): ) +@pytest.fixture +def stochastic_free_form_fins(calisto_free_form_fins): + """This fixture is used to create a StochasticFreeFormFins object for the + Calisto rocket. + + Parameters + ---------- + calisto_free_form_fins : FreeFormFins + This is another fixture. + + Returns + ------- + StochasticFreeFormFins + The stochastic free form fins object + """ + return StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, + shape_points=0.0005, + cant_angle=(0, 0.5), + ) + + @pytest.fixture def stochastic_tail(calisto_tail): """This fixture is used to create a StochasticTail object for the diff --git a/tests/integration/simulation/test_monte_carlo.py b/tests/integration/simulation/test_monte_carlo.py index bcfb59505..505c30f40 100644 --- a/tests/integration/simulation/test_monte_carlo.py +++ b/tests/integration/simulation/test_monte_carlo.py @@ -1,4 +1,5 @@ # pylint: disable=unused-argument +import json import os from unittest.mock import patch @@ -6,6 +7,9 @@ import numpy as np import pytest +from rocketpy.rocket.components import Components +from rocketpy.simulation import MonteCarlo + plt.rcParams.update({"figure.max_open_warning": 0}) @@ -263,3 +267,44 @@ def test_monte_carlo_simulate_convergence(monte_carlo_calisto): assert monte_carlo_calisto.num_of_loaded_sims <= 20 finally: _post_test_file_cleanup() + + +@pytest.mark.slow +def test_monte_carlo_simulate_free_form_fins( + stochastic_environment, + stochastic_calisto, + stochastic_free_form_fins, + stochastic_flight, + tmp_path, +): + """A free-form fin set must survive a whole Monte Carlo run: it has to be + sampled, flown, and written to the inputs file as an outline rather than as + a single point (see #953).""" + + stochastic_calisto.aerodynamic_surfaces = Components() + stochastic_calisto.add_free_form_fins( + stochastic_free_form_fins, position=(-1.04956, 0.001) + ) + + filename = str(tmp_path / "monte_carlo_free_form_fins") + monte_carlo = MonteCarlo( + filename=filename, + environment=stochastic_environment, + rocket=stochastic_calisto, + flight=stochastic_flight, + ) + monte_carlo.simulate(number_of_simulations=2, append=False) + + assert monte_carlo.num_of_loaded_sims == 2 + + nominal = np.asarray(stochastic_free_form_fins.obj.shape_points, dtype=float) + with open(filename + ".inputs.txt", encoding="utf-8") as file: + lines = file.read().splitlines() + assert len(lines) == 2 + for line in lines: + surfaces = json.loads(line)["aerodynamic_surfaces"] + outlines = [s["shape_points"] for s in surfaces if "shape_points" in s] + assert len(outlines) == 1 + sampled = np.asarray(outlines[0], dtype=float) + assert sampled.shape == nominal.shape + assert not np.allclose(sampled, nominal) diff --git a/tests/unit/stochastic/test_stochastic_aero_surfaces.py b/tests/unit/stochastic/test_stochastic_aero_surfaces.py index d63feb76c..ab979b699 100644 --- a/tests/unit/stochastic/test_stochastic_aero_surfaces.py +++ b/tests/unit/stochastic/test_stochastic_aero_surfaces.py @@ -1,4 +1,14 @@ -from rocketpy.rocket.aero_surface import NoseCone, RailButtons, Tail, TrapezoidalFins +import numpy as np +import pytest + +from rocketpy.rocket.aero_surface import ( + FreeFormFins, + NoseCone, + RailButtons, + Tail, + TrapezoidalFins, +) +from rocketpy.stochastic import StochasticFreeFormFins ## NOSE CONE @@ -46,6 +56,233 @@ class creates a StochasticTrapezoidalFins object from the randomly generated assert isinstance(obj, TrapezoidalFins) +## FREE FORM FINS + +NOMINAL_SHAPE = [(0, 0), (0.08, 0.1), (0.12, 0.1), (0.12, 0)] + + +def test_stochastic_free_form_fins_create_object(stochastic_free_form_fins): + """Test create object method of StochasticFreeFormFins class. + + This test checks if the create_object method of the StochasticFreeFormFins + class creates a FreeFormFins object from the randomly generated input + arguments. + + Parameters + ---------- + stochastic_free_form_fins : StochasticFreeFormFins + StochasticFreeFormFins object to be tested. + + Returns + ------- + None + """ + obj = stochastic_free_form_fins.create_object() + assert isinstance(obj, FreeFormFins) + + +def test_stochastic_free_form_fins_nominal_shape_is_preserved(calisto_free_form_fins): + """With nothing to randomize, the created fin set must keep the outline of + the object it was built from.""" + stochastic = StochasticFreeFormFins(free_form_fins=calisto_free_form_fins) + + created = stochastic.create_object() + + assert np.allclose( + np.asarray(created.shape_points, dtype=float), + np.asarray(calisto_free_form_fins.shape_points, dtype=float), + ) + + +@pytest.mark.parametrize( + "shape_points", + [ + 0.001, + (0.001, "normal"), + (NOMINAL_SHAPE, 0.001), + (NOMINAL_SHAPE, 0.001, "normal"), + ], + ids=["scalar", "std_and_dist", "outline_and_std", "outline_std_and_dist"], +) +def test_stochastic_free_form_fins_perturbs_the_whole_outline( + calisto_free_form_fins, shape_points +): + """A fin outline is only meaningful as a complete set of points, so every + accepted format must randomize all of them and keep the (n, 2) shape.""" + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, shape_points=shape_points + ) + stochastic._set_stochastic(42) + + created = stochastic.create_object() + + nominal = np.asarray(NOMINAL_SHAPE, dtype=float) + sampled = np.asarray(created.shape_points, dtype=float) + assert sampled.shape == nominal.shape + # Every coordinate is drawn on its own, so none of the four points is left + # exactly where it was, apart from the root's y (see the test below). + assert not np.allclose(sampled[:, 0], nominal[:, 0]) + assert not np.allclose(sampled[1:3, 1], nominal[1:3, 1]) + # A standard deviation of a millimetre must not turn into a new fin. + assert np.abs(sampled - nominal).max() < 0.01 + + +@pytest.mark.parametrize( + "shape_points", + [0.001, (0.001, "normal"), (NOMINAL_SHAPE, 0.001, "laplace")], + ids=["scalar", "std_and_dist", "outline_std_and_dist"], +) +def test_stochastic_free_form_fins_keeps_the_root_on_the_body_line( + calisto_free_form_fins, shape_points +): + """FreeFormFins measures the span from y = 0 and slices the chords over that + interval, so a perturbed root point must not drift off the body line: it + would put part of the fin inside the airframe and inflate the span the + chords are measured against. + """ + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, shape_points=shape_points + ) + stochastic._set_stochastic(3) + + for _ in range(50): + sampled = np.asarray(stochastic.create_object().shape_points, dtype=float) + # The first and last points of the nominal outline are on the body line. + assert sampled[0, 1] == 0 + assert sampled[-1, 1] == 0 + assert (sampled[:, 1] >= 0).all() + + +def test_stochastic_free_form_fins_bare_outline_is_a_single_candidate( + calisto_free_form_fins, +): + """A bare outline is a list, which the base class would otherwise read as a + list of candidate values and sample a single (x, y) point from.""" + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, shape_points=NOMINAL_SHAPE + ) + + created = stochastic.create_object() + + assert np.allclose( + np.asarray(created.shape_points, dtype=float), + np.asarray(NOMINAL_SHAPE, dtype=float), + ) + + +def test_stochastic_free_form_fins_chooses_between_outlines(calisto_free_form_fins): + """A list of outlines is a set of candidate shapes to choose from.""" + taller = [(0, 0), (0.06, 0.12), (0.12, 0.12), (0.12, 0)] + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, + shape_points=[NOMINAL_SHAPE, taller], + ) + stochastic._set_stochastic(42) + + spans = {round(stochastic.create_object().span, 4) for _ in range(50)} + + assert spans == {0.1, 0.12} + + +def test_stochastic_free_form_fins_chooses_between_outlines_of_different_lengths( + calisto_free_form_fins, +): + """Candidate outlines need not have the same number of points: choosing + between a three-point and a four-point fin is the plainest form of choosing + between shapes, and numpy raises on that ragged list if it is converted + whole instead of one candidate at a time. + """ + triangle = [(0, 0), (0.08, 0.1), (0.12, 0)] + quadrilateral = [(0, 0), (0.06, 0.12), (0.12, 0.12), (0.12, 0)] + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, + shape_points=[triangle, quadrilateral], + ) + stochastic._set_stochastic(42) + + point_counts = {len(stochastic.create_object().shape_points) for _ in range(50)} + + assert point_counts == {3, 4} + + +def test_stochastic_free_form_fins_accepts_an_array_outline(calisto_free_form_fins): + """A sampled outline comes back as an array, so feeding one straight back in + as the nominal outline must work.""" + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, + shape_points=np.asarray(NOMINAL_SHAPE, dtype=float), + ) + + created = stochastic.create_object() + + assert np.allclose( + np.asarray(created.shape_points, dtype=float), + np.asarray(NOMINAL_SHAPE, dtype=float), + ) + + +def test_stochastic_free_form_fins_empty_list_means_the_nominal_outline( + calisto_free_form_fins, +): + """An empty list means "take the nominal value and do not randomize" for + every other stochastic input, and this one is no different.""" + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, shape_points=[] + ) + + created = stochastic.create_object() + + assert np.allclose( + np.asarray(created.shape_points, dtype=float), + np.asarray(calisto_free_form_fins.shape_points, dtype=float), + ) + + +@pytest.mark.parametrize( + "shape_points", + [ + "not_an_outline", + [[(0, 0), (0.1, 0.1)]], + [(0, 0), (0.1, 0.1)], + [(0, 0, 0), (0.1, 0.1, 0), (0.1, 0, 0)], + [(0, 0), (1, 1, 1), (2, 0)], + [[("a", "b"), ("c", "d"), ("e", "f")]], + (0.001,), + (NOMINAL_SHAPE, 0.001, "normal", 1), + (NOMINAL_SHAPE, "normal"), + (0.001, 5), + (NOMINAL_SHAPE, 0.001, 7), + (0.001, "uniform"), + (NOMINAL_SHAPE, 0.001, "wald"), + ], + ids=[ + "string", + "too_few_points", + "bare_outline_too_few_points", + "three_dimensional_points", + "ragged_outline", + "non_numeric_points", + "tuple_too_short", + "tuple_too_long", + "outline_with_string_std", + "std_with_non_string_dist", + "outline_with_non_string_dist", + "bounded_distribution", + "shape_parameter_distribution", + ], +) +def test_stochastic_free_form_fins_rejects_invalid_shape_points( + calisto_free_form_fins, shape_points +): + """An outline that cannot mean a fin shape, or a distribution that cannot + mean a deviation around one, must fail during validation rather than + reaching FreeFormFins or the sampler.""" + with pytest.raises(AssertionError): + StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, shape_points=shape_points + ) + + ## TAIL diff --git a/tests/unit/stochastic/test_stochastic_model.py b/tests/unit/stochastic/test_stochastic_model.py index 9e35a5330..8bb360c48 100644 --- a/tests/unit/stochastic/test_stochastic_model.py +++ b/tests/unit/stochastic/test_stochastic_model.py @@ -1,5 +1,7 @@ import pytest +from rocketpy.stochastic import StochasticFreeFormFins + @pytest.mark.parametrize( "fixture_name", @@ -10,6 +12,7 @@ "stochastic_environment_custom_sampler", "stochastic_tail", "stochastic_calisto", + "stochastic_free_form_fins", ], ) def test_visualize_attributes(request, fixture_name): @@ -21,3 +24,29 @@ def test_visualize_attributes(request, fixture_name): report = fixture.visualize_attributes() assert isinstance(report, str) assert report + + +def test_list_choices_are_reproducible(calisto_free_form_fins): + """Choosing between the candidate values of a list input must come from the + model's own generator, so that the same seed replays the same choices. + + The interpreter-wide ``random.choice`` was used, which ``_set_stochastic`` + does not reseed: a fixed-seed run picked different values every time, and + Monte Carlo workers forked from one process walked a single shared stream + instead of sampling independently. + """ + taller = [(0, 0), (0.06, 0.12), (0.12, 0.12), (0.12, 0)] + stochastic = StochasticFreeFormFins( + free_form_fins=calisto_free_form_fins, + shape_points=[calisto_free_form_fins.shape_points, taller], + ) + + def spans(seed): + stochastic._set_stochastic(seed) + return [round(stochastic.create_object().span, 4) for _ in range(20)] + + assert spans(7) == spans(7) + assert spans(7) != spans(8) + # Both candidates must stay reachable, or the assertions above would also + # hold for a generator that always returned the same one. + assert set(spans(7)) == {0.1, 0.12} diff --git a/tests/unit/stochastic/test_stochastic_rocket.py b/tests/unit/stochastic/test_stochastic_rocket.py index c96122f04..d15eb0bb6 100644 --- a/tests/unit/stochastic/test_stochastic_rocket.py +++ b/tests/unit/stochastic/test_stochastic_rocket.py @@ -1,6 +1,15 @@ +import numpy as np +import pytest + +from rocketpy.rocket.aero_surface import FreeFormFins from rocketpy.rocket.parachute import Parachute from rocketpy.rocket.rocket import Rocket -from rocketpy.stochastic import StochasticParachute, StochasticRocket +from rocketpy.stochastic import ( + StochasticFreeFormFins, + StochasticParachute, + StochasticRocket, + StochasticTrapezoidalFins, +) def test_str(stochastic_calisto): @@ -123,3 +132,62 @@ def test_configured_geometry_survives_without_being_randomized(calisto_robust): flown = stochastic.create_object().parachutes[0] assert (flown.radius, flown.height, flown.porosity) == (2.0, 1.5, 0.05) + + +def test_a_deterministic_surface_is_wrapped_in_its_stochastic_model( + calisto_robust, calisto_trapezoidal_fins +): + """`_add_surfaces` used to wrap deterministic surfaces with a `component=` + keyword none of the stochastic classes accept, so passing any plain + aerodynamic surface raised a TypeError instead of being wrapped.""" + stochastic = StochasticRocket(rocket=calisto_robust) + + stochastic.add_trapezoidal_fins(calisto_trapezoidal_fins) + + added = stochastic.aerodynamic_surfaces.get_tuple_by_type(StochasticTrapezoidalFins) + assert len(added) == 1 + assert added[0].component.obj is calisto_trapezoidal_fins + + +def test_add_free_form_fins_reaches_the_created_rocket( + calisto_robust, stochastic_free_form_fins +): + """The fin set added to the stochastic rocket must be the one the created + rocket flies, with the outline randomized as a block.""" + stochastic = StochasticRocket(rocket=calisto_robust) + stochastic.add_free_form_fins(stochastic_free_form_fins, position=(-1.04956, 0.001)) + stochastic._set_stochastic(42) + + rocket = stochastic.create_object() + + fin_sets = rocket.aerodynamic_surfaces.get_tuple_by_type(FreeFormFins) + assert len(fin_sets) == 1 + flown = fin_sets[0].component + nominal = np.asarray(stochastic_free_form_fins.obj.shape_points, dtype=float) + sampled = np.asarray(flown.shape_points, dtype=float) + assert sampled.shape == nominal.shape + assert not np.allclose(sampled, nominal) + + +def test_add_free_form_fins_rejects_other_surfaces(calisto_robust, calisto_tail): + stochastic = StochasticRocket(rocket=calisto_robust) + + with pytest.raises(AssertionError): + stochastic.add_free_form_fins(calisto_tail) + + +def test_add_free_form_fins_wraps_a_deterministic_fin_set(calisto_robust): + """A plain FreeFormFins must be wrapped in its own stochastic model, the + same way the other surfaces are.""" + fins = calisto_robust.add_free_form_fins( + n=4, + shape_points=[(0, 0), (0.08, 0.1), (0.12, 0.1), (0.12, 0)], + position=-1.04956, + ) + stochastic = StochasticRocket(rocket=calisto_robust) + + stochastic.add_free_form_fins(fins) + + added = stochastic.aerodynamic_surfaces.get_tuple_by_type(StochasticFreeFormFins) + assert len(added) == 1 + assert added[0].component.obj is fins