diff --git a/aviary/api.py b/aviary/api.py index ecb3342f6e..90a163216e 100644 --- a/aviary/api.py +++ b/aviary/api.py @@ -112,12 +112,9 @@ ) from aviary.mission.two_dof.ode.two_dof_ode import TwoDOFODE from aviary.mission.two_dof.ode.accel_ode import AccelODE as TwoDOFAccelerationODE -from aviary.mission.two_dof.ode.ascent_ode import AscentODE as TwoDOFAscentODE from aviary.mission.two_dof.ode.breguet_cruise_ode import BreguetCruiseODE from aviary.mission.two_dof.ode.flight_ode import FlightODE as TwoDOFFlightODE -from aviary.mission.two_dof.ode.flight_path_ode import FlightPathODE as TwoDOFFlightPathODE -from aviary.mission.two_dof.ode.groundroll_ode import GroundrollODE as TwoDOFGroundrollODE -from aviary.mission.two_dof.ode.rotation_ode import RotationODE as TwoDOFRotationODE +from aviary.mission.two_dof.ode.takeoff_ode import TakeOffODE as TwoDOFTakeOffODE from aviary.mission.two_dof.ode.landing_ode import LandingSegment as TwoDOFSimplifiedLanding from aviary.mission.two_dof.ode.taxi_ode import TaxiSegment as AnalyticTaxi @@ -161,12 +158,8 @@ # Phase builders from aviary.mission.two_dof.phases.accel_phase import AccelPhase as TwoDOFAccelerationPhase -from aviary.mission.two_dof.phases.ascent_phase import AscentPhase as TwoDOFAscentPhase from aviary.mission.two_dof.phases.flight_phase import FlightPhase as TwoDOFFlightPhase -from aviary.mission.two_dof.phases.groundroll_phase import ( - GroundrollPhase as TwoDOFGroundrollPhase, -) -from aviary.mission.two_dof.phases.rotation_phase import RotationPhase as TwoDOFRotationPhase +from aviary.mission.two_dof.phases.takeoff_phase import TakeoffPhase as TwoDOFTakeoffPhase # Trajectory builders from aviary.mission.height_energy.phases.detailed_landing_phases import ( diff --git a/aviary/docs/getting_started/input_csv_phase_info.ipynb b/aviary/docs/getting_started/input_csv_phase_info.ipynb index 18324e9ed0..0a1ce68df0 100644 --- a/aviary/docs/getting_started/input_csv_phase_info.ipynb +++ b/aviary/docs/getting_started/input_csv_phase_info.ipynb @@ -137,13 +137,11 @@ "from aviary.models.missions.height_energy_default import phase_info as HE_phase_info\n", "from aviary.models.missions.two_dof_default import phase_info as TwoDOF_phase_info\n", "from aviary.mission.flight_phase_builder import FlightPhaseOptions\n", - "from aviary.mission.height_energy.phases.groundroll_phase import GroundrollPhaseOptions as FGopt\n", + "from aviary.mission.solved_two_dof.phases.groundroll_phase import GroundrollPhaseOptions as FGopt\n", "from aviary.mission.two_dof.phases.accel_phase import AccelPhaseOptions\n", - "from aviary.mission.two_dof.phases.ascent_phase import AscentPhaseOptions\n", "from aviary.mission.two_dof.phases.flight_phase import FlightPhaseOptions\n", - "from aviary.mission.two_dof.phases.groundroll_phase import GroundrollPhaseOptions as GGopt\n", - "from aviary.mission.two_dof.phases.rotation_phase import RotationPhaseOptions\n", "from aviary.mission.two_dof.phases.simple_cruise_phase import SimpleCruisePhaseOptions\n", + "from aviary.mission.two_dof.phases.takeoff_phase import TakeoffPhaseOptions\n", "from aviary.mission.solved_two_dof.phases.solved_twodof_phase import SolvedTwoDOFPhaseOptions\n", "from aviary.utils.doctape import glue_keys\n", "\n", @@ -166,14 +164,12 @@ "dummy_phase_info.update(TwoDOF_phase_info)\n", "\n", "dummy_phase_info.update(AccelPhaseOptions().items())\n", - "dummy_phase_info.update(AscentPhaseOptions().items())\n", "dummy_phase_info.update(FlightPhaseOptions().items())\n", "dummy_phase_info.update(FGopt().items())\n", "dummy_phase_info.update(FlightPhaseOptions().items())\n", - "dummy_phase_info.update(GGopt().items())\n", - "dummy_phase_info.update(RotationPhaseOptions().items())\n", "dummy_phase_info.update(SimpleCruisePhaseOptions().items())\n", "dummy_phase_info.update(SolvedTwoDOFPhaseOptions().items())\n", + "dummy_phase_info.update(TakeoffPhaseOptions().items())\n", "glue_keys(dummy_phase_info)" ] }, diff --git a/aviary/docs/user_guide/phase_info.ipynb b/aviary/docs/user_guide/phase_info.ipynb index 493d0e48ec..d21aa59488 100644 --- a/aviary/docs/user_guide/phase_info.ipynb +++ b/aviary/docs/user_guide/phase_info.ipynb @@ -12,25 +12,22 @@ "source": [ "# Testing Cell\n", "from aviary.mission.flight_phase_builder import FlightPhaseOptions\n", - "from aviary.mission.height_energy.phases.groundroll_phase import GroundrollPhaseOptions as FGopt\n", + "from aviary.mission.solved_two_dof.phases.groundroll_phase import GroundrollPhaseOptions as FGopt\n", "from aviary.mission.two_dof.phases.accel_phase import AccelPhaseOptions\n", - "from aviary.mission.two_dof.phases.ascent_phase import AscentPhaseOptions\n", - "from aviary.mission.two_dof.phases.flight_phase import FlightPhaseOptions\n", - "from aviary.mission.two_dof.phases.groundroll_phase import GroundrollPhaseOptions as GGopt\n", - "from aviary.mission.two_dof.phases.rotation_phase import RotationPhaseOptions\n", + "from aviary.mission.two_dof.phases.flight_phase import FlightPhaseOptions as GFlightPhaseOptions\n", "from aviary.mission.two_dof.phases.simple_cruise_phase import SimpleCruisePhaseOptions\n", + "from aviary.mission.two_dof.phases.takeoff_phase import TakeoffPhaseOptions\n", "from aviary.mission.solved_two_dof.phases.solved_twodof_phase import SolvedTwoDOFPhaseOptions\n", "from aviary.utils.doctape import glue_keys\n", "\n", "dummy_phase_info = {}\n", "dummy_phase_info.update(AccelPhaseOptions().items())\n", - "dummy_phase_info.update(AscentPhaseOptions().items())\n", "dummy_phase_info.update(FlightPhaseOptions().items())\n", "dummy_phase_info.update(FGopt().items())\n", "dummy_phase_info.update(FlightPhaseOptions().items())\n", - "dummy_phase_info.update(GGopt().items())\n", - "dummy_phase_info.update(RotationPhaseOptions().items())\n", + "dummy_phase_info.update(GFlightPhaseOptions().items())\n", "dummy_phase_info.update(SimpleCruisePhaseOptions().items())\n", + "dummy_phase_info.update(TakeoffPhaseOptions().items())\n", "dummy_phase_info.update(SolvedTwoDOFPhaseOptions().items())\n", "glue_keys(dummy_phase_info)" ] @@ -122,7 +119,7 @@ "\n", "- **{glue:md}`no_descent`**: Set to True to prevent the aircraft from descending during the phase. This can be used to prevent unexpected descent during a climb phase. This is only available when using the height_energy equations of motion.\n", "\n", - "- **{glue:md}'transcription'**: Set the Dymos transcription to use for the phase. Aviary supports 'Collocation' or 'PicardShooting'. [See Dymos documentation for more details](https://openmdao.github.io/dymos/getting_started/transcriptions.html)\n", + "- **{glue:md}`transcription`**: Set the Dymos transcription to use for the phase. Aviary supports 'Collocation' or 'PicardShooting'. [See Dymos documentation for more details](https://openmdao.github.io/dymos/getting_started/transcriptions.html)\n", "\n", "### Specialized phase settings\n", "\n", @@ -201,9 +198,9 @@ "\n", "- **{glue:md}`normal_ref0`**: Additive scale factor \"ref0\" for the normal force constraint. Only available in rotation phases (i.e, phases with \"rotation\" in their name).\n", "\n", - "- **{glue:md}`pitch_constraint_bounds`**: Lower and upper bounds of the pitch constraint. Only available in ascent phases (i.e, phases with \"ascent\" in their name).\n", + "- **{glue:md}`pitch_constraint_bounds`**: Lower and upper bounds of the pitch constraint. Only available in takeoff phases.\n", "\n", - "- **{glue:md}`pitch_constraint_ref`**: Multiplicative scale factor \"ref\" for the pitch constraint. Only available in ascent phases (i.e, phases with \"ascent\" in their name).\n", + "- **{glue:md}`pitch_constraint_ref`**: Multiplicative scale factor \"ref\" for the pitch constraint. Only available in takeoff phases.\n", "\n", "### Complete List of All Phase Info Options\n", "\n", diff --git a/aviary/docs/user_guide/phase_info_conversion.ipynb b/aviary/docs/user_guide/phase_info_conversion.ipynb index 07660b803c..df180ba6ce 100644 --- a/aviary/docs/user_guide/phase_info_conversion.ipynb +++ b/aviary/docs/user_guide/phase_info_conversion.ipynb @@ -12,25 +12,22 @@ "source": [ "# Testing Cell\n", "from aviary.mission.flight_phase_builder import FlightPhaseOptions\n", - "from aviary.mission.height_energy.phases.groundroll_phase import GroundrollPhaseOptions as FGopt\n", + "from aviary.mission.solved_two_dof.phases.groundroll_phase import GroundrollPhaseOptions as FGopt\n", "from aviary.mission.two_dof.phases.accel_phase import AccelPhaseOptions\n", - "from aviary.mission.two_dof.phases.ascent_phase import AscentPhaseOptions\n", - "from aviary.mission.two_dof.phases.flight_phase import FlightPhaseOptions\n", - "from aviary.mission.two_dof.phases.groundroll_phase import GroundrollPhaseOptions as GGopt\n", - "from aviary.mission.two_dof.phases.rotation_phase import RotationPhaseOptions\n", + "from aviary.mission.two_dof.phases.flight_phase import FlightPhaseOptions as GFlightPhaseOptions\n", "from aviary.mission.two_dof.phases.simple_cruise_phase import SimpleCruisePhaseOptions\n", + "from aviary.mission.two_dof.phases.takeoff_phase import TakeoffPhaseOptions\n", "from aviary.mission.solved_two_dof.phases.solved_twodof_phase import SolvedTwoDOFPhaseOptions\n", "from aviary.utils.doctape import glue_keys\n", "\n", "dummy_phase_info = {}\n", "dummy_phase_info.update(AccelPhaseOptions().items())\n", - "dummy_phase_info.update(AscentPhaseOptions().items())\n", "dummy_phase_info.update(FlightPhaseOptions().items())\n", "dummy_phase_info.update(FGopt().items())\n", "dummy_phase_info.update(FlightPhaseOptions().items())\n", - "dummy_phase_info.update(GGopt().items())\n", - "dummy_phase_info.update(RotationPhaseOptions().items())\n", + "dummy_phase_info.update(GFlightPhaseOptions().items())\n", "dummy_phase_info.update(SimpleCruisePhaseOptions().items())\n", + "dummy_phase_info.update(TakeoffPhaseOptions().items())\n", "dummy_phase_info.update(SolvedTwoDOFPhaseOptions().items())\n", "glue_keys(dummy_phase_info)" ] @@ -604,7 +601,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "This phase_info has been updated to the new format below:" + "This phase_info has been updated to the new format below. Note that the new \"phase_builder\" option allows different specific equations of motion to be selected. Formerly, these were selected based on the name that the user gave to each phase." ] }, { @@ -619,6 +616,8 @@ " 'groundroll': {\n", " 'subsystem_options': {'core_aerodynamics': {'method': 'low_speed'}},\n", " 'user_options': {\n", + " 'phase_builder': PhaseType.TWO_DOF_TAKEOFF,\n", + " 'ground_roll': True,\n", " 'num_segments': 1,\n", " 'order': 3,\n", " 'time_initial': (0.0, 's'),\n", @@ -645,6 +644,8 @@ " 'rotation': {\n", " 'subsystem_options': {'core_aerodynamics': {'method': 'low_speed'}},\n", " 'user_options': {\n", + " 'phase_builder': PhaseType.TWO_DOF_TAKEOFF,\n", + " 'rotation': True,\n", " 'num_segments': 1,\n", " 'order': 3,\n", " 'time_duration_bounds': ((1, 100), 's'),\n", @@ -674,6 +675,7 @@ " 'ascent': {\n", " 'subsystem_options': {'core_aerodynamics': {'method': 'low_speed'}},\n", " 'user_options': {\n", + " 'phase_builder': PhaseType.TWO_DOF_TAKEOFF,\n", " 'num_segments': 4,\n", " 'order': 3,\n", " 'velocity_bounds': ((0, 700), 'kn'),\n", @@ -712,6 +714,7 @@ " 'accel': {\n", " 'subsystem_options': {'core_aerodynamics': {'method': 'cruise'}},\n", " 'user_options': {\n", + " 'phase_builder': PhaseType.ACCEL,\n", " 'num_segments': 1,\n", " 'order': 3,\n", " 'alt': (500, 'ft'),\n", diff --git a/aviary/docs/user_guide/phase_info_detailed.ipynb b/aviary/docs/user_guide/phase_info_detailed.ipynb index 825fa13d4d..41493de753 100644 --- a/aviary/docs/user_guide/phase_info_detailed.ipynb +++ b/aviary/docs/user_guide/phase_info_detailed.ipynb @@ -65,7 +65,9 @@ "\n", "The TWO_DEGREES_OF_FREEDOM phases contain some options that are unique to each phase.\n", "\n", - "### AccelPhase\n" + "### AccelPhase\n", + "\n", + "This phase is only used for level accelerating flight." ] }, { @@ -85,41 +87,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "### AscentPhase\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "om.show_options_table('aviary.mission.two_dof.phases.ascent_phase.AscentPhaseOptions')" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### SimpleCruisePhase\n", - "\n", - "This is the best choice for a level cruise when using the two DOF equations. SimpleCruisePhase integrates mass using the specified transcription, but the distance is computed from the fixed velocity. This can be selected by setting the \"phase_type\" option to `PhaseType.SIMPLE_CRUISE`.\n" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "om.show_options_table('aviary.mission.two_dof.phases.simple_cruise_phase.SimpleCruisePhaseOptions')" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### BreguetCruisePhase (Optional)\n", + "### BreguetCruisePhase\n", "\n", "The `BreguetCruisePhase` computes distance analytically using a form of the Breguet Range equation. This will generally be a little faster than `SimpleCruisePhase`, but it comes with a limitation that no external subsystems can contain a state that requires integration. This can be selected by setting the \"phase_type\" option to `PhaseType.BREGUET_RANGE`.\n" ] @@ -157,7 +125,9 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "### GroundrollPhase\n" + "### SimpleCruisePhase\n", + "\n", + "This is the best choice for a level cruise when using the two DOF equations. SimpleCruisePhase integrates mass using the specified transcription, but the distance is computed from the fixed velocity. This can be selected by setting the \"phase_type\" option to `PhaseType.SIMPLE_CRUISE`.\n" ] }, { @@ -166,14 +136,16 @@ "metadata": {}, "outputs": [], "source": [ - "om.show_options_table('aviary.mission.two_dof.phases.groundroll_phase.GroundrollPhaseOptions')" + "om.show_options_table('aviary.mission.two_dof.phases.simple_cruise_phase.SimpleCruisePhaseOptions')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "### RotationPhase\n" + "### TakeoffPhase\n", + "\n", + "This phase can be used for all takeoff phases and includes dynamics for ground roll, rotation, and ascent." ] }, { @@ -182,8 +154,15 @@ "metadata": {}, "outputs": [], "source": [ - "om.show_options_table('aviary.mission.two_dof.phases.rotation_phase.RotationPhaseOptions')" + "om.show_options_table('aviary.mission.two_dof.phases.takeoff_phase.TakeoffPhaseOptions')" ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] } ], "metadata": { diff --git a/aviary/mission/two_dof/ode/groundroll_eom.py b/aviary/mission/solved_two_dof/ode/groundroll_eom.py similarity index 100% rename from aviary/mission/two_dof/ode/groundroll_eom.py rename to aviary/mission/solved_two_dof/ode/groundroll_eom.py diff --git a/aviary/mission/two_dof/ode/groundroll_ode.py b/aviary/mission/solved_two_dof/ode/groundroll_ode.py similarity index 98% rename from aviary/mission/two_dof/ode/groundroll_ode.py rename to aviary/mission/solved_two_dof/ode/groundroll_ode.py index 9768b4b7f7..27ff86fe69 100644 --- a/aviary/mission/two_dof/ode/groundroll_ode.py +++ b/aviary/mission/solved_two_dof/ode/groundroll_ode.py @@ -1,7 +1,7 @@ import numpy as np import openmdao.api as om -from aviary.mission.two_dof.ode.groundroll_eom import GroundrollEOM +from aviary.mission.solved_two_dof.ode.groundroll_eom import GroundrollEOM from aviary.mission.two_dof.ode.params import ParamPort from aviary.mission.two_dof.ode.two_dof_ode import TwoDOFODE from aviary.subsystems.aerodynamics.aerodynamics_builder import AerodynamicsBuilder diff --git a/aviary/mission/two_dof/ode/test/test_groundroll_eom.py b/aviary/mission/solved_two_dof/ode/test/test_groundroll_eom.py similarity index 94% rename from aviary/mission/two_dof/ode/test/test_groundroll_eom.py rename to aviary/mission/solved_two_dof/ode/test/test_groundroll_eom.py index f182934483..369e74b4a1 100644 --- a/aviary/mission/two_dof/ode/test/test_groundroll_eom.py +++ b/aviary/mission/solved_two_dof/ode/test/test_groundroll_eom.py @@ -4,7 +4,7 @@ import openmdao.api as om from openmdao.utils.assert_utils import assert_check_partials, assert_near_equal -from aviary.mission.two_dof.ode.groundroll_eom import GroundrollEOM +from aviary.mission.solved_two_dof.ode.groundroll_eom import GroundrollEOM from aviary.variable_info.variables import Aircraft, Dynamic @@ -59,12 +59,12 @@ class GroundrollEOMTestCase2(unittest.TestCase): """Test mass-weight conversion.""" def setUp(self): - import aviary.mission.two_dof.ode.groundroll_eom as gr + import aviary.mission.solved_two_dof.ode.groundroll_eom as gr gr.GRAV_ENGLISH_LBM = 1.1 def tearDown(self): - import aviary.mission.two_dof.ode.groundroll_eom as gr + import aviary.mission.solved_two_dof.ode.groundroll_eom as gr gr.GRAV_ENGLISH_LBM = 1.0 diff --git a/aviary/mission/two_dof/ode/test/test_groundroll_ode.py b/aviary/mission/solved_two_dof/ode/test/test_groundroll_ode.py similarity index 97% rename from aviary/mission/two_dof/ode/test/test_groundroll_ode.py rename to aviary/mission/solved_two_dof/ode/test/test_groundroll_ode.py index a006d169f9..7e39ea452b 100644 --- a/aviary/mission/two_dof/ode/test/test_groundroll_ode.py +++ b/aviary/mission/solved_two_dof/ode/test/test_groundroll_ode.py @@ -3,7 +3,7 @@ import openmdao.api as om from openmdao.utils.assert_utils import assert_check_partials -from aviary.mission.two_dof.ode.groundroll_ode import GroundrollODE +from aviary.mission.solved_two_dof.ode.groundroll_ode import GroundrollODE from aviary.mission.two_dof.ode.params import set_params_for_unit_tests from aviary.subsystems.propulsion.utils import build_engine_deck from aviary.utils.test_utils.default_subsystems import get_default_mission_subsystems diff --git a/aviary/mission/height_energy/phases/groundroll_phase.py b/aviary/mission/solved_two_dof/phases/groundroll_phase.py similarity index 96% rename from aviary/mission/height_energy/phases/groundroll_phase.py rename to aviary/mission/solved_two_dof/phases/groundroll_phase.py index 432104eb9c..72836eb348 100644 --- a/aviary/mission/height_energy/phases/groundroll_phase.py +++ b/aviary/mission/solved_two_dof/phases/groundroll_phase.py @@ -1,17 +1,17 @@ import dymos as dm -from aviary.mission.two_dof.ode.groundroll_ode import GroundrollODE from aviary.mission.initial_guess_builders import ( InitialGuessIntegrationVariable, InitialGuessPolynomialControl, InitialGuessState, ) from aviary.mission.phase_builder import PhaseBuilder, register +from aviary.mission.phase_utils import add_subsystem_variables_to_phase +from aviary.mission.solved_two_dof.ode.groundroll_ode import GroundrollODE from aviary.utils.aviary_options_dict import AviaryOptionsDictionary from aviary.utils.aviary_values import AviaryValues from aviary.variable_info.variable_meta_data import _MetaData from aviary.variable_info.variables import Dynamic -from aviary.mission.phase_utils import add_subsystem_variables_to_phase # Solved 2DOF uses this builder. # @@ -69,7 +69,11 @@ def declare_options(self): @register class GroundrollPhase(PhaseBuilder): - """A phase builder for a two degree of freedom (2DOF) phase.""" + """ + A phase builder for a two degree of freedom (2DOF) ground roll phase. + + This is used exclusively by the Solved 2DOF Phase Builder. + """ __slots__ = ('subsystems', 'meta_data') diff --git a/aviary/mission/solved_two_dof_problem_configurator.py b/aviary/mission/solved_two_dof_problem_configurator.py index 6340f8c33d..ffc6b69961 100644 --- a/aviary/mission/solved_two_dof_problem_configurator.py +++ b/aviary/mission/solved_two_dof_problem_configurator.py @@ -1,4 +1,4 @@ -from aviary.mission.height_energy.phases.groundroll_phase import ( +from aviary.mission.solved_two_dof.phases.groundroll_phase import ( GroundrollPhase as GroundrollPhaseVelocityIntegrated, ) from aviary.mission.solved_two_dof.phases.solved_twodof_phase import SolvedTwoDOFPhase diff --git a/aviary/mission/two_dof/ode/ascent_eom.py b/aviary/mission/two_dof/ode/ascent_eom.py deleted file mode 100644 index 3be401f7cf..0000000000 --- a/aviary/mission/two_dof/ode/ascent_eom.py +++ /dev/null @@ -1,346 +0,0 @@ -import numpy as np -import openmdao.api as om - -from aviary.constants import GRAV_ENGLISH_GASP, GRAV_ENGLISH_LBM, MU_TAKEOFF -from aviary.variable_info.functions import add_aviary_input -from aviary.variable_info.variables import Aircraft, Dynamic - - -class AscentEOM(om.ExplicitComponent): - """ - Ascent equation of motion: Compute the velocity rate, flight path angle rate, - altitude rate, distance rate, normal forces, fuselage pitch angle, load factor - and angle of attack rate. - """ - - def initialize(self): - self.options.declare('num_nodes', types=int) - - def setup(self): - nn = self.options['num_nodes'] - - self.add_input(Dynamic.Vehicle.MASS, val=np.ones(nn), desc='aircraft mass', units='lbm') - self.add_input( - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - val=np.ones(nn), - desc=Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - units='lbf', - ) - self.add_input( - Dynamic.Vehicle.LIFT, - val=np.ones(nn), - desc=Dynamic.Vehicle.LIFT, - units='lbf', - ) - self.add_input( - Dynamic.Vehicle.DRAG, - val=np.ones(nn), - desc=Dynamic.Vehicle.DRAG, - units='lbf', - ) - self.add_input(Dynamic.Mission.VELOCITY, val=np.ones(nn), desc='Velocity', units='ft/s') - self.add_input( - Dynamic.Mission.FLIGHT_PATH_ANGLE, - val=np.ones(nn), - desc='flight path angle', - units='rad', - ) - add_aviary_input(self, Aircraft.Wing.INCIDENCE, val=0, units='deg') - - self.add_input( - Dynamic.Vehicle.ANGLE_OF_ATTACK, - val=np.ones(nn), - desc='angle of attack', - units='deg', - ) - - self.add_output( - Dynamic.Mission.VELOCITY_RATE, - val=np.ones(nn), - desc='Velocity rate', - units='ft/s**2', - ) - self.add_output( - Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, - val=np.ones(nn), - desc='flight path angle rate', - units='rad/s', - ) - self.add_output( - Dynamic.Mission.ALTITUDE_RATE, - val=np.ones(nn), - desc='altitude rate', - units='ft/s', - ) - self.add_output( - Dynamic.Mission.DISTANCE_RATE, val=np.ones(nn), desc='distance rate', units='ft/s' - ) - self.add_output('normal_force', val=np.ones(nn), desc='normal forces', units='lbf') - self.add_output('fuselage_pitch', val=np.ones(nn), desc='fuselage pitch angle', units='deg') - self.add_output('load_factor', val=np.ones(nn), desc='load factor', units='unitless') - self.add_output( - 'angle_of_attack_rate', val=np.ones(nn), desc='angle of attack rate', units='deg/s' - ) - - def setup_partials(self): - arange = np.arange(self.options['num_nodes']) - - self.declare_partials( - Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, - [ - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - Dynamic.Vehicle.ANGLE_OF_ATTACK, - Dynamic.Vehicle.LIFT, - Dynamic.Vehicle.MASS, - Dynamic.Mission.FLIGHT_PATH_ANGLE, - Dynamic.Mission.VELOCITY, - ], - rows=arange, - cols=arange, - ) - self.declare_partials(Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, [Aircraft.Wing.INCIDENCE]) - self.declare_partials( - 'load_factor', - [ - Dynamic.Vehicle.LIFT, - Dynamic.Vehicle.MASS, - Dynamic.Mission.FLIGHT_PATH_ANGLE, - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - Dynamic.Vehicle.ANGLE_OF_ATTACK, - ], - rows=arange, - cols=arange, - ) - self.declare_partials('load_factor', [Aircraft.Wing.INCIDENCE]) - - self.declare_partials( - Dynamic.Mission.VELOCITY_RATE, - [ - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - Dynamic.Vehicle.ANGLE_OF_ATTACK, - Dynamic.Vehicle.DRAG, - Dynamic.Vehicle.MASS, - Dynamic.Mission.FLIGHT_PATH_ANGLE, - Dynamic.Vehicle.LIFT, - ], - rows=arange, - cols=arange, - ) - self.declare_partials(Dynamic.Mission.VELOCITY_RATE, [Aircraft.Wing.INCIDENCE]) - self.declare_partials( - Dynamic.Mission.ALTITUDE_RATE, - [Dynamic.Mission.VELOCITY, Dynamic.Mission.FLIGHT_PATH_ANGLE], - rows=arange, - cols=arange, - ) - self.declare_partials( - Dynamic.Mission.DISTANCE_RATE, - [Dynamic.Mission.VELOCITY, Dynamic.Mission.FLIGHT_PATH_ANGLE], - rows=arange, - cols=arange, - ) - self.declare_partials( - 'normal_force', - [ - Dynamic.Vehicle.MASS, - Dynamic.Vehicle.LIFT, - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - Dynamic.Vehicle.ANGLE_OF_ATTACK, - ], - rows=arange, - cols=arange, - ) - self.declare_partials('normal_force', [Aircraft.Wing.INCIDENCE]) - self.declare_partials( - 'fuselage_pitch', - Dynamic.Mission.FLIGHT_PATH_ANGLE, - rows=arange, - cols=arange, - val=180 / np.pi, - ) - self.declare_partials( - 'fuselage_pitch', - Dynamic.Vehicle.ANGLE_OF_ATTACK, - rows=arange, - cols=arange, - val=1, - ) - self.declare_partials('fuselage_pitch', Aircraft.Wing.INCIDENCE, val=-1) - - def compute(self, inputs, outputs): - weight = inputs[Dynamic.Vehicle.MASS] * GRAV_ENGLISH_LBM - thrust = inputs[Dynamic.Vehicle.Propulsion.THRUST_TOTAL] - incremented_lift = inputs[Dynamic.Vehicle.LIFT] - incremented_drag = inputs[Dynamic.Vehicle.DRAG] - TAS = inputs[Dynamic.Mission.VELOCITY] - gamma = inputs[Dynamic.Mission.FLIGHT_PATH_ANGLE] - i_wing = inputs[Aircraft.Wing.INCIDENCE] - alpha = inputs[Dynamic.Vehicle.ANGLE_OF_ATTACK] - - mu = 0.0 - - nn = self.options['num_nodes'] - - thrust_along_flightpath = thrust * np.cos((alpha - i_wing) * np.pi / 180) - thrust_across_flightpath = thrust * np.sin((alpha - i_wing) * np.pi / 180) - normal_force = weight - incremented_lift - thrust_across_flightpath - normal_force[normal_force < 0] = 0.0 - - outputs[Dynamic.Mission.VELOCITY_RATE] = ( - ( - thrust_along_flightpath - - incremented_drag - - weight * np.sin(gamma) - - mu * normal_force - ) - * GRAV_ENGLISH_GASP - / weight - ) - - outputs[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE] = ( - (thrust_across_flightpath + incremented_lift - weight * np.cos(gamma)) - * GRAV_ENGLISH_GASP - / (TAS * weight) - ) - - outputs[Dynamic.Mission.ALTITUDE_RATE] = TAS * np.sin(gamma) - outputs[Dynamic.Mission.DISTANCE_RATE] = TAS * np.cos(gamma) - outputs['normal_force'] = normal_force - outputs['fuselage_pitch'] = gamma * 180 / np.pi - i_wing + alpha - - outputs['angle_of_attack_rate'] = np.zeros(nn) - - load_factor = (incremented_lift + thrust_across_flightpath) / (weight * np.cos(gamma)) - - outputs['load_factor'] = load_factor - - def compute_partials(self, inputs, J): - nn = self.options['num_nodes'] - - weight = inputs[Dynamic.Vehicle.MASS] * GRAV_ENGLISH_LBM - thrust = inputs[Dynamic.Vehicle.Propulsion.THRUST_TOTAL] - incremented_lift = inputs[Dynamic.Vehicle.LIFT] - incremented_drag = inputs[Dynamic.Vehicle.DRAG] - TAS = inputs[Dynamic.Mission.VELOCITY] - gamma = inputs[Dynamic.Mission.FLIGHT_PATH_ANGLE] - i_wing = inputs[Aircraft.Wing.INCIDENCE] - alpha = inputs[Dynamic.Vehicle.ANGLE_OF_ATTACK] - - mu = MU_TAKEOFF - mu = 0.0 - - thrust_along_flightpath = thrust * np.cos((alpha - i_wing) * np.pi / 180) - thrust_across_flightpath = thrust * np.sin((alpha - i_wing) * np.pi / 180) - - dTAlF_dThrust = np.cos((alpha - i_wing) * np.pi / 180) - dTAlF_dAlpha = -thrust * np.sin((alpha - i_wing) * np.pi / 180) * np.pi / 180 - dTAlF_dIwing = thrust * np.sin((alpha - i_wing) * np.pi / 180) * np.pi / 180 - - dTAcF_dThrust = np.sin((alpha - i_wing) * np.pi / 180) - dTAcF_dAlpha = thrust * np.cos((alpha - i_wing) * np.pi / 180) * np.pi / 180 - dTAcF_dIwing = -thrust * np.cos((alpha - i_wing) * np.pi / 180) * np.pi / 180 - - J[ - Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - ] = dTAcF_dThrust * GRAV_ENGLISH_GASP / (TAS * weight) - J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Vehicle.ANGLE_OF_ATTACK] = ( - dTAcF_dAlpha * GRAV_ENGLISH_GASP / (TAS * weight) - ) - J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Aircraft.Wing.INCIDENCE] = ( - dTAcF_dIwing * GRAV_ENGLISH_GASP / (TAS * weight) - ) - J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Vehicle.LIFT] = GRAV_ENGLISH_GASP / ( - TAS * weight - ) - J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Vehicle.MASS] = ( - (GRAV_ENGLISH_GASP / TAS) - * GRAV_ENGLISH_LBM - * (-thrust_across_flightpath / weight**2 - incremented_lift / weight**2) - ) - J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = ( - weight * np.sin(gamma) * GRAV_ENGLISH_GASP / (TAS * weight) - ) - J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Mission.VELOCITY] = -( - (thrust_across_flightpath + incremented_lift - weight * np.cos(gamma)) - * GRAV_ENGLISH_GASP - / (TAS**2 * weight) - ) - - J['load_factor', Dynamic.Vehicle.LIFT] = 1 / (weight * np.cos(gamma)) - J['load_factor', Dynamic.Vehicle.MASS] = ( - -(incremented_lift + thrust_across_flightpath) - / (weight**2 * np.cos(gamma)) - * GRAV_ENGLISH_LBM - ) - J['load_factor', Dynamic.Mission.FLIGHT_PATH_ANGLE] = ( - -(incremented_lift + thrust_across_flightpath) - / (weight * (np.cos(gamma)) ** 2) - * (-np.sin(gamma)) - ) - J['load_factor', Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = dTAcF_dThrust / ( - weight * np.cos(gamma) - ) - J['load_factor', Dynamic.Vehicle.ANGLE_OF_ATTACK] = dTAcF_dAlpha / (weight * np.cos(gamma)) - J['load_factor', Aircraft.Wing.INCIDENCE] = dTAcF_dIwing / (weight * np.cos(gamma)) - - normal_force1 = weight - incremented_lift - thrust_across_flightpath - normal_force = np.where(normal_force1 < 0, np.zeros(nn), normal_force1) - - dNF_dWeight = np.ones(nn) - dNF_dWeight[normal_force1 < 0] = 0 - - dNF_dLift = -np.ones(nn) - dNF_dLift[normal_force1 < 0] = 0 - - dNF_dThrust = -np.ones(nn) * dTAcF_dThrust - dNF_dThrust[normal_force1 < 0] = 0 - - dNF_dAlpha = -np.ones(nn) * dTAcF_dAlpha - dNF_dAlpha[normal_force1 < 0] = 0 - - dNF_dIwing = -np.ones(nn) * dTAcF_dIwing - dNF_dIwing[normal_force1 < 0] = 0 - - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = ( - (dTAlF_dThrust - mu * dNF_dThrust) * GRAV_ENGLISH_GASP / weight - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.ANGLE_OF_ATTACK] = ( - (dTAlF_dAlpha - mu * dNF_dAlpha) * GRAV_ENGLISH_GASP / weight - ) - J[Dynamic.Mission.VELOCITY_RATE, Aircraft.Wing.INCIDENCE] = ( - (dTAlF_dIwing - mu * dNF_dIwing) * GRAV_ENGLISH_GASP / weight - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.DRAG] = -GRAV_ENGLISH_GASP / weight - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.MASS] = ( - GRAV_ENGLISH_GASP - * GRAV_ENGLISH_LBM - * ( - weight * (-np.sin(gamma) - mu * dNF_dWeight) - - ( - thrust_along_flightpath - - incremented_drag - - weight * np.sin(gamma) - - mu * normal_force - ) - ) - / weight**2 - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = ( - -np.cos(gamma) * GRAV_ENGLISH_GASP - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.LIFT] = ( - GRAV_ENGLISH_GASP * (-mu * dNF_dLift) / weight - ) - - J[Dynamic.Mission.ALTITUDE_RATE, Dynamic.Mission.VELOCITY] = np.sin(gamma) - J[Dynamic.Mission.ALTITUDE_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = TAS * np.cos(gamma) - - J[Dynamic.Mission.DISTANCE_RATE, Dynamic.Mission.VELOCITY] = np.cos(gamma) - J[Dynamic.Mission.DISTANCE_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = -TAS * np.sin(gamma) - - J['normal_force', Dynamic.Vehicle.MASS] = dNF_dWeight * GRAV_ENGLISH_LBM - J['normal_force', Dynamic.Vehicle.LIFT] = dNF_dLift - J['normal_force', Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = dNF_dThrust - J['normal_force', Dynamic.Vehicle.ANGLE_OF_ATTACK] = dNF_dAlpha - J['normal_force', Aircraft.Wing.INCIDENCE] = dNF_dIwing diff --git a/aviary/mission/two_dof/ode/ascent_ode.py b/aviary/mission/two_dof/ode/ascent_ode.py deleted file mode 100644 index c09cb2a353..0000000000 --- a/aviary/mission/two_dof/ode/ascent_ode.py +++ /dev/null @@ -1,88 +0,0 @@ -import numpy as np - -from aviary.mission.two_dof.ode.ascent_eom import AscentEOM -from aviary.mission.two_dof.ode.params import ParamPort -from aviary.mission.two_dof.ode.two_dof_ode import TwoDOFODE -from aviary.variable_info.enums import AlphaModes -from aviary.variable_info.variables import Aircraft, Dynamic - - -class AscentODE(TwoDOFODE): - """ODE for initial ascent. - - This phase is intended to model the portion of aircraft flight starting when the - aircraft detaches from the runway, and it is capable of retracting flaps and landing - gear as necessary. - """ - - def initialize(self): - super().initialize() - self.options.declare('alpha_mode', default=AlphaModes.DEFAULT, types=AlphaModes) - - def setup(self): - nn = self.options['num_nodes'] - alpha_mode = self.options['alpha_mode'] - - # TODO: paramport - ascent_params = ParamPort() - self.add_subsystem('params', ascent_params, promotes=['*']) - - self.add_atmosphere() - - kwargs = { - 'method': 'low_speed', - 'retract_gear': True, - 'retract_flaps': True, - } - self.options['subsystem_options'].setdefault('aerodynamics', {}).update(kwargs) - - self.add_subsystems() - - if alpha_mode is AlphaModes.DEFAULT: - # alpha as input - pass - else: - self.add_alpha_control(alpha_mode=alpha_mode, num_nodes=nn) - - self.add_subsystem( - 'ascent_eom', - AscentEOM(num_nodes=nn), - promotes_inputs=[ - Dynamic.Vehicle.MASS, - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - Dynamic.Vehicle.LIFT, - Dynamic.Vehicle.DRAG, - Dynamic.Mission.VELOCITY, - Dynamic.Mission.FLIGHT_PATH_ANGLE, - Dynamic.Vehicle.ANGLE_OF_ATTACK, - ] - + ['aircraft:*'], - promotes_outputs=[ - Dynamic.Mission.VELOCITY_RATE, - Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, - Dynamic.Mission.ALTITUDE_RATE, - Dynamic.Mission.DISTANCE_RATE, - 'angle_of_attack_rate', - 'normal_force', - 'fuselage_pitch', - 'load_factor', - ], - ) - - self.add_excess_rate_comps(nn) - - ParamPort.set_default_vals(self) - self.set_input_defaults('t_init_flaps', val=47.5) - self.set_input_defaults('t_init_gear', val=37.3) - self.set_input_defaults(Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(nn), units='deg') - self.set_input_defaults(Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(nn), units='deg') - self.set_input_defaults(Dynamic.Mission.ALTITUDE, val=np.zeros(nn), units='ft') - self.set_input_defaults(Dynamic.Mission.VELOCITY, val=np.zeros(nn), units='kn') - self.set_input_defaults('t_curr', val=np.zeros(nn), units='s') - self.set_input_defaults('aero_ramps.flap_factor:final_val', val=0.0) - self.set_input_defaults('aero_ramps.gear_factor:final_val', val=0.0) - self.set_input_defaults('aero_ramps.flap_factor:initial_val', val=1.0) - self.set_input_defaults('aero_ramps.gear_factor:initial_val', val=1.0) - self.set_input_defaults( - Dynamic.Vehicle.MASS, val=np.ones(nn), units='kg' - ) # val here is nominal diff --git a/aviary/mission/two_dof/ode/rotation_eom.py b/aviary/mission/two_dof/ode/rotation_eom.py deleted file mode 100644 index 863a827da1..0000000000 --- a/aviary/mission/two_dof/ode/rotation_eom.py +++ /dev/null @@ -1,269 +0,0 @@ -import numpy as np -import openmdao.api as om - -from aviary.constants import GRAV_ENGLISH_GASP, GRAV_ENGLISH_LBM, MU_TAKEOFF -from aviary.variable_info.functions import add_aviary_input -from aviary.variable_info.variables import Aircraft, Dynamic - - -class RotationEOM(om.ExplicitComponent): - """2-degree of freedom rotation EOM.""" - - def initialize(self): - self.options.declare('num_nodes', types=int) - - def setup(self): - nn = self.options['num_nodes'] - - self.add_input(Dynamic.Vehicle.MASS, val=np.ones(nn), desc='aircraft mass', units='lbm') - self.add_input( - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - val=np.ones(nn), - desc=Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - units='lbf', - ) - self.add_input( - Dynamic.Vehicle.LIFT, - val=np.ones(nn), - desc=Dynamic.Vehicle.LIFT, - units='lbf', - ) - self.add_input( - Dynamic.Vehicle.DRAG, - val=np.ones(nn), - desc=Dynamic.Vehicle.DRAG, - units='lbf', - ) - self.add_input( - Dynamic.Mission.VELOCITY, - val=np.ones(nn), - desc='true air speed', - units='ft/s', - ) - self.add_input( - Dynamic.Mission.FLIGHT_PATH_ANGLE, - val=np.ones(nn), - desc='flight path angle', - units='rad', - ) - - add_aviary_input(self, Aircraft.Wing.INCIDENCE, val=0.0, units='deg') - self.add_input( - Dynamic.Vehicle.ANGLE_OF_ATTACK, - val=np.ones(nn), - desc='angle of attack', - units='deg', - ) - - self.add_output( - Dynamic.Mission.VELOCITY_RATE, - val=np.ones(nn), - desc='TAS rate', - units='ft/s**2', - ) - self.add_output( - Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, - val=np.ones(nn), - desc='flight path angle rate', - units='rad/s', - ) - self.add_output( - Dynamic.Mission.ALTITUDE_RATE, - val=np.ones(nn), - desc='altitude rate', - units='ft/s', - ) - self.add_output( - Dynamic.Mission.DISTANCE_RATE, val=np.ones(nn), desc='distance rate', units='ft/s' - ) - self.add_output('normal_force', val=np.ones(nn), desc='normal forces', units='lbf') - self.add_output('fuselage_pitch', val=np.ones(nn), desc='fuselage pitch angle', units='deg') - self.add_output( - 'angle_of_attack_rate', val=np.ones(nn), desc='angle of attack rate', units='deg/s' - ) - - self.declare_partials('angle_of_attack_rate', ['*']) - - def setup_partials(self): - arange = np.arange(self.options['num_nodes']) - - self.declare_partials(Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, '*') - self.declare_partials( - Dynamic.Mission.VELOCITY_RATE, - [ - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - Dynamic.Vehicle.ANGLE_OF_ATTACK, - Dynamic.Vehicle.DRAG, - Dynamic.Vehicle.MASS, - Dynamic.Mission.FLIGHT_PATH_ANGLE, - Dynamic.Vehicle.LIFT, - ], - rows=arange, - cols=arange, - ) - self.declare_partials(Dynamic.Mission.VELOCITY_RATE, [Aircraft.Wing.INCIDENCE]) - self.declare_partials( - Dynamic.Mission.ALTITUDE_RATE, - [Dynamic.Mission.VELOCITY, Dynamic.Mission.FLIGHT_PATH_ANGLE], - rows=arange, - cols=arange, - ) - self.declare_partials( - Dynamic.Mission.DISTANCE_RATE, - [Dynamic.Mission.VELOCITY, Dynamic.Mission.FLIGHT_PATH_ANGLE], - rows=arange, - cols=arange, - ) - - self.declare_partials( - 'normal_force', - [ - Dynamic.Vehicle.MASS, - Dynamic.Vehicle.LIFT, - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, - Dynamic.Vehicle.ANGLE_OF_ATTACK, - ], - rows=arange, - cols=arange, - ) - self.declare_partials('normal_force', [Aircraft.Wing.INCIDENCE]) - self.declare_partials( - 'fuselage_pitch', - Dynamic.Mission.FLIGHT_PATH_ANGLE, - rows=arange, - cols=arange, - val=180 / np.pi, - ) - self.declare_partials( - 'fuselage_pitch', - Dynamic.Vehicle.ANGLE_OF_ATTACK, - rows=arange, - cols=arange, - val=1, - ) - self.declare_partials('fuselage_pitch', Aircraft.Wing.INCIDENCE, val=-1) - - def compute(self, inputs, outputs): - weight = inputs[Dynamic.Vehicle.MASS] * GRAV_ENGLISH_LBM - thrust = inputs[Dynamic.Vehicle.Propulsion.THRUST_TOTAL] - incremented_lift = inputs[Dynamic.Vehicle.LIFT] - incremented_drag = inputs[Dynamic.Vehicle.DRAG] - TAS = inputs[Dynamic.Mission.VELOCITY] - gamma = inputs[Dynamic.Mission.FLIGHT_PATH_ANGLE] - i_wing = inputs[Aircraft.Wing.INCIDENCE] - alpha = inputs[Dynamic.Vehicle.ANGLE_OF_ATTACK] - - mu = MU_TAKEOFF - - nn = self.options['num_nodes'] - - thrust_along_flightpath = thrust * np.cos((alpha - i_wing) * np.pi / 180) - thrust_across_flightpath = thrust * np.sin((alpha - i_wing) * np.pi / 180) - normal_force = np.clip( - weight - incremented_lift - thrust_across_flightpath, a_min=0.0, a_max=None - ) - - outputs[Dynamic.Mission.VELOCITY_RATE] = ( - ( - thrust_along_flightpath - - incremented_drag - - weight * np.sin(gamma) - - mu * normal_force - ) - * GRAV_ENGLISH_GASP - / weight - ) - outputs[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE] = np.zeros(nn) - - outputs[Dynamic.Mission.ALTITUDE_RATE] = TAS * np.sin(gamma) - outputs[Dynamic.Mission.DISTANCE_RATE] = TAS * np.cos(gamma) - outputs['normal_force'] = normal_force - outputs['fuselage_pitch'] = gamma * 180 / np.pi - i_wing + alpha - outputs['angle_of_attack_rate'] = np.full(nn, 3.33) - - def compute_partials(self, inputs, J): - mu = MU_TAKEOFF - - weight = inputs[Dynamic.Vehicle.MASS] * GRAV_ENGLISH_LBM - thrust = inputs[Dynamic.Vehicle.Propulsion.THRUST_TOTAL] - incremented_lift = inputs[Dynamic.Vehicle.LIFT] - incremented_drag = inputs[Dynamic.Vehicle.DRAG] - TAS = inputs[Dynamic.Mission.VELOCITY] - gamma = inputs[Dynamic.Mission.FLIGHT_PATH_ANGLE] - i_wing = inputs[Aircraft.Wing.INCIDENCE] - alpha = inputs[Dynamic.Vehicle.ANGLE_OF_ATTACK] - - nn = self.options['num_nodes'] - - thrust_along_flightpath = thrust * np.cos((alpha - i_wing) * np.pi / 180) - thrust_across_flightpath = thrust * np.sin((alpha - i_wing) * np.pi / 180) - - dTAlF_dThrust = np.cos((alpha - i_wing) * np.pi / 180) - dTAlF_dAlpha = -thrust * np.sin((alpha - i_wing) * np.pi / 180) * np.pi / 180 - dTAlF_dIwing = thrust * np.sin((alpha - i_wing) * np.pi / 180) * np.pi / 180 - - dTAcF_dThrust = np.sin((alpha - i_wing) * np.pi / 180) - dTAcF_dAlpha = thrust * np.cos((alpha - i_wing) * np.pi / 180) * np.pi / 180 - dTAcF_dIwing = -thrust * np.cos((alpha - i_wing) * np.pi / 180) * np.pi / 180 - - normal_force = np.clip( - weight - incremented_lift - thrust_across_flightpath, a_min=0.0, a_max=None - ) - - dNF_dWeight = np.ones(nn) - dNF_dWeight[normal_force < 0] = 0 - - dNF_dLift = -np.ones(nn) - dNF_dLift[normal_force < 0] = 0 - - dNF_dThrust = -np.ones(nn) * dTAcF_dThrust - dNF_dThrust[normal_force < 0] = 0 - - dNF_dAlpha = -np.ones(nn) * dTAcF_dAlpha - dNF_dAlpha[normal_force < 0] = 0 - - dNF_dIwing = -np.ones(nn) * dTAcF_dIwing - dNF_dIwing[normal_force < 0] = 0 - - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = ( - (dTAlF_dThrust - mu * dNF_dThrust) * GRAV_ENGLISH_GASP / weight - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.ANGLE_OF_ATTACK] = ( - (dTAlF_dAlpha - mu * dNF_dAlpha) * GRAV_ENGLISH_GASP / weight - ) - J[Dynamic.Mission.VELOCITY_RATE, Aircraft.Wing.INCIDENCE] = ( - (dTAlF_dIwing - mu * dNF_dIwing) * GRAV_ENGLISH_GASP / weight - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.DRAG] = -GRAV_ENGLISH_GASP / weight - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.MASS] = ( - GRAV_ENGLISH_GASP - * GRAV_ENGLISH_LBM - * ( - weight * (-np.sin(gamma) - mu * dNF_dWeight) - - ( - thrust_along_flightpath - - incremented_drag - - weight * np.sin(gamma) - - mu * normal_force - ) - ) - / weight**2 - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = ( - -np.cos(gamma) * GRAV_ENGLISH_GASP - ) - J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.LIFT] = ( - GRAV_ENGLISH_GASP * (-mu * dNF_dLift) / weight - ) - - J[Dynamic.Mission.ALTITUDE_RATE, Dynamic.Mission.VELOCITY] = np.sin(gamma) - J[Dynamic.Mission.ALTITUDE_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = TAS * np.cos(gamma) - - J[Dynamic.Mission.DISTANCE_RATE, Dynamic.Mission.VELOCITY] = np.cos(gamma) - J[Dynamic.Mission.DISTANCE_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = -TAS * np.sin(gamma) - - J['normal_force', Dynamic.Vehicle.MASS] = dNF_dWeight * GRAV_ENGLISH_LBM - J['normal_force', Dynamic.Vehicle.LIFT] = dNF_dLift - J['normal_force', Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = dNF_dThrust - J['normal_force', Dynamic.Vehicle.ANGLE_OF_ATTACK] = dNF_dAlpha - J['normal_force', Aircraft.Wing.INCIDENCE] = dNF_dIwing diff --git a/aviary/mission/two_dof/ode/rotation_ode.py b/aviary/mission/two_dof/ode/rotation_ode.py deleted file mode 100644 index 989e7ed72d..0000000000 --- a/aviary/mission/two_dof/ode/rotation_ode.py +++ /dev/null @@ -1,43 +0,0 @@ -import numpy as np -import openmdao.api as om - -from aviary.mission.two_dof.ode.params import ParamPort -from aviary.mission.two_dof.ode.rotation_eom import RotationEOM -from aviary.mission.two_dof.ode.two_dof_ode import TwoDOFODE -from aviary.variable_info.variables import Aircraft, Dynamic - - -class RotationODE(TwoDOFODE): - """ODE for takeoff rotation. - - This phase spans the time from when the aircraft is touching the runway but has - begun to rotate to liftoff. - """ - - def setup(self): - nn = self.options['num_nodes'] - - # TODO: paramport - self.add_subsystem('params', ParamPort(), promotes=['*']) - - self.add_atmosphere() - - kwargs = {'method': 'low_speed'} - self.options['subsystem_options'].setdefault('aerodynamics', {}).update(kwargs) - - self.add_subsystems() - - self.add_subsystem('rotation_eom', RotationEOM(num_nodes=nn), promotes=['*']) - - ParamPort.set_default_vals(self) - self.set_input_defaults('t_init_flaps', val=47.5, units='s') - self.set_input_defaults('t_init_gear', val=37.3, units='s') - self.set_input_defaults(Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.ones(nn), units='deg') - self.set_input_defaults(Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(nn), units='deg') - self.set_input_defaults(Dynamic.Mission.ALTITUDE, val=np.zeros(nn), units='ft') - self.set_input_defaults(Dynamic.Mission.VELOCITY, val=np.zeros(nn), units='kn') - self.set_input_defaults('t_curr', val=np.zeros(nn), units='s') - self.set_input_defaults('aero_ramps.flap_factor:final_val', val=1.0) - self.set_input_defaults('aero_ramps.gear_factor:final_val', val=1.0) - self.set_input_defaults('aero_ramps.flap_factor:initial_val', val=1.0) - self.set_input_defaults('aero_ramps.gear_factor:initial_val', val=1.0) diff --git a/aviary/mission/two_dof/ode/takeoff_eom.py b/aviary/mission/two_dof/ode/takeoff_eom.py new file mode 100644 index 0000000000..4f71181d8c --- /dev/null +++ b/aviary/mission/two_dof/ode/takeoff_eom.py @@ -0,0 +1,472 @@ +import numpy as np +import openmdao.api as om + +from aviary.constants import GRAV_ENGLISH_GASP, GRAV_ENGLISH_LBM, MU_TAKEOFF +from aviary.variable_info.functions import add_aviary_input +from aviary.variable_info.variables import Aircraft, Dynamic + +DEG2RAD = np.pi / 180.0 +RAD2DEG = 1.0 / DEG2RAD + + +class TakeoffEOM(om.ExplicitComponent): + """ + 2-degree of freedom EOM for takeoff phases. + + Compute the rates for the velocity, altitude, flight path angle, and angle of attack states. + This can be used for the groundroll, rotation, and ascent phases. The angle of attack rate + is fixed in this phase, and is set to the value in the "rotation_pitch_rate" option, which is + 3.33 degrees/second when "rotation" is True. Otherwise, the angle of attack is fixed for the + duration of the phase. + """ + + def __init__(self, **kwargs): + super().__init__(**kwargs) + + self._mu = MU_TAKEOFF + + def initialize(self): + self.options.declare('num_nodes', types=int) + + self.options.declare( + 'ground_roll', + types=bool, + default=False, + desc='True if the aircraft is confined to the ground. Angle of attack is fixed and ' + 'removed as an input.', + ) + + self.options.declare( + 'rotation', + types=bool, + default=False, + desc='True if the aircraft is pitching up, but the rear wheels are still on the ' + 'ground.', + ) + + # TODO: Make this a hiearchy variable with a default. + self.options.declare( + 'rotation_pitch_rate', + types=float, + default=3.33, + desc='Pitch rate during rotation in degrees/second.', + ) + + def setup(self): + nn = self.options['num_nodes'] + ground_roll = self.options['ground_roll'] + + self.add_input(Dynamic.Vehicle.MASS, val=np.ones(nn), desc='aircraft mass', units='lbm') + self.add_input( + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + val=np.ones(nn), + desc=Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + units='lbf', + ) + self.add_input( + Dynamic.Vehicle.LIFT, + val=np.ones(nn), + desc=Dynamic.Vehicle.LIFT, + units='lbf', + ) + self.add_input( + Dynamic.Vehicle.DRAG, + val=np.ones(nn), + desc=Dynamic.Vehicle.DRAG, + units='lbf', + ) + self.add_input( + Dynamic.Mission.VELOCITY, + val=np.ones(nn), + desc='true air speed', + units='ft/s', + ) + self.add_input( + Dynamic.Mission.FLIGHT_PATH_ANGLE, + val=np.ones(nn), + desc='flight path angle', + units='rad', + ) + + add_aviary_input(self, Aircraft.Wing.INCIDENCE, val=0) + + self.add_output( + Dynamic.Mission.VELOCITY_RATE, + val=np.ones(nn), + desc='TAS rate', + units='ft/s**2', + ) + + self.add_output( + Dynamic.Mission.ALTITUDE_RATE, + val=np.ones(nn), + desc='altitude rate', + units='ft/s', + ) + + self.add_output( + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, + val=np.ones(nn), + desc='flight path angle rate', + units='rad/s', + ) + + if not ground_roll: + self.add_input( + Dynamic.Vehicle.ANGLE_OF_ATTACK, + val=np.ones(nn), + desc='angle of attack', + units='deg', + ) + + self.add_output( + Dynamic.Mission.DISTANCE_RATE, + val=np.ones(nn), + desc='distance rate', + units='ft/s', + ) + + self.add_output('normal_force', val=np.ones(nn), desc='normal forces', units='lbf') + self.add_output('fuselage_pitch', val=np.ones(nn), desc='fuselage pitch angle', units='deg') + self.add_output('load_factor', val=np.ones(nn), desc='load factor', units='unitless') + self.add_output( + 'angle_of_attack_rate', val=np.ones(nn), desc='angle of attack rate', units='deg/s' + ) + + def setup_partials(self): + ground_roll = self.options['ground_roll'] + + arange = np.arange(self.options['num_nodes'], dtype=int) + + self.declare_partials( + 'load_factor', + [ + Dynamic.Vehicle.LIFT, + Dynamic.Vehicle.MASS, + Dynamic.Mission.FLIGHT_PATH_ANGLE, + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + ], + rows=arange, + cols=arange, + ) + self.declare_partials('load_factor', [Aircraft.Wing.INCIDENCE]) + + self.declare_partials( + Dynamic.Mission.VELOCITY_RATE, + [ + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + Dynamic.Vehicle.DRAG, + Dynamic.Vehicle.MASS, + Dynamic.Mission.FLIGHT_PATH_ANGLE, + Dynamic.Vehicle.LIFT, + ], + rows=arange, + cols=arange, + ) + self.declare_partials(Dynamic.Mission.VELOCITY_RATE, [Aircraft.Wing.INCIDENCE]) + + self.declare_partials( + Dynamic.Mission.ALTITUDE_RATE, + [Dynamic.Mission.VELOCITY, Dynamic.Mission.FLIGHT_PATH_ANGLE], + rows=arange, + cols=arange, + ) + + self.declare_partials( + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, + [ + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + Dynamic.Vehicle.LIFT, + Dynamic.Vehicle.MASS, + Dynamic.Mission.FLIGHT_PATH_ANGLE, + Dynamic.Mission.VELOCITY, + ], + rows=arange, + cols=arange, + ) + self.declare_partials(Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, [Aircraft.Wing.INCIDENCE]) + + if not ground_roll: + self.declare_partials( + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, + Dynamic.Vehicle.ANGLE_OF_ATTACK, + rows=arange, + cols=arange, + ) + + self.declare_partials( + 'normal_force', + Dynamic.Vehicle.ANGLE_OF_ATTACK, + rows=arange, + cols=arange, + ) + self.declare_partials( + 'fuselage_pitch', + Dynamic.Vehicle.ANGLE_OF_ATTACK, + rows=arange, + cols=arange, + ) + + self.declare_partials( + 'load_factor', + Dynamic.Vehicle.ANGLE_OF_ATTACK, + rows=arange, + cols=arange, + ) + self.declare_partials( + Dynamic.Mission.VELOCITY_RATE, + Dynamic.Vehicle.ANGLE_OF_ATTACK, + rows=arange, + cols=arange, + ) + + self.declare_partials( + Dynamic.Mission.DISTANCE_RATE, + [Dynamic.Mission.VELOCITY, Dynamic.Mission.FLIGHT_PATH_ANGLE], + rows=arange, + cols=arange, + ) + self.declare_partials( + 'normal_force', + [ + Dynamic.Vehicle.MASS, + Dynamic.Vehicle.LIFT, + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + ], + rows=arange, + cols=arange, + ) + self.declare_partials('normal_force', [Aircraft.Wing.INCIDENCE]) + + self.declare_partials( + 'fuselage_pitch', + Dynamic.Mission.FLIGHT_PATH_ANGLE, + rows=arange, + cols=arange, + val=RAD2DEG, + ) + self.declare_partials('fuselage_pitch', [Aircraft.Wing.INCIDENCE]) + + self.declare_partials('angle_of_attack_rate', ['*'], 0.0) + + def compute(self, inputs, outputs): + ground_roll = self.options['ground_roll'] + rotation = self.options['rotation'] + + weight = inputs[Dynamic.Vehicle.MASS] * GRAV_ENGLISH_LBM + thrust = inputs[Dynamic.Vehicle.Propulsion.THRUST_TOTAL] + incremented_lift = inputs[Dynamic.Vehicle.LIFT] + incremented_drag = inputs[Dynamic.Vehicle.DRAG] + TAS = inputs[Dynamic.Mission.VELOCITY] + gamma = inputs[Dynamic.Mission.FLIGHT_PATH_ANGLE] + i_wing = inputs[Aircraft.Wing.INCIDENCE] + + if ground_roll: + alpha = inputs[Aircraft.Wing.INCIDENCE] + else: + alpha = inputs[Dynamic.Vehicle.ANGLE_OF_ATTACK] + + thrust_along_flightpath = thrust * np.cos((alpha - i_wing) * DEG2RAD) + thrust_across_flightpath = thrust * np.sin((alpha - i_wing) * DEG2RAD) + + if ground_roll or rotation: + mu = MU_TAKEOFF + normal_force = weight - incremented_lift - thrust_across_flightpath + normal_force[normal_force < 0] = 0.0 + outputs['normal_force'] = normal_force + + else: + mu = 0.0 + normal_force = 0.0 + outputs['normal_force'][:] = normal_force + + cos_gamma = np.cos(gamma) + sin_gamma = np.sin(gamma) + + outputs[Dynamic.Mission.VELOCITY_RATE] = ( + (thrust_along_flightpath - incremented_drag - weight * sin_gamma - mu * normal_force) + * GRAV_ENGLISH_GASP + / weight + ) + + if ground_roll or rotation: + outputs[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE][:] = 0.0 + else: + outputs[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE] = ( + (thrust_across_flightpath + incremented_lift - weight * cos_gamma) + * GRAV_ENGLISH_GASP + / (TAS * weight) + ) + + outputs[Dynamic.Mission.ALTITUDE_RATE] = TAS * sin_gamma + outputs[Dynamic.Mission.DISTANCE_RATE] = TAS * cos_gamma + + outputs['fuselage_pitch'] = gamma * RAD2DEG - i_wing + alpha + + load_factor = (incremented_lift + thrust_across_flightpath) / (weight * cos_gamma) + outputs['load_factor'] = load_factor + + if rotation: + outputs['angle_of_attack_rate'][:] = self.options['rotation_pitch_rate'] + else: + outputs['angle_of_attack_rate'][:] = 0.0 + + def compute_partials(self, inputs, J): + ground_roll = self.options['ground_roll'] + rotation = self.options['rotation'] + nn = self.options['num_nodes'] + + mu = MU_TAKEOFF if (ground_roll or rotation) else 0.0 + + weight = inputs[Dynamic.Vehicle.MASS] * GRAV_ENGLISH_LBM + thrust = inputs[Dynamic.Vehicle.Propulsion.THRUST_TOTAL] + incremented_lift = inputs[Dynamic.Vehicle.LIFT] + incremented_drag = inputs[Dynamic.Vehicle.DRAG] + TAS = inputs[Dynamic.Mission.VELOCITY] + gamma = inputs[Dynamic.Mission.FLIGHT_PATH_ANGLE] + i_wing = inputs[Aircraft.Wing.INCIDENCE] + + if ground_roll: + alpha = i_wing + else: + alpha = inputs[Dynamic.Vehicle.ANGLE_OF_ATTACK] + + cos_alpha = np.cos((alpha - i_wing) * DEG2RAD) + sin_alpha = np.sin((alpha - i_wing) * DEG2RAD) + cos_gamma = np.cos(gamma) + sin_gamma = np.sin(gamma) + + thrust_along_flightpath = thrust * cos_alpha + thrust_across_flightpath = thrust * sin_alpha + + dTAlF_dThrust = cos_alpha + dTAlF_dAlpha = -thrust * sin_alpha * DEG2RAD + dTAlF_dIwing = thrust * sin_alpha * DEG2RAD + + dTAcF_dThrust = sin_alpha + dTAcF_dAlpha = thrust * cos_alpha * DEG2RAD + dTAcF_dIwing = -thrust * cos_alpha * DEG2RAD + + J['load_factor', Dynamic.Vehicle.LIFT] = 1 / (weight * cos_gamma) + J['load_factor', Dynamic.Vehicle.MASS] = ( + -(incremented_lift + thrust_across_flightpath) + / (weight**2 * cos_gamma) + * GRAV_ENGLISH_LBM + ) + J['load_factor', Dynamic.Mission.FLIGHT_PATH_ANGLE] = ( + -(incremented_lift + thrust_across_flightpath) + / (weight * (cos_gamma) ** 2) + * (-sin_gamma) + ) + J['load_factor', Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = dTAcF_dThrust / ( + weight * cos_gamma + ) + + if ground_roll or rotation: + normal_force = weight - incremented_lift - thrust_across_flightpath + else: + normal_force = np.zeros(nn, dtype=TAS.dtype) + + idx = np.where(normal_force < 0) + normal_force[idx] = 0.0 + + dNF_dWeight = np.ones(nn, dtype=TAS.dtype) + dNF_dWeight[idx] = 0 + + dNF_dLift = -np.ones(nn, dtype=TAS.dtype) + dNF_dLift[idx] = 0 + + dNF_dThrust = -np.ones(nn, dtype=TAS.dtype) * dTAcF_dThrust + dNF_dThrust[idx] = 0 + + dNF_dIwing = -np.ones(nn, dtype=TAS.dtype) * dTAcF_dIwing + dNF_dIwing[idx] = 0 + + J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = ( + (dTAlF_dThrust - mu * dNF_dThrust) * GRAV_ENGLISH_GASP / weight + ) + + J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.DRAG] = -GRAV_ENGLISH_GASP / weight + J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.MASS] = ( + GRAV_ENGLISH_GASP + * GRAV_ENGLISH_LBM + * ( + weight * (-sin_gamma - mu * dNF_dWeight) + - ( + thrust_along_flightpath + - incremented_drag + - weight * sin_gamma + - mu * normal_force + ) + ) + / weight**2 + ) + J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = ( + -cos_gamma * GRAV_ENGLISH_GASP + ) + J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.LIFT] = ( + GRAV_ENGLISH_GASP * (-mu * dNF_dLift) / weight + ) + + J[Dynamic.Mission.ALTITUDE_RATE, Dynamic.Mission.VELOCITY] = sin_gamma + J[Dynamic.Mission.ALTITUDE_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = TAS * cos_gamma + + if not (ground_roll or rotation): + # OF flight path angle rate + J[ + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + ] = dTAcF_dThrust * GRAV_ENGLISH_GASP / (TAS * weight) + J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Vehicle.ANGLE_OF_ATTACK] = ( + dTAcF_dAlpha * GRAV_ENGLISH_GASP / (TAS * weight) + ) + J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Aircraft.Wing.INCIDENCE] = ( + dTAcF_dIwing * GRAV_ENGLISH_GASP / (TAS * weight) + ) + J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Vehicle.LIFT] = GRAV_ENGLISH_GASP / ( + TAS * weight + ) + J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Vehicle.MASS] = ( + (GRAV_ENGLISH_GASP / TAS) + * GRAV_ENGLISH_LBM + * (-thrust_across_flightpath / weight**2 - incremented_lift / weight**2) + ) + J[ + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, + Dynamic.Mission.FLIGHT_PATH_ANGLE, + ] = weight * sin_gamma * GRAV_ENGLISH_GASP / (TAS * weight) + J[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, Dynamic.Mission.VELOCITY] = -( + (thrust_across_flightpath + incremented_lift - weight * cos_gamma) + * GRAV_ENGLISH_GASP + / (TAS**2 * weight) + ) + + if not ground_roll: + # WRT incidence angle + J[Dynamic.Mission.VELOCITY_RATE, Aircraft.Wing.INCIDENCE] = ( + (dTAlF_dIwing - mu * dNF_dIwing) * GRAV_ENGLISH_GASP / weight + ) + if rotation: + J['normal_force', Aircraft.Wing.INCIDENCE] = dNF_dIwing + J['fuselage_pitch', Aircraft.Wing.INCIDENCE] = -1 + J['load_factor', Aircraft.Wing.INCIDENCE] = dTAcF_dIwing / (weight * cos_gamma) + + # WRT angle of attack + dNF_dAlpha = -np.ones(nn, dtype=TAS.dtype) * dTAcF_dAlpha + dNF_dAlpha[idx] = 0 + + J[Dynamic.Mission.VELOCITY_RATE, Dynamic.Vehicle.ANGLE_OF_ATTACK] = ( + (dTAlF_dAlpha - mu * dNF_dAlpha) * GRAV_ENGLISH_GASP / weight + ) + if rotation: + J['normal_force', Dynamic.Vehicle.ANGLE_OF_ATTACK] = dNF_dAlpha + J['fuselage_pitch', Dynamic.Vehicle.ANGLE_OF_ATTACK] = 1 + J['load_factor', Dynamic.Vehicle.ANGLE_OF_ATTACK] = dTAcF_dAlpha / (weight * cos_gamma) + + J[Dynamic.Mission.DISTANCE_RATE, Dynamic.Mission.VELOCITY] = cos_gamma + J[Dynamic.Mission.DISTANCE_RATE, Dynamic.Mission.FLIGHT_PATH_ANGLE] = -TAS * sin_gamma + + if ground_roll or rotation: + J['normal_force', Dynamic.Vehicle.MASS] = dNF_dWeight * GRAV_ENGLISH_LBM + J['normal_force', Dynamic.Vehicle.LIFT] = dNF_dLift + J['normal_force', Dynamic.Vehicle.Propulsion.THRUST_TOTAL] = dNF_dThrust diff --git a/aviary/mission/two_dof/ode/takeoff_ode.py b/aviary/mission/two_dof/ode/takeoff_ode.py new file mode 100644 index 0000000000..078e889a6a --- /dev/null +++ b/aviary/mission/two_dof/ode/takeoff_ode.py @@ -0,0 +1,270 @@ +import numpy as np +import openmdao.api as om + +from aviary.mission.two_dof.ode.params import ParamPort +from aviary.mission.two_dof.ode.takeoff_eom import TakeoffEOM +from aviary.mission.two_dof.ode.two_dof_ode import TwoDOFODE +from aviary.subsystems.aerodynamics.aerodynamics_builder import AerodynamicsBuilder +from aviary.subsystems.mass.mass_to_weight import MassToWeight +from aviary.subsystems.propulsion.propulsion_builder import PropulsionBuilder +from aviary.variable_info.enums import AlphaModes, SpeedType +from aviary.variable_info.variables import Aircraft, Dynamic, Mission + + +class TakeOffODE(TwoDOFODE): + """ODE using 2D aircraft equations of motion with states distance, alt, TAS, and gamma. + + Control is managed via angle-of-attack (alpha). + This ODE is used for two-dof takeoff phases, and supports ground roll, rotation, and + ascending flight. + """ + + def initialize(self): + super().initialize() + + self.options.declare('alpha_mode', default=AlphaModes.DEFAULT, types=AlphaModes) + + self.options.declare( + 'input_speed_type', + default=SpeedType.TAS, + types=SpeedType, + desc='Whether the speed is given as a equivalent airspeed, true airspeed, or Mach number', + ) + + self.options.declare( + 'ground_roll', + types=bool, + default=False, + desc='True if the aircraft is confined to the ground. Angle of attack is fixed and ' + 'removed as an input.', + ) + + self.options.declare( + 'rotation', + types=bool, + default=False, + desc='True if the aircraft is pitching up, but the rear wheels are still on the ' + 'ground.', + ) + + self.options.declare( + 'clean', + types=bool, + default=False, + desc='If true then no flaps or gear are included. Useful for high-speed flight phases.', + ) + + def setup(self): + nn = self.options['num_nodes'] + aviary_options = self.options['aviary_options'] + alpha_mode = self.options['alpha_mode'] + input_speed_type = self.options['input_speed_type'] + ground_roll = self.options['ground_roll'] + rotation = self.options['rotation'] + + EOM_inputs = [ + Dynamic.Vehicle.MASS, + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, + Dynamic.Vehicle.LIFT, + Dynamic.Vehicle.DRAG, + Dynamic.Mission.VELOCITY, + Dynamic.Mission.FLIGHT_PATH_ANGLE, + ] + ['aircraft:*'] + if not self.options['ground_roll']: + EOM_inputs.append(Dynamic.Vehicle.ANGLE_OF_ATTACK) + + subsystems = self.options['subsystems'] + subsystem_options = self.options['subsystem_options'] + + # TODO: paramport + self.add_subsystem('params', ParamPort(), promotes=['*']) + + self.add_atmosphere(input_speed_type=input_speed_type) + + if ground_roll: + # Angle of attack equals the incidence angle. + self.add_subsystem( + 'init_alpha', + om.ExecComp( + 'alpha = i_wing', + i_wing={'units': 'deg', 'val': 1.1}, + alpha={'units': 'deg', 'val': 1.1 * np.ones(nn)}, + ), + promotes=[ + ('i_wing', Aircraft.Wing.INCIDENCE), + ('alpha', Dynamic.Vehicle.ANGLE_OF_ATTACK), + ], + ) + + if alpha_mode is AlphaModes.DEFAULT: + # alpha as input + pass + else: + if alpha_mode is AlphaModes.REQUIRED_LIFT: + self.add_subsystem( + 'calc_weight', + MassToWeight(num_nodes=nn), + promotes_inputs=[('mass', Dynamic.Vehicle.MASS)], + promotes_outputs=['weight'], + ) + self.add_subsystem( + 'calc_lift', + om.ExecComp( + 'required_lift = weight*cos(alpha + gamma) - thrust*sin(i_wing)', + required_lift={'val': 0, 'units': 'lbf'}, + weight={'val': 0, 'units': 'lbf'}, + thrust={'val': 0, 'units': 'lbf'}, + alpha={'val': 0, 'units': 'rad'}, + gamma={'val': 0, 'units': 'rad'}, + i_wing={'val': 0, 'units': 'rad'}, + ), + promotes_inputs=[ + 'weight', + ('thrust', Dynamic.Vehicle.Propulsion.THRUST_TOTAL), + ('alpha', Dynamic.Vehicle.ANGLE_OF_ATTACK), + ('gamma', Dynamic.Mission.FLIGHT_PATH_ANGLE), + ('i_wing', Aircraft.Wing.INCIDENCE), + ], + promotes_outputs=['required_lift'], + ) + self.add_alpha_control( + alpha_mode=alpha_mode, + target_load_factor=1, + atol=1e-6, + rtol=1e-12, + num_nodes=nn, + print_level=0, + ) + + for subsystem in subsystems: + name = subsystem.name + kwargs = {} + + if isinstance(subsystem, AerodynamicsBuilder): + kwargs = {'method': 'low_speed'} + if self.options['clean']: + kwargs['method'] = 'cruise' + kwargs['output_alpha'] = False + + if not (ground_roll or rotation): + kwargs['retract_gear'] = True + kwargs['retract_flaps'] = True + + if name in subsystem_options: + kwargs.update(subsystem_options[name]) + + system = subsystem.build_mission(num_nodes=nn, aviary_inputs=aviary_options, **kwargs) + + if system is not None: + if isinstance(subsystem, PropulsionBuilder): + self.add_subsystem( + name, + system, + promotes_inputs=subsystem.mission_inputs(**kwargs), + promotes_outputs=subsystem.mission_outputs(**kwargs), + ) + else: + self.add_subsystem( + name, + system, + promotes_inputs=subsystem.mission_inputs(**kwargs), + promotes_outputs=subsystem.mission_outputs(**kwargs), + ) + + if ground_roll: + eom_name = 'groundroll_eom' + elif rotation: + eom_name = 'rotation_eom' + else: + eom_name = 'ascent_eom' + + self.add_subsystem( + eom_name, + TakeoffEOM( + num_nodes=nn, + ground_roll=ground_roll, + rotation=rotation, + ), + promotes_inputs=EOM_inputs, + promotes_outputs=[ + Dynamic.Mission.VELOCITY_RATE, + Dynamic.Mission.DISTANCE_RATE, + Dynamic.Mission.ALTITUDE_RATE, + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, + 'normal_force', + 'fuselage_pitch', + 'load_factor', + 'angle_of_attack_rate', + ], + ) + + if ground_roll: + self.add_subsystem( + 'exec', + om.ExecComp( + 'over_a = velocity / velocity_rate', + velocity_rate={'units': 'kn/s', 'val': np.ones(nn)}, + velocity={'units': 'kn', 'val': np.ones(nn)}, + over_a={'units': 's', 'val': np.ones(nn)}, + has_diag_partials=True, + ), + promotes=['*'], + ) + + self.add_subsystem( + 'exec2', + om.ExecComp( + 'dt_dv = 1 / velocity_rate', + velocity_rate={'units': 'kn/s', 'val': np.ones(nn)}, + dt_dv={'units': 's/kn', 'val': np.ones(nn)}, + has_diag_partials=True, + ), + promotes=['*'], + ) + + self.add_subsystem( + 'exec3', + om.ExecComp( + 'dmass_dv = mass_rate * dt_dv', + mass_rate={'units': 'lbm/s', 'val': np.ones(nn)}, + dt_dv={'units': 's/kn', 'val': np.ones(nn)}, + dmass_dv={'units': 'lbm/kn', 'val': np.ones(nn)}, + has_diag_partials=True, + ), + promotes_outputs=[ + 'dmass_dv', + ], + promotes_inputs=[ + ('mass_rate', Dynamic.Vehicle.Propulsion.FUEL_FLOW_RATE_NEGATIVE_TOTAL), + 'dt_dv', + ], + ) + + if not (ground_roll or rotation): + self.add_excess_rate_comps(nn) + + ParamPort.set_default_vals(self) + + self.set_input_defaults(Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(nn), units='rad') + self.set_input_defaults(Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(nn), units='deg') + self.set_input_defaults(Dynamic.Mission.VELOCITY, val=np.zeros(nn), units='kn') + self.set_input_defaults(Dynamic.Mission.ALTITUDE, val=np.zeros(nn), units='ft') + self.set_input_defaults(Dynamic.Vehicle.MASS, val=np.zeros(nn), units='lbm') + + # TODO: Some of these are backdoor defaults. + if not self.options['clean']: + self.set_input_defaults('t_init_flaps', val=47.5, units='s') + self.set_input_defaults('t_init_gear', val=37.3, units='s') + self.set_input_defaults('t_curr', val=np.zeros(nn), units='s') + if ground_roll or rotation: + self.set_input_defaults('aero_ramps.flap_factor:final_val', val=1.0) + self.set_input_defaults('aero_ramps.gear_factor:final_val', val=1.0) + else: + self.set_input_defaults('aero_ramps.flap_factor:final_val', val=0.0) + self.set_input_defaults('aero_ramps.gear_factor:final_val', val=0.0) + self.set_input_defaults('aero_ramps.flap_factor:initial_val', val=1.0) + self.set_input_defaults('aero_ramps.gear_factor:initial_val', val=1.0) + + if ground_roll: + self.set_input_defaults(Dynamic.Mission.VELOCITY_RATE, val=np.zeros(nn), units='kn/s') + self.set_input_defaults(Aircraft.Wing.INCIDENCE, val=1.0, units='deg') diff --git a/aviary/mission/two_dof/ode/test/test_ascent_eom.py b/aviary/mission/two_dof/ode/test/test_ascent_eom.py deleted file mode 100644 index 8a290f83c9..0000000000 --- a/aviary/mission/two_dof/ode/test/test_ascent_eom.py +++ /dev/null @@ -1,91 +0,0 @@ -import unittest - -import numpy as np -import openmdao.api as om -from openmdao.utils.assert_utils import assert_check_partials, assert_near_equal - -from aviary.mission.two_dof.ode.ascent_eom import AscentEOM -from aviary.variable_info.variables import Aircraft, Dynamic - - -class AscentEOMTestCase(unittest.TestCase): - def setUp(self): - self.prob = om.Problem() - self.prob.model.add_subsystem('group', AscentEOM(num_nodes=2), promotes=['*']) - self.prob.model.set_input_defaults( - Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm' - ) - self.prob.model.set_input_defaults( - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' - ) - self.prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') - self.prob.model.set_input_defaults( - Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf' - ) - self.prob.model.set_input_defaults( - Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s' - ) - self.prob.model.set_input_defaults( - Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' - ) - self.prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') - self.prob.model.set_input_defaults( - Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg' - ) - - self.prob.setup(check=False, force_alloc_complex=True) - - def test_case1(self): - tol = 1e-6 - self.prob.run_model() - - assert_near_equal( - self.prob[Dynamic.Mission.VELOCITY_RATE], - np.array([2.202965, 2.202965]), - tol, - ) - assert_near_equal( - self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], - np.array([-3.216328, -3.216328]), - tol, - ) - - partial_data = self.prob.check_partials(out_stream=None, method='cs') - assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) - - -class AscentEOMTestCase2(unittest.TestCase): - """Test mass-weight conversion.""" - - def setUp(self): - import aviary.mission.two_dof.ode.ascent_eom as ascent - - ascent.GRAV_ENGLISH_LBM = 1.1 - - def tearDown(self): - import aviary.mission.two_dof.ode.ascent_eom as ascent - - ascent.GRAV_ENGLISH_LBM = 1.0 - - def test_case1(self): - prob = om.Problem() - prob.model.add_subsystem('group', AscentEOM(num_nodes=2), promotes=['*']) - prob.model.set_input_defaults(Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm') - prob.model.set_input_defaults( - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' - ) - prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') - prob.model.set_input_defaults(Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf') - prob.model.set_input_defaults(Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s') - prob.model.set_input_defaults( - Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' - ) - prob.model.set_input_defaults(Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg') - prob.setup(check=False, force_alloc_complex=True) - - partial_data = prob.check_partials(out_stream=None, method='cs') - assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) - - -if __name__ == '__main__': - unittest.main() diff --git a/aviary/mission/two_dof/ode/test/test_ascent_ode.py b/aviary/mission/two_dof/ode/test/test_ascent_ode.py deleted file mode 100644 index e0fd34f916..0000000000 --- a/aviary/mission/two_dof/ode/test/test_ascent_ode.py +++ /dev/null @@ -1,78 +0,0 @@ -import unittest - -import numpy as np -import openmdao.api as om -from openmdao.utils.assert_utils import assert_check_partials, assert_near_equal - -from aviary.mission.two_dof.ode.ascent_ode import AscentODE -from aviary.mission.two_dof.ode.params import set_params_for_unit_tests -from aviary.subsystems.propulsion.utils import build_engine_deck -from aviary.utils.aviary_values import AviaryValues -from aviary.utils.test_utils.default_subsystems import get_default_mission_subsystems -from aviary.variable_info.functions import setup_model_options -from aviary.variable_info.options import get_option_defaults -from aviary.variable_info.variables import Aircraft, Dynamic - - -class AscentODETestCase(unittest.TestCase): - def setUp(self): - self.prob = om.Problem() - - aviary_options = get_option_defaults() - aviary_options.set_val(Aircraft.Engine.GLOBAL_THROTTLE, True) - default_mission_subsystems = get_default_mission_subsystems( - 'GASP', [build_engine_deck(aviary_options)] - ) - - self.prob.model = AscentODE( - num_nodes=2, aviary_options=aviary_options, subsystems=default_mission_subsystems - ) - - setup_model_options( - self.prob, AviaryValues({Aircraft.Engine.NUM_ENGINES: ([2], 'unitless')}) - ) - - def test_ascent_partials(self): - # Test partial derivatives - self.prob.setup(check=False, force_alloc_complex=True) - - self.prob.set_val(Dynamic.Mission.VELOCITY, [100, 100], units='kn') - self.prob.set_val('t_curr', [1, 2], units='s') - self.prob.set_val('interference_independent_of_shielded_area', 1.89927266) - self.prob.set_val('drag_loss_due_to_shielded_wing_area', 68.02065834) - self.prob.set_val(Aircraft.Wing.FORM_FACTOR, 1.25) - self.prob.set_val(Aircraft.VerticalTail.FORM_FACTOR, 1.25) - self.prob.set_val(Aircraft.HorizontalTail.FORM_FACTOR, 1.25) - - set_params_for_unit_tests(self.prob) - - self.prob.run_model() - - tol = tol = 1e-6 - assert_near_equal( - self.prob[Dynamic.Mission.VELOCITY_RATE], - np.array([642156.99315828, 642156.99315828]), - tol, - ) - assert_near_equal( - self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], - np.array([2260.37849562, 2260.37849562]), - tol, - ) - assert_near_equal(self.prob[Dynamic.Mission.ALTITUDE_RATE], np.array([0.0, 0.0]), tol) - assert_near_equal( - self.prob[Dynamic.Mission.DISTANCE_RATE], np.array([168.781, 168.781]), tol - ) - assert_near_equal(self.prob['angle_of_attack_rate'], np.array([0.0, 0.0]), tol) - assert_near_equal(self.prob['normal_force'], np.array([0.0, 0.0]), tol) - assert_near_equal(self.prob['fuselage_pitch'], np.array([0.0, 0.0]), tol) - assert_near_equal(self.prob['load_factor'], np.array([11849.10281268, 11849.10281268]), tol) - - partial_data = self.prob.check_partials( - out_stream=None, method='cs', excludes=['*params*', '*aero*'] - ) - assert_check_partials(partial_data, atol=1e-8, rtol=1e-8) - - -if __name__ == '__main__': - unittest.main() diff --git a/aviary/mission/two_dof/ode/test/test_rotation_eom.py b/aviary/mission/two_dof/ode/test/test_rotation_eom.py deleted file mode 100644 index a4d6ff36df..0000000000 --- a/aviary/mission/two_dof/ode/test/test_rotation_eom.py +++ /dev/null @@ -1,94 +0,0 @@ -import unittest - -import numpy as np -import openmdao.api as om -from openmdao.utils.assert_utils import assert_check_partials, assert_near_equal - -from aviary.mission.two_dof.ode.rotation_eom import RotationEOM -from aviary.variable_info.variables import Aircraft, Dynamic - - -class RotationEOMTestCase(unittest.TestCase): - def setUp(self): - self.prob = om.Problem() - self.prob.model.add_subsystem('group', RotationEOM(num_nodes=2), promotes=['*']) - self.prob.model.set_input_defaults( - Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm' - ) - self.prob.model.set_input_defaults( - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' - ) - self.prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') - self.prob.model.set_input_defaults( - Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf' - ) - self.prob.model.set_input_defaults( - Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s' - ) - self.prob.model.set_input_defaults( - Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' - ) - self.prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') - self.prob.model.set_input_defaults( - Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg' - ) - - self.prob.setup(check=False, force_alloc_complex=True) - - def test_case1(self): - tol = 1e-6 - self.prob.run_model() - - assert_near_equal( - self.prob[Dynamic.Mission.VELOCITY_RATE], - np.array([1.5597, 1.5597]), - tol, - ) - assert_near_equal( - self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], np.array([0.0, 0.0]), tol - ) - assert_near_equal(self.prob[Dynamic.Mission.ALTITUDE_RATE], np.array([0.0, 0.0]), tol) - assert_near_equal(self.prob[Dynamic.Mission.DISTANCE_RATE], np.array([10.0, 10.0]), tol) - assert_near_equal(self.prob['normal_force'], np.array([175200.0, 175200.0]), tol) - assert_near_equal(self.prob['fuselage_pitch'], np.array([0.0, 0.0]), tol) - - partial_data = self.prob.check_partials(out_stream=None, method='cs') - assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) - - -class RotationEOMTestCase2(unittest.TestCase): - """Test mass-weight conversion.""" - - def setUp(self): - import aviary.mission.two_dof.ode.rotation_eom as rotation - - rotation.GRAV_ENGLISH_LBM = 1.1 - - def tearDown(self): - import aviary.mission.two_dof.ode.rotation_eom as rotation - - rotation.GRAV_ENGLISH_LBM = 1.0 - - def test_case1(self): - prob = om.Problem() - prob.model.add_subsystem('group', RotationEOM(num_nodes=2), promotes=['*']) - prob.model.set_input_defaults(Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm') - prob.model.set_input_defaults( - Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' - ) - prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') - prob.model.set_input_defaults(Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf') - prob.model.set_input_defaults(Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s') - prob.model.set_input_defaults( - Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' - ) - prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') - prob.model.set_input_defaults(Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg') - prob.setup(check=False, force_alloc_complex=True) - - partial_data = prob.check_partials(out_stream=None, method='cs') - assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) - - -if __name__ == '__main__': - unittest.main() diff --git a/aviary/mission/two_dof/ode/test/test_rotation_ode.py b/aviary/mission/two_dof/ode/test/test_rotation_ode.py deleted file mode 100644 index 01fbcc76cf..0000000000 --- a/aviary/mission/two_dof/ode/test/test_rotation_ode.py +++ /dev/null @@ -1,79 +0,0 @@ -import unittest - -import numpy as np -import openmdao.api as om -from openmdao.utils.assert_utils import assert_check_partials, assert_near_equal - -from aviary.mission.two_dof.ode.params import set_params_for_unit_tests -from aviary.mission.two_dof.ode.rotation_ode import RotationODE -from aviary.subsystems.propulsion.utils import build_engine_deck -from aviary.utils.test_utils.default_subsystems import get_default_mission_subsystems -from aviary.variable_info.functions import setup_model_options -from aviary.variable_info.options import get_option_defaults -from aviary.variable_info.variables import Aircraft, Dynamic - - -class RotationODETestCase(unittest.TestCase): - """Test 2-degree of freedom rotation ODE.""" - - def setUp(self): - self.prob = om.Problem() - - aviary_options = get_option_defaults() - aviary_options.set_val(Aircraft.Engine.GLOBAL_THROTTLE, True) - default_mission_subsystems = get_default_mission_subsystems( - 'GASP', [build_engine_deck(aviary_options)] - ) - - self.prob.model = RotationODE( - num_nodes=2, - aviary_options=get_option_defaults(), - subsystems=default_mission_subsystems, - ) - setup_model_options(self.prob, aviary_options) - - def test_rotation_partials(self): - # Check partial derivatives - self.prob.setup(check=False, force_alloc_complex=True) - - self.prob.set_val(Aircraft.Wing.INCIDENCE, 1.5, units='deg') - self.prob.set_val(Dynamic.Vehicle.MASS, [100000, 100000], units='lbm') - self.prob.set_val(Dynamic.Vehicle.ANGLE_OF_ATTACK, [1.5, 1.5], units='deg') - self.prob.set_val(Dynamic.Mission.VELOCITY, [100, 100], units='kn') - self.prob.set_val('t_curr', [1, 2], units='s') - self.prob.set_val('interference_independent_of_shielded_area', 1.89927266) - self.prob.set_val('drag_loss_due_to_shielded_wing_area', 68.02065834) - self.prob.set_val(Aircraft.Wing.FORM_FACTOR, 1.25) - self.prob.set_val(Aircraft.VerticalTail.FORM_FACTOR, 1.25) - self.prob.set_val(Aircraft.HorizontalTail.FORM_FACTOR, 1.25) - - set_params_for_unit_tests(self.prob) - - self.prob.run_model() - - tol = 1e-6 - assert_near_equal( - self.prob[Dynamic.Mission.VELOCITY_RATE], - np.array([13.68875852, 13.68875852]), - tol, - ) - assert_near_equal( - self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], np.array([0.0, 0.0]), tol - ) - assert_near_equal(self.prob[Dynamic.Mission.ALTITUDE_RATE], np.array([0.0, 0.0]), tol) - assert_near_equal( - self.prob[Dynamic.Mission.DISTANCE_RATE], np.array([168.781, 168.781]), tol - ) - assert_near_equal( - self.prob['normal_force'], np.array([66936.59676831, 66936.59676831]), tol - ) - assert_near_equal(self.prob['fuselage_pitch'], np.array([0.0, 0.0]), tol) - - partial_data = self.prob.check_partials( - out_stream=None, method='cs', excludes=['*params*', '*aero*'] - ) - assert_check_partials(partial_data, atol=1e-8, rtol=1e-8) - - -if __name__ == '__main__': - unittest.main() diff --git a/aviary/mission/two_dof/ode/test/test_takeoff_eom.py b/aviary/mission/two_dof/ode/test/test_takeoff_eom.py new file mode 100644 index 0000000000..d0bd280da1 --- /dev/null +++ b/aviary/mission/two_dof/ode/test/test_takeoff_eom.py @@ -0,0 +1,255 @@ +import unittest + +import numpy as np +import openmdao.api as om +from openmdao.utils.assert_utils import assert_check_partials, assert_near_equal + +from aviary.mission.two_dof.ode.takeoff_eom import TakeoffEOM +from aviary.variable_info.variables import Aircraft, Dynamic + + +class GroundrollEOMTestCase(unittest.TestCase): + def setUp(self): + self.prob = om.Problem() + self.prob.model.add_subsystem( + 'group', TakeoffEOM(num_nodes=2, ground_roll=True), promotes=['*'] + ) + self.prob.model.set_input_defaults( + Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm' + ) + self.prob.model.set_input_defaults( + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' + ) + self.prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') + self.prob.model.set_input_defaults( + Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf' + ) + self.prob.model.set_input_defaults( + Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s' + ) + self.prob.model.set_input_defaults( + Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' + ) + self.prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') + + self.prob.setup(check=False, force_alloc_complex=True) + + def test_case1(self): + tol = 1e-6 + self.prob.run_model() + + assert_near_equal( + self.prob[Dynamic.Mission.VELOCITY_RATE], + np.array([1.5597, 1.5597]), + tol, + ) + assert_near_equal( + self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], np.array([0.0, 0.0]), tol + ) + assert_near_equal(self.prob[Dynamic.Mission.ALTITUDE_RATE], np.array([0.0, 0.0]), tol) + assert_near_equal(self.prob[Dynamic.Mission.DISTANCE_RATE], np.array([10.0, 10.0]), tol) + assert_near_equal(self.prob['normal_force'], np.array([175200.0, 175200.0]), tol) + assert_near_equal(self.prob['fuselage_pitch'], np.array([0.0, 0.0]), tol) + partial_data = self.prob.check_partials(out_stream=None, method='cs') + assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) + + +class GroundrollEOMTestCase2(unittest.TestCase): + """Test mass-weight conversion.""" + + def setUp(self): + import aviary.mission.two_dof.ode.takeoff_eom as gr + + gr.GRAV_ENGLISH_LBM = 1.1 + + def tearDown(self): + import aviary.mission.two_dof.ode.takeoff_eom as gr + + gr.GRAV_ENGLISH_LBM = 1.0 + + def test_case1(self): + prob = om.Problem() + prob.model.add_subsystem('group', TakeoffEOM(num_nodes=2, ground_roll=True), promotes=['*']) + prob.model.set_input_defaults(Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm') + prob.model.set_input_defaults( + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' + ) + prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') + prob.model.set_input_defaults(Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf') + prob.model.set_input_defaults(Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s') + prob.model.set_input_defaults( + Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' + ) + prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') + + prob.setup(check=False, force_alloc_complex=True) + + partial_data = prob.check_partials(out_stream=None, method='cs') + assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) + + +class RotationEOMTestCase(unittest.TestCase): + def setUp(self): + self.prob = om.Problem() + self.prob.model.add_subsystem( + 'group', TakeoffEOM(num_nodes=2, rotation=True), promotes=['*'] + ) + self.prob.model.set_input_defaults( + Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm' + ) + self.prob.model.set_input_defaults( + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' + ) + self.prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') + self.prob.model.set_input_defaults( + Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf' + ) + self.prob.model.set_input_defaults( + Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s' + ) + self.prob.model.set_input_defaults( + Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' + ) + self.prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') + self.prob.model.set_input_defaults( + Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg' + ) + + self.prob.setup(check=False, force_alloc_complex=True) + + def test_case1(self): + tol = 1e-6 + self.prob.run_model() + + assert_near_equal( + self.prob[Dynamic.Mission.VELOCITY_RATE], + np.array([1.5597, 1.5597]), + tol, + ) + assert_near_equal( + self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], np.array([0.0, 0.0]), tol + ) + assert_near_equal(self.prob[Dynamic.Mission.ALTITUDE_RATE], np.array([0.0, 0.0]), tol) + assert_near_equal(self.prob[Dynamic.Mission.DISTANCE_RATE], np.array([10.0, 10.0]), tol) + assert_near_equal(self.prob['normal_force'], np.array([175200.0, 175200.0]), tol) + assert_near_equal(self.prob['fuselage_pitch'], np.array([0.0, 0.0]), tol) + + partial_data = self.prob.check_partials(out_stream=None, method='cs') + assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) + + +class RotationEOMTestCase2(unittest.TestCase): + """Test mass-weight conversion.""" + + def setUp(self): + import aviary.mission.two_dof.ode.takeoff_eom as rotation + + rotation.GRAV_ENGLISH_LBM = 1.1 + + def tearDown(self): + import aviary.mission.two_dof.ode.takeoff_eom as rotation + + rotation.GRAV_ENGLISH_LBM = 1.0 + + def test_case1(self): + prob = om.Problem() + prob.model.add_subsystem('group', TakeoffEOM(num_nodes=2, rotation=True), promotes=['*']) + prob.model.set_input_defaults(Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm') + prob.model.set_input_defaults( + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' + ) + prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') + prob.model.set_input_defaults(Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf') + prob.model.set_input_defaults(Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s') + prob.model.set_input_defaults( + Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' + ) + prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') + prob.model.set_input_defaults(Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg') + prob.setup(check=False, force_alloc_complex=True) + + partial_data = prob.check_partials(out_stream=None, method='cs') + assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) + + +class AscentEOMTestCase(unittest.TestCase): + def setUp(self): + self.prob = om.Problem() + self.prob.model.add_subsystem('group', TakeoffEOM(num_nodes=2), promotes=['*']) + self.prob.model.set_input_defaults( + Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm' + ) + self.prob.model.set_input_defaults( + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' + ) + self.prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') + self.prob.model.set_input_defaults( + Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf' + ) + self.prob.model.set_input_defaults( + Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s' + ) + self.prob.model.set_input_defaults( + Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' + ) + self.prob.model.set_input_defaults(Aircraft.Wing.INCIDENCE, val=0, units='deg') + self.prob.model.set_input_defaults( + Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg' + ) + + self.prob.setup(check=False, force_alloc_complex=True) + + def test_case1(self): + tol = 1e-6 + self.prob.run_model() + + assert_near_equal( + self.prob[Dynamic.Mission.VELOCITY_RATE], + np.array([2.202965, 2.202965]), + tol, + ) + assert_near_equal( + self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], + np.array([-3.216328, -3.216328]), + tol, + ) + + partial_data = self.prob.check_partials(out_stream=None, method='cs') + assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) + + +class AscentEOMTestCase2(unittest.TestCase): + """Test mass-weight conversion.""" + + def setUp(self): + import aviary.mission.two_dof.ode.takeoff_eom as ascent + + ascent.GRAV_ENGLISH_LBM = 1.1 + + def tearDown(self): + import aviary.mission.two_dof.ode.takeoff_eom as ascent + + ascent.GRAV_ENGLISH_LBM = 1.0 + + def test_case1(self): + prob = om.Problem() + prob.model.add_subsystem('group', TakeoffEOM(num_nodes=2), promotes=['*']) + prob.model.set_input_defaults(Dynamic.Vehicle.MASS, val=175400 * np.ones(2), units='lbm') + prob.model.set_input_defaults( + Dynamic.Vehicle.Propulsion.THRUST_TOTAL, val=22000 * np.ones(2), units='lbf' + ) + prob.model.set_input_defaults(Dynamic.Vehicle.LIFT, val=200 * np.ones(2), units='lbf') + prob.model.set_input_defaults(Dynamic.Vehicle.DRAG, val=10000 * np.ones(2), units='lbf') + prob.model.set_input_defaults(Dynamic.Mission.VELOCITY, val=10 * np.ones(2), units='ft/s') + prob.model.set_input_defaults( + Dynamic.Mission.FLIGHT_PATH_ANGLE, val=np.zeros(2), units='rad' + ) + prob.model.set_input_defaults(Dynamic.Vehicle.ANGLE_OF_ATTACK, val=np.zeros(2), units='deg') + prob.setup(check=False, force_alloc_complex=True) + + partial_data = prob.check_partials(out_stream=None, method='cs') + assert_check_partials(partial_data, atol=1e-12, rtol=1e-12) + + +if __name__ == '__main__': + unittest.main() diff --git a/aviary/mission/two_dof/ode/test/test_takeoff_ode.py b/aviary/mission/two_dof/ode/test/test_takeoff_ode.py new file mode 100644 index 0000000000..d88564b672 --- /dev/null +++ b/aviary/mission/two_dof/ode/test/test_takeoff_ode.py @@ -0,0 +1,201 @@ +import unittest + +import numpy as np +import openmdao.api as om +from openmdao.utils.assert_utils import assert_check_partials, assert_near_equal + +from aviary.mission.two_dof.ode.takeoff_ode import TakeOffODE +from aviary.mission.two_dof.ode.params import set_params_for_unit_tests +from aviary.subsystems.propulsion.utils import build_engine_deck +from aviary.utils.aviary_values import AviaryValues +from aviary.utils.test_utils.default_subsystems import get_default_mission_subsystems +from aviary.utils.test_utils.IO_test_util import check_prob_outputs +from aviary.variable_info.functions import setup_model_options +from aviary.variable_info.options import get_option_defaults +from aviary.variable_info.variables import Aircraft, Dynamic + + +class GroundrollODETestCase(unittest.TestCase): + """Test groundroll ODE.""" + + def setUp(self): + self.prob = om.Problem() + + aviary_options = get_option_defaults() + aviary_options.set_val(Aircraft.Engine.GLOBAL_THROTTLE, True) + default_mission_subsystems = get_default_mission_subsystems( + 'GASP', [build_engine_deck(aviary_options)] + ) + + self.prob.model = TakeOffODE( + num_nodes=2, + ground_roll=True, + aviary_options=get_option_defaults(), + subsystems=default_mission_subsystems, + ) + + setup_model_options(self.prob, aviary_options) + + def test_groundroll_partials(self): + # Check partial derivatives + self.prob.setup(check=False, force_alloc_complex=True) + + set_params_for_unit_tests(self.prob) + + self.prob.set_val('t_curr', [1, 2], units='s') + self.prob.set_val('aircraft:wing:incidence', 0, units='deg') + self.prob.set_val('interference_independent_of_shielded_area', 1.89927266) + self.prob.set_val('drag_loss_due_to_shielded_wing_area', 68.02065834) + self.prob.set_val(Aircraft.Wing.FORM_FACTOR, 1.25) + self.prob.set_val(Aircraft.VerticalTail.FORM_FACTOR, 1.25) + self.prob.set_val(Aircraft.HorizontalTail.FORM_FACTOR, 1.25) + self.prob.set_val(Dynamic.Mission.VELOCITY, [75, 150], units='kn') + self.prob.set_val(Dynamic.Vehicle.MASS, [100000, 100000], units='lbm') + + self.prob.run_model() + + testvals = { + Dynamic.Mission.VELOCITY_RATE: [14.58304081, 11.87430892], + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE: [0.0, 0.0], + Dynamic.Mission.ALTITUDE_RATE: [0.0, 0.0], + Dynamic.Mission.DISTANCE_RATE: [126.58573928, 253.17147857], + 'normal_force': [85313.25425063, 41138.11842255], + 'fuselage_pitch': [0.0, 0.0], + 'dmass_dv': [-0.4852005, -0.60896963], + } + check_prob_outputs(self.prob, testvals, rtol=1e-6) + + partial_data = self.prob.check_partials( + out_stream=None, method='cs', excludes=['*params*', '*aero*'] + ) + assert_check_partials(partial_data, atol=1e-8, rtol=1e-8) + + +class RotationODETestCase(unittest.TestCase): + """Test 2-degree of freedom rotation ODE.""" + + def setUp(self): + self.prob = om.Problem() + + aviary_options = get_option_defaults() + aviary_options.set_val(Aircraft.Engine.GLOBAL_THROTTLE, True) + default_mission_subsystems = get_default_mission_subsystems( + 'GASP', [build_engine_deck(aviary_options)] + ) + + self.prob.model = TakeOffODE( + num_nodes=2, + rotation=True, + aviary_options=get_option_defaults(), + subsystems=default_mission_subsystems, + ) + setup_model_options(self.prob, aviary_options) + + def test_rotation_partials(self): + # Check partial derivatives + self.prob.setup(check=False, force_alloc_complex=True) + + self.prob.set_val(Aircraft.Wing.INCIDENCE, 1.5, units='deg') + self.prob.set_val(Dynamic.Vehicle.MASS, [100000, 100000], units='lbm') + self.prob.set_val(Dynamic.Vehicle.ANGLE_OF_ATTACK, [1.5, 1.5], units='deg') + self.prob.set_val(Dynamic.Mission.VELOCITY, [100, 100], units='kn') + self.prob.set_val('t_curr', [1, 2], units='s') + self.prob.set_val('interference_independent_of_shielded_area', 1.89927266) + self.prob.set_val('drag_loss_due_to_shielded_wing_area', 68.02065834) + self.prob.set_val(Aircraft.Wing.FORM_FACTOR, 1.25) + self.prob.set_val(Aircraft.VerticalTail.FORM_FACTOR, 1.25) + self.prob.set_val(Aircraft.HorizontalTail.FORM_FACTOR, 1.25) + + set_params_for_unit_tests(self.prob) + + self.prob.run_model() + + tol = 1e-6 + assert_near_equal( + self.prob[Dynamic.Mission.VELOCITY_RATE], + np.array([13.68875852, 13.68875852]), + tol, + ) + assert_near_equal( + self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], np.array([0.0, 0.0]), tol + ) + assert_near_equal(self.prob[Dynamic.Mission.ALTITUDE_RATE], np.array([0.0, 0.0]), tol) + assert_near_equal( + self.prob[Dynamic.Mission.DISTANCE_RATE], np.array([168.781, 168.781]), tol + ) + assert_near_equal( + self.prob['normal_force'], np.array([66936.59676831, 66936.59676831]), tol + ) + assert_near_equal(self.prob['fuselage_pitch'], np.array([0.0, 0.0]), tol) + + partial_data = self.prob.check_partials( + out_stream=None, method='cs', excludes=['*params*', '*aero*'] + ) + assert_check_partials(partial_data, atol=1e-8, rtol=1e-8) + + +class AscentODETestCase(unittest.TestCase): + def setUp(self): + self.prob = om.Problem() + + aviary_options = get_option_defaults() + aviary_options.set_val(Aircraft.Engine.GLOBAL_THROTTLE, True) + default_mission_subsystems = get_default_mission_subsystems( + 'GASP', [build_engine_deck(aviary_options)] + ) + + self.prob.model = TakeOffODE( + num_nodes=2, aviary_options=aviary_options, subsystems=default_mission_subsystems + ) + + setup_model_options( + self.prob, AviaryValues({Aircraft.Engine.NUM_ENGINES: ([2], 'unitless')}) + ) + + def test_ascent_partials(self): + # Test partial derivatives + self.prob.setup(check=False, force_alloc_complex=True) + + # TODO: These values are kind of hokey, but were in the previous test. + self.prob.set_val(Dynamic.Mission.VELOCITY, [100, 100], units='kn') + self.prob.set_val(Dynamic.Vehicle.MASS, [1, 1], units='kg') + + self.prob.set_val('t_curr', [1, 2], units='s') + self.prob.set_val('interference_independent_of_shielded_area', 1.89927266) + self.prob.set_val('drag_loss_due_to_shielded_wing_area', 68.02065834) + self.prob.set_val(Aircraft.Wing.FORM_FACTOR, 1.25) + self.prob.set_val(Aircraft.VerticalTail.FORM_FACTOR, 1.25) + self.prob.set_val(Aircraft.HorizontalTail.FORM_FACTOR, 1.25) + + set_params_for_unit_tests(self.prob) + + self.prob.run_model() + + tol = tol = 1e-6 + assert_near_equal( + self.prob[Dynamic.Mission.VELOCITY_RATE], + np.array([642156.99315828, 642156.99315828]), + tol, + ) + assert_near_equal( + self.prob[Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE], + np.array([2260.37849562, 2260.37849562]), + tol, + ) + assert_near_equal(self.prob[Dynamic.Mission.ALTITUDE_RATE], np.array([0.0, 0.0]), tol) + assert_near_equal( + self.prob[Dynamic.Mission.DISTANCE_RATE], np.array([168.781, 168.781]), tol + ) + assert_near_equal(self.prob['angle_of_attack_rate'], np.array([0.0, 0.0]), tol) + assert_near_equal(self.prob['normal_force'], np.array([0.0, 0.0]), tol) + assert_near_equal(self.prob['fuselage_pitch'], np.array([0.0, 0.0]), tol) + assert_near_equal(self.prob['load_factor'], np.array([11849.10281268, 11849.10281268]), tol) + + partial_data = self.prob.check_partials( + out_stream=None, method='cs', excludes=['*params*', '*aero*'] + ) + assert_check_partials(partial_data, atol=1e-8, rtol=1e-8) + + +if __name__ == '__main__': + unittest.main() diff --git a/aviary/mission/two_dof/phases/groundroll_phase.py b/aviary/mission/two_dof/phases/groundroll_phase.py deleted file mode 100644 index c7c27c1688..0000000000 --- a/aviary/mission/two_dof/phases/groundroll_phase.py +++ /dev/null @@ -1,132 +0,0 @@ -from aviary.mission.two_dof.ode.groundroll_ode import GroundrollODE -from aviary.mission.initial_guess_builders import ( - InitialGuessControl, - InitialGuessIntegrationVariable, - InitialGuessState, -) -from aviary.mission.phase_builder import PhaseBuilder -from aviary.utils.aviary_options_dict import AviaryOptionsDictionary -from aviary.utils.aviary_values import AviaryValues -from aviary.variable_info.variables import Dynamic -from aviary.mission.phase_utils import add_subsystem_variables_to_phase - - -class GroundrollPhaseOptions(AviaryOptionsDictionary): - def declare_options(self): - self.declare( - name='num_segments', - types=int, - default=1, - desc='The number of segments in transcription creation in Dymos. ' - 'The default value is 1.', - ) - - self.declare( - name='order', - types=int, - default=3, - desc='The order of polynomials for interpolation in the transcription ' - 'created in Dymos. The default value is 3.', - ) - - defaults = { - 'mass_bounds': (0.0, 200_000.0), - 'mass_ref': 100_000.0, - 'mass_defect_ref': 100.0, - } - self.add_state_options('mass', units='lbm', defaults=defaults) - - defaults = { - 'time_duration_bounds': (1.0, 100.0), - } - self.add_time_options(units='s', defaults=defaults) - - # NOTE: All GASP phases before accel are in 'ft'. - defaults = { - 'distance_bounds': (0.0, 4000.0), - 'distance_ref': 3000.0, - } - self.add_state_options('distance', units='ft', defaults=defaults) - - defaults = { - 'velocity_bounds': (0.0, 1000.0), - 'velocity_ref': 100.0, - } - self.add_state_options('velocity', units='kn', defaults=defaults) - - -class GroundrollPhase(PhaseBuilder): - """ - A phase builder for a groundroll phase in a mission simulation. - - This class extends the PhaseBuilder class, providing specific implementations for - the groundroll phase of a 2-degree of freedom flight mission. - - Attributes - ---------- - Inherits all attributes from PhaseBuilder. - - Methods - ------- - Inherits all methods from PhaseBuilder. - Additional method overrides and new methods specific to the groundroll phase are included. - """ - - default_name = 'groundroll_phase' - default_ode_class = GroundrollODE - default_options_class = GroundrollPhaseOptions - - _initial_guesses_meta_data_ = {} - - def build_phase(self, aviary_options: AviaryValues = None): - phase = self.phase = super().build_phase(aviary_options) - - # Add states - self.add_state('velocity', Dynamic.Mission.VELOCITY, Dynamic.Mission.VELOCITY_RATE) - self.add_state( - 'mass', - Dynamic.Vehicle.MASS, - Dynamic.Vehicle.Propulsion.FUEL_FLOW_RATE_NEGATIVE_TOTAL, - ) - self.add_state('distance', Dynamic.Mission.DISTANCE, Dynamic.Mission.DISTANCE_RATE) - - add_subsystem_variables_to_phase(phase, self.name, self.subsystems) - - phase.add_parameter('t_init_gear', units='s', static_target=True, opt=False, val=100) - phase.add_parameter('t_init_flaps', units='s', static_target=True, opt=False, val=100) - - # boundary/path constraints + controls - # the final TAS is constrained externally to define the transition to the rotation - # phase - - phase.add_timeseries_output('time', units='s', output_name='time') - phase.add_timeseries_output(Dynamic.Vehicle.Propulsion.THRUST_TOTAL, units='lbf') - - phase.add_timeseries_output('normal_force') - phase.add_timeseries_output(Dynamic.Atmosphere.MACH) - phase.add_timeseries_output('EAS', units='kn') - - phase.add_timeseries_output(Dynamic.Vehicle.LIFT) - phase.add_timeseries_output('CL') - phase.add_timeseries_output('CD') - phase.add_timeseries_output('fuselage_pitch', output_name='theta', units='deg') - - return phase - - -# Adding initial guess metadata -GroundrollPhase._add_initial_guess_meta_data( - InitialGuessIntegrationVariable(), desc='initial guess for time options' -) -GroundrollPhase._add_initial_guess_meta_data( - InitialGuessState('velocity'), desc='initial guess for true airspeed state' -) -GroundrollPhase._add_initial_guess_meta_data( - InitialGuessState('mass'), desc='initial guess for mass state' -) -GroundrollPhase._add_initial_guess_meta_data( - InitialGuessState('distance'), desc='initial guess for distance state' -) -GroundrollPhase._add_initial_guess_meta_data( - InitialGuessControl('throttle'), desc='initial guess for throttle' -) diff --git a/aviary/mission/two_dof/phases/rotation_phase.py b/aviary/mission/two_dof/phases/rotation_phase.py deleted file mode 100644 index 30e89a2f46..0000000000 --- a/aviary/mission/two_dof/phases/rotation_phase.py +++ /dev/null @@ -1,186 +0,0 @@ -import numpy as np - -from aviary.mission.two_dof.ode.rotation_ode import RotationODE -from aviary.mission.initial_guess_builders import ( - InitialGuessControl, - InitialGuessIntegrationVariable, - InitialGuessState, -) -from aviary.mission.phase_builder import PhaseBuilder -from aviary.utils.aviary_options_dict import AviaryOptionsDictionary -from aviary.utils.aviary_values import AviaryValues -from aviary.variable_info.variables import Dynamic -from aviary.mission.phase_utils import add_subsystem_variables_to_phase - - -class RotationPhaseOptions(AviaryOptionsDictionary): - def declare_options(self): - self.declare( - name='num_segments', - types=int, - default=None, - desc='The number of segments in transcription creation in Dymos. ', - ) - - self.declare( - name='order', - types=int, - default=None, - desc='The order of polynomials for interpolation in the transcription ' - 'created in Dymos.', - ) - - defaults = { - 'mass_ref': 100_000.0, - 'mass_bounds': (0.0, 190_000.0), - } - self.add_state_options('mass', units='lbm', defaults=defaults) - - # NOTE: All GASP phases before accel are in 'ft'. - defaults = { - 'distance_ref': 3000.0, - 'distance_bounds': (0.0, 10.0e3), - } - self.add_state_options('distance', units='ft', defaults=defaults) - - defaults = { - 'velocity_ref': 100.0, - 'velocity_bounds': (0.0, 1000.0), - } - self.add_state_options('velocity', units='kn', defaults=defaults) - - defaults = { - 'angle_of_attack_defect_ref': 0.01, - 'angle_of_attack_bounds': (0.0, 25 * np.pi / 180), - } - self.add_state_options('angle_of_attack', units='rad', defaults=defaults) - - defaults = { - 'time_duration_bounds': (1.0, 100.0), - } - self.add_time_options(units='s', defaults=defaults) - - self.declare( - 'reserve', - types=bool, - default=False, - desc='Designate this phase as a reserve phase and contributes its fuel burn ' - 'towards the reserve mission fuel requirements. Reserve phases should be ' - 'be placed after all non-reserve phases in the phase_info.', - ) - - self.declare( - name='target_distance', - default=None, - units='m', - desc='The total distance traveled by the aircraft from takeoff to landing ' - 'for the primary mission, not including reserve missions. This value must ' - 'be positive.', - ) - - # The options below have not yet been revamped. - - self.declare( - name='normal_ref', - default=1.0, - units='lbf', - desc='Scale factor ref for the normal force constraint.', - ) - - self.declare( - name='normal_ref0', - default=0.0, - units='lbf', - desc='Scale factor ref0 for the normal force constraint.', - ) - - -class RotationPhase(PhaseBuilder): - """ - A phase builder for a rotation phase in a 2-degree of freedom mission simulation. - - This class extends the PhaseBuilder class, providing specific implementations for - the rotation phase of a flight mission. - - Attributes - ---------- - Inherits all attributes from PhaseBuilder. - - Methods - ------- - Inherits all methods from PhaseBuilder. - Additional method overrides and new methods specific to the rotation phase are included. - """ - - default_name = 'rotation_phase' - default_ode_class = RotationODE - default_options_class = RotationPhaseOptions - - _initial_guesses_meta_data_ = {} - - def build_phase(self, aviary_options: AviaryValues = None): - phase = self.phase = super().build_phase(aviary_options) - - # Retrieve user options values - user_options = self.user_options - normal_ref = user_options.get_val('normal_ref', units='lbf') - normal_ref0 = user_options.get_val('normal_ref0', units='lbf') - - # Add states - self.add_state('angle_of_attack', Dynamic.Vehicle.ANGLE_OF_ATTACK, 'angle_of_attack_rate') - self.add_state('velocity', Dynamic.Mission.VELOCITY, Dynamic.Mission.VELOCITY_RATE) - self.add_state( - 'mass', - Dynamic.Vehicle.MASS, - Dynamic.Vehicle.Propulsion.FUEL_FLOW_RATE_NEGATIVE_TOTAL, - ) - self.add_state('distance', Dynamic.Mission.DISTANCE, Dynamic.Mission.DISTANCE_RATE) - - add_subsystem_variables_to_phase(phase, self.name, self.subsystems) - - # Add parameters - phase.add_parameter('t_init_gear', units='s', static_target=True, opt=False, val=100) - phase.add_parameter('t_init_flaps', units='s', static_target=True, opt=False, val=100) - - # Add boundary constraints - phase.add_boundary_constraint( - 'normal_force', - loc='final', - equals=0, - units='lbf', - ref=normal_ref, - ref0=normal_ref0, - ) - - # Add timeseries outputs - phase.add_timeseries_output(Dynamic.Vehicle.Propulsion.THRUST_TOTAL, units='lbf') - phase.add_timeseries_output('normal_force') - phase.add_timeseries_output(Dynamic.Atmosphere.MACH) - phase.add_timeseries_output('EAS', units='kn') - phase.add_timeseries_output(Dynamic.Vehicle.LIFT) - phase.add_timeseries_output('CL') - phase.add_timeseries_output('CD') - phase.add_timeseries_output('fuselage_pitch', output_name='theta', units='deg') - - return phase - - -# Adding initial guess metadata -RotationPhase._add_initial_guess_meta_data( - InitialGuessIntegrationVariable(), desc='initial guess for time options' -) -RotationPhase._add_initial_guess_meta_data( - InitialGuessState('angle_of_attack'), desc='initial guess for angle of attack state' -) -RotationPhase._add_initial_guess_meta_data( - InitialGuessState('velocity'), desc='initial guess for true airspeed state' -) -RotationPhase._add_initial_guess_meta_data( - InitialGuessState('mass'), desc='initial guess for mass state' -) -RotationPhase._add_initial_guess_meta_data( - InitialGuessState('distance'), desc='initial guess for distance state' -) -RotationPhase._add_initial_guess_meta_data( - InitialGuessControl('throttle'), desc='initial guess for throttle' -) diff --git a/aviary/mission/two_dof/phases/ascent_phase.py b/aviary/mission/two_dof/phases/takeoff_phase.py similarity index 51% rename from aviary/mission/two_dof/phases/ascent_phase.py rename to aviary/mission/two_dof/phases/takeoff_phase.py index 322bbe3174..6063fd11b8 100644 --- a/aviary/mission/two_dof/phases/ascent_phase.py +++ b/aviary/mission/two_dof/phases/takeoff_phase.py @@ -1,19 +1,19 @@ import numpy as np -from aviary.mission.two_dof.ode.ascent_ode import AscentODE from aviary.mission.initial_guess_builders import ( InitialGuessControl, InitialGuessIntegrationVariable, InitialGuessState, ) from aviary.mission.phase_builder import PhaseBuilder +from aviary.mission.phase_utils import add_subsystem_variables_to_phase +from aviary.mission.two_dof.ode.takeoff_ode import TakeOffODE from aviary.utils.aviary_options_dict import AviaryOptionsDictionary from aviary.utils.aviary_values import AviaryValues from aviary.variable_info.variables import Dynamic -from aviary.mission.phase_utils import add_subsystem_variables_to_phase -class AscentPhaseOptions(AviaryOptionsDictionary): +class TakeoffPhaseOptions(AviaryOptionsDictionary): def declare_options(self): self.declare( name='num_segments', @@ -26,28 +26,33 @@ def declare_options(self): self.declare( name='order', types=int, - default=None, + default=3, desc='The order of polynomials for interpolation in the transcription ' - 'created in Dymos.', + 'created in Dymos. The default value is 3.', ) defaults = { + 'mass_bounds': (0.0, 200_000.0), 'mass_ref': 100_000.0, - 'mass_defect_ref': 1.0e2, - 'mass_bounds': (0.0, 190_000.0), + 'mass_defect_ref': 100.0, } self.add_state_options('mass', units='lbm', defaults=defaults) + defaults = { + 'time_duration_bounds': (1.0, 100.0), + } + self.add_time_options(units='s', defaults=defaults) + # NOTE: All GASP phases before accel are in 'ft'. defaults = { + 'distance_bounds': (0.0, 10000.0), 'distance_ref': 3000.0, - 'distance_bounds': (0.0, 10.0e3), } self.add_state_options('distance', units='ft', defaults=defaults) defaults = { - 'velocity_ref': 1.0e2, 'velocity_bounds': (0.0, 1000.0), + 'velocity_ref': 100.0, } self.add_state_options('velocity', units='kn', defaults=defaults) @@ -69,9 +74,26 @@ def declare_options(self): 'angle_of_attack_bounds': (np.deg2rad(-30), np.deg2rad(30)), 'angle_of_attack_optimize': True, } + # NOTE: Sometimes alpha is a control, and sometimes a state. This actually works + # for making sure both sets of options are in there. + self.add_state_options('angle_of_attack', units='rad', defaults=defaults) self.add_control_options('angle_of_attack', units='rad', defaults=defaults) - self.add_time_options(units='s') + self.declare( + 'ground_roll', + types=bool, + default=False, + desc='True if the aircraft is confined to the ground. Angle of attack is fixed and ' + 'removed as an input.', + ) + + self.declare( + 'rotation', + types=bool, + default=False, + desc='True if the aircraft is pitching up, but the rear wheels are still on the ' + 'ground.', + ) self.declare( 'reserve', @@ -93,6 +115,20 @@ def declare_options(self): # The options below have not yet been revamped. + self.declare( + name='normal_ref', + default=1.0, + units='lbf', + desc='Scale factor ref for the normal force constraint.', + ) + + self.declare( + name='normal_ref0', + default=0.0, + units='lbf', + desc='Scale factor ref0 for the normal force constraint.', + ) + self.declare( name='pitch_constraint_bounds', default=(0.0, 15.0), @@ -111,12 +147,11 @@ def declare_options(self): ) -class AscentPhase(PhaseBuilder): +class TakeoffPhase(PhaseBuilder): """ - A phase builder for an ascent phase in a 2-degree of freedom mission simulation. + A phase builder for a two DOF takeoff phase. - This class extends the PhaseBuilder class, providing specific implementations for - the ascent phase of a flight mission. + This can be used to build takeoff phases and includes support for ground roll and rotation. Attributes ---------- @@ -125,12 +160,11 @@ class AscentPhase(PhaseBuilder): Methods ------- Inherits all methods from PhaseBuilder. - Additional method overrides and new methods specific to the ascent phase are included. """ - default_name = 'ascent_phase' - default_ode_class = AscentODE - default_options_class = AscentPhaseOptions + default_name = 'takeoff_phase' + default_ode_class = TakeOffODE + default_options_class = TakeoffPhaseOptions _initial_guesses_meta_data_ = {} @@ -139,16 +173,27 @@ def build_phase(self, aviary_options: AviaryValues = None): # Retrieve user options values user_options = self.user_options - + ground_roll = user_options.get_val('ground_roll') + rotation = user_options.get_val('rotation') pitch_constraint_bounds = user_options.get_val('pitch_constraint_bounds', units='deg') pitch_constraint_ref = user_options.get_val('pitch_constraint_ref', units='deg') - self.add_state( - 'flight_path_angle', - Dynamic.Mission.FLIGHT_PATH_ANGLE, - Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, - ) - self.add_state('altitude', Dynamic.Mission.ALTITUDE, Dynamic.Mission.ALTITUDE_RATE) + # Add states + if rotation: + self.add_state( + 'angle_of_attack', + Dynamic.Vehicle.ANGLE_OF_ATTACK, + 'angle_of_attack_rate', + ) + + if not (ground_roll or rotation): + self.add_state( + 'flight_path_angle', + Dynamic.Mission.FLIGHT_PATH_ANGLE, + Dynamic.Mission.FLIGHT_PATH_ANGLE_RATE, + ) + self.add_state('altitude', Dynamic.Mission.ALTITUDE, Dynamic.Mission.ALTITUDE_RATE) + self.add_state('velocity', Dynamic.Mission.VELOCITY, Dynamic.Mission.VELOCITY_RATE) self.add_state( 'mass', @@ -159,75 +204,117 @@ def build_phase(self, aviary_options: AviaryValues = None): add_subsystem_variables_to_phase(phase, self.name, self.subsystems) - self.add_control( - 'angle_of_attack', - Dynamic.Vehicle.ANGLE_OF_ATTACK, + # Add controls + if not (ground_roll or rotation): + self.add_control( + 'angle_of_attack', + Dynamic.Vehicle.ANGLE_OF_ATTACK, + ) + + # Add parameters + # TODO: These are backdoor defaults. + # Although, I can see how this is useful for keeping flaps/gear from engaging before + # ascent. + if ground_roll or rotation: + t_init_gear = 100.0 + t_init_flaps = 100.0 + else: + t_init_gear = 38.25 + t_init_flaps = 48.21 + + phase.add_parameter( + 't_init_gear', + units='s', + static_target=True, + opt=False, + val=t_init_gear, ) - - phase.add_path_constraint('load_factor', upper=1.10, lower=0.0) - - phase.add_path_constraint( - 'fuselage_pitch', - 'theta', - lower=pitch_constraint_bounds[0], - upper=pitch_constraint_bounds[1], - units='deg', - ref=pitch_constraint_ref, + phase.add_parameter( + 't_init_flaps', + units='s', + static_target=True, + opt=False, + val=t_init_flaps, ) - phase.add_parameter('t_init_gear', units='s', static_target=True, opt=False, val=38.25) - - phase.add_parameter('t_init_flaps', units='s', static_target=True, opt=False, val=48.21) - + # Add boundary constraints + if rotation: + normal_ref = user_options.get_val('normal_ref', units='lbf') + normal_ref0 = user_options.get_val('normal_ref0', units='lbf') + + phase.add_boundary_constraint( + 'normal_force', + loc='final', + equals=0, + units='lbf', + ref=normal_ref, + ref0=normal_ref0, + ) + + if not (ground_roll or rotation): + phase.add_path_constraint('load_factor', upper=1.10, lower=0.0) + + phase.add_path_constraint( + 'fuselage_pitch', + 'theta', + lower=pitch_constraint_bounds[0], + upper=pitch_constraint_bounds[1], + units='deg', + ref=pitch_constraint_ref, + ) + + # Add timeseries outputs + phase.add_timeseries_output('time', units='s', output_name='time') phase.add_timeseries_output(Dynamic.Vehicle.Propulsion.THRUST_TOTAL, units='lbf') + phase.add_timeseries_output('normal_force') phase.add_timeseries_output(Dynamic.Atmosphere.MACH) phase.add_timeseries_output('EAS', units='kn') + phase.add_timeseries_output(Dynamic.Vehicle.LIFT) phase.add_timeseries_output('CL') phase.add_timeseries_output('CD') + phase.add_timeseries_output('fuselage_pitch', output_name='theta', units='deg') return phase + def _extra_ode_init_kwargs(self): + """Return extra kwargs required for initializing the ODE.""" + # TODO: support external_subsystems and meta_data in the base class + return { + 'ground_roll': self.user_options.get_val('ground_roll'), + 'rotation': self.user_options.get_val('rotation'), + } + # Adding initial guess metadata -AscentPhase._add_initial_guess_meta_data( +TakeoffPhase._add_initial_guess_meta_data( InitialGuessIntegrationVariable(), desc='initial guess for time options' ) - -AscentPhase._add_initial_guess_meta_data( - InitialGuessState('flight_path_angle'), desc='initial guess for flight path angle state' -) - -AscentPhase._add_initial_guess_meta_data( +TakeoffPhase._add_initial_guess_meta_data( InitialGuessState('altitude'), desc='initial guess for altitude state' ) - -AscentPhase._add_initial_guess_meta_data( +TakeoffPhase._add_initial_guess_meta_data( InitialGuessState('velocity'), desc='initial guess for true airspeed state' ) - -AscentPhase._add_initial_guess_meta_data( +TakeoffPhase._add_initial_guess_meta_data( InitialGuessState('mass'), desc='initial guess for mass state' ) - -AscentPhase._add_initial_guess_meta_data( +TakeoffPhase._add_initial_guess_meta_data( InitialGuessState('distance'), desc='initial guess for distance state' ) - -AscentPhase._add_initial_guess_meta_data( - InitialGuessControl('angle_of_attack'), - desc='initial guess for angle of attack control', +TakeoffPhase._add_initial_guess_meta_data( + InitialGuessControl('throttle'), desc='initial guess for throttle' ) - -AscentPhase._add_initial_guess_meta_data( +TakeoffPhase._add_initial_guess_meta_data( + InitialGuessState('flight_path_angle'), desc='initial guess for flight path angle state' +) +TakeoffPhase._add_initial_guess_meta_data( + InitialGuessState('angle_of_attack'), desc='initial guess for angle of attack state' +) +TakeoffPhase._add_initial_guess_meta_data( InitialGuessControl('tau_gear'), desc='when the gear is retracted' ) - -AscentPhase._add_initial_guess_meta_data( +TakeoffPhase._add_initial_guess_meta_data( InitialGuessControl('tau_flaps'), desc='when the flaps are retracted' ) - -AscentPhase._add_initial_guess_meta_data( - InitialGuessControl('throttle'), desc='initial guess for throttle' -) diff --git a/aviary/mission/two_dof_problem_configurator.py b/aviary/mission/two_dof_problem_configurator.py index 441246d18c..67ff54c00f 100644 --- a/aviary/mission/two_dof_problem_configurator.py +++ b/aviary/mission/two_dof_problem_configurator.py @@ -5,15 +5,13 @@ from aviary.mission.two_dof.ode.params import ParamPort from aviary.mission.two_dof.ode.taxi_ode import TaxiSegment from aviary.mission.two_dof.phases.accel_phase import AccelPhase -from aviary.mission.two_dof.phases.ascent_phase import AscentPhase from aviary.mission.two_dof.phases.breguet_cruise_phase import ( BreguetCruisePhase, ElectricCruisePhase, ) from aviary.mission.two_dof.phases.flight_phase import FlightPhase -from aviary.mission.two_dof.phases.groundroll_phase import GroundrollPhase -from aviary.mission.two_dof.phases.rotation_phase import RotationPhase from aviary.mission.two_dof.phases.simple_cruise_phase import SimpleCruisePhase +from aviary.mission.two_dof.phases.takeoff_phase import TakeoffPhase from aviary.mission.two_dof.polynomial_fit import PolynomialFit from aviary.mission.problem_configurator import ProblemConfiguratorBase from aviary.mission.utils import process_guess_var @@ -219,24 +217,25 @@ def get_phase_builder(self, aviary_group, phase_name, phase_options): builder = phase_options['user_options']['phase_builder'] if builder is PhaseType.BREGUET_RANGE: + # TODO: This is the only one currently accessed by name string, to preserve + # some legacy tests. if 'electric_cruise' in phase_name: phase_builder = ElectricCruisePhase else: phase_builder = BreguetCruisePhase - return phase_builder elif builder is PhaseType.SIMPLE_CRUISE: phase_builder = SimpleCruisePhase - return phase_builder - - if 'groundroll' in phase_name: - phase_builder = GroundrollPhase - elif 'rotation' in phase_name: - phase_builder = RotationPhase - elif 'accel' in phase_name: - phase_builder = AccelPhase - elif 'ascent' in phase_name: - phase_builder = AscentPhase + + elif builder is PhaseType.TWO_DOF_TAKEOFF: + phase_builder = TakeoffPhase + + elif builder is PhaseType.ACCEL: + phase_builder = AccelPhase + + else: + phase_builder = FlightPhase + else: # All other phases are flight phases. phase_builder = FlightPhase @@ -264,7 +263,8 @@ def set_phase_options(self, aviary_group, phase_name, phase_idx, phase, user_opt comm : MPI.Comm or MPI Communicator from OpenMDAO problem. """ - if user_options['phase_builder'] is PhaseType.BREGUET_RANGE: + phase_builder = user_options['phase_builder'] + if phase_builder is PhaseType.BREGUET_RANGE: # The Breguet Range Cruise phase integrates over mass instead of time. # We rely on mass being monotonically non-increasing across the phase. phase.set_time_options( @@ -291,63 +291,47 @@ def set_phase_options(self, aviary_group, phase_name, phase_idx, phase, user_opt fix_duration = duration is not None fix_initial = initial is not None - if 'ascent' in phase_name: - phase.set_time_options( - units='s', - targets='t_curr', - input_initial=True, - input_duration=True, - ) - - elif 'descent' in phase_name: - phase.set_time_options( - duration_bounds=duration_bounds, - fix_initial=fix_initial, - input_initial=input_initial, - units='s', - duration_ref=duration_ref, - ) - + if initial_bounds[0] is not None and initial_bounds[1]: + # Upper bound is good for a ref. + time_initial_ref = initial_bounds[1] else: - if initial_bounds[0] is not None and initial_bounds[1]: - # Upper bound is good for a ref. - time_initial_ref = initial_bounds[1] - else: - time_initial_ref = 600.0 + time_initial_ref = 600.0 - if duration_ref is None: - # Why are we overriding this? - duration_ref = 0.5 * (duration_bounds[0] + duration_bounds[1]) + if duration_ref is None: + # Why are we overriding this? + duration_ref = 0.5 * (duration_bounds[0] + duration_bounds[1]) - input_initial = phase_idx > 0 - extra_args = {} - - if fix_initial or input_initial: - if comm.size > 1: - # Phases are disconnected to run in parallel, so initial ref is - # valid. - initial_ref = time_initial_ref - else: - # Redundant on a fixed input; raises a warning if specified. - initial_ref = None + input_initial = phase_idx > 0 + extra_args = {} + if fix_initial or input_initial: + if comm.size > 1: + # Phases are disconnected to run in parallel, so initial ref is + # valid. + initial_ref = time_initial_ref else: - extra_args['initial_bounds'] = initial_bounds - - if user_options['phase_builder'] is PhaseType.SIMPLE_CRUISE: - extra_args['targets'] = 'time' + # Redundant on a fixed input; raises a warning if specified. + initial_ref = None - phase.set_time_options( - fix_initial=fix_initial, - fix_duration=fix_duration, - units=time_units, - duration_bounds=duration_bounds, - duration_ref=duration_ref, - initial_ref=initial_ref, - **extra_args, - ) + else: + extra_args['initial_bounds'] = initial_bounds + + if phase_builder is PhaseType.SIMPLE_CRUISE: + extra_args['targets'] = 'time' + elif phase_builder is PhaseType.TWO_DOF_TAKEOFF: + extra_args['targets'] = 't_curr' + + phase.set_time_options( + fix_initial=fix_initial, + fix_duration=fix_duration, + units=time_units, + duration_bounds=duration_bounds, + duration_ref=duration_ref, + initial_ref=initial_ref, + **extra_args, + ) - if user_options['phase_builder'] is not PhaseType.SIMPLE_CRUISE: + if phase_builder is not PhaseType.SIMPLE_CRUISE: phase.add_control( Dynamic.Vehicle.Propulsion.THROTTLE, targets=Dynamic.Vehicle.Propulsion.THROTTLE, @@ -359,7 +343,7 @@ def set_phase_options(self, aviary_group, phase_name, phase_idx, phase, user_opt # The issue is that aero methods are hardcoded for GASP mission phases # instead of being defaulted somewhere, so they don't use phase_info # aviary_group.mission_info[phase_name]['phase_type'] = phase_name - if phase_name in ['ascent', 'groundroll', 'rotation']: + if phase_builder is PhaseType.TWO_DOF_TAKEOFF: # safely add in default method in way that doesn't overwrite existing method # and create nested structure if it doesn't already exist aviary_group.mission_info[phase_name].setdefault('subsystem_options', {}).setdefault( @@ -470,17 +454,18 @@ def link_phases(self, aviary_group, phases, connect_directly=True): (phase1 in aviary_group.reserve_phases) == (phase2 in aviary_group.reserve_phases) ) - and not ({'groundroll', 'rotation'} & {phase1, phase2}) - and not ('accel', 'climb1') == (phase1, phase2) + and not info1['user_options'].get('ground_roll') + and phase_builder1 is not PhaseType.TWO_DOF_TAKEOFF ): # required for convergence of FwGm states_to_link[Dynamic.Mission.ALTITUDE] = connect_directly # if either phase is rotation, we need to connect velocity # ascent to accel also requires velocity - if 'rotation' in (phase1, phase2) or ('ascent', 'accel') == (phase1, phase2): + # TODO: This may need refined now that trajectories are more flexible. + if phase_builder1 is PhaseType.TWO_DOF_TAKEOFF: states_to_link[Dynamic.Mission.VELOCITY] = connect_directly # if the first phase is rotation, we also need alpha - if phase1 == 'rotation': + if info1['user_options'].get('rotation'): states_to_link[Dynamic.Vehicle.ANGLE_OF_ATTACK] = False for state, connected in states_to_link.items(): @@ -572,8 +557,8 @@ def link_phases(self, aviary_group, phases, connect_directly=True): aviary_group.promotes('landing', inputs=param_list) aviary_group.connect('taxi.mass', 'vrot.mass') - if 'ascent' in aviary_group.mission_info: - self._add_groundroll_eq_constraint(aviary_group) + if len(phases) > 1: + self._add_groundroll_eq_constraint(aviary_group, phases) def check_trajectory(self, aviary_group): """ @@ -586,11 +571,20 @@ def check_trajectory(self, aviary_group): """ pass - def _add_groundroll_eq_constraint(self, aviary_group): + def _add_groundroll_eq_constraint(self, aviary_group, phases): """ Add an equality constraint to the problem to ensure that the TAS at the end of the groundroll phase is equal to the rotation velocity at the start of the rotation phase. """ + phase1 = phases[0] + info1 = aviary_group.mission_info[phase1]['user_options'] + groundroll = info1.get('ground_roll') + + if not groundroll: + # Ground roll should be the first phase. If it isn't, then the aircraft is starting + # in flight. + return + aviary_group.add_subsystem( 'groundroll_boundary', om.EQConstraintComp( @@ -602,7 +596,7 @@ def _add_groundroll_eq_constraint(self, aviary_group): ) aviary_group.connect(Mission.Takeoff.ROTATION_VELOCITY, 'groundroll_boundary.rhs:velocity') aviary_group.connect( - 'traj.groundroll.states:velocity', + f'traj.{phase1}.states:velocity', 'groundroll_boundary.lhs:velocity', src_indices=[-1], flat_src_indices=True, @@ -794,9 +788,6 @@ def set_phase_initial_guesses( units=units, ) - if guess_key in control_keys: - pass - # Set initial guess for state variables elif guess_key in state_keys: target_prob.set_val( diff --git a/aviary/models/aircraft/blended_wing_body/generic_BWB_2dof_phase_info.py b/aviary/models/aircraft/blended_wing_body/generic_BWB_2dof_phase_info.py index e1ef9a191c..36a7a91188 100644 --- a/aviary/models/aircraft/blended_wing_body/generic_BWB_2dof_phase_info.py +++ b/aviary/models/aircraft/blended_wing_body/generic_BWB_2dof_phase_info.py @@ -5,6 +5,8 @@ 'groundroll': { 'subsystem_options': {'core_aerodynamics': {'method': 'low_speed'}}, 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'ground_roll': True, 'num_segments': 1, 'order': 3, 'time_initial': (0.0, 's'), @@ -31,6 +33,8 @@ 'rotation': { 'subsystem_options': {'core_aerodynamics': {'method': 'low_speed'}}, 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'rotation': True, 'num_segments': 1, 'order': 3, 'time_duration_bounds': ((1, 50), 's'), @@ -59,10 +63,12 @@ 'ascent': { 'subsystem_options': {'core_aerodynamics': {'method': 'low_speed'}}, 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, 'num_segments': 4, 'order': 3, 'velocity_bounds': ((0, 300), 'kn'), 'velocity_ref': (200, 'kn'), + 'time_duration_ref': (10, 's'), 'mass_bounds': ((0, None), 'lbm'), 'mass_ref': (150_000, 'lbm'), 'mass_defect_ref': (150_000, 'lbm'), @@ -97,6 +103,7 @@ 'accel': { 'subsystem_options': {'core_aerodynamics': {'method': 'cruise'}}, 'user_options': { + 'phase_builder': PhaseType.ACCEL, 'num_segments': 1, 'order': 3, 'alt': (500, 'ft'), diff --git a/aviary/models/aircraft/large_turboprop_freighter/electrified_phase_info.py b/aviary/models/aircraft/large_turboprop_freighter/electrified_phase_info.py index ec9e6dd0ba..7b2faa8fb1 100644 --- a/aviary/models/aircraft/large_turboprop_freighter/electrified_phase_info.py +++ b/aviary/models/aircraft/large_turboprop_freighter/electrified_phase_info.py @@ -1,4 +1,4 @@ -from aviary.variable_info.enums import SpeedType, ThrottleAllocation +from aviary.variable_info.enums import SpeedType, PhaseType, ThrottleAllocation # Energy method energy_phase_info = { @@ -80,6 +80,8 @@ two_dof_phase_info = { 'groundroll': { 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'ground_roll': True, 'num_segments': 1, 'order': 3, 'time_duration_ref': (50.0, 's'), @@ -104,6 +106,8 @@ }, 'rotation': { 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'rotation': True, 'num_segments': 1, 'order': 3, 'time_duration_bounds': ((1, 100), 's'), @@ -133,11 +137,13 @@ }, 'ascent': { 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, 'num_segments': 4, 'order': 3, 'velocity_bounds': ((0, 700), 'kn'), 'velocity_ref': (200, 'kn'), 'velocity_ref0': (0, 'kn'), + 'time_duration_ref': (10, 's'), 'mass_bounds': ((0, None), 'lbm'), 'mass_ref': (150_000, 'lbm'), 'mass_defect_ref': (150_000, 'lbm'), @@ -170,6 +176,7 @@ }, 'accel': { 'user_options': { + 'phase_builder': PhaseType.ACCEL, 'num_segments': 1, 'order': 3, 'alt': (500, 'ft'), @@ -248,6 +255,7 @@ }, 'electric_cruise': { 'user_options': { + 'phase_builder': PhaseType.BREGUET_RANGE, 'alt_cruise': (21_000, 'ft'), 'mach_cruise': 0.475, }, diff --git a/aviary/models/aircraft/large_turboprop_freighter/phase_info.py b/aviary/models/aircraft/large_turboprop_freighter/phase_info.py index 9be4677d11..e4effdb48c 100644 --- a/aviary/models/aircraft/large_turboprop_freighter/phase_info.py +++ b/aviary/models/aircraft/large_turboprop_freighter/phase_info.py @@ -77,6 +77,8 @@ two_dof_phase_info = { 'groundroll': { 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'ground_roll': True, 'num_segments': 1, 'order': 3, 'time_duration_ref': (50.0, 's'), @@ -101,6 +103,8 @@ }, 'rotation': { 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'rotation': True, 'num_segments': 1, 'order': 3, 'time_duration_bounds': ((1, 100), 's'), @@ -130,11 +134,13 @@ }, 'ascent': { 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, 'num_segments': 4, 'order': 3, 'velocity_bounds': ((0, 700), 'kn'), 'velocity_ref': (200, 'kn'), 'velocity_ref0': (0, 'kn'), + 'time_duration_ref': (10, 's'), 'mass_bounds': ((0, None), 'lbm'), 'mass_ref': (150_000, 'lbm'), 'mass_defect_ref': (150_000, 'lbm'), @@ -167,6 +173,7 @@ }, 'accel': { 'user_options': { + 'phase_builder': PhaseType.ACCEL, 'num_segments': 1, 'order': 3, 'alt': (500, 'ft'), diff --git a/aviary/models/missions/two_dof_default.py b/aviary/models/missions/two_dof_default.py index 19a1c9ac76..e96f7e96f2 100644 --- a/aviary/models/missions/two_dof_default.py +++ b/aviary/models/missions/two_dof_default.py @@ -8,6 +8,8 @@ 'groundroll': { 'subsystem_options': {'aerodynamics': {'method': 'low_speed'}}, 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'ground_roll': True, 'num_segments': 1, 'order': 3, 'time_initial': (0.0, 's'), @@ -34,6 +36,8 @@ 'rotation': { 'subsystem_options': {'aerodynamics': {'method': 'low_speed'}}, 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, + 'rotation': True, 'num_segments': 1, 'order': 3, 'time_duration_bounds': ((1, 100), 's'), @@ -63,10 +67,12 @@ 'ascent': { 'subsystem_options': {'aerodynamics': {'method': 'low_speed'}}, 'user_options': { + 'phase_builder': PhaseType.TWO_DOF_TAKEOFF, 'num_segments': 4, 'order': 3, 'velocity_bounds': ((0, 700), 'kn'), 'velocity_ref': (200, 'kn'), + 'time_duration_ref': (10, 's'), 'mass_bounds': ((0, None), 'lbm'), 'mass_ref': (150_000, 'lbm'), 'mass_defect_ref': (150_000, 'lbm'), @@ -101,6 +107,7 @@ 'accel': { 'subsystem_options': {'aerodynamics': {'method': 'cruise'}}, 'user_options': { + 'phase_builder': PhaseType.ACCEL, 'num_segments': 1, 'order': 3, 'alt': (500, 'ft'), diff --git a/aviary/validation_cases/benchmark_tests/test_bench_off_design.py b/aviary/validation_cases/benchmark_tests/test_bench_off_design.py index 547d5c5cf6..983f39f38b 100644 --- a/aviary/validation_cases/benchmark_tests/test_bench_off_design.py +++ b/aviary/validation_cases/benchmark_tests/test_bench_off_design.py @@ -108,32 +108,32 @@ def test_fallout_mission_changed(self): assert_near_equal( prob_fallout.get_val(Mission.Summary.FUEL_MASS, 'lbm'), 46164.07769361, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Mission.Summary.TOTAL_FUEL_MASS, 'lbm'), 29031.53317628, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Mission.Summary.OPERATING_MASS, 'lbm'), 97743.46682372, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Aircraft.CrewPayload.CARGO_MASS, 'lbm'), 5000, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Aircraft.CrewPayload.TOTAL_PAYLOAD_MASS, 'lbm'), 23225, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Aircraft.CrewPayload.PASSENGER_PAYLOAD_MASS, 'lbm'), 18225, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Aircraft.Design.EMPTY_MASS, 'lbm'), @@ -198,17 +198,17 @@ def test_alternate_mission_changed(self): assert_near_equal( prob_alternate.get_val(Mission.Summary.FUEL_MASS, 'lbm'), 33139.01658754, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_alternate.get_val(Mission.Summary.TOTAL_FUEL_MASS, 'lbm'), 23663.00633808, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_alternate.get_val(Mission.Summary.OPERATING_MASS, 'lbm'), 97743.52792979, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_alternate.get_val(Aircraft.CrewPayload.CARGO_MASS, 'lbm'), @@ -218,12 +218,12 @@ def test_alternate_mission_changed(self): assert_near_equal( prob_alternate.get_val(Aircraft.CrewPayload.TOTAL_PAYLOAD_MASS, 'lbm'), 36250, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_alternate.get_val(Aircraft.CrewPayload.PASSENGER_PAYLOAD_MASS, 'lbm'), 33750, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_alternate.get_val(Aircraft.Design.EMPTY_MASS, 'lbm'), @@ -238,7 +238,7 @@ def test_alternate_mission_changed(self): assert_near_equal( prob_alternate.get_val(Mission.Summary.GROSS_MASS, 'lbm'), 157656.53426787, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_alternate.aviary_inputs.get_val(Aircraft.CrewPayload.NUM_FIRST_CLASS), @@ -318,7 +318,7 @@ def compare_results(self, comparison_prob): ] for var in prob_var_list: - assert_near_equal(comparison_prob.get_val(var), self.prob.get_val(var), tolerance=1e-6) + assert_near_equal(comparison_prob.get_val(var), self.prob.get_val(var), tolerance=1e-5) for var in inputs_var_list: assert_near_equal( @@ -354,28 +354,28 @@ def test_fallout_mission_changed(self): ) assert_near_equal( prob_fallout.get_val(Mission.Summary.FUEL_MASS, 'lbm'), - 40505.09154366, - tolerance=1e-6, + 40505.91552026, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Mission.Summary.TOTAL_FUEL_MASS, 'lbm'), - 39910.0046378, - tolerance=1e-6, + 39909.74193096, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Mission.Summary.OPERATING_MASS, 'lbm'), - 95089.9953622, - tolerance=1e-6, + 95090.25806904, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Aircraft.CrewPayload.CARGO_MASS, 'lbm'), 5000, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Aircraft.CrewPayload.TOTAL_PAYLOAD_MASS, 'lbm'), 20000, - tolerance=1e-6, + tolerance=1e-5, ) assert_near_equal( prob_fallout.get_val(Aircraft.CrewPayload.PASSENGER_PAYLOAD_MASS, 'lbm'), @@ -432,17 +432,17 @@ def test_alternate_mission_changed(self): assert_near_equal(prob_alternate.get_val(Mission.Summary.RANGE), 1800, tolerance=1e-6) assert_near_equal( prob_alternate.get_val(Mission.Summary.FUEL_MASS, 'lbm'), - 40505.07149946, + 40505.91552026, tolerance=1e-6, ) assert_near_equal( prob_alternate.get_val(Mission.Summary.TOTAL_FUEL_MASS, 'lbm'), - 21484.2904701, + 21484.97566914, tolerance=1e-6, ) assert_near_equal( prob_alternate.get_val(Mission.Summary.OPERATING_MASS, 'lbm'), - 95089.98897155, + 95090.25806904, tolerance=1e-6, ) assert_near_equal( @@ -473,7 +473,7 @@ def test_alternate_mission_changed(self): ) assert_near_equal( prob_alternate.get_val(Mission.Summary.GROSS_MASS, 'lbm'), - 148674.27944164, + 148675.23373818, tolerance=1e-6, ) assert_near_equal( @@ -558,7 +558,7 @@ def test_payload_range(self): # unittest.main() test = Test2DOFOffDesign() test.setUp() - test.test_fallout_mission_changed() + test.test_alternate_mission_changed() # test = PayloadRangeTest() # test.test_payload_range() diff --git a/aviary/variable_info/enums.py b/aviary/variable_info/enums.py index 3c0b60b708..0b9deb38a8 100644 --- a/aviary/variable_info/enums.py +++ b/aviary/variable_info/enums.py @@ -170,16 +170,24 @@ class PhaseType(Enum): PhaseType is used for replacing a default phase and its equations of motion with a different one. - DEFAULT: Use the default phase builder for this EquationsOfMotion. + ACCEL: Use a phase_builder for accelerating level flight. BREGUET_RANGE: Use a phase builder that implements the Breguet Range equations. + DEFAULT: Use the default phase builder for this EquationsOfMotion. + SIMPLE_CRUISE: Use a phase builder that implements a single DOF (mass) cruise. + + TWO_DOF_TAKEOFF: Use a phase builder that implements two DOF equations suitable for takeoff + phases, including ground_roll and rotation phases. Angle of attack is an additional control + unless ground_roll=True. """ - DEFAULT = 'default' + ACCEL = 'accel' BREGUET_RANGE = 'breguet_range' + DEFAULT = 'default' SIMPLE_CRUISE = 'simple_cruise' + TWO_DOF_TAKEOFF = 'two_dof_takeoff' class ProblemType(Enum):