I tried to solve the simple implicit equation using OpenMDAO. The equations are shown below,
x*z + z - 4 = 0
y = x + 2*z
The solution is z = 2.666667, y = 5.833333 for x = 0.5.
For this case, I have used the code which is shown below,
from __future__ import print_function
from openmdao.api import Component, Group, Problem, Newton, ScipyGMRES
class SimpleEquationSystem(Component):
"""Solve the Equation
x*z + z - 4 = 0
y = x + 2*z
Solution: z = 2.666667, y = 5.833333 for x = 0.5
"""
def __init__(self):
super(SimpleEquationSystem, self).__init__()
self.add_param('x', 0.5)
self.add_state('y', 0.0)
self.add_state('z', 0.0)
self.iter=0
def solve_nonlinear(self, params, unknowns, resids):
"""This component does no calculation on its own. It mainly holds the
initial value of the state. An OpenMDAO solver outside of this
component varies it to drive the residual to zero."""
pass
def apply_nonlinear(self, params, unknowns, resids):
""" Report the residual """
self.iter+=1
x=params['x']
y = unknowns['y']
z = unknowns['z']
resids['y'] = x*z + z - 4
resids['z'] = x + 2*z - y
print('y_%d' % self.iter,'=%f' %resids['y'], 'z_%d' % self.iter, '=%f' %resids['z'])
print('x' ,'=%f' %x, 'y', '=%f' %y, 'z', '=%f' %z)
top = Problem()
root = top.root = Group()
root.add('comp', SimpleEquationSystem())
# Tell these components to finite difference
root.comp.deriv_options['type'] = 'fd'
root.comp.deriv_options['form'] = 'central'
root.comp.deriv_options['step_size'] = 1.0e-4
root.nl_solver = Newton()
root.ln_solver = ScipyGMRES()
top.setup()
top.print_all_convergence(level=1, depth=2)
top.run()
print('Solution x=%.2f, y=%.2f, z=%.2f' % (top['comp.x'], top['comp.y'], top['comp.z']))
I wrote a code based on Solving an Implicit Relationship with a Newton Solver method. To ran this code, I got a solution like this,
##############################################
Setup: Checking root problem for potential issues...
No recorders have been specified, so no data will be saved.
The following parameters have no associated unknowns:
comp.x
The following components have no connections:
comp
Setup: Check of root problem complete.
##############################################
y_1 =-4.000000 z_1 =0.500000
x =0.500000 y =0.000000 z =0.000000
y_2 =-4.000000 z_2 =0.499900
x =0.500000 y =0.000100 z =0.000000
y_3 =-4.000000 z_3 =0.500100
x =0.500000 y =-0.000100 z =0.000000
y_4 =-3.999850 z_4 =0.500200
x =0.500000 y =0.000000 z =0.000100
y_5 =-4.000150 z_5 =0.499800
x =0.500000 y =0.000000 z =-0.000100
[root] LN: GMRES 1 | 0 0
[root] LN: GMRES 1 | Converged in 1 iterations
y_6 =-inf z_6 =-inf
x =0.500000 y =inf z =-inf
y_7 =-inf z_7 =-inf
x =0.500000 y =inf z =-inf
y_8 =-inf z_8 =-inf
x =0.500000 y =inf z =-inf
y_9 =-inf z_9 =-inf
x =0.500000 y =inf z =-inf
y_10 =-inf z_10 =-inf
x =0.500000 y =inf z =-inf
[root] LN: GMRES 1000 | nan nan
[root] LN: GMRES 1000 | Converged in 1000 iterations
y_11 =nan z_11 =nan
x =0.500000 y =nan z =nan
[root] NL: NEWTON 2 | nan nan (nan)
[root] NL: NEWTON 2 | FAILED to converge after 2 iterations
Solution x=0.50, y=nan, z=nan
C:\ProgramData\Anaconda3\Anaconda3\lib\site-packages\openmdao\core\system.py:750: RuntimeWarning: invalid value encountered in subtract
resultvec.vec[:] -= cache2
Could you know, How to solve this problem with respective iteration?. And tell me, How to construct the solve_nonlinear & apply_nonlinear for this case.
Looks like the problem clears up when you add an Indepvarcomp for the independent variable:
root.add('p1', IndepVarComp('x', 0.5))
root.add('comp', SimpleEquationSystem())
root.connect('p1.x', 'comp.x')
Then it converges quickly.
[root] NL: NEWTON 1 | 2.73692411e-11 6.7894731e-12 (6.41396046825)
[root] NL: NEWTON 1 | Converged in 1 iterations
I agree that you shouldn't have to add that for this to work, and it is probably a bug. I have verified that this problem has been fixed in our latest development version of OpenMDAO.