Skip to content

Commit

Permalink
Removed default_(non)linear_solver options.
Browse files Browse the repository at this point in the history
Added auto_solvers options to Trajectory and Phases, allowing users to disable the behavior if desired.
Added parallel_phases option to Trajectory, which toggles whether the top level phase container is a parallel group.
  • Loading branch information
robfalck committed Nov 28, 2023
1 parent cba8c29 commit fd9f601
Show file tree
Hide file tree
Showing 7 changed files with 113 additions and 149 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
def make_traj(transcription='gauss-lobatto', transcription_order=3, compressed=False,
connected=False, default_nonlinear_solver=None, default_linear_solver=None):
"""
Build a traejctory for the finite burn orbit raise problem.
Build a trajectory for the finite burn orbit raise problem.
Parameters
----------
Expand All @@ -31,8 +31,13 @@ def make_traj(transcription='gauss-lobatto', transcription_order=3, compressed=F
t = {'gauss-lobatto': dm.GaussLobatto(num_segments=5, order=transcription_order, compressed=compressed),
'radau': dm.Radau(num_segments=5, order=transcription_order, compressed=compressed)}

traj = dm.Trajectory(default_nonlinear_solver=default_nonlinear_solver,
default_linear_solver=default_linear_solver)
traj = dm.Trajectory()

if default_nonlinear_solver is not None:
traj.phases.nonlinear_solver = default_nonlinear_solver

if default_linear_solver is not None:
traj.phases.linear_solver = default_linear_solver

traj.add_parameter('c', opt=False, val=1.5, units='DU/TU',
targets={'burn1': ['c'], 'burn2': ['c']})
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -114,7 +114,6 @@ def test_ex_two_burn_orbit_raise_connected(self):
show_output=False, restart='dymos_solution.db',
connected=True, solution_record_file='dymos_solution2.db',
simulation_record_file='dymos_simulation2.db')
#

case2 = om.CaseReader('dymos_solution2.db').get_case('final')
sim_case2 = om.CaseReader('dymos_simulation2.db').get_case('final')
Expand All @@ -127,41 +126,30 @@ def test_ex_two_burn_orbit_raise_connected(self):
def test_restart_from_solution_radau_to_connected(self):
optimizer = 'IPOPT'

expected_warnings = \
[(om.OpenMDAOWarning,
"'traj' <class Trajectory>: Setting phases.nonlinear_solver to `om.NonlinearBlockJac(iprint=0)`.\n"
"Connected phases in parallel require a non-default nonlinear solver.\n"
"Use traj.options[\'default_nonlinear_solver\'] to explicitly set the solver."),
(om.OpenMDAOWarning,
"'traj' <class Trajectory>: Setting phases.linear_solver to `om.PETScKrylov()`.\n"
"Connected phases in parallel require a non-default linear solver.\n"
"Use traj.options[\'default_linear_solver\'] to explicitly set the solver.")]

with assert_warnings(expected_warnings):
p = two_burn_orbit_raise_problem(transcription='radau', transcription_order=3,
compressed=False, optimizer=optimizer,
show_output=False, connected=True)
p = two_burn_orbit_raise_problem(transcription='radau', transcription_order=3,
compressed=False, optimizer=optimizer,
show_output=False, connected=True)

if p.model.traj.phases.burn2 in p.model.traj.phases._subsystems_myproc:
assert_near_equal(p.get_val('traj.burn2.states:deltav')[0], 0.3995,
tolerance=4.0E-3)
if p.model.traj.phases.burn2 in p.model.traj.phases._subsystems_myproc:
assert_near_equal(p.get_val('traj.burn2.states:deltav')[0], 0.3995,
tolerance=4.0E-3)

case1 = om.CaseReader('dymos_solution.db').get_case('final')
sim_case1 = om.CaseReader('dymos_simulation.db').get_case('final')
case1 = om.CaseReader('dymos_solution.db').get_case('final')
sim_case1 = om.CaseReader('dymos_simulation.db').get_case('final')

# Run again without an actual optimizer
p = two_burn_orbit_raise_problem(transcription='radau', transcription_order=3,
compressed=False, optimizer=optimizer, run_driver=False,
show_output=False, restart='dymos_solution.db',
connected=True, solution_record_file='dymos_solution2.db',
simulation_record_file='dymos_simulation2.db')
# Run again without an actual optimizer
p = two_burn_orbit_raise_problem(transcription='radau', transcription_order=3,
compressed=False, optimizer=optimizer, run_driver=False,
show_output=False, restart='dymos_solution.db',
connected=True, solution_record_file='dymos_solution2.db',
simulation_record_file='dymos_simulation2.db')

case2 = om.CaseReader('dymos_solution2.db').get_case('final')
sim_case2 = om.CaseReader('dymos_simulation2.db').get_case('final')
case2 = om.CaseReader('dymos_solution2.db').get_case('final')
sim_case2 = om.CaseReader('dymos_simulation2.db').get_case('final')

# Verify that the second case has the same inputs and outputs
assert_cases_equal(case1, case2, tol=1.0E-8)
assert_cases_equal(sim_case1, sim_case2, tol=1.0E-8)
# Verify that the second case has the same inputs and outputs
assert_cases_equal(case1, case2, tol=1.0E-8)
assert_cases_equal(sim_case1, sim_case2, tol=1.0E-8)


if __name__ == '__main__': # pragma: no cover
Expand Down
11 changes: 2 additions & 9 deletions dymos/phase/phase.py
Original file line number Diff line number Diff line change
Expand Up @@ -114,15 +114,8 @@ def initialize(self):
desc='Keyword arguments provided when initializing the ODE System')
self.options.declare('transcription', types=TranscriptionBase,
desc='Transcription technique of the optimal control problem.')
self.options.declare('default_nonlinear_solver',
default=None, allow_none=True, recordable=False,
desc='A nonlinear solver to be used when the phase has implicit behavior due'
'to the use of `solve_segments` or `input_initial` in one or more states,'
'or the use of `set_duration_balance`.')
self.options.declare('default_linear_solver', default=None, allow_none=True, recordable=False,
desc='A linear solver to be used when the phase has implicit behavior due'
'to the use of `solve_segments` or `input_initial` in one or more states,'
'or the use of `set_duration_balance`.')
self.options.declare('auto_solvers', types=bool, default=True,
desc='If True, attempt to automatically assign solvers if necessary.')

def add_state(self, name, units=_unspecified, shape=_unspecified,
rate_source=_unspecified, targets=_unspecified,
Expand Down
5 changes: 1 addition & 4 deletions dymos/trajectory/test/test_trajectory.py
Original file line number Diff line number Diff line change
Expand Up @@ -1779,7 +1779,7 @@ def compute(self, inputs, outputs):
p.model.add_design_var('radius', lower=0.01, upper=0.10,
ref0=0.01, ref=0.10, units='m')

traj = p.model.add_subsystem('traj', dm.Trajectory())
traj = p.model.add_subsystem('traj', dm.Trajectory(parallel_phases=True))

transcription = dm.Radau(num_segments=5, order=3, compressed=True)

Expand Down Expand Up @@ -1856,9 +1856,6 @@ def compute(self, inputs, outputs):
p.model.connect('size_comp.mass', 'traj.parameters:m')
p.model.connect('size_comp.S', 'traj.parameters:S')

# A linear solver at the top level can improve performance.
p.model.linear_solver = om.DirectSolver()

# Finish Problem Setup
p.setup()

Expand Down
55 changes: 28 additions & 27 deletions dymos/trajectory/trajectory.py
Original file line number Diff line number Diff line change
Expand Up @@ -61,23 +61,20 @@ def __init__(self, **kwargs):
self._has_connected_phases = False

self.parameter_options = {}
self.phases = om.ParallelGroup()
self.phases = om.ParallelGroup() if self.options['parallel_phases'] else om.Group()

def initialize(self):
"""
Declare any options for Trajectory.
"""
self.options.declare('sim_mode', types=bool, default=False,
desc='Used internally by Dymos when invoking simulate on a trajectory')
self.options.declare('default_nonlinear_solver',
default=None, allow_none=True, recordable=False,
desc='A nonlinear solver to be used when Phases are connected but being '
'run in parallel. If not specified, Dymos will automatically use '
'NonlinearBlockJac in this situation.')
self.options.declare('default_linear_solver', default=None, allow_none=True, recordable=False,
desc='A linear solver to be used when Phases are connected but being '
'run in parallel. If not specified, Dymos will automatically use '
'PETScKrylov in this situation.')
self.options.declare('parallel_phases',
default=True, types=bool,
desc='If True, the top-level container of all phases will be a ParallelGroup,'
'otherwise it will be a standard OpenMDAO Group.')
self.options.declare('auto_solvers', types=bool, default=True,
desc='If True, attempt to automatically assign solvers if necessary.')

def add_phase(self, name, phase, **kwargs):
"""
Expand Down Expand Up @@ -1016,26 +1013,30 @@ def _configure_solvers(self):
These solvers can be changed through the
'default_nonlinear_solver' and 'default_linear_solver' options.
"""
if self._has_connected_phases and isinstance(self.phases, om.ParallelGroup):
if not self.options['auto_solvers']:
return

warn = False
if self._has_connected_phases and isinstance(self.phases, om.ParallelGroup) and MPI and self.comm.size > 1:
msg = (f'{self.msginfo}\n'
f' Non-default solvers are required\n'
f' connected phases in parallel: True\n')
if isinstance(self.phases.nonlinear_solver, om.NonlinearRunOnce):
if self.options['default_nonlinear_solver'] is None:
msg = f'{self.msginfo}: Setting phases.nonlinear_solver to `om.NonlinearBlockJac(iprint=0)`.\n' \
f'Connected phases in parallel require a non-default nonlinear solver.\n' \
f'Use {self.pathname}.options[\'default_nonlinear_solver\'] to explicitly set the solver.'
om.issue_warning(msg)
self.phases.nonlinear_solver = om.NonlinearBlockJac(iprint=0)
else:
self.phases.nonlinear_solver = self.options['default_nonlinear_solver']
warn = True
msg += (f' Setting \n'
f' {self.pathname}.phases.nonlinear_solver = om.NonlinearBlockJac(iprint=0)\n'
f' Explicitly set {self.pathname}.phases.nonlinear_solver to override.\n')
self.phases.nonlinear_solver = om.NonlinearBlockJac(iprint=0)

if isinstance(self.phases.linear_solver, om.LinearRunOnce):
if self.options['default_linear_solver'] is None:
msg = f'{self.msginfo}: Setting phases.linear_solver to `om.PETScKrylov()`.\n' \
f'Connected phases in parallel require a non-default linear solver.\n' \
f'Use {self.pathname}.options[\'default_linear_solver\'] to explicitly set the solver.'
om.issue_warning(msg)
self.phases.linear_solver = om.PETScKrylov()
else:
self.phases.linear_solver = self.options['default_linear_solver']
warn = True
msg += (f' Setting\n'
f' {self.pathname}.phases.linear_solver = om.PETScKrylov(iprint=0)\n'
f' Explicitly set {self.pathname}.phases.linear_solver to override.\n')
self.phases.linear_solver = om.PETScKrylov(iprint=0)

if warn:
om.issue_warning(msg)

def configure(self):
"""
Expand Down
50 changes: 16 additions & 34 deletions dymos/transcriptions/pseudospectral/pseudospectral_base.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,12 @@ class PseudospectralBase(TranscriptionBase):
**kwargs : dict
Dictionary of optional arguments.
"""
def __init__(self, **kwargs):
super().__init__(**kwargs)

self.any_solved_segs = False
self.any_connected_opt_segs = False

def initialize(self):
"""
Declare transcription options.
Expand Down Expand Up @@ -80,8 +86,6 @@ def setup_states(self, phase):
"""
grid_data = self.grid_data

self.any_solved_segs = False
self.any_connected_opt_segs = False
for options in phase.state_options.values():
# Transcription solve_segments overrides state solve_segments if its not set
if options['solve_segments'] is None:
Expand Down Expand Up @@ -494,45 +498,23 @@ def setup_solvers(self, phase):
"""
pass

def configure_solvers(self, phase):
def configure_solvers(self, phase, requires_solvers=None):
"""
Configure the solvers.
Parameters
----------
phase : dymos.Phase
The phase object to which this transcription instance applies.
"""
super().configure_solvers(phase)
if self.any_solved_segs or self.any_connected_opt_segs:
if isinstance(phase.nonlinear_solver, om.NonlinearRunOnce):
if phase.options['default_nonlinear_solver'] is None:
msg = f'{phase.msginfo}: Setting {phase.pathname}.nonlinear_solver to ' \
f'`om.NewtonSolver()`.\n' \
f'A phase requires a non-default nonlinear solver when a state utilizes solve_segments or ' \
f'input_initial, or when implicit duration is used.\n' \
f'Use `phase.options[\'default_nonlinear_solver\']` to explicitly set the solver.'
om.issue_warning(msg)
phase.nonlinear_solver = om.NewtonSolver(iprint=0)
phase.nonlinear_solver.options['solve_subsystems'] = True
phase.nonlinear_solver.options['maxiter'] = 1000
phase.nonlinear_solver.options['iprint'] = 2
phase.nonlinear_solver.options['stall_limit'] = 5
phase.nonlinear_solver.linesearch = om.ArmijoGoldsteinLS()
else:
phase.nonlinear_solver = phase.options['default_nonlinear_solver']

if isinstance(phase.linear_solver, om.LinearRunOnce):
if phase.options['default_linear_solver'] is None:
msg = f'{phase.msginfo}: Setting {phase.pathname}.linear_solver to ' \
f'`om.DirectSolver()`.\n' \
f'A phase requires a non-default linear solver when a state utilizes solve_segments or ' \
f'input_initial, or when implicit duration is used.\n' \
f'Use `phase.options[\'default_linear_solver\']` to explicitly set the solver.'
om.issue_warning(msg)
phase.linear_solver = om.DirectSolver()
else:
phase.linear_solver = self.options['default_linear_solver']
requires_solvers : dict[str: bool]
A dictionary mapping a string descriptor of a reason why a solver is required,
and whether a solver is required.
"""
req_solvers = {'solved segments': self.any_solved_segs,
'input initial': self.any_connected_opt_segs}
if requires_solvers is not None:
req_solvers.update(requires_solvers)
super().configure_solvers(phase, requires_solvers=req_solvers)

def setup_timeseries_outputs(self, phase):
"""
Expand Down
Loading

0 comments on commit fd9f601

Please sign in to comment.