固定次数 RHS 评估的最优 ODE 方法

计算科学 数值分析
2021-11-27 23:57:25

在实践中,数值求解 IVP 的运行时间

x˙(t)=f(t,x(t)) for t[t0,t1]
x(t0)=x0
通常受评估右侧 (RHS) 的持续时间支配f. 因此,让我们假设所有其他操作都是即时的(即没有计算成本)。如果求解 IVP 的总运行时间是有限的,那么这相当于限制了f对一些NN.

我们只对最终值感兴趣x(t1).

我正在寻找可以帮助我在这种环境中选择最佳 ODE 方法的理论和实践结果。

如果,例如,N=2那么我们可以使用两个显式的欧拉宽度步长来求解 IVP(t1t0)/2或一步宽度t1t0使用中点法。我还不清楚哪个更可取。对于较大的N,当然也可以考虑多步法、迭代龙格-库塔方案等。

我正在寻找的是类似于现有结果的结果,例如,对于正交规则:我们可以选择n权重{wi}及相关点{xi}使得正交规则i=1nwig(xi)对所有多项式都是精确的g这样deg(g)2n1.

因此,鉴于允许的 RHS 评估数量有限,我正在寻找 ODE 方法的全局准确性的上限或下限f. 如果边界仅适用于某些类别的 RHS 或对解决方案构成额外限制,则可以x(就像求积法则的结果一样,它只适用于一定程度的多项式)。

编辑:一些背景信息:这是用于硬实时应用程序,即结果x(t1)必须在已知的截止日期之前可用。因此限制了 RHS 评估的数量N作为主要的成本因素。通常我们的问题是僵硬的并且相对较小。

EDIT2:不幸的是,我没有精确的时间要求,但可以安全地假设N会相当小(肯定<100,可能更接近10)。考虑到实时性要求,我们必须在模型的准确性之间进行权衡(更好的模型会导致 RHS 的执行时间更长,从而导致更低的N) 和 ODE 方法的准确性(更好的方法需要更高的N)。

4个回答

我认为回答您问题的关键参考是Hosea 和 Shampine 的这篇论文现在我将介绍一些背景。

通常,在对 IVP 进行数值积分时可以使用的步长可能会受到稳定性或准确性的限制。如果要在稳定性方面选择最佳求解器,则需要考虑绝对稳定性区域对于一步法,这是

S={zC:|P(z)|1}.

这里P(z)是方法的稳定性函数;参见例如 Hairer et 的文本。人。稳定的必要条件是λhS在哪里λ范围超过 的 jacobian 的特征值fh是步长。这并不总是非线性问题的充分条件,但它通常是一个很好的经验法则并在实践中使用。

有关寻找(显式)方法允许大稳定步长的问题的广泛处理,请参阅我的这篇关于稳定性多项式的论文这篇关于可压缩流体模拟的龙格库塔方法优化的论文

如果您发现最大的稳定步长已经为您提供了足够的精度,那么稳定性是相关的问题。另一方面,步长可能会受到您的精度要求的限制。通常所做的是本地错误控制。解决方案是使用两种方法计算的,它们的差异用作对不太准确的方法中的误差的估计。自适应地选择步长以尽可能接近规定的公差。

两个理论度量是预测准确率效率的关键。第一个是该方法的准确度顺序,它描述了当步长减小时误差接近零的速率。第二个是准确率效率指数(参见上面第一句话中链接的 Hosea 和 Shampine 的论文),它考虑了误差项中出现的常数,并允许在相同顺序的方法之间进行比较。

可以使用NodePy以简单且自动化的方式计算各种方法的准确性和稳定性效率(免责声明:NodePy 由我开发)。

在这个方向上没有很多结果,因为它比仅仅修复精度更困难,因为稳定性考虑通常需要您选择小于所需精度所需的时间步长。因此,结果分为刚性和非刚性案例。前一种情况下,时间步长和 RHS 评估要求通常不受精度控制,而在后一种情况下,它们是。

我将专注于显式方法,因为隐式情况远不那么明显,您需要使用多少 RHS 评估..这完全取决于您决定如何解决结果系统。

对于非刚性系统:

显式 Runge-Kutta 方法存在阶段限制,即需要多少阶段(RHS 评估)才能达到一定的准确度。在四阶之后,阶段的数量超过了精度的数量级,并且差距继续扩大。Butcher 的 ODE 大书: http: //books.google.com/books/about/Numerical_Methods_for_Ordinary_Different.html ?id=opd2NkBmMxsC

很好地解释了其中一些“不存在”的证据。

您的正交规则示例导致了多步类型的方法,例如 Adams-bashforth,或者导致了现在称为频谱延迟校正的方法。对于 adams-bashforth,您每一步只需要一个 RHS 评估,但由于这些方法的稳定性区域通常非常小,您通常最终在 RHS 评估方面做与具有相同的 Runge-Kutta 方法相同的工作量命令。

这是一篇关于光谱延迟校正的论文:

https://www.google.com/search?q=spectral+deferred+correction&aq=f&oq=spectral+deferred+correction&aqs=chrome.0.57j0l2j62.3336j0&sourceid=chrome&ie=UTF-8

我不确定这些集成方法如何与标准显式方法进行对比,它们通常需要更多内存来保存正交节点处的解决方案状态,因此我自己从未使用过它们。

对于刚性系统:

有“优化”的时间步进器,但不幸的是,关于这些时间步进器的精确理论结果仅限于一些简单的情况(甚至那些被证明不是微不足道的工作)。三个标准结果表明,对于具有S阶段:它可以包含在其稳定区域中的最负实轴是一个长度区间2S2,它所能包含的最虚轴是一个长度区间S1,并且它可以包含的与假想轴相切的最大圆具有半径S(所有这些也是相互排斥的)。

Reid 和 David 已经回答了技术问题,但无论如何我还是想提供一些背景信息。今天,对于您能想到的大多数 ODE,您可以达到基本上完美的精度(例如,1014) 使用大多数优秀的 ODE 求解器包,如果评估的话,基本上没有时间f(x)免费。

当然也有例外(非常大的系统,非常僵硬的系统),但社区中的一种普遍看法是,为“标准”系统设计 ODE 求解器的问题已得到解决。因此,我认为您提出的问题不是一个非常有趣的问题——它解决了 ODE 求解器设计中不再重要的一个组件。这也可以解释缺乏有关该主题的文献。

我的推理如下。由于您的问题很僵硬且不大,因此很可能您应该使用一些隐式方法。如果您使用隐式方法,那么您在每一步中至少求解一个线性系统。这意味着如果雅可比不是结构化的,那么每一步都将采取O(dim3)或者O(dim2)除了评估 RHS之外,每一步都会失败(取决于你 LU 分解雅可比行列式的频率) 。

所以第一点是确定你的 RHS 是否真的比基础线性代数更昂贵。

第二点:从文献中得知,基于“昂贵”方法(即显式 RK 方法)的求解器有时比“更便宜”的求解器(显式多步方法)执行得更快。

总而言之,我认为您不应该只考虑 RHS 评估计数。