我有一些离散数据，在，中不等距。$x$$y=f(x)$

我想在某个时候使用数值有限差分法来计算的二阶导数。$y$

我正在使用 Fornberg 方法，这里和这里都有很好的描述，并且在 Fortran 中工作。

Program fornberg

  implicit none
  integer, parameter :: dp = selected_real_kind(32,307)
  integer, parameter :: Nrows =1d4 
  integer (kind=dp) :: j, counts,m,k, N, i, nn, mn
  real(kind=dp), dimension(Nrows,2) :: SomeData
  real(kind=dp), dimension(:,:), allocatable ::carray
  real(kind=dp), dimension(:), allocatable :: xdata,ydata
  real(kind=dp) :: c1,c2,c3,c4,c5, z


  !Load the data from an unformatted binary file


  open(unit = 10, file='example.dat', form = 'unformatted')
  read(10) SomeData
  close(10)


  !Count the number of non-zero rows
  do j=1,Nrows
  if (SomeData(j,1) .EQ. 0) then
  counts = j-1
  EXIT
  endif
  enddo


  !Allocate the x and y arrays to hold the data with indexing starting at 0
  ALLOCATE(xdata(0:counts-1)) !specific zero indexing
  ALLOCATE(ydata(0:counts-1))


  !Populate the arrays
  do j = 1,counts
  xdata(j-1) = SomeData(j,1)
  ydata(j-1)= SomeData(j,2)
  enddo


  !Define the point at which we want the derivative evaluated
  z = xdata(0)

  !Length of data array
  N = counts
  nn = N - 1

  !Define the maximum order of the derivative
  m = 2


  !Set up zeroes array
  ALLOCATE(carray(0:N-1, 0:m)) !zero indexing

  !Determine the weights via the Fornberg algorithm

  c1 = 1.0_dp
  c4 = xdata(0) - z
  carray = 0.0_dp
  carray(0,0) = 1.0_dp

  do i = 1,nn

  mn = min(i,m)

  c2 = 1.0_dp
  c5 = c4
  c4 = xdata(i) - z




  do j = 0, i-1


      c3 = xdata(i) - xdata(j)
      c2 = c2*c3


      if ( j .EQ. i-1) then

          do k=mn,1,-1
              carray(i,k) = c1*(k*carray(i-1,k-1) - c5*carray(i-1,k))/c2
          enddo
          carray(i,0) = -c1*c5*carray(i-1,0)/c2
      endif

      do k = mn,1,-1
          carray(j,k) = (c4*carray(j,k) - k*carray(j,k-1))/c3
      enddo
      carray(j,0) = c4*carray(j,0)/c3
  enddo

  c1 = c2



  enddo

end program fornberg

我的问题是二阶导数的权重很快变得很大。

查看代码，这是c2 = c2*c3命令的直接结果。对于c3 > 1大量迭代（数据集约为 300 行），我对权重如何“合理”感到困惑。

任何指导将不胜感激。如有必要，我还可以提供数据集。