Skip to content

Commit

Permalink
Simplify ssto (#534)
Browse files Browse the repository at this point in the history
* Modified SSTO problem to be a single component rather than a group

* Updated problem to use wildcard for derivative declaration

* Removed anayltical derivatives. Removed cb as an option and replaced it with options for rho_ref and h_scale. Removed targets and units where not longer necessary

* fixed test_upgrade_guide problem to match changes to LaunchVehicleODE

* added tolerance to coloring

* Modified scaling on problem

* updated some docs for kaushiks SSTO case

* formatting fix for pep8

Co-authored-by: Kaushik Ponnapalli <kaushik.s.ponnapalli@nasa.gov>
Co-authored-by: Rob Falck <rfalck@nasa.gov>
  • Loading branch information
3 people authored Feb 8, 2021
1 parent 9ee496b commit 41c23c5
Show file tree
Hide file tree
Showing 8 changed files with 146 additions and 277 deletions.
31 changes: 14 additions & 17 deletions dymos/examples/ssto/doc/test_doc_ssto_earth.py
Original file line number Diff line number Diff line change
@@ -1,8 +1,6 @@
import unittest

import matplotlib.pyplot as plt
plt.switch_backend('Agg')
plt.style.use('ggplot')

from openmdao.utils.assert_utils import assert_near_equal
from openmdao.utils.testing_utils import use_tempdirs
Expand All @@ -23,7 +21,7 @@ def test_doc_ssto_earth(self):
#
p = om.Problem(model=om.Group())
p.driver = om.pyOptSparseDriver()
p.driver.declare_coloring()
p.driver.declare_coloring(tol=1.0E-12)

from dymos.examples.ssto.launch_vehicle_ode import LaunchVehicleODE

Expand All @@ -33,7 +31,6 @@ def test_doc_ssto_earth(self):
traj = dm.Trajectory()

phase = dm.Phase(ode_class=LaunchVehicleODE,
ode_init_kwargs={'central_body': 'earth'},
transcription=dm.GaussLobatto(num_segments=12, order=3, compressed=False))

traj.add_phase('phase0', phase)
Expand All @@ -44,19 +41,19 @@ def test_doc_ssto_earth(self):
#
phase.set_time_options(initial_bounds=(0, 0), duration_bounds=(10, 500))

phase.add_state('x', fix_initial=True, ref=1.0E5, defect_ref=1.0,
rate_source='eom.xdot', units='m')
phase.add_state('y', fix_initial=True, ref=1.0E5, defect_ref=1.0,
rate_source='eom.ydot', targets=['atmos.y'], units='m')
phase.add_state('vx', fix_initial=True, ref=1.0E3, defect_ref=1.0,
rate_source='eom.vxdot', targets=['eom.vx'], units='m/s')
phase.add_state('vy', fix_initial=True, ref=1.0E3, defect_ref=1.0,
rate_source='eom.vydot', targets=['eom.vy'], units='m/s')
phase.add_state('m', fix_initial=True, ref=1.0E3, defect_ref=1.0,
rate_source='eom.mdot', targets=['eom.m'], units='kg')

phase.add_control('theta', units='rad', lower=-1.57, upper=1.57, targets=['eom.theta'])
phase.add_parameter('thrust', units='N', opt=False, val=2100000.0, targets=['eom.thrust'])
phase.add_state('x', fix_initial=True, ref=1.0E5, defect_ref=10000.0,
rate_source='xdot')
phase.add_state('y', fix_initial=True, ref=1.0E5, defect_ref=10000.0,
rate_source='ydot')
phase.add_state('vx', fix_initial=True, ref=1.0E3, defect_ref=1000.0,
rate_source='vxdot')
phase.add_state('vy', fix_initial=True, ref=1.0E3, defect_ref=1000.0,
rate_source='vydot')
phase.add_state('m', fix_initial=True, ref=1.0E3, defect_ref=100.0,
rate_source='mdot')

phase.add_control('theta', units='rad', lower=-1.57, upper=1.57, targets=['theta'])
phase.add_parameter('thrust', units='N', opt=False, val=2100000.0, targets=['thrust'])

#
# Set the options for our constraints and objective
Expand Down
10 changes: 0 additions & 10 deletions dymos/examples/ssto/doc/test_doc_ssto_linear_tangent_guidance.py
Original file line number Diff line number Diff line change
Expand Up @@ -202,16 +202,6 @@ def compute_partials(self, inputs, jacobian):
jacobian['theta', 'b_ctrl'] = 1.0 / denom
jacobian['theta', 'time_phase'] = a / denom

# @dm.declare_time(units='s', targets=['guidance.time_phase'])
# @dm.declare_state('x', rate_source='eom.xdot', units='m')
# @dm.declare_state('y', rate_source='eom.ydot', units='m')
# @dm.declare_state('vx', rate_source='eom.vxdot', targets=['eom.vx'], units='m/s')
# @dm.declare_state('vy', rate_source='eom.vydot', targets=['eom.vy'], units='m/s')
# @dm.declare_state('m', rate_source='eom.mdot', targets=['eom.m'], units='kg')
# @dm.declare_parameter('thrust', targets=['eom.thrust'], units='N')
# @dm.declare_parameter('a_ctrl', targets=['guidance.a_ctrl'], units='1/s')
# @dm.declare_parameter('b_ctrl', targets=['guidance.b_ctrl'], units=None)
# @dm.declare_parameter('Isp', targets=['eom.Isp'], units='s')
class LaunchVehicleLinearTangentODE(om.Group):
"""
The LaunchVehicleLinearTangentODE for this case consists of a guidance component and
Expand Down
155 changes: 0 additions & 155 deletions dymos/examples/ssto/launch_vehicle_2d_eom_comp.py

This file was deleted.

125 changes: 106 additions & 19 deletions dymos/examples/ssto/launch_vehicle_ode.py
Original file line number Diff line number Diff line change
@@ -1,34 +1,121 @@
import openmdao.api as om
import numpy as np

from .log_atmosphere_comp import LogAtmosphereComp
from .launch_vehicle_2d_eom_comp import LaunchVehicle2DEOM


class LaunchVehicleODE(om.Group):
class LaunchVehicleODE(om.ExplicitComponent):

def initialize(self):
self.options.declare('num_nodes', types=int,
desc='Number of nodes to be evaluated in the RHS')

self.options.declare('central_body', values=['earth', 'moon'], default='earth',
desc='The central gravitational body for the launch vehicle.')
self.options.declare('g', types=float, default=9.80665,
desc='Gravitational acceleration, m/s**2')

self.options.declare('rho_ref', types=float, default=1.225,
desc='Reference atmospheric density, kg/m**3')

self.options.declare('h_scale', types=float, default=8.44E3,
desc='Reference altitude, m')

self.options.declare('CD', types=float, default=0.5,
desc='coefficient of drag')

self.options.declare('S', types=float, default=7.069,
desc='aerodynamic reference area (m**2)')

def setup(self):
nn = self.options['num_nodes']
cb = self.options['central_body']

if cb == 'earth':
rho_ref = 1.225
h_scale = 8.44E3
elif cb == 'moon':
rho_ref = 0.0
h_scale = 1.0
else:
raise RuntimeError('Unrecognized value for central_body: {0}'.format(cb))
self.add_input('y',
val=np.zeros(nn),
desc='altitude',
units='m')

self.add_input('vx',
val=np.zeros(nn),
desc='x velocity',
units='m/s')

self.add_input('vy',
val=np.zeros(nn),
desc='y velocity',
units='m/s')

self.add_input('m',
val=np.zeros(nn),
desc='mass',
units='kg')

self.add_input('theta',
val=np.zeros(nn),
desc='pitch angle',
units='rad')

self.add_input('thrust',
val=2100000 * np.ones(nn),
desc='thrust',
units='N')

self.add_input('Isp',
val=265.2 * np.ones(nn),
desc='specific impulse',
units='s')
# Outputs
self.add_output('xdot',
val=np.zeros(nn),
desc='velocity component in x',
units='m/s')

self.add_output('ydot',
val=np.zeros(nn),
desc='velocity component in y',
units='m/s')

self.add_output('vxdot',
val=np.zeros(nn),
desc='x acceleration magnitude',
units='m/s**2')

self.add_output('vydot',
val=np.zeros(nn),
desc='y acceleration magnitude',
units='m/s**2')

self.add_output('mdot',
val=np.zeros(nn),
desc='mass rate of change',
units='kg/s')

self.add_output('rho',
val=np.zeros(nn),
desc='density',
units='kg/m**3')

# Setup partials
# Complex-step derivatives
self.declare_coloring(wrt='*', method='cs', show_sparsity=True)

def compute(self, inputs, outputs):

theta = inputs['theta']
cos_theta = np.cos(theta)
sin_theta = np.sin(theta)
vx = inputs['vx']
vy = inputs['vy']
m = inputs['m']
F_T = inputs['thrust']
Isp = inputs['Isp']
y = inputs['y']

self.add_subsystem('atmos',
LogAtmosphereComp(num_nodes=nn, rho_ref=rho_ref, h_scale=h_scale))
g = self.options['g']
rho_ref = self.options['rho_ref']
h_scale = self.options['h_scale']

self.add_subsystem('eom', LaunchVehicle2DEOM(num_nodes=nn, central_body=cb))
CDA = self.options['CD'] * self.options['S']

self.connect('atmos.rho', 'eom.rho')
outputs['rho'] = rho_ref * np.exp(-y / h_scale)
outputs['xdot'] = vx
outputs['ydot'] = vy
outputs['vxdot'] = (F_T * cos_theta - 0.5 * CDA * outputs['rho'] * vx**2) / m
outputs['vydot'] = (F_T * sin_theta - 0.5 * CDA * outputs['rho'] * vy**2) / m - g
outputs['mdot'] = -F_T / (g * Isp)
Loading

0 comments on commit 41c23c5

Please sign in to comment.