Finite Difference Approximation of 1st Order Error Bi-Harmonic Operator

Hello,

I am working on a summer project at my university which involves simulating a pair of reaction-diffusion equations. I have chosen to approximate the differentials with the finite difference method. I have successfully applied the method to the laplacian operator (shown below), however, I appear to be having issue with the bi-harmonic operator (shown below).

Laplacian_2D_1stOrder

function out = Laplacian2D(in, res)

h = res;
k = res;

f_xx = (circshift(in,[ 0, 1]) - 2 .* circshift(in,[ 0, 0]) + circshift(in,[ 0, -1])).*h;
f_yy = (circshift(in,[ 1, 0]) - 2 .* circshift(in,[ 0, 0]) + circshift(in,[ -1, 0])).*k;

out = f_xx + f_yy;

end

BiHarmonic_2D_1stOrder

function out = Biharmonic2D(in, res)

h = res;
k = res;

f_xxxx = (circshift(in,[ 0, -2]) - 4.*circshift(in,[ 0, -1]) + 6.*circshift(in,[ 0, 0]) - 4.*circshift(in,[ 0, 1]) + circshift(in,[ 0, 2])).*res;
f_yyyy = (circshift(in,[ -2, 0]) - 4.*circshift(in,[ -1, 0]) + 6.*circshift(in,[ 0, 0]) - 4.*circshift(in,[ 1, 0]) + circshift(in,[ 2, 0])).*res;
f_xxyy = (circshift(in,[ -1, -1]) - 2.*circshift(in,[ 0, -1]) + circshift(in,[ 1, -1]) - 2.*circshift(in,[ -1, 0]) + 4.*circshift(in,[ 0, 0]) - 2.*circshift(in,[ 1, 0]) + circshift(in,[ -1, 1]) - 2.*circshift(in,[ 0, 1]) + circshift(in,[ 1, 1])).*res;

out = f_xxxx + f_xxyy + f_xxyy + f_yyyy;

end

Quoting from Wikipedia:

In mathematics, finite-difference methods (FDM) are numerical methods for solving differential equations by approximating them with difference equations, in which finite differences approximate the derivatives. FDMs are thus discretization methods.

I derived the approximations by following the work of David Eberly "Derivative Approximation by Finite Difference" Link: http://read.pudn.com/downloads203/sourcecode/game/955182/3D%20Geometry%20Tuts/FiniteDifferences.pdf

The x and y resolutions are the same. I found that different resolutions led to spatial stretching. (I think this is obvious but would be happy to learn if there is a reason why or a situation where one would have different x and y spatial resolutions?)

When I call my functions, they appear to work, but the Laplacian appears far better behaved than the bi-harmonic operator. I tested both on the MATLAB Peaks function and compared them to MATHEMATICA's built in laplacian and hiharmonic operator functions and they returned the same results (roughly, I assume the difference is between my approximation and MATHEMATICA's more accurate differentiation).

My questions to the community are:

1. Have I implemented these correctly, specifically the resolutions, but also the addition at the end of the functions to return the overall approximation?
2. Is there a more efficient way to implement them within MATLAB? I know that MATLAB has a laplacian operator but I was unable to find a Bi-Harmonic function, hence deriving it myself.
3. When applying these to PDEs, is the error order acceptable? The PDEs will be applied to a pair of Cahn-Hilliard equations:

Should they therefore be of a higher order for a stability or a convergence reason?

Thank you for your assisstance. Regards, Matthew

Answer 1

[Elizabeth Reese]

I believe you have these implemented correctly and that both have second-order accuracy. I am not sure how you are defining res, but make sure that there is a correct translation from resolution (generally points/distance) to grid spacing ( h = distance between points). Then, the Laplacian should be over h^2 and the biharmonic should be over h^4.

As far as efficiency, you can still use the built-in del2 function to compute the discrete Laplacian in MATLAB. You can apply it twice to your data to get the BiHarmonic.

As far as the x and y resolutions are concerned, this will depend on the mesh that you are working with but they will be the same if you have a uniform mesh.

If you are interested in improving the accuracy, I recommend taking a look at the following page and just expanding your finite differences to include more points.

https://en.wikipedia.org/wiki/Finite_difference_coefficient

[Matthew Noble]

Dear Liz,

Thank you very much, really helpful and great food-for-thought!

Regards, Matthew