0

I have a query on numerical differentiation that's beyond the language used. Suppose I have an array of n points x and f(x) and I want to take the first derivative of f(x). Every method will consume points making the derivative array shorter than the function, so how can a "lengthen" the array in a smart way. For example I want to take the derivate using the Five-point stencil, that is

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

where f(n) is the function evaluated at the n-th point. So with this method the f' array is shorter by 4 points. How can I lengthen this array in a smart way and if it's possible in a way that produces an error similar as this 5 point stencil method?

Vladimir F Героям слава
  • 57,977
  • 4
  • 76
  • 119
Daniel
  • 105
  • 1
  • 9
  • At endpoints you can use non-centered formulas. https://en.wikipedia.org/wiki/Finite_difference_coefficient#Forward_and_backward_finite_difference It is probably also possible to get the formulas for the intermediate cases, but that might be splitting hairs. –  Oct 29 '16 at 15:10

1 Answers1

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

The middle of the vector is easy, just the ends needs something special.

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

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

Of course it should not take long to realise that you could get rid of Temp entirely and do the ends manually. It depends on how long your vectors are, and whether you need some alignment. In some parallel world you may want the temp array as a function, For a simple serial implementation the temp may be conceptually easier to grasp. Also you mentioned lengthening the array, which is really accordioning the vector by grabbing each end and moving your paws further apart. Extending it means that you have chin-scratch to keep track of the index, where as in the above implimentation the X, X', and Temp are all aligned in index value. Since fortran can go from any#:AnyOther# this the perfect example of when you may want to do that.

Holmz
  • 714
  • 7
  • 14