Hi PF!(adsbygoogle = window.adsbygoogle || []).push({});

I'm trying to finite difference a FTCS PDE $$h_t+\partial_x\left[h^3(h_{xxx}+(Kf'(h)-G)h_x\right] = 0$$ where ##K## and ##G## are constants. ##f'(h) = -(n(h^*/h)^n-m(h^*/h)^m)/h##. Boundary conditions are ##h_x=h_{xxx}=0## at both ends of the ##x## domain (however long you want to make it) and some initial height, call it ##h_0(1+\epsilon \cos(2 \pi x/\lambda)##. If we index such that ##h_1## is the first gridpoint, then the two BC's imply ##h_{-1}=h_3## and ##h_0 = h_2##, and similarly at the other end of the domain. Can you tell me if the following, which is my approach, is correct?

I should say, since ##h_0## and ##h_{-1}## don't exist in our indexing, to satisfy the boundary conditions I concatenate ##h## by 4 components, 2 on each side. These are ghost points that I use to satisfy the BC. Code wrote in MATLAB.

However, a rough stability check implies a parameter ##dt/dx^4## which is very large for the specified ##dx## value. Does anyone have any ideas on how to handle this problem, perhaps a different way?

Code (Text):

%% Parameters

hstar = 0.01;

n = 3;

m = 2;

theta = 50*pi/180;

K = 2*(1-cos(theta))/hstar;

G = 1;

eps = 10^(-3);

lambda = 2.66;

%% Numerics

tf = 200;% final time

dt = 0.01;% time step

t = 0:dt:tf;% time vector

xmin = 0;

xmax = 2.5;

dx = 0.01;% spacial step

x = xmin:dx:xmax;% x-vector

xc = xmin-2*dx:dx:xmax+2*dx;% concatenate h by 4 fictitious values (two at each

% domain endpoint) to accomodate BC. Then h actually starts at i=3 and

% ends at i = N-2.

%% IC

h0 = 0.1;% IC

hi = h0*(1+eps*cos(2*pi*xc/lambda));% IC

h = hi;

% preallocate

hnew = h;

h2x = zeros(size(x));

hxx = zeros(size(x));

hxxxx = zeros(size(x));

fh = zeros(size(x));

fhh = zeros(size(x));

N = length(h);

%% FTCS difference equations for eq. (4)

for j = 1:size(t,2)

fh = (m*(hstar./h).^m-n*(hstar./h).^n)./h;% f'(h)

fhh = ((hstar./h).^n*n*(n+1)-(hstar./h).^m*m*(m+1))./h.^2;% f''(h)

% BC

h(1) = h(5);

h(2) = h(4);

h(N) = h(N-4);

h(N-1) = h(N-3);

% FTCS

for i = 3:N-2

h2x(i) = (h(i+1)^2-h(i-1)^2)/(2*dx);% (h^2)'(x)

hxx(i) = (h(i-1)-2*h(i)+h(i+1))/dx^2;% h''(x)

hxxxx(i) = (h(i-2)-4*h(i-1)+6*h(i)-4*h(i+1)+h(i+2))/dx^4;% h''''(x)

hnew(i) = h(i) - dt*h(i)^3*(hxxxx(i)+(K*fh(i)-G)*hxx(i)+K*h2x(i)*fhh(i));

end% for i

h = hnew;% update time

end% for j

h = h(3:N-2);% delete fictitious points used for BC

**Physics Forums - The Fusion of Science and Community**

Join Physics Forums Today!

The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

# A Finite Difference

Tags:

Have something to add?

Draft saved
Draft deleted

Loading...

Similar Threads - Finite Difference | Date |
---|---|

A Convergence order of central finite difference scheme | Nov 8, 2017 |

A Better way to find Finite Difference | Oct 16, 2017 |

A Finite difference of fourth order partial differential | Aug 30, 2017 |

A Finite difference formulation ideas (journal verification) | May 15, 2017 |

**Physics Forums - The Fusion of Science and Community**