Setting Nonlinear and Linear Solvers#
A nonlinear solver, like NonlinearBlockGS or Newton, is used to converge the nonlinear analysis. A nonlinear solver is needed whenever there is a cyclic dependency between components in your model. It might also be needed if you have an ImplicitComponent in your model that expects the framework to handle its convergence.
Whenever you use a nonlinear solver on a Group or Component, if you’re going to be working with analytic derivatives, you will also need a linear solver. A linear solver, like LinearBlockGS or DirectSolver, is used to solve the linear system that provides total derivatives across the model.
You can add nonlinear and linear solvers at any level of the model hierarchy, letting you build a hierarchical solver setup to efficiently converge your model and solve for total derivatives across it.
Solvers for the Sellar Problem#
The Sellar Problem has two components with a cyclic dependency, so the appropriate nonlinear solver is necessary. We’ll use the Newton nonlinear solver, which requires derivatives, so we’ll also use the DirectSolver linear solver.
SellarDerivatives
class definition
class SellarDerivatives(om.Group):
"""
Group containing the Sellar MDA. This version uses the disciplines with derivatives.
"""
def setup(self):
self.add_subsystem('d1', SellarDis1withDerivatives(), promotes=['x', 'z', 'y1', 'y2'])
self.add_subsystem('d2', SellarDis2withDerivatives(), promotes=['z', 'y1', 'y2'])
self.add_subsystem('obj_cmp', om.ExecComp('obj = x**2 + z[1] + y1 + exp(y2)', obj=0.0,
x=0.0, z=np.array([0.0, 0.0]), y1=0.0, y2=0.0),
promotes=['obj', 'x', 'z', 'y1', 'y2'])
self.add_subsystem('con_cmp1', om.ExecComp('con1 = 3.16  y1', con1=0.0, y1=0.0),
promotes=['con1', 'y1'])
self.add_subsystem('con_cmp2', om.ExecComp('con2 = y2  24.0', con2=0.0, y2=0.0),
promotes=['con2', 'y2'])
self.set_input_defaults('x', 1.0)
self.set_input_defaults('z', np.array([5.0, 2.0]))
import openmdao.api as om
from openmdao.test_suite.components.sellar_feature import SellarDerivatives
prob = om.Problem()
model = prob.model = SellarDerivatives()
model.nonlinear_solver = newton = om.NewtonSolver(solve_subsystems=False)
# using a different linear solver for Newton with a looser tolerance
newton.linear_solver = om.ScipyKrylov(atol=1e4)
# used for analytic derivatives
model.linear_solver = om.DirectSolver()
prob.setup()
prob.run_model()
print(prob.get_val('y1'))
print(prob.get_val('y2'))
NL: Newton Converged in 3 iterations
[25.58830237]
[12.05848815]
Some models have more complex coupling. There could be toplevel cycles between groups as well as lowerlevel groups that have cycles of their own. openmdao.test_suite.components.double_sellar.DoubleSellar is a simple example of this kind of model structure. In these problems, you might want to specify a more complex hierarchical solver structure for both nonlinear and linear solvers.
DoubleSellar
class definition
class DoubleSellar(om.Group):
def __init__(self, units=None, scaling=None, **kwargs):
super().__init__(**kwargs)
self.add_subsystem('g1', SubSellar(units=units, scaling=scaling))
self.add_subsystem('g2', SubSellar(units=units, scaling=scaling))
self.connect('g1.y2', 'g2.x')
self.connect('g2.y2', 'g1.x')
# Converge the outer loop with Gauss Seidel, with a looser tolerance.
self.nonlinear_solver = om.NewtonSolver(solve_subsystems=False)
self.linear_solver = om.DirectSolver()
from openmdao.test_suite.components.double_sellar import DoubleSellar
prob = om.Problem()
model = prob.model = DoubleSellar()
# each SubSellar group converges itself
g1 = model.g1
g1.nonlinear_solver = om.NewtonSolver(solve_subsystems=False)
g1.linear_solver = om.DirectSolver() # used for derivatives
g2 = model.g2
g2.nonlinear_solver = om.NewtonSolver(solve_subsystems=False)
g2.linear_solver = om.DirectSolver()
# Converge the outer loop with Gauss Seidel, with a looser tolerance.
model.nonlinear_solver = om.NonlinearBlockGS(rtol=1.0e5)
model.linear_solver = om.ScipyKrylov()
model.linear_solver.precon = om.LinearBlockGS()
prob.setup()
prob.run_model()
print(prob.get_val('g1.y1'))
print(prob.get_val('g1.y2'))
print(prob.get_val('g2.y1'))
print(prob.get_val('g2.y2'))

 ==
 g1
 ==
 NL: Newton Converged in 3 iterations

 ==
 g2
 ==
 NL: Newton Converged in 3 iterations

 ==
 g1
 ==
 NL: Newton Converged in 3 iterations

 ==
 g2
 ==
 NL: Newton Converged in 3 iterations

 ==
 g1
 ==
 NL: Newton Converged in 3 iterations

 ==
 g2
 ==
 NL: Newton Converged in 2 iterations

 ==
 g1
 ==
 NL: Newton Converged in 2 iterations

 ==
 g2
 ==
 NL: Newton Converged in 2 iterations

 ==
 g1
 ==
 NL: Newton Converged in 2 iterations

 ==
 g2
 ==
 NL: Newton Converged in 2 iterations

 ==
 g1
 ==
 NL: Newton Converged in 2 iterations

 ==
 g2
 ==
 NL: Newton Converged in 2 iterations

 ==
 g1
 ==
 NL: Newton Converged in 2 iterations

 ==
 g2
 ==
 NL: Newton Converged in 2 iterations

 ==
 g1
 ==
 NL: Newton Converged in 2 iterations

 ==
 g2
 ==
 NL: Newton Converged in 2 iterations

 ==
 g1
 ==
 NL: Newton Converged in 2 iterations

 ==
 g2
 ==
 NL: Newton Converged in 1 iterations

 ==
 g1
 ==
 NL: Newton Converged in 1 iterations

 ==
 g2
 ==
 NL: Newton Converged in 1 iterations
NL: NLBGSSolver 'NL: NLBGS' on system '' failed to converge in 10 iterations.
[0.64000398]
[0.80000249]
[0.64000221]
[0.80000138]