2017-09-08 13 views
-1

私はOpenMDAOを使って以下に示す暗黙の方程式を解こうとしました。 「:double_scalarsで遭遇無効値RuntimeWarning」方程式は、以下この場合のOpenMDAOの配列構造として方程式系を解く方法は?

x[i]*z[i] + z[i] - 4 = 0, 
y[i] = x[i] + 2*z[i] 

The solution is (for i=2) z = [2.666667, 2], y = [5.833333, 5] for x = [0.5, 1.0]. 

に示されているが、私はメッセージを以下に示すコード、

from __future__ import print_function 
from openmdao.api import Component, Group, Problem, Newton, ScipyGMRES, IndepVarComp 

class SimpleEquationSystem(Component): 
    """Solve the Equation 
      x[i]*z[i] + z[i] - 4 = 0 
      y[i] = x[i] + 2*z[i] 

     Solution: z = [2.666667, 2], y = [5.833333, 5] for x = [0.5, 1.0] 
    """ 

    def __init__(self): 
     super(SimpleEquationSystem, self).__init__() 


     self.add_param('x', [0.5, 1.0]) 
     self.add_state('y', [0.0,0.0]) 
     self.add_state('z', [0.0,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 
     for i in range(2): 
      x=params['x'][i] 
      y = unknowns['y'][i] 
      z = unknowns['z'][i] 


      resids['y'][i] = x*z + z - 4 
      resids['z'][i] = x + 2*z - y 

     print('y_%d' % self.iter,'=%s' %resids['y'], 'z_%d' % self.iter, '=%s' %resids['z']) 
     print('x' ,'=%s' %x, 'y', '=%s' %y, 'z', '=%s' %z) 

top = Problem() 
root = top.root = Group() 
root.add('comp', SimpleEquationSystem()) 
root.add('p1', IndepVarComp('x', [0.5, 1.0])) 
root.connect('p1.x', 'comp.x') 
# 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.nl_solver.options['maxiter']=int(200) 
root.ln_solver = ScipyGMRES() 

top.setup() 
top.print_all_convergence(level=1, depth=2) 
top.run() 
print('Solution x=%s, y=%s, z=%s' % (top['comp.x'], top['comp.y'], top['comp.z'])) 

このコードエラーアウトを使用しています。上記のコードを使用している間、OpenMDAOでこのエラーが発生するはずです。

以下に示すように、ベクトルの代わりに残差方程式をスカラーとして実装すると、正常に動作します。

resids['z1']= params['x1'] + 2*unknowns['z1'] - unknowns['y1'] 
resids['z2']= params['x2'] + 2*unknowns['z2'] - unknowns['y2'] 
resids['y1']= params['x1']*unknowns['z1'] + unknowns['z1'] - 4 
resids['y2']= params['x2']*unknowns['z2'] + unknowns['z2'] - 4 

しかし、forループで処理する方が簡単なように、残差をベクトルとして使用したいと考えています。この問題を解決するために私を助けてもらえますか?

+0

あなたは、最小完全、かつ検証例HTTPSを生成する必要があります。 //stackoverflow.com/help/mcve。投稿が長すぎます。 – DyZ

答えて

1

変数を宣言するときは、それを宣言する変数型です。あなたはそれをOpenMDAOがリストとして解釈するリストを与えました。リストは分化をサポートするデータ型ではないので、ニュートンは何もすることができません。また

import numpy as np 

self.add_param('x', np.array([0.5, 1.0])) 
self.add_state('y', np.array([0.0,0.0])) 
self.add_state('z', np.array([0.0,0.0])) 

:あなたは、以下の変更を加えることで、それらnumpyの配列にする必要があり

root.add('p1', IndepVarComp('x', np.array([0.5, 1.0]))) 
+0

はい、動作しています。ありがとう、@ケネス・ムーア –

関連する問題