如何将离散时间估计集成到 OpenMDAO 组件中？

Question

我正在尝试熟悉 OpenMDAO。我发现很难理解的一件事是 integrated/state 变量如何在 OpenMDAO 的单独组件中工作。我认为这很重要，我想了解基础知识。

假设我有一个推力恒定且质量不断变化的火箭，我想用离散时间和 0.1 秒的时间步长对前 10 秒的飞行进行建模。假设不断变化的质量也是时间的函数，但也假装它是一个依赖于许多事物的困难变量，所以我们希望在不同的组件中计算它并简单地作为该组件的输入。

-如何确保质量在速度的离散时间步长计算中更新，同时处于不同的组件中？

下面的代码是我试图解决这个示例问题的平庸尝试，但我现在相当确定 Mass 是一个简单的静态输入，而它应该随时间以某种任意关系发生变化。（假设推力不变）

from openmdao.api import ExplicitComponent

class VelocityComponent(ExplicitComponent):

    def setup(self):
        self.add_input('T', desc='Propulsion Thrust')
        self.add_input('M', desc='Instanteneous Mass')\

        self.add_output('v', desc='Satellite Velocity')

        self.declare_partials('*','*')

    def compute(self, inputs, outputs)

        v = 10                          #some initial velocity value
        t = 0                           #initial time value
        tstep = 0.1
        tend = 10

        for i in range(0,tend/tstep):
            a = inputs['T']/inputs['M'] #calculate acceleration  
            v += a                      #update velocity
            t += tstep                  #next time step 

        outputs['v'] = v

速度 v 应与随时间变化的加速度 a 相结合，而不是恒定加速度。

PS：由于我对所有这一切都比较陌生，但愿意学习，因此非常感谢任何可以帮助像我这样的初学者使用 OpenMDAO 的资源提示。

PSPS: 我已阅读 OpenMDAO 文档的初学者和高级用户指南，但找不到包含集成变量的示例。旧文档有一个 engine and transmission system 的示例，它的引擎组件确实包含状态变量和一些离散的集成步骤，但它使用旧的 OpenMDAO 版本，我不知道它在新版本中如何工作（如果我什至正确理解旧的）

Answer 1

这是关于 ODE 积分的讨论，您选择了一个相当复杂的主题，因为 1) 有很多不同的方法可以完成积分（例如显式 Euler、RK4、BDF ...） 2）通过时间积分携带解析导数是非常棘手的。

对于#2，困难的根源正是您所确定的问题。如果您还想将一组不同的组件组织成一个组来构造 ODE，则不能在单个组件内使用简单的 for 循环结构。

好消息是已经编写了一个库来为您处理所有时间积分：Dymos。截至 2019 年 4 月，该库由 OpenMDAO 团队自己积极开发，并且仍在进行某些模式 API 修订和功能添加。尽管 APIs 有点流畅，我还是建议您看一下那里的示例。对于复杂的问题，这是最好的选择。

但是，您无需额外的库即可实现简单的时间步长欧拉积分。诀窍是您为每个时间步消除 ODE 的一个实例，并将新状态从一个实例传递到下一个实例。这是一个在恒定重力下坠落物体的简单示例：

from openmdao.api import IndepVarComp, Problem, ExplicitComponent


class Cannonball(ExplicitComponent): 

    def initialize(self): 

        self.options.declare('delta_t', default=0.1)

    def setup(self): 

        self.add_input('Yi', units='m',   desc='position at the start of the time-step')
        self.add_input('Vi', units='m/s', desc='velocity at the start of the time-step')

        self.add_output('Ye', units='m',   desc='position at the end of the time-step')
        self.add_output('Ve', units='m/s', desc='velocity at the end of the time-step')

        self.declare_partials(of='*', wrt='*', method='cs')


    def compute(self, inputs, outputs): 

        dt = self.options['delta_t']

        outputs['Ve'] = 9.81 * dt + inputs['Vi']
        outputs['Ye'] = 0.5 * 9.81 * dt**2 + inputs['Vi'] * dt + inputs['Yi']


if __name__ == "__main__": 
    import numpy as np
    import matplotlib.pylab as plt

    N_TIMES = 10

    p = Problem()

    ivc = p.model.add_subsystem('init_conditions', IndepVarComp(), promotes=['*'])
    ivc.add_output('Y0', 100., units='m')
    ivc.add_output('V0', 0, units='m/s')

    p.model.connect('Y0', 't_0.Yi')
    p.model.connect('V0', 't_0.Vi')

    for i in range(N_TIMES): 
        p.model.add_subsystem(f't_{i}', Cannonball())


    for i in range(N_TIMES-1): 
        p.model.connect(f't_{i}.Ye', f't_{i+1}.Yi')
        p.model.connect(f't_{i}.Ve', f't_{i+1}.Vi')

    p.setup()

    p.run_model()


    # collect the data into an array for plotting
    Y = [p['Y0'],]
    V = [p['V0'],]
    for i in range(N_TIMES): 
        Y.append(p[f't_{i}.Ye'])
        V.append(p[f't_{i}.Ve'])

    times = np.arange(N_TIMES+1) * .01 # delta_t

    fig, ax = plt.subplots()
    ax.plot(times, Y)
    ax.set_ylabel('velocity (m/s')
    ax.set_xlabel('time (s)')

    plt.show()

这将为您的模型提供时间结构的前馈（使用 OpenMDAO's built-in N2 viewer).

你可以看到你得到了预期的时间二次方位置。

您可以通过将不同的 ODE 方案编码到此组件（例如 RK4）中来进行更复杂的集成。您还可以编写一个由多个组件组成的更复杂的组作为您的时间步长，然后剔除该组的多个副本。

我想强调的是，虽然上面的示例有助于加深理解，但这并不是在 OpenMDAO 中实际进行时间积分的非常有效的方法。 Dymos 在内部做事非常不同，使用及时向量化的组件以提高效率。 None 更少，如果你真的对 OpenMDAO 中最简单的时间积分方案感兴趣......就是这样。

如何将离散时间估计集成到 OpenMDAO 组件中？

How to integrate discrete time estimation into OpenMDAO Components?

openmdao