solving system of non linear PDEs

Hello

I provide a solution to Allen cahn equation using spectral methods and Euler's formula.

The Allen Cahn Equation is with BC , and BC

Now my question is if I have two systems of nonlinear partial differential equations (PDEs)., like with BC , . How can I use the spectral method and Euler's formula to solve these problems?

Thank you. I truly appreciate any help

This script outlines the solution of the Allen-Cahn equation using spectral methods and Euler's formula:

% Allen-Cahn eq. u_t = u_xx + u - u^3, u(-1)=-1, u(1)=1
% (compare p6.m and p32.m).The challenge idea here is that the BC satisfy
% the IC. so Zeroing out the first and last rows of D2, This ensures that
% when you comput D2*(v - x) ,the first and last rows of this expression are
% always zero. So the time derivative at the boundary points is zero, i.e., those boundary values do not 
% change during time stepping, (  v + dt*(eps*D2*(v-x) + v - v.^3) = 0 at u(1) and u(end))


% Differentiation matrix and initial data:
N = 20; [D,x] = cheb(N); D2 = D^2; % use full-size matrix
D2([1 N+1],:) = zeros(2,N+1); % You're assigning values to two specific rows of the matrix D2 — the first and last.
                              % D2([1 N+1],:) means: "access rows 1 and N+1, and all columns. zeros(2,N+1) You're setting those rows
                              % to a matrix of zeros with size 2 by N+1


eps = 0.01; dt = min([.01,50*N^(-4)/eps]);
t = 0; v = .53*x + .47*sin(-1.5*pi*x);

% Solve PDE by Euler formula and plot results:
tmax = 100; tplot = 2; nplots = round(tmax/tplot);
plotgap = round(tplot/dt); dt = tplot/plotgap;
xx = -1:.025:1; vv = polyval(polyfit(x,v,N),xx);
plotdata = [vv; zeros(nplots,length(xx))]; tdata = t;

for i = 1:nplots
for n = 1:plotgap
t = t+dt; v = v + dt*(eps*D2*(v-x) + v - v.^3); % D2*(v-x): the function u(x,t)−x vanishes at the boundaries
end
vv = polyval(polyfit(x,v,N),xx);
plotdata(i+1,:) = vv; tdata = [tdata; t];
end
clf, subplot('position',[.1 .4 .8 .5])
mesh(xx,tdata,plotdata), grid on, axis([-1 1 0 tmax -1 1]),
view(-60,55), colormap([0 0 0]), xlabel x, ylabel t, zlabel u

[Torsten]

Your second equation does not need boundary conditions - it's an ODE in disguise. Only an initial condition for s at t=0 is required.

[Erm]

I assumed the second equation because the main task is to understand how to solve system of partial differential equations (PDEs) using spectral methods with Euler. In general, if we have a system of equations like the one mentioned, we need four boundary conditions since we have two variables involved like and , and we need two initial conditions as well since we have and .

[Torsten]

I get substantially different results if I use the MATLAB standard solver "pdepe" for your problem. Maybe the cheb.m I added to your code is incorrect. If not, first take care that your code is able to solve only one PDE correctly.

% Allen-Cahn eq. u_t = u_xx + u - u^3, u(-1)=-1, u(1)=1
% (compare p6.m and p32.m).The challenge idea here is that the BC satisfy
% the IC. so Zeroing out the first and last rows of D2, This ensures that
% when you comput D2*(v - x) ,the first and last rows of this expression are
% always zero. So the time derivative at the boundary points is zero, i.e., those boundary values do not 
% change during time stepping, (  v + dt*(eps*D2*(v-x) + v - v.^3) = 0 at u(1) and u(end))


% Differentiation matrix and initial data:
N = 20; [D,x] = cheb(N); D2 = D^2; % use full-size matrix
D2([1 N+1],:) = zeros(2,N+1); % You're assigning values to two specific rows of the matrix D2 — the first and last.
                              % D2([1 N+1],:) means: "access rows 1 and N+1, and all columns. zeros(2,N+1) You're setting those rows
                              % to a matrix of zeros with size 2 by N+1


eps = 0.01; dt = min([.01,50*N^(-4)/eps]);
t = 0; v = .53*x + .47*sin(-1.5*pi*x);

% Solve PDE by Euler formula and plot results:
tmax = 100; tplot = 2; nplots = round(tmax/tplot);
plotgap = round(tplot/dt); dt = tplot/plotgap;
xx = -1:.025:1; vv = polyval(polyfit(x,v,N),xx);
plotdata = [vv; zeros(nplots,length(xx))]; tdata = t;

for i = 1:nplots
for n = 1:plotgap
t = t+dt; v = v + dt*(eps*D2*(v-x) + v - v.^3); % D2*(v-x): the function u(x,t)−x vanishes at the boundaries
end
vv = polyval(polyfit(x,v,N),xx);
plotdata(i+1,:) = vv; tdata = [tdata; t];
end
%clf, subplot('position',[.1 .4 .8 .5])
%mesh(xx,tdata,plotdata), grid on, axis([-1 1 0 tmax -1 1]),
%view(-60,55), colormap([0 0 0]), xlabel x, ylabel t, zlabel u
pdefun = @(x,t,u,dudx)deal(1,dudx,u-u^3);
pdeic = @(x).53*x + .47*sin(-1.5*pi*x);
pdebc = @(xl,ul,xr,ur,t) deal(ul+1,0,ur-1,0);
m = 0;
sol = pdepe(m,pdefun,pdeic,pdebc,xx,tdata);
u = sol(:,:,1);
figure(1)
plot(xx,[plotdata(1,:);plotdata(2,:);plotdata(3,:);plotdata(end,:)])

figure(2)
plot(xx,[u(1,:);u(2,:);u(3,:);u(end,:)])

●

function [D,x] = cheb(N)
% CHEB compute D = differentiation matrix, x = Chebyshev grid, from
% Trefethen

if N==0, D=0; x=1; return, end
x = cos(pi*(0:N)/N)';
c = [2; ones(N-1,1); 2].*(-1).^(0:N)';
X = repmat(x,1,N+1);
dX = X-X';
D = (c*(1./c)')./(dX+eye(N+1));
D = D-diag(sum(D'));
end

[Erm]

Cheb.m is correct.

Can I use the Euler method for a system of nonlinear partial differential equations (PDEs)?

Could you recommend a book or reference?

Thank you so much for your help,I really appreciate it

[Torsten]

You can replace the spatial derivatives by finite-difference approximations and then use the Euler method as time integrator (also for systems of PDEs) .

In case of the above PDE:

uxx(x_i ) ~ (u(x_(i+1))-2*u(x_i)+ u(x_(i-1)))/dx^2

for 2 <= i <= N-1,

thus you have to solve for x_new(x_i):

u_new(x_1) = -1

u_new(x_N) = 1

(u_new(x_i)-u_old(x_i))/dt = (u_old(x_(i+1))-2*u_old(x_i)+ u_old(x_(i-1)))/dx^2 + u_old(x_i) - u_old(x_i)^3

for 2 <= i <= N-1.

But usually the system above is stiff, and a very small time step dt is required for the Explicit Euler method in time to keep the solution stable.

Look up "Method-of-lines" for more information about the solution method. The advantage in using this method is that you can use a stiff MATLAB integrator like ODE15S for the time integration and that you only need to care about the spatial discretization.

[Erm]

thank you. I appreciate your help.

I will read about MOF.