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 Direct linear solver.

from openmdao.api import Problem, NewtonSolver, ScipyKrylov, DirectSolver
from openmdao.test_suite.components.sellar import SellarDerivatives

prob = Problem()
model = prob.model = SellarDerivatives()

model.nonlinear_solver = newton = NewtonSolver()

# using a different linear solver for Newton with a looser tolerance
newton.linear_solver = ScipyKrylov(atol=1e-4)

# used for analytic derivatives
model.linear_solver = DirectSolver()

prob.setup()
prob.run_model()
NL: NLBGS Converged in 7 iterations
print(prob['y1'])
[ 25.58830237]
print(prob['y2'])
[ 12.05848815]

Some models have more complex coupling. There could be top-level cycles between groups as well as lower-level groups that have cycles of their own. The openmdao.test_suite.components.double_sellar.DoubleSellar (TODO: Link to problem page) 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.

from openmdao.api import Problem, NewtonSolver, ScipyKrylov, DirectSolver, NonlinearBlockGS, LinearBlockGS
from openmdao.test_suite.components.double_sellar import DoubleSellar

prob = Problem()
model = prob.model = DoubleSellar()

# each SubSellar group converges itself
g1 = model.g1
g1.nonlinear_solver = NewtonSolver()
g1.linear_solver = DirectSolver()  # used for derivatives

g2 = model.g2
g2.nonlinear_solver = NewtonSolver()
g2.linear_solver = DirectSolver()

# Converge the outer loop with Gauss Seidel, with a looser tolerance.
model.nonlinear_solver = NonlinearBlockGS(rtol=1.0e-5)
model.linear_solver = ScipyKrylov()
model.linear_solver.precon = LinearBlockGS()

prob.setup()
prob.run_model()
|  
|  ==
|  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: NLBGS Failed to Converge in 10 iterations
print(prob['g1.y1'])
[ 0.64000398]
print(prob['g1.y2'])
[ 0.80000249]
print(prob['g2.y1'])
[ 0.64000221]
print(prob['g2.y2'])
[ 0.80000138]

Note

Preconditioning for iterative linear solvers is a complex topic. The structure of the preconditioner should follow the model hierarchy itself, but developing an effective and efficient preconditioner is not trivial. If you’re having trouble converging the linear solves with an iterative solver, you should try using the Direct solver instead. Before changing solvers, first verify that all your partials derivatives are correct with the check_partials method.


You can also specify solvers as part of the initialization of a Group

class DoubleSellar(Group):

    def __init__(self, units=None, scaling=None, **kwargs):
        super(DoubleSellar, self).__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 = NewtonSolver()
        self.linear_solver = DirectSolver()

Tags

Solver