A Discretization for a fourth-order PDE (and solution)

1. Jul 1, 2017

maistral

Hi. I have this PDE that governs an L x L plate (similar to the Poisson equation, it seems)

with boundary conditions F = 0 and F" = 0 along the edges. I have successfully solved the problem by setting up an equality W = ∇2F then I solved the two PDEs simultaneously:
W = ∇2F (boundary condition F = 0 along the edges)
φ = ∇2W (boundary condition W = 0 along the edges since W = F").

Simple discretization using second order derivatives followed, then poof. Answer.

What I would like to know is how to solve this without using substitutions (I mean, directly attacking the differentials using fourth-order finite differences). What I did is I took the second order derivative of the second-order finite differences with respect to the same variables (for the non-mixed derivative) and the other variable (for the mixed derivative).

My problem is coming from the boundary conditions. The boundary F = 0 is a no-brainer. The problem is how should I implement the second order boundary condition F" = 0. Does anyone have an idea on how to execute this properly?

What I did is that I tried using the second order three-point formula for the second derivative using double the steps then equated it to zero, then applied the resulting formula to the ghost points formed outside the boundary. ie.
(F(x+2Δx,y) - 2F(x,y) + F(x-2Δx,y)) / (2Δx)2 = 0 say, for points beside the right border then I solved for the ghost point F(x+2Δx,y). I applied the same to the points beside the top, bottom, and left borders.

And, this happened:

How should I implement the second order condition? I hope for any kind and generous response. Thanks!

2. Jul 3, 2017

Telemachus

Perhaps the algorithm is unstable. You could try with the Von Neumann stability analysis to see if the problem is in the scheme. However, the instability could be due to the boundary conditions. I think what you did is fine, but the problem is in the algorithm (of course this is under the assumption that you have made no coding mistakes). You should look for a stable algorithm to solve this equation, if that is actually the problem that makes everything blow up.

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

3. Jul 6, 2017

maistral

Yeah, apparently the algorithm is unstable. A simple adjustment to Gauss-Seidel SOR made it work.

Thanks!