复杂步骤在 OpenMDAO 组件中未按预期工作
Complex-step not working as expected in OpenMDAO component
我正在尝试使用 complex-step 来获取 OpenMDAO 中一个简单组件的导数。虽然我有解析导数,但我想用复步比较性能。这是一个更大的设计问题的一部分。
这是一个最小的例子:
import numpy as np
from openmdao.api import IndepVarComp, Component, Problem, Group
class SpatialBeamDisp(Component):
def __init__(self, ny):
super(SpatialBeamDisp, self).__init__()
self.ny = ny
self.add_param('disp_aug', val=np.zeros(((self.ny+1)*6), dtype='complex'))
self.add_output('disp', val=np.zeros((self.ny, 6), dtype='complex'))
# Comment out this line to use analytic derivatives
self.deriv_options['type'] = 'cs'
def solve_nonlinear(self, params, unknowns, resids):
# Obtain the relevant portions of disp_aug and store the reshaped
# displacements in disp
unknowns['disp'] = params['disp_aug'][:-6].reshape((self.ny, 6))
def linearize(self, params, unknowns, resids):
jac = self.alloc_jacobian()
n = self.ny * 6
jac['disp', 'disp_aug'] = np.hstack((np.eye((n)), np.zeros((n, 6))))
return jac
top = Problem()
root = top.root = Group()
n = 5
disp_aug = np.random.random(((n+1) * 6))
root.add('disp_input', IndepVarComp('disp_aug', disp_aug), promotes=['*'])
root.add('disp_', SpatialBeamDisp(n), promotes=['*'])
top.setup()
top.run_once()
top.check_partial_derivatives(compact_print=True)
运行 此代码按原样生成不正确的雅可比行列式,同时注释掉 self.deriv_options['type' = 'cs'
行并使用解析表达式生成正确的雅可比行列式。
我在 Ubuntu 上使用 OpenMDAO 1.7.3 和 numpy 1.10.2。
我是否错误地设置了复阶导数?如果是这样,我应该如何编写此组件才能支持复杂步骤?
该问题是 memory/pointer 相关的问题。本质上,当您执行切片和重塑时,您最终会看到原始数组的视图,然后将其分配给 unknowns['disp']
。这个赋值打破了 OpenMDAO 指向它用于进行复杂步骤计算的原始数据的指针。您可以通过更改分配以请求将数据副本放入当前内存位置来明确修复它:
def solve_nonlinear(self, params, unknowns, resids):
# Obtain the relevant portions of disp_aug and store the reshaped
# displacements in disp
unknowns['disp'][:] = params['disp_aug'][:-6].reshape((self.ny, 6))
OpenMDAO 的 CS 代码中一定存在潜在的错误,这使得此分配保护成为必要。该问题不会出现在 fd 上,只会出现在 cs 上。
我正在尝试使用 complex-step 来获取 OpenMDAO 中一个简单组件的导数。虽然我有解析导数,但我想用复步比较性能。这是一个更大的设计问题的一部分。
这是一个最小的例子:
import numpy as np
from openmdao.api import IndepVarComp, Component, Problem, Group
class SpatialBeamDisp(Component):
def __init__(self, ny):
super(SpatialBeamDisp, self).__init__()
self.ny = ny
self.add_param('disp_aug', val=np.zeros(((self.ny+1)*6), dtype='complex'))
self.add_output('disp', val=np.zeros((self.ny, 6), dtype='complex'))
# Comment out this line to use analytic derivatives
self.deriv_options['type'] = 'cs'
def solve_nonlinear(self, params, unknowns, resids):
# Obtain the relevant portions of disp_aug and store the reshaped
# displacements in disp
unknowns['disp'] = params['disp_aug'][:-6].reshape((self.ny, 6))
def linearize(self, params, unknowns, resids):
jac = self.alloc_jacobian()
n = self.ny * 6
jac['disp', 'disp_aug'] = np.hstack((np.eye((n)), np.zeros((n, 6))))
return jac
top = Problem()
root = top.root = Group()
n = 5
disp_aug = np.random.random(((n+1) * 6))
root.add('disp_input', IndepVarComp('disp_aug', disp_aug), promotes=['*'])
root.add('disp_', SpatialBeamDisp(n), promotes=['*'])
top.setup()
top.run_once()
top.check_partial_derivatives(compact_print=True)
运行 此代码按原样生成不正确的雅可比行列式,同时注释掉 self.deriv_options['type' = 'cs'
行并使用解析表达式生成正确的雅可比行列式。
我在 Ubuntu 上使用 OpenMDAO 1.7.3 和 numpy 1.10.2。
我是否错误地设置了复阶导数?如果是这样,我应该如何编写此组件才能支持复杂步骤?
该问题是 memory/pointer 相关的问题。本质上,当您执行切片和重塑时,您最终会看到原始数组的视图,然后将其分配给 unknowns['disp']
。这个赋值打破了 OpenMDAO 指向它用于进行复杂步骤计算的原始数据的指针。您可以通过更改分配以请求将数据副本放入当前内存位置来明确修复它:
def solve_nonlinear(self, params, unknowns, resids):
# Obtain the relevant portions of disp_aug and store the reshaped
# displacements in disp
unknowns['disp'][:] = params['disp_aug'][:-6].reshape((self.ny, 6))
OpenMDAO 的 CS 代码中一定存在潜在的错误,这使得此分配保护成为必要。该问题不会出现在 fd 上,只会出现在 cs 上。