I'm working on a finite difference scheme to simulate my model equations. As part of it I have to discretise the linear inhomogeneous diffusion operator:
$$\frac{\partial}{\partial x} \left( D(x) \frac{\partial u}{\partial x}\right)$$
Using a central difference on the outer derivative I obtain: $$\frac{\partial}{\partial x} \left( D(x) \frac{\partial u}{\partial x}\right) \bigg|_r^n \approx \frac{1}{2h} \left[ \left(D\frac{\partial u}{\partial x}\right)\bigg|_{r+1}^n - \left(D\frac{\partial u}{\partial x}\right)\bigg|_{r-1}^n\right]$$ Subsequently applying a forward and backward difference to the two derivatives at $x_{r-1}^n$ and $x_{r+1}^n$ gives me the following approximation: \begin{align} \frac{\partial}{\partial x} \left( D(x) \frac{\partial u}{\partial x}\right) \bigg|_r^n &\approx \frac{1}{2h^2} \left[ D_{r+1}^n\left(u_{r+1}^n - u_r^n\right) - D_{r-1}^n\left(u_r^n - u_{r-1}^n\right) \right]\\ &=\frac{1}{2h^2} \left[D_{r-1}^n u_{r-1}^n - \left(D_{r-1}^n + D_{r+1}^n \right) u_r^n +D_{r+1}^n u_{r+1}^n \right] \end{align}
However, if I now let $D=1$ everywhere so that I go back to homogeneous linear diffusion this gives: $$\frac{\partial^2 u}{\partial x^2} \bigg|_r^n \approx \frac{1}{2h^2} \left[u_{r-1}^n - 2 u_r^n + u_{r+1}^n \right],$$ which is 1/2 times what one would expect.
I suspect there's something wrong in the way I mix central, forward and backward differences, but I don't see what. How can the factor of 1/2 be explained?
Thanks very much for your help!
Maxi
You used a first order (one-sided) approximation of the first derivative in the second step; apparently this is not enough. How about
$ \partial_x (D \partial_x) u |_r \approx $
$ \frac1h \left[ (D \partial_x u)_{r+1/2} - (D \partial_x u)_{r-1/2} \right] \approx $
$ \frac1h \left[ \frac1h D_{r+1/2}(u_{r+1} - u_r) - \frac1h D_{r-1/2} (u_r - u_{r-1}) \right] $
?
That gives the usual central finite difference when $D$ is constant.
Btw., your operator is a linear inhomogeneous diffusion operator in divergence form.