如何解决MDO 测试服：Golinski 的使用 OpenMDAO 的减速器问题
最近，我开始在多学科设计优化领域工作。我在 MDO 测试服中使用 OpenMDAO 框架对 Golinski 的减速器进行重量优化。我想为这个问题应用 MDF 架构。我指的是 Tedford 和 Martins 的论文“Benchmarking Multidisciplinary Design Optimization Algorithms”用于问题的表述和分解。他们将这个问题分解为三个学科及其各自的限制。
在编码时，我从 OpenMDAO 文档中提到了卖家问题。我做了三个学科和一个组（speed_mda()）来实现多学科分析。我在 speed_mda() 组中添加了目标函数和约束作为子系统。我已经将它们与学科输出（耦合变量）建立了联系。但是我并没有将这些限制应用于个别学科（实际上，我不知道该怎么做）。所以我已经将所有这些都应用到了顶级组中。我违反了一些限制，得到了输出 2713.678。
这是我的代码：
# Discipline 1
class speed_1(om.ExplicitComponent):
def setup(self):
self.add_input('z1',val = 0)
self.add_input('z2',val = 0)
self.add_output('y1',val = 1)
def setup_partials(self):
# Finite difference all partials.
self.declare_partials('*','*',method='fd')
def compute(self,inputs,outputs):
outputs['y1'] = max(27/(inputs['z1']**2*(inputs['z2'])),397.5/(inputs['z1']**2*inputs['z2']**2),5*inputs['z1'],2.6)
# Discipline 2
class speed_2(om.ExplicitComponent):
def setup(self):
self.add_input('z1',val = 0)
self.add_input('x21',val = 0)
self.add_output('y2',val = 0)
def setup_partials(self):
self.declare_partials('*',outputs):
outputs['y2'] = max((1.93*inputs['x21']**3/(inputs['z1']*inputs['z2']))**0.25,1/(0.5*(((1.69*10**7)**2)*inputs['x21']**2/(inputs['z1']**2*inputs['z2']**2) + 745)**0.5)**(0.3333),2.9)
# Discipline 3
class speed_3(om.ExplicitComponent):
def setup(self):
self.add_input('z1',val = 0)
self.add_input('x31',val = 0)
self.add_output('y3',outputs):
outputs['y3'] = max((1.93*inputs['x31']**3/(inputs['z1']*inputs['z2']))**0.25,1/(85*(((1.69*10**7)**2)*inputs['x31']**2/(inputs['z1']**2*inputs['z2']**2) + 1.575*(10**8))**0.5)**(0.3333),5)
class speed_mda(om.Group):
def setup(self):
# Adding all discipline to MDA
cycle = self.add_subsystem('cycle',om.Group(),promotes_inputs=['z1','z2','x21','x31'])
cycle.add_subsystem('d1',speed_1(),promotes_inputs = ['z1','z2'])
cycle.add_subsystem('d2',speed_2(),'x21'])
cycle.add_subsystem('d3',speed_3(),'x31'])
# No need of connections for the discipline
cycle.set_input_defaults('x21',7.8)
cycle.set_input_defaults('x31',8.3)
cycle.set_input_defaults('z1',0.75)
cycle.set_input_defaults('z2',22.0)
# Add solver to MDA: Nonlinear Block Gauss Seidel is a gradient free solver
cycle.nonlinear_solver = om.NonlinearBlockGS()
# Adding obj. function and constraints as a subsystem
self.add_subsystem('obj_fun',om.ExecComp('obj = (0.7854*y1*z1**2)*(3.3333*z2**2+14.933*z2-43.0934) - 1.5079*y1*(y2**2+y3**2)+7.477*(y2**3+y3**3)+0.7854*(x21*y2**2+x31*y3**2)',z1=0.0,z2=0.0,x21=0.0,x31=0.0),promotes=['x21','x31','z1','obj'])
self.add_subsystem('con1',om.ExecComp('c1 = z1*z2 - 40.0'),promotes=['c1']) # Global
self.add_subsystem('con10',om.ExecComp('c10 = y1 - 12.0*z1'),promotes=['c10']) # 1
self.add_subsystem('con11',om.ExecComp('c11 = y1 - 3.6'),promotes=['c11'])
self.add_subsystem('con12',om.ExecComp('c12 = y2 - 3.9'),promotes=['c12']) # 2
self.add_subsystem('con13',om.ExecComp('c13 = 2.85*y2 - x21'),promotes=['c13'])
self.add_subsystem('con14',om.ExecComp('c14 = y3 - 5.5'),promotes=['c14']) # 3
self.add_subsystem('con15',om.ExecComp('c15 = 2.09*y3 - x31'),promotes=['c15'])
# Connect outputs from MDA (coupled variables) to obj. function and constraints
self.connect('cycle.d1.y1',['obj_fun.y1','con10.y1','con11.y1'])
self.connect('cycle.d2.y2',['obj_fun.y2','con12.y2','con13.y2'])
self.connect('cycle.d3.y3',['obj_fun.y3','con14.y3','con15.y3'])
# Form topmost group (Problem) and add above MDF model to it
prob = om.Problem()
model = prob.model = speed_mda()
prob.driver = om.ScipyOptimizeDriver()
prob.driver.options['optimizer'] = 'SLSQP'
prob.driver.options['tol'] = 1e-9
prob.driver.options['disp'] = True
model.add_design_var('x21',lower=7.3,upper=8.3)
model.add_design_var('x31',upper=8.3)
model.add_design_var('z1',lower=0.7,upper=0.8)
model.add_design_var('z2',lower=17.0,upper=28.0)
model.add_objective('obj')
model.add_constraint('c1',upper = 0)
model.add_constraint('c10',upper = 0)
model.add_constraint('c11',upper = 0)
model.add_constraint('c12',upper = 0)
model.add_constraint('c13',upper = 0)
model.add_constraint('c14',upper = 0)
model.add_constraint('c15',upper = 0)
prob.model.approx_totals()
prob.setup()
prob.set_solver_print(level=0)
prob.set_val('x21',7.8)
prob.set_val('x31',8.3)
prob.set_val('z1',0.75)
prob.set_val('z2',22.0)
prob.run_model()
prob.run_driver()
print('minimum found at')
print((prob.get_val('z1')[0],prob.get_val('z2')[0],prob.get_val('x21')[0],prob.get_val('x31')[0],prob.get_val('cycle.d1.y1')[0],prob.get_val('cycle.d2.y2')[0],prob.get_val('cycle.d3.y3')[0]))
print('minumum objective')
print(prob.get_val('obj')[0])
我得到以下输出：
Positive directional derivative for linesearch (Exit mode 8)
Current function value: [2713.67806668]
Iterations: 8
Function evaluations: 4
Gradient evaluations: 4
Optimization FAILED.
Positive directional derivative for linesearch
-----------------------------------
minimum found at
(0.7,17.00007920093152,7.300007915120889,7.300015825246984,3.5,2.9,5.0)
minumum objective
2713.6780666813797
这里，目标值远小于实际值，表示优化失败。我寻找了上述错误，但我最初的猜测是在界限内。我的所有问题约束都不满足。我也尝试对个别学科施加限制，但我做不到。我不知道实际问题是什么，可能我犯了一些基本的概念错误。可以请任何人帮助我。
解决方法
那篇论文中提出减速器问题的方式有些奇怪。 max 函数的使用在技术上根本无法区分，这使得在基于梯度的优化中实现问题的方法不太理想
此外，该公式看起来与我在其他论文中看到的对减速器问题的其他描述非常不同。这个公式可以追溯到 MDO 测试问题套件的原始工作，其中不是真正多学科的问题被分解为单独的块，并添加了额外的约束以确保兼容性。就本文而言，我认为问题表述的变化导致了一些不太理想的问题结构。我建议您寻找更合适的公式，例如 this one
无论如何，当我从他们论文的最佳状态设置给定的输入时，我没有得到他们报告的最佳值。我得到了一个较低的数字，所以你的代码肯定有一些微妙的不同。不知何故，您的代码返回了较低的值，因此请仔细检查您的方程式。
我在 OpenMDAO V3.8 上运行您的代码并得到以下结果：
/Users/jsgray/work/packages/OpenMDAO/openmdao/core/total_jac.py:1713: UserWarning:Constraints or objectives ['con1.c1','con12.c12','con13.c13','con14.c14','con15.c15'] cannot be impacted by the design variables of the problem.
Positive directional derivative for linesearch (Exit mode 8)
Current function value: [2713.66402046]
Iterations: 9
Function evaluations: 5
Gradient evaluations: 5
Optimization FAILED.
Positive directional derivative for linesearch
-----------------------------------
minimum found at
(0.7,17.000000001268262,7.3,3.5,2.9,5.0)
minumum objective
2713.6640204584155
所以我看到了与您相同的值，但我还收到了一个有用的附加警告（在 V3.8 中新增）关于不受任何设计变量影响的约束。当我注释掉这些约束时，结果变成了
Optimization terminated successfully (Exit mode 0)
Current function value: [2713.66402024]
Iterations: 10
Function evaluations: 6
Gradient evaluations: 6
Optimization Complete
-----------------------------------
minimum found at
(0.7,17.0,5.0)
minumum objective
2713.6640202393
这是与以前相同的答案，但没有来自优化器的可怕警告。因此，您看到的错误是由于存在大量约束，虽然本质上满足这些约束，但实际上优化器无法控制这些约束。这会导致全零行出现在雅可比矩阵中，从而使优化问题变得奇异。虽然 SLSQP 能够解决奇点问题，但它引起了足够多的数值问题，导致它发出警告。
版权声明：本文内容由互联网用户自发贡献，该文观点与技术仅代表作者本人。本站仅提供信息存储空间服务，不拥有所有权，不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容， 请发送邮件至 dio@foxmail.com 举报，一经查实，本站将立刻删除。
