2016-09-15 59 views
0

是否可以向OpenMDAO問題添加約束?在下面的例子中,我想限制目標函數低於-3.16。我從另一個文件sellar_backend.py導入了sellar問題。我可以在不修改sellar_backend.py的情況下添加這個約束嗎?在驅動程序/概率級別向OpenMDAO添加約束

sellar_backend.py

import numpy as np 
from openmdao.api import Problem, ScipyOptimizer, Group, ExecComp, IndepVarComp, Component 
from openmdao.api import Newton, ScipyGMRES 
class SellarDis1(Component): 
    """Component containing Discipline 1.""" 

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

     # Global Design Variable 
     self.add_param('z', val=np.zeros(2)) 

     # Local Design Variable 
     self.add_param('x', val=0.) 

     # Coupling parameter 
     self.add_param('y2', val=1.0) 

     # Coupling output 
     self.add_output('y1', val=1.0) 

    def solve_nonlinear(self, params, unknowns, resids): 
     """Evaluates the equation 
     y1 = z1**2 + z2 + x1 - 0.2*y2""" 

     z1 = params['z'][0] 
     z2 = params['z'][1] 
     x1 = params['x'] 
     y2 = params['y2'] 

     unknowns['y1'] = z1**2 + z2 + x1 - 0.2*y2 

    def linearize(self, params, unknowns, resids): 
     """ Jacobian for Sellar discipline 1.""" 
     J = {} 

     J['y1','y2'] = -0.2 
     J['y1','z'] = np.array([[2*params['z'][0], 1.0]]) 
     J['y1','x'] = 1.0 

     return J 


class SellarDis2(Component): 
    """Component containing Discipline 2.""" 

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

     # Global Design Variable 
     self.add_param('z', val=np.zeros(2)) 

     # Coupling parameter 
     self.add_param('y1', val=1.0) 

     # Coupling output 
     self.add_output('y2', val=1.0) 

    def solve_nonlinear(self, params, unknowns, resids): 
     """Evaluates the equation 
     y2 = y1**(.5) + z1 + z2""" 

     z1 = params['z'][0] 
     z2 = params['z'][1] 
     y1 = params['y1'] 

     # Note: this may cause some issues. However, y1 is constrained to be 
     # above 3.16, so lets just let it converge, and the optimizer will 
     # throw it out 
     y1 = abs(y1) 

     unknowns['y2'] = y1**.5 + z1 + z2 

    def linearize(self, params, unknowns, resids): 
     """ Jacobian for Sellar discipline 2.""" 
     J = {} 

     J['y2', 'y1'] = .5*params['y1']**-.5 

     #Extra set of brackets below ensure we have a 2D array instead of a 1D array 
     # for the Jacobian; Note that Jacobian is 2D (num outputs x num inputs). 
     J['y2', 'z'] = np.array([[1.0, 1.0]]) 

     return J 

class StateConnection(Component): 
    """ Define connection with an explicit equation""" 

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

     # Inputs 
     self.add_param('y2_actual', 1.0) 

     # States 
     self.add_state('y2_command', val=1.0) 

    def apply_nonlinear(self, params, unknowns, resids): 
     """ Don't solve; just calculate the residual.""" 

     y2_actual = params['y2_actual'] 
     y2_command = unknowns['y2_command'] 

     resids['y2_command'] = y2_actual - y2_command 

    def solve_nonlinear(self, params, unknowns, resids): 
     """ This is a dummy comp that doesn't modify its state.""" 
     pass 

    def linearize(self, params, unknowns, resids): 
     """Analytical derivatives.""" 

     J = {} 

     # State equation 
     J[('y2_command', 'y2_command')] = -1.0 
     J[('y2_command', 'y2_actual')] = 1.0 

     return J 

class SellarStateConnection(Group): 
    """ Group containing the Sellar MDA. This version uses the disciplines 
    with derivatives.""" 

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

     self.add('px', IndepVarComp('x', 1.0), promotes=['x']) 
     self.add('pz', IndepVarComp('z', np.array([5.0, 2.0])), promotes=['z']) 

     self.add('state_eq', StateConnection()) 
     self.add('d1', SellarDis1(), promotes=['x', 'z', 'y1']) 
     self.add('d2', SellarDis2(), promotes=['z', 'y1']) 
     self.connect('state_eq.y2_command', 'd1.y2') 
     self.connect('d2.y2', 'state_eq.y2_actual') 

     self.add('obj_cmp', ExecComp('obj = x**2 + z[1] + y1 + exp(-y2)', 
            z=np.array([0.0, 0.0]), x=0.0, y1=0.0, y2=0.0), 
        promotes=['x', 'z', 'y1', 'obj']) 
     self.connect('d2.y2', 'obj_cmp.y2') 

     self.add('con_cmp1', ExecComp('con1 = 3.16 - y1'), promotes=['con1', 'y1']) 
     self.add('con_cmp2', ExecComp('con2 = y2 - 24.0'), promotes=['con2']) 
     self.connect('d2.y2', 'con_cmp2.y2') 

     self.nl_solver = Newton() 
     self.ln_solver = ScipyGMRES() 

example.py

from sellar_backend import * 
top = Problem() 
top.root = SellarStateConnection() 

top.driver = ScipyOptimizer() 
top.driver.options['optimizer'] = 'SLSQP' 
top.driver.options['tol'] = 1.0e-8 

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

# This is my best attempt so far at adding a constraint at this level  
top.add('new_constraint', ExecComp('new_con = -3.16 - obj'), promotes=['*']) 
top.driver.add_constraint('new_constraint', upper=0.0) 

top.driver.add_objective('obj') 
top.driver.add_constraint('con1', upper=0.0) 
top.driver.add_constraint('con2', upper=0.0) 

top.setup() 
top.run() 

print("\n") 
print("Minimum found at (%f, %f, %f)" % (top['z'][0], \ 
             top['z'][1], \ 
             top['x'])) 
print("Coupling vars: %f, %f" % (top['y1'], top['d2.y2'])) 
print("Minimum objective: ", top['obj']) 

這種失敗AttributeError: 'Problem' object has no attribute 'add'。在問題層面添加這個新約束將會非常方便。

回答

1

您近距離了。 add方法在top.root上,而不是頂部。您將組件添加到組中,在這種情況下是問題的根組。

# This is my best attempt so far at adding a constraint at this level 
top.root.add('new_constraint', ExecComp('new_con = -3.16 - obj'), promotes=['*']) 
top.driver.add_constraint('new_con', upper=0.0) 

而且,由於你在第一線推廣「*」,要約束的數量被稱爲「new_con」。

+0

謝謝!我什麼時候不想使用'promote = ['*']'?總是使用它並且永遠不會再考慮它是很誘人的。 – kilojoules

+0

'promote = ['*']'有點像在Python中使用'from foobar import *'。在某些情況下,這很好,但是如果你的團隊有很多組件,那麼你意想不到的後果的機會就會增加。在某些情況下,您確實想要將所有組件參數和未知數展示給父組邊界,只要注意不要發生意外的名稱衝突。 –

+0

添加到@RobFalck所說的內容,而'promote = ['*]'非常方便,可能會導致一些問題。首先,你可以得到意想不到的連接,但如果你小心這是可以避免的。更有問題的是,當你與別人分享你的模型。在這種情況下,他們無法通過查看文件輕鬆查看與什麼相關的內容。後一個問題可以通過使用我們的新模型查看器稍微緩解,但它仍然不是很好的做法:http://openmdao.readthedocs.io/en/latest/usr-guide/tutorials/visualizing-model-connections.html –