Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Added automatic parameter generation to dymos #17

Open
wants to merge 11 commits into
base: master
Choose a base branch
from
Open
33 changes: 22 additions & 11 deletions dymos/examples/balanced_field/balanced_field_ode.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,17 +22,28 @@ def setup(self):
nn = self.options['num_nodes']

# Scalar (constant) inputs
self.add_input('rho', val=1.225, desc='atmospheric density at runway', units='kg/m**3')
self.add_input('S', val=124.7, desc='aerodynamic reference area', units='m**2')
self.add_input('CD0', val=0.03, desc='zero-lift drag coefficient', units=None)
self.add_input('CL0', val=0.5, desc='zero-alpha lift coefficient', units=None)
self.add_input('CL_max', val=2.0, desc='maximum lift coefficient for linear fit', units=None)
self.add_input('alpha_max', val=np.radians(10), desc='angle of attack at CL_max', units='rad')
self.add_input('h_w', val=1.0, desc='height of the wing above the CG', units='m')
self.add_input('AR', val=9.45, desc='wing aspect ratio', units=None)
self.add_input('e', val=0.801, desc='Oswald span efficiency factor', units=None)
self.add_input('span', val=35.7, desc='Wingspan', units='m')
self.add_input('T', val=1.0, desc='thrust', units='N')
self.add_input('rho', val=1.225, desc='atmospheric density at runway', units='kg/m**3',
tags=['dymos.static_target'])
self.add_input('S', val=124.7, desc='aerodynamic reference area', units='m**2',
tags=['dymos.static_target'])
self.add_input('CD0', val=0.03, desc='zero-lift drag coefficient', units=None,
tags=['dymos.static_target'])
self.add_input('CL0', val=0.5, desc='zero-alpha lift coefficient', units=None,
tags=['dymos.static_target'])
self.add_input('CL_max', val=2.0, desc='maximum lift coefficient for linear fit', units=None,
tags=['dymos.static_target'])
self.add_input('alpha_max', val=np.radians(10), desc='angle of attack at CL_max', units='rad',
tags=['dymos.static_target'])
self.add_input('h_w', val=1.0, desc='height of the wing above the CG', units='m',
tags=['dymos.static_target'])
self.add_input('AR', val=9.45, desc='wing aspect ratio', units=None,
tags=['dymos.static_target'])
self.add_input('e', val=0.801, desc='Oswald span efficiency factor', units=None,
tags=['dymos.static_target'])
self.add_input('span', val=35.7, desc='Wingspan', units='m',
tags=['dymos.static_target'])
self.add_input('T', val=1.0, desc='thrust', units='N',
tags=['dymos.static_target'])

# Dynamic inputs (can assume a different value at every node)
self.add_input('m', shape=(nn,), desc='aircraft mass', units='kg')
Expand Down
39 changes: 37 additions & 2 deletions dymos/phase/phase.py
Original file line number Diff line number Diff line change
Expand Up @@ -70,6 +70,9 @@ def __init__(self, from_phase=None, **kwargs):
self.timeseries_ec_vars = {}
self.timeseries_options = PhaseTimeseriesOptionsDictionary()

# A dictionary that keeps track of all targets to which we've connected in the ODE.
self._ode_connections = {}

# Dictionaries of variable options that are set by the user via the API
# These will be applied over any defaults specified by decorators on the ODE
if from_phase is None:
Expand Down Expand Up @@ -1839,8 +1842,7 @@ def setup(self):
if self.polynomial_control_options:
transcription.setup_polynomial_controls(self)

if self.parameter_options:
transcription.setup_parameters(self)
transcription.setup_parameters(self)

transcription.setup_states(self)
self._check_ode()
Expand Down Expand Up @@ -1889,6 +1891,8 @@ def configure(self):

transcription.configure_ode(self)

transcription.configure_automatic_parameters(self)

transcription.configure_defects(self)

_configure_constraint_introspection(self)
Expand Down Expand Up @@ -2024,6 +2028,37 @@ def _check_parameter_options(self):
f"phase '{self.name}': {', '.join(invalid_options)}",
RuntimeWarning)

def _connect_to_ode(self, src_name, ode_tgt_name,
src_indices=None, flat_src_indices=None):
"""
Connect source src_name to target tgt_name in the ODE and cache the name of the ODE target.

Parameters
----------
src_name : str
Name of the source variable to connect.
ode_tgt_name : str or [str, ... ] or (str, ...)
Name of the target variable(s) to connect relative to the ODE system.
src_indices : int or list of ints or tuple of ints or int ndarray or Iterable or None
The global indices of the source variable to transfer data from.
The shapes of the target and src_indices must match, and form of the
entries within is determined by the value of 'flat_src_indices'.
flat_src_indices : bool
If True, each entry of src_indices is assumed to be an index into the
flattened source. Otherwise it must be a tuple or list of size equal
to the number of dimensions of the source.
"""
ode_paths = self.options['transcription']._ode_paths

for ode_path, default_src_idxs in ode_paths.items():
if isinstance(ode_tgt_name, str):
ode_tgt_name = [ode_tgt_name]
for tgt in ode_tgt_name:
src_idxs = None if src_indices is None else src_indices[ode_path]
super().connect(src_name=src_name, tgt_name=f'{ode_path}.{tgt}',
src_indices=src_idxs, flat_src_indices=flat_src_indices)
self._ode_connections[tgt] = src_name

def interpolate(self, xs=None, ys=None, nodes='all', kind='linear', axis=0):
"""
Return an array of values on interpolated to the given node subset of the phase.
Expand Down
12 changes: 8 additions & 4 deletions dymos/phase/test/test_sized_input_parameters.py
Original file line number Diff line number Diff line change
Expand Up @@ -55,7 +55,8 @@ def setup(self):
phase.set_time_options(initial_bounds=(0.0, 100.0), duration_bounds=(0., 100.), units='s')

phase.add_state('h', fix_initial=True, fix_final=True, lower=0.0, units='m', rate_source='eom.h_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s', rate_source='eom.v_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s',
targets=['eom.v'], rate_source='eom.v_dot')

phase.add_parameter('m', val=[[1, 2], [3, 4]], units='kg', targets='sum.m')

Expand Down Expand Up @@ -117,7 +118,8 @@ def setup(self):
phase.set_time_options(initial_bounds=(0.0, 100.0), duration_bounds=(0., 100.))

phase.add_state('h', fix_initial=True, fix_final=True, lower=0.0, units='m', rate_source='eom.h_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s', rate_source='eom.v_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s',
targets=['eom.v'], rate_source='eom.v_dot')

phase.add_parameter('m', val=[[1, 2], [3, 4]], units='kg', targets='sum.m', static_target=True)

Expand Down Expand Up @@ -178,7 +180,8 @@ def setup(self):
phase.set_time_options(initial_bounds=(0.0, 100.0), duration_bounds=(0., 100.), units='s')

phase.add_state('h', fix_initial=True, fix_final=True, lower=0.0, units='m', rate_source='eom.h_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s', rate_source='eom.v_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s',
targets=['eom.v'], rate_source='eom.v_dot')

phase.add_parameter('m', val=[[1, 2], [3, 4]], units='kg', targets='sum.m')

Expand Down Expand Up @@ -246,7 +249,8 @@ def setup(self):
phase.set_time_options(initial_bounds=(0.0, 100.0), duration_bounds=(0., 100.), units='s')

phase.add_state('h', fix_initial=True, fix_final=True, lower=0.0, units='m', rate_source='eom.h_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s', rate_source='eom.v_dot')
phase.add_state('v', fix_initial=True, fix_final=False, units='m/s',
targets=['eom.v'], rate_source='eom.v_dot')

phase.add_parameter('m', val=[[1, 2], [3, 4]], units='kg', targets='sum.m', static_target=True)

Expand Down
2 changes: 1 addition & 1 deletion dymos/trajectory/trajectory.py
Original file line number Diff line number Diff line change
Expand Up @@ -599,7 +599,7 @@ def _update_linkage_options_configure(self, linkage_options):
units[i] = phases[i].parameter_options[vars[i]]['units']
shapes[i] = phases[i].parameter_options[vars[i]]['shape']
else:
rhs_source = phases[i].options['transcription']._rhs_source
rhs_source = phases[i].options['transcription']._ode_paths[0]
sources[i] = f'{rhs_source}.{vars[i]}'
try:
meta = get_source_metadata(phases[i]._get_subsystem(rhs_source), vars[i], user_units=units[i],
Expand Down
55 changes: 24 additions & 31 deletions dymos/transcriptions/analytic/analytic.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,13 +2,11 @@
import openmdao.api as om

from ..transcription_base import TranscriptionBase
from ...utils.introspection import configure_analytic_states_introspection, get_promoted_vars, get_targets, \
get_source_metadata, configure_analytic_states_discovery
from ...utils.introspection import configure_analytic_states_introspection, get_promoted_vars, get_source_metadata, configure_analytic_states_discovery
from ...utils.indexing import get_src_indices_by_row
from ..grid_data import GridData
from .analytic_timeseries_output_comp import AnalyticTimeseriesOutputComp
from ..common import TimeComp, TimeseriesOutputGroup
from ..._options import options as dymos_options


class Analytic(TranscriptionBase):
Expand All @@ -25,7 +23,7 @@ class Analytic(TranscriptionBase):
"""
def __init__(self, **kwargs):
super(Analytic, self).__init__(**kwargs)
self._rhs_source = 'rhs'
self._ode_paths = {'rhs': self.grid_data.subset_node_indices['all']}

def init_grid(self):
"""
Expand Down Expand Up @@ -56,7 +54,10 @@ def setup_time(self, phase):
initial_val=phase.time_options['initial_val'],
duration_val=phase.time_options['duration_val'])

phase.add_subsystem('time', time_comp, promotes_inputs=['*'], promotes_outputs=['*'])
phase.add_subsystem('time', time_comp)

phase.connect('t_initial_val', 'time.t_initial')
phase.connect('t_duration_val', 'time.t_duration')

def configure_time(self, phase):
"""
Expand All @@ -75,16 +76,16 @@ def configure_time(self, phase):
phase.time.configure_io()

options = phase.time_options
ode = phase._get_subsystem(self._rhs_source)
ode = phase._get_subsystem(list(self._ode_paths.keys())[0])
ode_inputs = get_promoted_vars(ode, iotypes='input')

# The tuples here are (name, user_specified_targets, dynamic)
for name, targets, dynamic in [('t', options['targets'], True),
('t_phase', options['time_phase_targets'], True)]:
for name, targets, dynamic in [('time.t', options['targets'], True),
('time.t_phase', options['time_phase_targets'], True)]:
if targets:
src_idxs = self.grid_data.subset_node_indices['all'] if dynamic else None
phase.connect(name, [f'rhs.{t}' for t in targets], src_indices=src_idxs,
flat_src_indices=True if dynamic else None)
phase._connect_to_ode(name, targets, src_indices={'rhs': src_idxs},
flat_src_indices=True if dynamic else None)

for name, targets in [('t_initial', options['t_initial_targets']),
('t_duration', options['t_duration_targets'])]:
Expand All @@ -94,18 +95,12 @@ def configure_time(self, phase):
if shape == (1,):
src_idxs = None
flat_src_idxs = None
src_shape = None
else:
src_idxs = np.zeros(self.grid_data.subset_num_nodes['all'])
flat_src_idxs = True
src_shape = (1,)

phase.promotes('rhs', inputs=[(t, name)], src_indices=src_idxs,
flat_src_indices=flat_src_idxs, src_shape=src_shape)
if targets:
phase.set_input_defaults(name=name,
val=np.ones((1,)),
units=options['units'])
phase._connect_to_ode(name, t, src_indices={'rhs': src_idxs},
flat_src_indices=flat_src_idxs)

def setup_controls(self, phase):
"""
Expand Down Expand Up @@ -323,7 +318,7 @@ def setup_timeseries_outputs(self, phase):
timeseries_group = TimeseriesOutputGroup(has_expr=has_expr, timeseries_output_comp=timeseries_comp)
phase.add_subsystem(name, subsys=timeseries_group)

phase.connect('dt_dstau', f'{name}.dt_dstau', flat_src_indices=True)
phase.connect('time.dt_dstau', f'{name}.dt_dstau', flat_src_indices=True)

def configure_defects(self, phase):
"""
Expand Down Expand Up @@ -372,11 +367,11 @@ def _get_timeseries_var_source(self, var, output_name, phase):

# Determine the path to the variable
if var_type == 't':
path = 't'
path = 'time.t'
src_units = time_units
src_shape = (1,)
elif var_type == 't_phase':
path = 't_phase'
path = 'time.t_phase'
src_units = time_units
src_shape = (1,)
elif var_type == 'parameter':
Expand Down Expand Up @@ -435,8 +430,7 @@ def get_parameter_connections(self, name, phase):
src_idxs = get_src_indices_by_row(src_idxs_raw, options['shape'])
src_idxs = np.squeeze(src_idxs, axis=0)

rhs_tgts = [f'rhs.{t}' for t in options['targets']]
connection_info.append((rhs_tgts, (src_idxs,)))
connection_info.append((options['targets'], {'rhs': src_idxs}))

return connection_info

Expand Down Expand Up @@ -511,21 +505,20 @@ def _get_objective_src(self, var, loc, phase, ode_outputs=None):
var_type = phase.classify_var(var)

if ode_outputs is None:
ode_outputs = get_promoted_vars(phase._get_subsystem(self._rhs_source), 'output')
ode_outputs = get_promoted_vars(phase._get_subsystem(list(self._ode_paths.keys())[0]), 'output')

if var_type == 't':
shape = (1,)
units = time_units
linear = True
constraint_path = 't'
obj_path = 't'
elif var_type == 't_phase':
shape = (1,)
units = time_units
linear = True
constraint_path = 't_phase'
obj_path = 't_phase'
elif var_type == 'state':
constraint_path = f'{self._rhs_source}.{var}'
src_path = phase.state_options[var]['source']
obj_path = f'rhs.{var}'
meta = get_source_metadata(ode_outputs, var, user_units=None, user_shape=None)
shape = meta['shape']
units = meta['units']
Expand All @@ -534,13 +527,13 @@ def _get_objective_src(self, var, loc, phase, ode_outputs=None):
shape = phase.parameter_options[var]['shape']
units = phase.parameter_options[var]['units']
linear = True
constraint_path = f'parameter_vals:{var}'
obj_path = f'parameter_vals:{var}'
else:
# Failed to find variable, assume it is in the ODE. This requires introspection.
constraint_path = f'{self._rhs_source}.{var}'
obj_path = f'{self._ode_paths[0]}.{var}'
meta = get_source_metadata(ode_outputs, var, user_units=None, user_shape=None)
shape = meta['shape']
units = meta['units']
linear = False

return constraint_path, shape, units, linear
return obj_path, shape, units, linear
Loading