2016-04-19 14 views
0

MDO Test suiteの問題を解決するために、OpenMDAOを検証手段として使用しています。私が直面している問題は、特定の問題に対するオプティマイザ/アルゴリズムの選択です。複数の問題を抱えるオプティマイザの選択

以下は、IDFメソッドによって策定されたスイートのPropane Combustion Problemのコードです。私はpyOptSparseDriverで利用可能なSLSQP & NSGA2を使って解決しようとしました。しかし、両方のオプティマイザが最適以下のソリューションを提供しています。

from __future__ import print_function, division 
    import numpy as np 
    import time 
    from openmdao.api import * 

    class PropaneDis1(Component): 
     def __init__(self): 
      super(PropaneDis1, self).__init__() 

      self.add_param('z1', 1.0) 
      self.add_param('y3', np.ones(3)) 

      self.add_output('y1', val=np.ones(2)) 

     def solve_nonlinear(self, params, unknowns, resids): 

      z1 = params['z1'] 
      y31 = params['y3'][0] 

      y12 = 3.0 - z1 
      y11 = z1*y31/y12 

      unknowns['y1'] = np.array([y11, y12]) 

    class PropaneDis2(Component): 
     def __init__(self): 
      super(PropaneDis2, self).__init__() 

      self.add_param('z1', val= 0.0) 
      self.add_param('y1', val= np.ones(2)) 
      self.add_param('y3', val= np.ones(3)) 

      self.add_output('y2', val=np.ones(2)) 

     def solve_nonlinear(self, params, unknowns, resids): 

      z1 = params['z1'] 
      y12 = params['y1'][1] 
      y33 = params['y3'][2] 

      y21 = 0.1*z1*y33/(40.0*y12) 
      y22 = 0.1*(z1**2)*y33/(40.0*(y12**2)) 

      unknowns['y2'] = np.array([y21, y22]) 

    class PropaneDis3(Component): 
     def __init__(self): 
      super(PropaneDis3, self).__init__() 

      self.add_param('z1', val=0.0) 
      self.add_param('y1', val=np.ones(2)) 
      self.add_param('y2', val=np.ones(2)) 
      self.add_param('x3', val=np.ones(3)) 

      self.add_output('y3', val=np.ones(3)) 

     def solve_nonlinear(self, params, unknowns, resids): 

      z1 = params['z1'] 
      y11 = params['y1'][0] 
      x31 = params['x3'][0] 
      y12 = params['y1'][1] 
      x32 = params['x3'][1] 
      x33 = params['x3'][2] 
      y21 = params['y2'][0] 
      y22 = params['y2'][1] 

      y31 = (8.0 - 2*y11 - x32 - x33)/2.0 
      y32 = 4*10.0 - 2*x31 
      y33 = z1 + y11 + x31 + y12 + x32 + x33 + y21 + y22 + y31 + y32 

      unknowns['y3'] = np.array([y31, y32, y33]) 


    class PropaneIDF(Group): 

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

      self.add('pz1', IndepVarComp('z1', val=1.0), promotes=['*']) 
      self.add('px3', IndepVarComp('x3', val=np.array([20.0,1.5,1.0])), promotes=['*']) 
      self.add('py1t', IndepVarComp('y1t', val=np.array([1.0,1.0])), promotes=['*']) 
      self.add('py2t', IndepVarComp('y2t', val=np.array([1.0,1.0])), promotes=['*']) 
      self.add('py3t', IndepVarComp('y3t', val=np.array([1.0,1.0,1.0])), promotes=['*']) 

      self.add('d1', PropaneDis1(), promotes=['z1', 'y1']) 
      self.add('d2', PropaneDis2(), promotes=['z1', 'y2']) 
      self.add('d3', PropaneDis3(), promotes=['z1', 'x3', 'y3']) 

      self.connect('y1t', 'd2.y1') 
      self.connect('y1t', 'd3.y1') 

      self.connect('y2t', 'd3.y2') 

      self.connect('y3t', 'd1.y3') 
      self.connect('y3t', 'd2.y3') 

      self.add('con_cmp1', ExecComp('f2=2*z1 + y1[0] + y1[1] + x3[2] + y2[0] + y3[1] + 2*y2[1] - 10.0' 
        , z1=0.0, y1=np.zeros(2), y2=np.zeros(2), y3=np.zeros(3), x3=np.zeros(3)) 
        , promotes=['*']) 
      self.add('con_cmp2', ExecComp('f6=abs(y1[0]*y1[1])**0.5 - x3[1]*abs(40.0*z1/y3[2])**0.5' 
        , z1=0.0, y1=np.zeros(2), y3=np.ones(3), x3=np.zeros(3)) 
        , promotes=['*']) 
      self.add('con_cmp3', ExecComp('f7=abs(z1*y1[0])**0.5 - x3[2]*abs(40.0*y1[1]/y3[2])**0.5' 
        , z1=0.0, y1=np.zeros(2), y3=np.ones(3), x3=np.zeros(3)) 
        , promotes=['*']) 
      self.add('con_cmp4', ExecComp('f9=z1*(x3[0]**0.5) - y1[1]*y3[1]*(40.0/y3[2])**0.5' 
        , z1=0.0, y1=np.zeros(2), y3=np.zeros(3), x3=np.zeros(3)) 
        , promotes=['*']) 

      self.add('con_cmp51', ExecComp('con51=y1t[0] - y1[0]', y1t=np.zeros(2), y1=np.zeros(2)), promotes=['*']) 
      self.add('con_cmp52', ExecComp('con52=y1t[1] - y1[1]', y1t=np.zeros(2), y1=np.zeros(2)), promotes=['*']) 

      self.add('con_cmp61', ExecComp('con61=y2t[0] - y2[0]', y2t=np.zeros(2), y2=np.zeros(2)), promotes=['*']) 
      self.add('con_cmp62', ExecComp('con62=y2t[1] - y2[1]', y2t=np.zeros(2), y2=np.zeros(2)), promotes=['*']) 

      self.add('con_cmp71', ExecComp('con71=y3t[0] - y3[0]', y3t=np.zeros(3), y3=np.zeros(3)), promotes=['*']) 
      self.add('con_cmp72', ExecComp('con72=y3t[1] - y3[1]', y3t=np.zeros(3), y3=np.zeros(3)), promotes=['*']) 
      self.add('con_cmp73', ExecComp('con73=y3t[2] - y3[2]', y3t=np.zeros(3), y3=np.zeros(3)), promotes=['*']) 

      self.add('obj_cmp', ExecComp('obj=2*z1 + y1[0] + y1[1] + x3[2] + y2[0] + y3[1] + 2*y2[1] - 10.0 +\ 
             abs(y1[0]*y1[1])**0.5 - x3[1]*abs(40.0*z1/y3[2])**0.5+\ 
             abs(z1*y1[0])**0.5 - x3[2]*abs(40.0*y1[1]/y3[2])**0.5+\ 
             z1*(x3[0]**0.5) - y1[1]*y3[1]*(40.0/y3[2])**0.5 ' 
        , z1=0.0, y1=np.zeros(2), y2=np.zeros(2), y3=np.zeros(3), x3=np.zeros(3)) 
        , promotes=['*']) 

      self.fd_options['force_fd']=True 
      self.fd_options['step_size']=1.0e-12 

    top = Problem() 
    top.root = PropaneIDF() 

    top.driver = pyOptSparseDriver() 
    top.driver.options['optimizer']='SLSQP' 

    top.driver.add_desvar('z1', lower=0.0, upper=20.0) 
    top.driver.add_desvar('x3', lower=np.array([0.0, 0.0, 0.0]), upper=np.array([20.0, 20.0, 20.0])) 
    top.driver.add_desvar('y1t', lower=np.array([0.0, 0.0]), upper=np.array([20.0, 20.0])) 
    top.driver.add_desvar('y2t', lower=np.array([0.0, 0.0]), upper=np.array([20.0, 20.0])) 
    top.driver.add_desvar('y3t', lower=np.array([0.0, 0.0, 0.0]), upper=np.array([20.0, 20.0, 20.0])) 

    top.driver.add_objective('obj') 

    top.driver.add_constraint('con51', equals=0.0) 
    top.driver.add_constraint('con52', equals=0.0) 
    top.driver.add_constraint('con61', equals=0.0) 
    top.driver.add_constraint('con62', equals=0.0) 
    top.driver.add_constraint('con71', equals=0.0) 
    top.driver.add_constraint('con72', equals=0.0) 
    top.driver.add_constraint('con73', equals=0.0) 

    top.driver.add_constraint('f2', lower=0.0) 
    top.driver.add_constraint('f6', lower=0.0) 
    top.driver.add_constraint('f7', lower=0.0) 
    top.driver.add_constraint('f9', lower=0.0) 
    # top.driver.add_constraint('f2', equals=0.0) 
    # top.driver.add_constraint('f6', equals=0.0) 
    # top.driver.add_constraint('f7', equals=0.0) 
    # top.driver.add_constraint('f9', equals=0.0) 

    top.setup() 
    tt = time.time() 

    top.run() 

SLSQPのためのソリューション - NSGA2ため

Optimization using pyOpt_sparse 
    ================================================================================ 

      Objective Function: _objfunc 

     Solution: 
    -------------------------------------------------------------------------------- 
     Total Time:     0.1535 
      User Objective Time :  0.0183 
      User Sensitivity Time :  0.0806 
      Interface Time :   0.0534 
      Opt Solver Time:   0.0013 
     Calls to Objective Function :  10 
     Calls to Sens Function :   7 

     Objectives: 
      Name  Value  Optimum 
      obj  14.3997    0 

     Variables (c - continuous, i - integer, d - discrete): 
       Name  Type  Value  Lower Bound Upper Bound 
      z1_0  c  2.130827  0.00e+00  2.00e+01 
      x3_0  c  20.000000  0.00e+00  2.00e+01 
      x3_1  c  0.000000  0.00e+00  2.00e+01 
      x3_2  c  -0.000000  0.00e+00  2.00e+01 
      y1t_0  c  3.645641  0.00e+00  2.00e+01 
      y1t_1  c  0.869179  0.00e+00  2.00e+01 
      y2t_0  c  0.180723  0.00e+00  2.00e+01 
      y2t_1  c  0.679352  0.00e+00  2.00e+01 
      y3t_0  c  1.691004  0.00e+00  2.00e+01 
      y3t_1  c  -0.000000  0.00e+00  2.00e+01 
      y3t_2  c  20.000000  0.00e+00  2.00e+01 

     Constraints (i - inequality, e - equality): 
      Name Type     Bounds 
      con51  i  0.00e+00 <= -0.499951 <= 0.00e+00 
      con52  i  0.00e+00 <= 0.000006 <= 0.00e+00 
      con61  i  0.00e+00 <= 0.058146 <= 0.00e+00 
      con62  i  0.00e+00 <= 0.378850 <= 0.00e+00 
      con71  i  0.00e+00 <= 1.336646 <= 0.00e+00 
      con72  i  0.00e+00 <= -0.000000 <= 0.00e+00 
      con73  i  0.00e+00 <= -7.860081 <= 0.00e+00 
      f2  i  0.00e+00 <= -0.000000 <= 1.00e+20 
      f6  i  0.00e+00 <= 1.898220 <= 1.00e+20 
      f7  i  0.00e+00 <= 2.972127 <= 1.00e+20 
      f9  i  0.00e+00 <= 9.529347 <= 1.00e+20 

    -------------------------------------------------------------------------------- 

ソリューション -

Optimization using pyOpt_sparse 
    ================================================================================ 

      Objective Function: _objfunc 

     Solution: 
    -------------------------------------------------------------------------------- 
     Total Time:     107.5302 
      User Objective Time :  62.5869 
      User Sensitivity Time :  0.0000 
      Interface Time :   41.6987 
      Opt Solver Time:   3.2446 
     Calls to Objective Function : 99981 
     Calls to Sens Function :   0 

     Objectives: 
      Name  Value  Optimum 
      obj  7.67208    0 

     Variables (c - continuous, i - integer, d - discrete): 
       Name  Type  Value  Lower Bound Upper Bound 
      z1_0  c  2.184051  0.00e+00  2.00e+01 
      x3_0  c  20.000000  0.00e+00  2.00e+01 
      x3_1  c  0.000000  0.00e+00  2.00e+01 
      x3_2  c  0.000000  0.00e+00  2.00e+01 
      y1t_0  c  3.757993  0.00e+00  2.00e+01 
      y1t_1  c  0.815949  0.00e+00  2.00e+01 
      y2t_0  c  0.193249  0.00e+00  2.00e+01 
      y2t_1  c  0.746901  0.00e+00  2.00e+01 
      y3t_0  c  1.592243  0.00e+00  2.00e+01 
      y3t_1  c  -0.000000  0.00e+00  2.00e+01 
      y3t_2  c  20.000000  0.00e+00  2.00e+01 

     Constraints (i - inequality, e - equality): 
      Name Type     Bounds 
      con51  i  0.00e+00 <= 0.173682 <= 0.00e+00 
      con52  i  0.00e+00 <= -0.520903 <= 0.00e+00 
      con61  i  0.00e+00 <= -0.122582 <= 0.00e+00 
      con62  i  0.00e+00 <= -0.292454 <= 0.00e+00 
      con71  i  0.00e+00 <= -0.120713 <= 0.00e+00 
      con72  i  0.00e+00 <= 0.000372 <= 0.00e+00 
      con73  i  0.00e+00 <= -7.729249 <= 0.00e+00 
      f2  i  0.00e+00 <= 0.143730 <= 1.00e+20 
      f6  i  0.00e+00 <= 1.641686 <= 1.00e+20 
      f7  i  0.00e+00 <= 1.957046 <= 1.00e+20 
      f9  i  0.00e+00 <= 3.929613 <= 1.00e+20 

    ---------------------------------------------------------------------- 

私は制約f2の問題を実行しようとしたような問題が定義していないことを知っています、 f6、f7 & f9下限の代わりに等価制約(= 0)として、私は最適を得ています。

他のオプティマイザを調べる必要がありますか?あるいは、オプティマイザの問題について私が混乱している他の問題もあります。

ありがとうございます。

+0

実装の詳細はどこにありますか?あなたが参照したページの方程式と一致するようには見えません。また、既知の解についても言及していません。 –

答えて

0

まず、SLSQPの出力を見ると、オプティマイザは実行可能なソリューションをまったく提供できませんでした。具体的には、con51〜con73の制約には多数の違反があります。これは、規律コンポーネント間の循環依存性を置き換える互換性の式です。したがって、オプティマイザは、多分野の結合を収束できなかったようです。

MDF設定に変換し、非線形ソルバーをニュートンに切り替えて、MDAサイクルを1つに収束させる方法を見てみましょう(私もオプティマイザを取り出しました) 。他の設計変数が初期条件に泊まりましたと仮定すると、私はそれは素敵な1E-10許容範囲に収束するようになった:

y1 [ 1. 2.] 
y2 [ 0.05886036 0.02943018] 
y3 [ 2.   38.   47.08829054] 

お知らせY3のための収束値はy3tを使用するときに指定したことを上限外の方法であることIDFの下で設計変数として:私は彼らがプロパンの問題で指定され表示されませんでしたよう

top.driver.add_desvar('y3t', lower=np.array([0.0, 0.0, 0.0]), 
          upper=np.array([10.0, 10.0, 10.0])) 

私は、あなたが上限を持っていることを確認していません。さらに重要なことは、あなたのコードに見られるものとMDOテストスイートのリンクで指定されているものを一致させることができなかったことです。私はいくつかのエラーがあるかもしれないと思うが、何かを確認することができませんでした。

関連する問題