数组的有限差分导数

Finite difference derivative of an array

我正在尝试对数组求导,但遇到了问题。该数组是二维的,xy 方向。我想使用中心差分离散化沿着 xy 求导数。该数组具有数字的随机值,没有值是 NaN。我将提供下面代码的基本部分来说明我的观点(假设数组 u 已定义并且已经输入了一些初始值)

integer :: i,j
integer, parameter :: nx=10, ny=10
real, dimension(-nx:nx, -ny:ny) :: u,v,w
real, parameter :: h

do i=-nx,nx
  do j=-ny,ny

    v = (u(i+1,j)-u(i-1,j))/(2*h)
    w = (u(i,j+1)-u(i,j-1))/(2*h)

  end do 
end do

注意,假设数组 u 在我找到 vw 之前已定义并填充。 vw 应该分别是数组 u 沿 x 和沿 y 的导数。这是获取数组导数的正确方法吗?

我在您的代码中发现了几个问题。

1.You 一定要小心左手边的东西。

v = (u(i+1,j)-u(i-1,j))/(2*h)

意味着整个数组 v 将被设置为相同的数字。你不希望这个循环。在一个循环中你想一次只设置一个点

v(i,j) = (u(i+1,j)-u(i-1,j)) / (2*h)

和 2) 您正在越界访问数组。您可以保留简单循环,但必须使用存储边界值的边界点 "ghost points"。如果我假设点 -nx,nx,-nyandny` 位于边界上,那么您只能使用域内的中心差来计算导数:

do i=-nx+1,nx-1
  do j=-ny+1,ny-1

    v(i,j) = (u(i+1,j)-u(i-1,j)) / (2*h)
    w(i,j) = (u(i,j+1)-u(i,j-1)) / (2*h)

  end do 
end do

如果你需要边界上的导数,你必须使用像

这样的on-sided difference
  do j=-ny+1,ny-1

    v(nx,j) = (u(nx,j)-u(nx-1,j)) / h
    w(nx,j) = (u(nx,j+1)-u(nx,j-1)) / h

  end do