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

Disabled reports by default for subproblems under simulate and ExplicitShooting #741

Merged
merged 8 commits into from
Apr 22, 2022
Original file line number Diff line number Diff line change
@@ -0,0 +1,230 @@
import pathlib
import unittest

import openmdao.api as om
import dymos as dm

from openmdao.utils.testing_utils import use_tempdirs
from dymos.examples.brachistochrone.brachistochrone_ode import BrachistochroneODE


@use_tempdirs
class TestSimulateReportsToggle(unittest.TestCase):

def test_no_sim_reports(self):
p = om.Problem(model=om.Group())

p.driver = om.pyOptSparseDriver()
p.driver.options['optimizer'] = 'SLSQP'
p.driver.declare_coloring(tol=1.0E-12)

t = dm.Radau(num_segments=10, order=3)

traj = dm.Trajectory()
phase = dm.Phase(ode_class=BrachistochroneODE, transcription=t)
p.model.add_subsystem('traj0', traj)
traj.add_phase('phase0', phase)

# p.model.add_subsystem('traj0', traj)

phase.set_time_options(fix_initial=True, duration_bounds=(.5, 10))

phase.add_state('x', fix_initial=True, fix_final=False)
phase.add_state('y', fix_initial=True, fix_final=False)

# Note that by omitting the targets here Dymos will automatically attempt to connect
# to a top-level input named 'v' in the ODE, and connect to nothing if it's not found.
phase.add_state('v', fix_initial=True, fix_final=False)

phase.add_control('theta',
continuity=True, rate_continuity=True,
units='deg', lower=0.01, upper=179.9)

phase.add_parameter('g', targets=['g'], units='m/s**2')

phase.add_boundary_constraint('x', loc='final', equals=10)
phase.add_boundary_constraint('y', loc='final', equals=5)
# Minimize time at the end of the phase
phase.add_objective('time_phase', loc='final', scaler=10)

p.setup()

phase.set_simulate_options(method='RK23')

p['traj0.phase0.t_initial'] = 0.0
p['traj0.phase0.t_duration'] = 2.0

p['traj0.phase0.states:x'] = phase.interp('x', [0, 10])
p['traj0.phase0.states:y'] = phase.interp('y', [10, 5])
p['traj0.phase0.states:v'] = phase.interp('v', [0, 9.9])
p['traj0.phase0.controls:theta'] = phase.interp('theta', [5, 100])
p['traj0.phase0.parameters:g'] = 9.80665

dm.run_problem(p, run_driver=True, simulate=True)

report_dir = pathlib.Path.cwd() / 'reports'
report_subdirs = [e for e in report_dir.iterdir() if e.is_dir()]

# Test that a report subdir was made
self.assertEqual(len(report_subdirs), 1)

def test_make_sim_reports(self):
p = om.Problem(model=om.Group())

p.driver = om.pyOptSparseDriver()
p.driver.options['optimizer'] = 'SLSQP'
p.driver.declare_coloring(tol=1.0E-12)

t = dm.Radau(num_segments=10, order=3)

traj = dm.Trajectory()
phase = dm.Phase(ode_class=BrachistochroneODE, transcription=t)
p.model.add_subsystem('traj0', traj)
traj.add_phase('phase0', phase)

# p.model.add_subsystem('traj0', traj)

phase.set_time_options(fix_initial=True, duration_bounds=(.5, 10))

phase.add_state('x', fix_initial=True, fix_final=False)
phase.add_state('y', fix_initial=True, fix_final=False)

# Note that by omitting the targets here Dymos will automatically attempt to connect
# to a top-level input named 'v' in the ODE, and connect to nothing if it's not found.
phase.add_state('v', fix_initial=True, fix_final=False)

phase.add_control('theta',
continuity=True, rate_continuity=True,
units='deg', lower=0.01, upper=179.9)

phase.add_parameter('g', targets=['g'], units='m/s**2')

phase.add_boundary_constraint('x', loc='final', equals=10)
phase.add_boundary_constraint('y', loc='final', equals=5)
# Minimize time at the end of the phase
phase.add_objective('time_phase', loc='final', scaler=10)

p.setup()

phase.set_simulate_options(method='RK23')

p['traj0.phase0.t_initial'] = 0.0
p['traj0.phase0.t_duration'] = 2.0

p['traj0.phase0.states:x'] = phase.interp('x', [0, 10])
p['traj0.phase0.states:y'] = phase.interp('y', [10, 5])
p['traj0.phase0.states:v'] = phase.interp('v', [0, 9.9])
p['traj0.phase0.controls:theta'] = phase.interp('theta', [5, 100])
p['traj0.phase0.parameters:g'] = 9.80665

dm.run_problem(p, run_driver=True, simulate=True, simulate_kwargs={'reports': True})

reports_dir = pathlib.Path.cwd() / 'reports'
report_subdirs = [e for e in reports_dir.iterdir() if e.is_dir()]

# Test that a report subdir was made
# There is the nominal problem, the simulation problem, and a subproblem for each segment in the simulation.
self.assertEqual(len(report_subdirs), 12)


@use_tempdirs
class TestExplicitShootingReportsToggle(unittest.TestCase):

def test_no_subprob_reports(self):
prob = om.Problem()

prob.driver = om.pyOptSparseDriver(optimizer='SLSQP')

tx = dm.ExplicitShooting(num_segments=3, grid='gauss-lobatto',
method='rk4', order=5,
num_steps_per_segment=5,
compressed=False)

phase = dm.Phase(ode_class=BrachistochroneODE, transcription=tx)

phase.set_time_options(units='s', fix_initial=True, duration_bounds=(1.0, 10.0))

# automatically discover states
phase.set_state_options('x', fix_initial=True)
phase.set_state_options('y', fix_initial=True)
phase.set_state_options('v', fix_initial=True)

phase.add_parameter('g', val=1.0, units='m/s**2', opt=True, lower=1, upper=9.80665)
phase.add_control('theta', val=45.0, units='deg', opt=True, lower=1.0E-6, upper=179.9,
ref=90., rate2_continuity=True)

phase.add_boundary_constraint('x', loc='final', equals=10.0)
phase.add_boundary_constraint('y', loc='final', equals=5.0)

prob.model.add_subsystem('phase0', phase)

phase.add_objective('time', loc='final')

prob.setup(force_alloc_complex=True)

prob.set_val('phase0.t_initial', 0.0)
prob.set_val('phase0.t_duration', 2)
prob.set_val('phase0.states:x', 0.0)
prob.set_val('phase0.states:y', 10.0)
prob.set_val('phase0.states:v', 1.0E-6)
prob.set_val('phase0.parameters:g', 1.0, units='m/s**2')
prob.set_val('phase0.controls:theta', phase.interp('theta', ys=[0.01, 90]), units='deg')

dm.run_problem(prob, run_driver=True, simulate=False)

reports_dir = pathlib.Path.cwd() / 'reports'
report_subdirs = [e for e in reports_dir.iterdir() if e.is_dir()]

# Test that a report subdir was made
# There is the nominal problem, the simulation problem, and a subproblem for each segment in the simulation.
self.assertEqual(len(report_subdirs), 1)

def test_make_subprob_reports(self):
prob = om.Problem()

prob.driver = om.pyOptSparseDriver(optimizer='SLSQP')

tx = dm.ExplicitShooting(num_segments=3, grid='gauss-lobatto',
method='rk4', order=5,
num_steps_per_segment=5,
compressed=False,
subprob_reports=True)

phase = dm.Phase(ode_class=BrachistochroneODE, transcription=tx)

phase.set_time_options(units='s', fix_initial=True, duration_bounds=(1.0, 10.0))

# automatically discover states
phase.set_state_options('x', fix_initial=True)
phase.set_state_options('y', fix_initial=True)
phase.set_state_options('v', fix_initial=True)

phase.add_parameter('g', val=1.0, units='m/s**2', opt=True, lower=1, upper=9.80665)
phase.add_control('theta', val=45.0, units='deg', opt=True, lower=1.0E-6, upper=179.9,
ref=90., rate2_continuity=True)

phase.add_boundary_constraint('x', loc='final', equals=10.0)
phase.add_boundary_constraint('y', loc='final', equals=5.0)

prob.model.add_subsystem('phase0', phase)

phase.add_objective('time', loc='final')

prob.setup(force_alloc_complex=True)

prob.set_val('phase0.t_initial', 0.0)
prob.set_val('phase0.t_duration', 2)
prob.set_val('phase0.states:x', 0.0)
prob.set_val('phase0.states:y', 10.0)
prob.set_val('phase0.states:v', 1.0E-6)
prob.set_val('phase0.parameters:g', 1.0, units='m/s**2')
prob.set_val('phase0.controls:theta', phase.interp('theta', ys=[0.01, 90]), units='deg')

dm.run_problem(prob, run_driver=True, simulate=False)

reports_dir = pathlib.Path.cwd() / 'reports'
report_subdirs = [e for e in reports_dir.iterdir() if e.is_dir()]

# Test that a report subdir was made
# There is the nominal problem, a subproblem for integration, and a subproblem for the derivatives.
self.assertEqual(len(report_subdirs), 3)
8 changes: 6 additions & 2 deletions dymos/phase/phase.py
Original file line number Diff line number Diff line change
Expand Up @@ -1852,7 +1852,8 @@ def interp(self, name=None, ys=None, xs=None, nodes=None, kind='linear', axis=0)
return res

def get_simulation_phase(self, times_per_seg=None, method=_unspecified, atol=_unspecified,
rtol=_unspecified, first_step=_unspecified, max_step=_unspecified):
rtol=_unspecified, first_step=_unspecified, max_step=_unspecified,
reports=False):
"""
Return a SolveIVPPhase instance.

Expand All @@ -1877,6 +1878,8 @@ def get_simulation_phase(self, times_per_seg=None, method=_unspecified, atol=_un
Initial step size for the integration.
max_step : float or _unspecified
Maximum step size for the integration.
reports : bool or None or str or Sequence
The reports setting for the subproblem run under each simulation segment.

Returns
-------
Expand All @@ -1890,7 +1893,8 @@ def get_simulation_phase(self, times_per_seg=None, method=_unspecified, atol=_un

sim_phase = dm.Phase(from_phase=self,
transcription=SolveIVP(grid_data=t.grid_data,
output_nodes_per_seg=times_per_seg))
output_nodes_per_seg=times_per_seg,
reports=reports))

# Copy over any simulation options from the simulate call. The fallback will be to
# phase.simulate_options, which are copied from the original phase.
Expand Down
8 changes: 5 additions & 3 deletions dymos/trajectory/trajectory.py
Original file line number Diff line number Diff line change
Expand Up @@ -1072,7 +1072,7 @@ def _print_constraints(phase, outstream):

def simulate(self, times_per_seg=10, method=_unspecified, atol=_unspecified, rtol=_unspecified,
first_step=_unspecified, max_step=_unspecified, record_file=None, case_prefix=None,
reset_iter_counts=True):
reset_iter_counts=True, reports=False):
"""
Simulate the Trajectory using scipy.integrate.solve_ivp.

Expand All @@ -1098,6 +1098,8 @@ def simulate(self, times_per_seg=10, method=_unspecified, atol=_unspecified, rto
Prefix to prepend to coordinates when recording.
reset_iter_counts : bool
If True and model has been run previously, reset all iteration counters.
reports : bool or None or str or Sequence
Reports setting for the subproblems run under simualate.

Returns
-------
Expand All @@ -1111,12 +1113,12 @@ def simulate(self, times_per_seg=10, method=_unspecified, atol=_unspecified, rto
for name, phs in self._phases.items():
sim_phs = phs.get_simulation_phase(times_per_seg=times_per_seg, method=method,
atol=atol, rtol=rtol, first_step=first_step,
max_step=max_step)
max_step=max_step, reports=reports)
sim_traj.add_phase(name, sim_phs)

sim_traj.parameter_options.update(self.parameter_options)

sim_prob = om.Problem(model=om.Group())
sim_prob = om.Problem(model=om.Group(), reports=reports)

traj_name = self.name if self.name else 'sim_traj'
sim_prob.model.add_subsystem(traj_name, sim_traj)
Expand Down
5 changes: 4 additions & 1 deletion dymos/transcriptions/explicit_shooting/explicit_shooting.py
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,8 @@ def initialize(self):
str(list(rk_methods.keys())))
self.options.declare('num_steps_per_segment', types=int,
default=10, desc='Number of integration steps in each segment')
self.options.declare('subprob_reports', default=False,
desc='Controls the reports made when running the subproblems for ExplicitShooting')

def init_grid(self):
"""
Expand Down Expand Up @@ -166,7 +168,8 @@ def setup_ode(self, phase):
num_steps_per_segment=self.options['num_steps_per_segment'],
grid_data=self.grid_data,
ode_init_kwargs=phase.options['ode_init_kwargs'],
standalone_mode=False)
standalone_mode=False,
reports=self.options['subprob_reports'])

phase.add_subsystem(name='integrator', subsys=integrator_comp, promotes_inputs=['*'])

Expand Down
9 changes: 6 additions & 3 deletions dymos/transcriptions/explicit_shooting/rk_integration_comp.py
Original file line number Diff line number Diff line change
Expand Up @@ -94,6 +94,8 @@ class RKIntegrationComp(om.ExplicitComponent):
standalone_mode : bool
When True, this component will perform its configuration during setup. This is useful
for unittesting this component when not embedded in a larger system.
reports : bool or None or str or Sequence
Controls the reports generated by the subproblems used during the integration.
**kwargs : dict
Additional keyword arguments passed to Group.

Expand All @@ -105,7 +107,7 @@ class RKIntegrationComp(om.ExplicitComponent):
def __init__(self, ode_class, time_options=None,
state_options=None, parameter_options=None, control_options=None,
polynomial_control_options=None, timeseries_options=None,
grid_data=None, standalone_mode=True, **kwargs):
grid_data=None, standalone_mode=True, reports=False, **kwargs):
super().__init__(**kwargs)
self.ode_class = ode_class
self.time_options = time_options
Expand All @@ -118,6 +120,7 @@ def __init__(self, ode_class, time_options=None,
self._deriv_subprob = None
self._grid_data = grid_data
self._DTYPE = float
self._reports = reports

self._inputs_cache = None

Expand Down Expand Up @@ -153,7 +156,7 @@ def _setup_subprob(self):
rk = rk_methods[self.options['method']]
num_stages = len(rk['b'])

self._eval_subprob = p = om.Problem(comm=self.comm)
self._eval_subprob = p = om.Problem(comm=self.comm, reports=self._reports)
p.model.add_subsystem('ode_eval',
ODEEvaluationGroup(self.ode_class, self.time_options,
self.state_options,
Expand All @@ -168,7 +171,7 @@ def _setup_subprob(self):
p.setup(force_alloc_complex=True if self._DTYPE is complex else False)
p.final_setup()

self._deriv_subprob = p = om.Problem(comm=self.comm)
self._deriv_subprob = p = om.Problem(comm=self.comm, reports=self._reports)
p.model.add_subsystem('ode_eval',
ODEEvaluationGroup(self.ode_class, self.time_options,
self.state_options,
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -30,9 +30,12 @@ class ODEIntegrationInterface(object):
The parameter options for the phase being simulated.
ode_init_kwargs : dict
Keyword argument dictionary passed to the ODE at initialization.
reports : bool or None or str or Sequence
The reports argument to be passed to the subproblems. By default, no subproblem reports are generated.
"""
def __init__(self, ode_class, time_options, state_options, control_options,
polynomial_control_options, parameter_options, ode_init_kwargs=None):
polynomial_control_options, parameter_options, ode_init_kwargs=None,
reports=False):

# Get the state vector. This isn't necessarily ordered
# so just pick the default ordering and go with it.
Expand Down Expand Up @@ -67,7 +70,8 @@ def __init__(self, ode_class, time_options, state_options, control_options,
control_options=control_options,
polynomial_control_options=polynomial_control_options,
parameter_options=parameter_options,
ode_init_kwargs=ode_init_kwargs))
ode_init_kwargs=ode_init_kwargs),
reports=reports)

def _unpack_state_vec(self, x):
"""
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -72,6 +72,8 @@ def initialize(self):
'segment. If an int (n) then results are provided at n '
'equally distributed points in time within each segment.')

self.options.declare('reports', default=False, desc='Reports setting for the subproblem.')

self.recording_options['options_excludes'] = ['ode_integration_interface']

def configure_io(self):
Expand Down Expand Up @@ -106,7 +108,8 @@ def configure_io(self):
control_options=self.options['control_options'],
polynomial_control_options=self.options['polynomial_control_options'],
parameter_options=self.options['parameter_options'],
ode_init_kwargs=self.options['ode_init_kwargs'])
ode_init_kwargs=self.options['ode_init_kwargs'],
reports=self.options['reports'])

self.add_input(name='time', val=np.ones(nnps_i),
units=self.options['time_options']['units'],
Expand Down