这些只是一些快速而肮脏的建议 . 希望有人会发现它们有用!
1. Do you have a symbolic function or a set of points?
如果您有符号功能,您可以分析计算衍生物 . (很可能,如果那么容易的话,你会做到这一点,你不会在这里寻找替代品 . )
如果您具有符号函数且无法通过分析计算导数,则可以始终在一组点上评估函数,并使用此页面上列出的其他方法来评估导数 .
在大多数情况下,你有一组点(xi,fi),并且必须使用以下方法之一....
2. Is your grid evenly or unevenly spaced?
如果网格间隔均匀,您可能需要使用有限差分格式(请参阅维基百科文章here或here),除非您使用的是周期性边界条件(见下文) . Here是在解决网格上的常微分方程的背景下对有限差分方法的一个不错的介绍(特别参见幻灯片9-14) . 这些方法通常在计算上有效,易于实现,并且该方法的误差可以简单地估计为用于导出它的泰勒展开的截断误差 .
如果您的网格间距不均匀,您仍然可以使用有限差分格式,但表达式更加困难,并且精度随着网格的统一程度而变化很大 . 如果您的网格非常不均匀,您可能需要使用较大的模板尺寸(更多相邻点)来计算给定的导数点 . 人们经常构造插值多项式(通常是Lagrange polynomial)并区分该多项式来计算导数 . 例如,参见this StackExchange问题 . 使用这些方法通常很难估计错误(尽管有些人试图这样做:here和here) . Fornberg's method在这些情况下通常非常有用....
必须注意您域名的边界,因为模板通常涉及域外的点 . 有些人引入"ghost points"或将边界条件与不同阶数的导数组合以消除这些并简化模板 . 另一种方法是使用右侧或左侧有限差分方法 .
Here's一种优秀的有限差分方法,包括低阶的中心,右侧和左侧方案 . 我在工作站附近打印出这个,因为我发现它非常有用 .
3. Is your domain periodic? Can you assume periodic boundary conditions?
如果您的域是周期性的,您可以使用傅里叶谱方法计算导数到非常高的阶数 . 这种技术在某种程度上牺牲了性能以获得高精度 . 实际上,如果您使用的是N点,那么您对导数的估计大约是第N阶的准确度 . 有关更多信息,请参阅(例如)this WikiBook .
傅立叶方法通常使用快速傅立叶变换(FFT)算法来实现大致O(N log(N))性能,而不是天然实现的离散傅立叶变换(DFT)可能采用的O(N ^ 2)算法 .
如果您的函数和域不是周期性的,则不应使用傅里叶谱方法 . 如果您尝试将其与非周期性函数一起使用,则会出现大的错误和不良的"ringing"现象 .
计算任何阶数的导数需要1)从网格空间到光谱空间的变换(O(N log(N))),2)傅里叶系数乘以其光谱波数(O(N)),和2)从光谱空间到网格空间的逆变换(再次为O(N log(N))) .
在将傅里叶系数乘以其光谱波数时必须小心 . FFT算法的每个实现似乎都有自己的频谱模式和归一化参数的排序 . 例如,请参阅Math StackExchange上this question的答案,以获取有关在MATLAB中执行此操作的说明 .
4. What level of accuracy are you looking for? Do you need to compute the derivatives within a given tolerance?
出于许多目的,一阶或二阶有限差分方案可能就足够了 . 为了获得更高的精度,您可以使用更高阶的泰勒展开式,从而降低高阶项 .
如果需要计算给定容差范围内的导数,可能需要查看具有所需错误的高阶方案 .
通常,减少误差的最佳方法是在有限差分方案中减少网格间距,但这并非总是可行的 .
请注意,高阶有限差分格式几乎总是需要更大的模板尺寸(更多的相邻点) . 这可能会导致边界问题 . (参见上面关于鬼点的讨论 . )
5. Does it matter to you that your derivative is evaluated on the same points as your function is defined?
MATLAB提供 diff 函数来计算相邻数组元素之间的差异 . 这可以用于通过一阶前向差分(或前向有限差分)方案来计算近似导数,但估计是低阶估计 . 如MATLAB的 diff (link)文档中所述,如果输入长度为N的数组,它将返回一个长度为N-1的数组 . 当您在N点上使用此方法估计导数时,您只能估算出N-1点的导数 . (请注意,如果它们按升序排序,则可以在不均匀网格上使用 . )
在大多数情况下,我们希望在所有点评估衍生物,这意味着我们想要使用除 diff 方法之外的其他东西 .
6. Do you need to calculate multiple orders of derivatives?
可以 Build 一个方程组,其中网格点函数值和这些点上的一阶和二阶导数都相互依赖 . 这可以通过像往常一样在相邻点处组合泰勒展开来找到,但是保持导数项而不是将它们抵消,并将它们与相邻点的那些连接在一起 . 这些方程式可以通过线性代数求解,不仅可以得到一阶导数,也可以得到第二个导数(或者更高阶,如果设置得当) . 我相信这些被称为组合有限差分方案,它们通常与紧致有限差分方案结合使用,这将在下面讨论 .
紧凑有限差分格式(link) . 在这些方案中, Build 一个设计矩阵并通过矩阵求解同时计算所有点的导数 . 它们被称为"compact"因为它们通常被设计为需要比具有可比精度的普通有限差分方案更少的模板点 . 因为它们涉及将所有点连接在一起的矩阵方程,确定据说紧凑有限差分方案具有"spectral-like resolution"(例如Lele's 1992 paper - 优秀!),这意味着它们通过依赖于所有节点值来模拟频谱方案,并且因此,它们在所有长度尺度上保持准确性 . 相反,典型的有限差分方法仅局部精确(例如,点#13处的导数通常不依赖于点#200处的函数值) .
目前的研究领域是如何在紧凑的模板中最好地解决多种衍生物 . 这种研究,组合,紧凑的有限差分方法的结果是强大的和广泛适用的,尽管许多研究人员倾向于针对特定需求(性能,准确性,稳定性或诸如流体动力学的特定研究领域)调整它们 .
Ready-to-Go Routines
如上所述,可以使用 diff 函数(文档link)来计算相邻数组元素之间的粗略导数 .
MATLAB的 gradient 例程(link到文档)是一个很好的选择,用于许多目的 . 它实现了二阶中心差分方案 . 它具有计算多维导数并支持任意网格间距的优点 . (感谢@thewaywewalk指出这个明显的遗漏!)
我用Fornberg的方法(见上文)开发了一个小程序( nderiv_fornberg )来计算任意网格间距的一维的有限差分 . 我觉得它很容易使用 . 它在边界处使用6个点的侧面模板,在内部使用中心的5点模板 . 它可以在MATLAB File Exchange here中找到 .
Conclusion
数值微分领域非常多样化 . 对于上面列出的每种方法,有许多变体具有各自的优点和缺点 . 这篇文章几乎不是数值微分的完整处理 .
每个应用程序都不同 . 希望这篇文章能够为感兴趣的读者提供有组织的考虑因素和资源清单,以便选择适合自己需要的方法 .
可以使用特定于MATLAB的代码片段和示例来改进此社区wiki .