Merge pull request #203 from OpenMDAO/develop

Develop

benchmark/benchmark_brachistochrone.py

@@ -1,20 +1,22 @@
 
 import unittest
 
-from openmdao.api import Problem, Group, pyOptSparseDriver, DirectSolver
+import openmdao.api as om
 from openmdao.utils.assert_utils import assert_rel_error
 import dymos as dm
 from dymos.examples.brachistochrone import BrachistochroneODE
 
 
 def brachistochrone_min_time(transcription='gauss-lobatto', num_segments=8, transcription_order=3,
                              compressed=True, optimizer='SLSQP', simul_derivs=True):
-    p = Problem(model=Group())
+    p = om.Problem(model=om.Group())
 
     # if optimizer == 'SNOPT':
-    p.driver = pyOptSparseDriver()
+    p.driver = om.pyOptSparseDriver()
     p.driver.options['optimizer'] = optimizer
-    p.driver.options['dynamic_simul_derivs'] = simul_derivs
+
+    if simul_derivs:
+        p.driver.declare_coloring()
 
     if transcription == 'gauss-lobatto':
         t = dm.GaussLobatto(num_segments=num_segments,
@@ -49,7 +51,7 @@ def brachistochrone_min_time(transcription='gauss-lobatto', num_segments=8, tran
     # Minimize time at the end of the phase
     phase.add_objective('time_phase', loc='final', scaler=10)
 
-    p.model.linear_solver = DirectSolver()
+    p.model.linear_solver = om.DirectSolver()
     p.setup(check=True)
 
     p['phase0.t_initial'] = 0.0

dymos/docs/examples/brachistochrone.rst

@@ -57,7 +57,7 @@ Finally, note that we are specifying rows and columns when declaring the partial
 Since our inputs and outputs are scalars *at each point in time*, and the value at an input at
 one time only directly impacts the values of an output at the same point in time, the partial
 derivative jacobian will be diagonal.  Specifying the partial derivatives as being sparse
-greatly improves the performance of Dymos when the driver option ``dynamic_simul_derivs`` is used.
+greatly improves the performance of Dymos when the driver function ``declare_coloring`` is used.
 Using a sparse optimizer like SNOPT or IPOPT can provide significant addition improvements in
 performance.
 

dymos/docs/examples/figures/ssto_linear_tangent_xdsm.py

@@ -2,12 +2,12 @@
 
 from pyxdsm.XDSM import XDSM
 
-from openmdao.api import Problem
+import openmdao.api as om
 from dymos.examples.ssto.launch_vehicle_linear_tangent_ode import LaunchVehicleLinearTangentODE
 
 
 def main():  # pragma: no cover
-    p = Problem()
+    p = om.Problem()
     p.model = LaunchVehicleLinearTangentODE(num_nodes=1)
     p.setup()
 

dymos/docs/examples/figures/ssto_xdsm.py

@@ -2,12 +2,12 @@
 
 from pyxdsm.XDSM import XDSM
 
-from openmdao.api import Problem
+import openmdao.api as om
 from dymos.examples.ssto.launch_vehicle_ode import LaunchVehicleODE
 
 
 def main():  # pragma: no cover
-    p = Problem()
+    p = om.Problem()
     p.model = LaunchVehicleODE(num_nodes=1)
     p.setup()
 

dymos/docs/examples/multibranch_trajectory.rst

@@ -25,7 +25,7 @@ phases are connected to the output of a single phase.  This way you can simulate
 trajectory paths in the same model. For this example, we will start with a single phase ("phase0")
 that simulates the model for one hour. Three follow-on phases will be linked to the output of the
 first phase: "phase1" will run as normal, "phase1_bfail" will fail one of the battery cells, and
-"phase1_bfail" will fail a motor. All three of these phases start where "phase0" leaves off, so
+"phase1_mfail" will fail a motor. All three of these phases start where "phase0" leaves off, so
 they share the same initial time and state of charge.
 
 

dymos/docs/feature_reference/feature_reference.rst

@@ -18,4 +18,5 @@ problems.
     phases/phases
     trajectories
     timeseries
+    tandem_phases
     simultaneous_derivs

dymos/docs/feature_reference/simultaneous_derivs.rst

@@ -38,11 +38,30 @@ Step 1: Using OpenMDAO's Simul-Coloring Capability
 
 OpenMDAO supports dynamic simul-coloring, meaning it can automatically run the Jacobian coloring
 algorithm before handing the problem to the optimizer.  To enable this capability, simply
-add the following lines to the driver.
+add the following line to the driver.
 
 .. code-block:: python
 
-    driver.options['dynamic_simul_derivs'] = True
+    driver.declare_coloring()
+
+By default the coloring algorithm will attempt to determine the sparsity pattern of the total jacobian
+by filling the partial jacobian matrices with random noise and searching for nonzeros in the resulting
+total jacobian.  At times this might report that it failed to converge on a number of nonzero entries.
+This is due to the introduction of noise during the matrix inversion by the system's linear solver.
+This can be remedied by using a different linear solver, such as PETScKrylov, or by telling the
+coloring algorithm to accept a given tolerance on the nonzero elements rather than trying to determine
+it automatically.  This can be accomplished with the following options to declare coloring:
+
+.. code-block:: python
+
+    driver.declare_coloring(tol=1.0E-12, orders=None)
+
+Setting `orders` to None prevents the automatic tolerance detection.  The value of `tol` is up to
+the user.  If it is too large, then some nonzero values will erroneously be determined to be zeros
+and the total derivative will be incorrect.  If `tol` is too small then the sparsity pattern may be
+overly conservative and degrade performance somewhat.  We recommend letting the coloring algorithm
+detect the sparsity automatically and only resorting to a fixed tolerance if necessary.
+
 
 The simul_coloring script outputs the following information about our problem:
 

dymos/docs/feature_reference/tandem_phases.rst

@@ -0,0 +1,52 @@
+==========================================================
+Tandem Phases:  Using different ODE simultaneously in time
+==========================================================
+
+Complex models sometimes encounter state variables which are best simulated on different time
+scales, with some state variables changing quickly (fast variables) and some evolving slowly (slow variables).
+For instance, and aircraft trajectory optimization which includes vehicle component temperatures might
+see relatively gradual changes in altitude over the course of a two hour flight while temperatures of
+some components seem to exhibit step-function-like behavior on the same scale.
+
+To accommodate both fast and slow variables in the same ODE, one would typically need to use a _dense_
+grid (with many segments/higher order segments).  This can be unnecessarily burdensome when there are
+many slow variables or evaluating their rates is particularly expensive.
+
+As a solution, Dymos allows the user to run two phases over the same range of times, where one
+phase may have a more sparse grid to accommodate the slow variables, and one has a more dense grid
+for the fast variables.
+
+To connect the two phases, state variable values are passed from the first (slow) phase to the second
+(fast) phase as non-optimal dynamic control variables.  These values are then used to evaluate the
+rates of the fast variables.  Since outputs from the first phase in generally will not fall on the
+appropriate grid points to be used by the second phase, interpolation is necessary.  This is one
+application of the interpolating timeseries component.
+
+In the following example, we solve the brachistochrone problem but do so to minimize the arclength
+of the resulting wire instead of the time required for the bead to travel along the wire.
+This is a trivial solution which should find a straight line from the starting point to the ending point.
+There are two phases involved, the first utilizes the standard ODE for the brachistochrone problem.
+The second integrates the arclength (:math:`S`) of the wire using the equation:
+
+.. math::
+
+    S = \int v \sin \theta  \sqrt{1 + \frac{1}{\tan^2 \theta}} \, dt
+
+.. embed-code::
+    dymos.examples.brachistochrone.doc.test_doc_brachistochrone_tandem_phases.BrachistochroneArclengthODE
+    :layout: code
+
+The trick is that the bead velocity (:math:`v`) is a state variable solved for in the first phase,
+and the wire angle (:math:`\theta`) is a control variable "owned" by the first phase.  In the
+second phase they are used as control variables with option ``opt=False`` so that their values are
+expected as inputs for the second phase.  We need to connect their values from the first phase
+to the second phase, at the :code:`control_input` node subset of the second phase.
+
+In the following example, we instantiate two phases and add an interpolating timeseries to the first phase
+which provides outputs at the :code:`control_input` nodes of the second phase.  Those values are
+then connected and the entire problem run. The result is that the position and velocity variables
+are solved on a relatively coarse grid while the arclength of the wire is solved on a much denser grid.
+
+.. embed-code::
+    dymos.examples.brachistochrone.doc.test_doc_brachistochrone_tandem_phases.TestDocTandemPhases.test_tandem_phases_for_docs
+    :layout: code, plot

dymos/docs/feature_reference/timeseries.rst

@@ -37,3 +37,14 @@ using the ``add_timeseries_output`` method on Phase.  These outputs are availabl
 .. automethod:: dymos.Phase.add_timeseries_output
     :noindex:
 
+Interpolated Timeseries Outputs
+===============================
+
+Sometimes a user may want to interpolate the results of a phase onto a different grid.  This is particularly
+useful in the context of tandem phases.  Additional timeseries may be added to a phase using the
+``add_timeseries`` method.  By default all timeseries will provide times, states, controls, and
+parameters on the specified output grid.  Adding other variables is accomplished using the
+``timeseries`` argument in the ``add_timeseries_output`` method.
+
+.. automethod:: dymos.Phase.add_timeseries
+    :noindex:

dymos/docs/feature_reference/trajectories.rst

@@ -61,9 +61,20 @@ at the trajectory level which maybe be connected to some output external to the
 When using Trajectory Design and Input parameters, their values are connected to each phase as an
 Input Parameter within the Phase.  Because ODEs in different phases may have different names
 for parameters (e.g. 'mass', 'm', 'm_total', etc) Dymos allows the user to specify the targeted
-ODE parameters on a phase-by-phase basis.
+ODE parameters on a phase-by-phase basis using the `custom_targets` option.  It can take on the
+following values.
 
-<EXAMPLE OF TARGETS>
+*  Left unspecified (`custom_targets = None`), a trajectory design or input parameter will be connected to a decorated ODE parameter of the same name in each phase.
+
+*  Otherwise, `custom_targets` is specified as a dictionary keyed by phase name.
+
+    * If the name of a phase is omitted from `custom_targets`, the trajectory design or input parameter will be connected to a decorated ODE parameter of the same name in that phase.
+
+    * If the phase is specified with a corresponding value of `None`, **the trajectory parameter will not be passed to the phase**.
+
+    * If the phase is specified and the corresponding value is a string, assume the given value is a decorated ODE parameter to which the trajectory parameter should be connected.
+
+    * If the phase is specified and the corresponding parameter value is a list, assume the list gives ODE-level targets to which the parameter should be connected in that phase.
 
 If a phase exists within the Trajectory that doesn't utilize the trajectory
 design/input parameters, it is simply ignored for the purposes of that phase.

dymos/docs/feature_reference/transcriptions/figures/gauss-lobatto/plot_high_order_gauss_lobatto.py

@@ -6,13 +6,13 @@
 matplotlib.use('Agg')
 import matplotlib.pyplot as plt
 
-from openmdao.api import Problem, Group
-from dymos import Phase, GaussLobatto
+import openmdao.api as om
+import dymos as dm
 from dymos.examples.brachistochrone.brachistochrone_ode import BrachistochroneODE
 
-p = Problem(model=Group())
-phase = Phase(ode_class=BrachistochroneODE,
-              transcription=GaussLobatto(num_segments=4, order=[3, 5, 3, 5]))
+p = om.Problem(model=om.Group())
+phase = dm.Phase(ode_class=BrachistochroneODE,
+                 transcription=dm.GaussLobatto(num_segments=4, order=[3, 5, 3, 5]))
 p.model.add_subsystem('phase0', phase)
 
 p.setup()

dymos/docs/feature_reference/transcriptions/figures/radau-pseudospectral/plot_radau-pseudospectral.py

@@ -6,13 +6,13 @@
 matplotlib.use('Agg')
 import matplotlib.pyplot as plt
 
-from openmdao.api import Problem, Group
-from dymos import Phase, Radau
+import openmdao.api as om
+import dymos as dm
 from dymos.examples.brachistochrone.brachistochrone_ode import BrachistochroneODE
 
-p = Problem(model=Group())
-phase = Phase(ode_class=BrachistochroneODE,
-              transcription=Radau(num_segments=4, order=[3, 5, 3, 5]))
+p = om.Problem(model=om.Group())
+phase = dm.Phase(ode_class=BrachistochroneODE,
+                 transcription=dm.Radau(num_segments=4, order=[3, 5, 3, 5]))
 p.model.add_subsystem('phase0', phase)
 
 p.setup()

dymos/examples/aircraft_steady_flight/aero/aero_coef_comp.py

@@ -1,13 +1,12 @@
 from __future__ import print_function, division, absolute_import
 
 import numpy as np
-
-from openmdao.api import MetaModelStructuredComp
+import openmdao.api as om
 
 from .crm_data import h_bp, alpha_bp, mach_bp, eta_bp, CL_data, CD_data, CM_data
 
 
-class AeroCoefComp(MetaModelStructuredComp):
+class AeroCoefComp(om.MetaModelStructuredComp):
     """ Interpolates aerodynamic coefficients for the NASA Common Research Model. """
 
     def setup(self):

dymos/examples/aircraft_steady_flight/aero/aero_forces_comp.py

@@ -2,10 +2,10 @@
 
 import numpy as np
 
-from openmdao.api import ExplicitComponent
+import openmdao.api as om
 
 
-class AeroForcesComp(ExplicitComponent):
+class AeroForcesComp(om.ExplicitComponent):
 
     def initialize(self):
         self.options.declare('num_nodes', types=int)

dymos/examples/aircraft_steady_flight/aero/aerodynamics_group.py

@@ -1,7 +1,6 @@
 from __future__ import print_function, division, absolute_import
 
-from openmdao.api import Group
-
+import openmdao.api as om
 from .aero_coef_comp import AeroCoefComp
 from .aero_forces_comp import AeroForcesComp
 from .mbi_aero_coef_comp import MBIAeroCoeffComp, setup_surrogates_all
@@ -12,7 +11,7 @@
     MBI = None
 
 
-class AerodynamicsGroup(Group):
+class AerodynamicsGroup(om.Group):
     """
     The purpose of the Aerodynamics is to compute the lift and
     drag forces on the aircraft.

dymos/examples/aircraft_steady_flight/aero/mbi_aero_coef_comp.py

@@ -5,8 +5,7 @@
 import os.path
 
 import numpy as np
-from openmdao.api import ExplicitComponent
-
+import openmdao.api as om
 try:
     import MBI
 except ImportError:
@@ -60,7 +59,7 @@ def setup_surrogates_all(model_name='CRM'):
     return [CL_arr, CD_arr, CM_arr, nums]
 
 
-class MBIAeroCoeffComp(ExplicitComponent):
+class MBIAeroCoeffComp(om.ExplicitComponent):
     """ Compute the lift, drag, and moment coefficients of the aircraft """
     def initialize(self):
         self.options.declare('vec_size', types=int)

dymos/examples/aircraft_steady_flight/aero/test/test_aero_coef_comp.py

@@ -5,7 +5,7 @@
 import numpy as np
 from numpy import array
 
-from openmdao.api import Problem, Group, IndepVarComp
+import openmdao.api as om
 from openmdao.utils.assert_utils import assert_rel_error
 
 from dymos.examples.aircraft_steady_flight.aero.aero_coef_comp import AeroCoefComp
@@ -17,11 +17,11 @@ def test_aero_coefs(self):
 
         nn = 100
 
-        prob = Problem(model=Group())
+        prob = om.Problem(model=om.Group())
 
         prob.model.add_subsystem(name='aero', subsys=AeroCoefComp(vec_size=nn, method='quintic'))
 
-        ivc = prob.model.add_subsystem(name='ivc', subsys=IndepVarComp(), promotes_outputs=['*'])
+        ivc = prob.model.add_subsystem(name='ivc', subsys=om.IndepVarComp(), promotes_outputs=['*'])
 
         ivc.add_output('mach', val=np.zeros(nn), units=None)
         ivc.add_output('alpha', val=np.zeros(nn), units='rad')

dymos/examples/aircraft_steady_flight/aero/test/test_mbi_aero_coef_comp.py

@@ -2,7 +2,7 @@
 
 import numpy as np
 
-from openmdao.api import Problem, Group, IndepVarComp
+import openmdao.api as om
 from openmdao.utils.assert_utils import assert_rel_error
 
 from dymos.examples.aircraft_steady_flight.aero.mbi_aero_coef_comp import setup_surrogates_all, \
@@ -22,7 +22,7 @@ def test_aero_coefs(self):
         NUM_NODES = 100
         MODEL = 'CRM'
 
-        prob = Problem(root=Group())
+        prob = om.Problem(model=om.Group())
 
         mbi_CL, mbi_CD, mbi_CM, mbi_num = setup_surrogates_all(MODEL)
 
@@ -32,16 +32,16 @@ def test_aero_coefs(self):
                                                          mbi_num=mbi_num))
 
         prob.model.add_subsystem(name='M_ivc',
-                                 subsys=IndepVarComp('M', val=np.zeros(NUM_NODES), units=None),
+                                 subsys=om.IndepVarComp('M', val=np.zeros(NUM_NODES), units=None),
                                  promotes=['M'])
         prob.model.add_subsystem(name='alpha_ivc',
-                                 subsys=IndepVarComp('alpha', val=np.zeros(NUM_NODES), units='rad'),
+                                 subsys=om.IndepVarComp('alpha', val=np.zeros(NUM_NODES), units='rad'),
                                  promotes=['alpha'])
         prob.model.add_subsystem(name='eta_ivc',
-                                 subsys=IndepVarComp('eta', val=np.zeros(NUM_NODES), units='rad'),
+                                 subsys=om.IndepVarComp('eta', val=np.zeros(NUM_NODES), units='rad'),
                                  promotes=['eta'])
         prob.model.add_subsystem(name='h_ivc',
-                                 subsys=IndepVarComp('h', val=np.zeros(NUM_NODES), units='km'),
+                                 subsys=om.IndepVarComp('h', val=np.zeros(NUM_NODES), units='km'),
                                  promotes=['h'])
 
         prob.model.connect('M', 'aero.M')