数值微分中如何处理边界点

Question

我有一个关于数值微分的查询，超出了所使用的语言。假设我有一个包含n个点x和f(x)的数组我想对f(x)求导。每个方法都会消耗点，使得导数数组比函数短，那么如何以一种智能的方式"延长"数组呢?例如，我想使用五点模板取衍生物，即

    f'(0) = 1/12 h (-f(-2) + 8 f(-1)- 8 f(1) + f(2))

，其中f(n)为函数在第n点处求值。因此，用这种方法，f'阵列缩短了4个点。我怎样才能以一种聪明的方式延长这个数组，如果它有可能以一种方式产生类似于这个5点模板方法的错误?

Answer 1

SUBROUTINE Diff(X, N, XPrime)
INTEGER           , INTENT(IN   ) :: N
REAL, DIMENSION(N), INTENT(IN   ) :: X
REAL, DIMENSION(N), INTENT(INOUT) :: XPrime
enter code here
REAL, DIMENSION(-1:N+2)           :: Temp
INTEGER                           :: I
!Use temp for X
Temp(1:N) = X
!... Temp(O)   = X(1) - (X(2) - X(1)  )
!... Temp(N+1) = X(N) + (X(N) - X(N-1))
!Your code here
!output XPrime from 1:N
END SUBROUTINE Diff

中间很容易，只是末端需要一些特殊的东西。

For X' maybe Temp(0:N+1).

For X " maybe Temp(-1:N+2).

当然，你很快就会意识到，你可以完全摆脱Temp，手动完成任务。这取决于向量的长度，以及是否需要对齐。在一些并行世界中，您可能希望temp数组作为一个函数，对于一个简单的串行实现，temp可能在概念上更容易掌握。你还提到延长数组，这是真正通过抓住每一个末端和移动你的爪子进一步分开的矢量。扩展它意味着您可以轻松地跟踪索引，在上面的实现中，X、X'和Temp都在索引值中对齐。由于fortran可以从任意#:AnyOther#开始，这是您可能想要这样做的完美示例。