问题描述
我正在尝试使用 SLSQP 来优化翼型的攻角，以将停滞点放置在所需的位置。这纯粹是一个测试用例，用于检查我计算停滞位置部分的方法是否有效。
使用 COBYLA 运行时，优化在 47 次迭代后收敛到正确的 alpha (6.04144912)。当使用 SLSQP 运行时，它完成一次迭代，然后挂起 很长时间（10、20 分钟或更长时间，我没有准确计时），并以错误的值退出。输出为：
Driver debug print for iter coord: rank0:ScipyOptimize_SLSQP|0
--------------------------------------------------------------
Design Vars
{'alpha': array([0.5])}
Nonlinear constraints
None
Linear constraints
None
Objectives
{'obj_cmp.obj': array([0.00023868])}
Driver debug print for iter coord: rank0:ScipyOptimize_SLSQP|1
--------------------------------------------------------------
Design Vars
{'alpha': array([0.5])}
Nonlinear constraints
None
Linear constraints
None
Objectives
{'obj_cmp.obj': array([0.00023868])}
Optimization terminated successfully. (Exit mode 0)
Current function value: 0.0002386835700364719
Iterations: 1
Function evaluations: 1
Gradient evaluations: 1
Optimization Complete
-----------------------------------
Finished optimisation
为什么 SLSQP 会出现这样的行为不端？据我所知，当我查看 check_partials() 时，没有不正确的分析导数。
代码很长，所以我把它放在Pastebin这里：
- 核心：https://pastebin.com/fKJpnWHp
- 无粘性：https://pastebin.com/7Cmac5GF
- 翼型坐标 (NACA64-012)：https://pastebin.com/UZHXEsr6
解决方法
你问了两个问题，他们的答案最终彼此无关：
- 为什么模型在使用 SLSQP 时很慢，而在使用 COBYLA 时却很快
- 为什么 SLSQP 在一次迭代后停止？
1) 为什么 SLSQP 这么慢？
COBYLA 是一种无梯度方法。 SLSQP 使用梯度。因此，可靠的赌注是当 SLSQP 要求衍生品（COBYLA 从未这样做过）时，放缓发生了。
那是我先去看看的地方。计算导数分两步进行：a) 计算每个分量的偏分量和 b) 用这些偏分量求解线性系统以计算总数。减速必须在这两个步骤之一中。
由于您可以毫不费力地运行 check_partials，因此步骤 (a) 不太可能是罪魁祸首。所以这意味着步骤 (b) 可能是我们需要加快速度的地方。
我在您的模型上运行了汇总实用程序 (openmdao summary core.py) 并看到了这一点：
============== Problem Summary ============
Groups: 9
Components: 36
Max tree depth: 4
Design variables: 1 Total size: 1
Nonlinear Constraints: 0 Total size: 0
equality: 0 0
inequality: 0 0
Linear Constraints: 0 Total size: 0
equality: 0 0
inequality: 0 0
Objectives: 1 Total size: 1
Input variables: 87 Total size: 1661820
Output variables: 44 Total size: 1169614
Total connections: 87 Total transfer data size: 1661820
然后我生成了你的模型的 N2 并看到了这个：
所以我们有一个长度为 1169614 个元素的输出向量，这意味着您的线性系统是一个大约为 1e6x1e6 的矩阵。这非常大，您正在使用 DirectSolver 尝试计算/存储它的分解。这就是减速的根源。使用 DirectSolvers 非常适合较小的模型（经验法则是输出向量应小于 10000 个元素）。对于较大的求解器，您需要更加小心并使用更高级的线性求解器。
在您的情况下，我们可以从 N2 中看到模型中的任何地方都没有耦合（N2 的下三角形中没有任何耦合）。像这样的纯前馈模型可以使用更简单、更快的 LinearRunOnce 求解器（如果您不设置任何其他内容，这是默认设置）。所以我关闭了你模型中的所有 DirectSolvers，导数变得有效即时。让你的 N2 看起来像这样：
最佳线性求解器的选择非常依赖于模型。要考虑的一个因素是计算成本，另一个是数值鲁棒性。 Section 5.3 of the OpenMDAO paper 中详细介绍了此问题，我不会在此处介绍所有内容。但这里非常简要地总结了关键考虑因素。
刚开始使用 OpenMDAO 时，使用 DirectSolver 既是最简单的，也是通常最快的选择。这很简单，因为它不需要考虑您的模型结构，而且速度很快，因为对于小模型，OpenMDAO 可以将雅可比矩阵组装成密集或稀疏矩阵，并提供直接分解。但是，对于较大的模型（或具有非常大的输出向量的模型），计算因式分解的成本高得令人望而却步。在这种情况下，您需要更有意地分解求解器结构，并使用其他线性求解器（有时与直接求解器结合使用——see Section 5.3 of OpenMDAO paper 和 this OpenMDAO doc）。
您表示您想使用 DirectSolver 来利用稀疏的 Jacobian 存储。这是一种很好的直觉，但 OpenMDAO 的结构方式无论如何都不是问题。我们现在已经陷入困境，但既然你问了，我会给出一个简短的总结解释。从 OpenMDAO 3.7 开始，只有 DirectSolver 需要一个组装的雅可比行列式（实际上，线性求解器本身决定了它所连接的任何系统）。所有其他 LinearSolvers 使用 DictionaryJacobian（它存储每个 sub-jac 键控到 [of-var,wrt-var] 对）。每个 sub-jac 可以存储为密集或稀疏（取决于您如何声明该特定偏导数）。字典雅可比矩阵实际上是稀疏矩阵的一种形式，尽管不是传统矩阵。这里的关键要点是，如果您使用 LinearRunOnce（或任何其他求解器），那么无论如何您都会获得内存高效的数据存储。只有 DirectSolver 可以转换为实际矩阵对象的更传统的组合。
关于内存分配的问题。我从 openmdao 文档中借用了这张图片
2) 为什么 SLSQP 在一次迭代后就停止了？
基于梯度的优化对缩放非常敏感。我在你允许的设计空间内绘制了你的目标函数并得到了这个：
所以我们可以看到最小值在 6 度左右，但目标值是 TINY（大约 1e-4）。
作为一般经验法则，将您的目标设置为 1 级左右是一个好主意（我们有一个 scaling report 功能可以帮助实现这一点）。我添加了一个关于您的目标数量级的参考：
p.model.add_objective('obj',ref=1e-4)
然后我得到了一个很好的结果：
Optimization terminated successfully (Exit mode 0)
Current function value: [3.02197589e-11]
Iterations: 7
Function evaluations: 9
Gradient evaluations: 7
Optimization Complete
-----------------------------------
Finished optimization
alpha = [6.04143334]
time: 2.1188600063323975 seconds
不幸的是，基于梯度的优化很难缩放。首先将目标/约束扩展到 1 阶是一个不错的经验法则，但通常需要针对更复杂的问题进行调整。