Monday, February 25, 2013

Improved Finite Difference Formula: Regular Grid

I am came across this excellent paper "Calculation of Weights in Finite Difference Formulas" by Bengt Fornberg, for deriving numerical differentiation formulae quickly.

If you don't have access to the website, check out this article on Scholarpedia (by Fornberg himself). You won't miss much.

Here, I will consider a very special case considered above.

Problem: Say you have a regular grid of equispaced points \(x_0, x_1, ..., x_n\) of \(n+1\) points, such that \(x_j = x_0 + j h\). Let \(f_j\) denote the function value at node \(x_j\). Let us say, we want to find the best approximation to the m-th derivative at a point \(x = x_0 + s h\), \(f^{(m)(x)}\), where \(x_0 \leq x \leq x_n\).

In the figure above, for example, n = 4 (5 points), m =2 (we are interested in an approximation to the second derivative) at x defined by s = 1.5.


We want to find the set of n+1 coefficients \(c_i\)s so that the explicit approximation is as good as possible
\[ f^{(m)}(x_0 + sh) = \sum_{i=0}^n c_i  f_i\]
Solution: The solution is a one-liner in Mathematica:

CoefficientList[Normal[Series[(x^s*Log[x]^m),{x,1,n}]/h^m],x]

Why does this work?

Set \(f(x) = e^{i \omega x}\) in the approximate formula above. Note that \(f^{(m)}(x_0 + sh) = (i \omega)^m e^{i \omega x_0} e^{i \omega s h}\), and \(f_j = e^{i \omega x_0} e^{i \omega j h}\).
\[(i \omega)^m e^{i \omega x_0} e^{i \omega s h} = \sum c_j e^{i \omega x_0} e^{i \omega j h} \]
Using the substitution \(e^{i \omega h} = \zeta\), which implies \( \ln \zeta = i \omega h \), we get
\[ \left( \frac{\ln \zeta}{h} \right)^m \zeta^s = \sum c_j \zeta^j \]
We want the approximation above to be as accurate as possible at small "h" or when \(\zeta = 1\). Pulling the "h" term aside, and recognizing that the RHS looks like a truncated Taylor series, the coefficients \(c_j\) can be obtained.

We can post facto throw in the factor of \(h^m\) in the denominator to fix things up.

No comments: