【问题标题】:Add constraint to OpenMDAO at the driver/prob level在驱动程序/问题级别向 OpenMDAO 添加约束
【发布时间】:2016-09-15 20:53:59
【问题描述】:

是否可以为 OpenMDAO 问题添加约束?在下面的示例中,我想将目标函数限制在 -3.16 以下。我从另一个文件sellar_backend.py 导入了卖方问题。我是否可以在不修改 sellar_backend.py 的情况下添加此约束?

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' 失败。在问题级别添加这个新约束会非常方便。

【问题讨论】:

    标签: python optimization constraints openmdao


    【解决方案1】:

    你很亲密。 add 方法位于 top.root,而不是 top。您将组件添加到组中,在本例中是问题的根组。

    # 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”。

    【讨论】:

    • 谢谢!我什么时候不想使用promotes=['*']?总是使用它而不再考虑它是很诱人的。
    • promotes=['*'] 有点像在 Python 中使用 from foobar import *。在某些情况下它很好,但如果你的团队有很多组件,那么你产生意外后果的机会就会增加。在某些情况下,您确实希望将所有组件参数和未知数暴露给父组边界,请注意这样做时不会出现意外的名称冲突。
    • 补充@RobFalck 所说的内容,虽然promote=['*] 非常方便,但它可能会导致一些问题。首先,您可能会得到意想不到的连接,但如果您小心,这是可以避免的。但是,当您与其他人共享您的模型时,问题会更大。在这种情况下,他们无法通过查看您的文件轻松查看与什么相关的内容。使用我们的新模型查看器可以在一定程度上缓解后一个问题,但它仍然不是很好的做法:openmdao.readthedocs.io/en/latest/usr-guide/tutorials/…
    猜你喜欢
    • 2017-11-24
    • 2022-01-13
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2017-07-01
    • 2020-04-20
    • 1970-01-01
    • 2016-01-17
    相关资源
    最近更新 更多