0

我对超出所使用语言的数值微分有疑问。假设我有一个包含 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 点模板方法的错误？

4

1 回答 1

0 
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

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

对于 X' 可能是 Temp(0 :N+1)。

对于 X'' 可能是 Temp(-1:N+2)。

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

于 2016-10-29T04:13:50.153 回答