# what is wrong in my pdepe code ?

2 views (last 30 days)
Kazi Haq on 21 Apr 2018
Answered: Kazi Haq on 21 Apr 2018
I was trying to solve a very simple pde with matlab pdepe solver . du/dt=D(d2u/dr2+1/r(du/dr)+0.01, IC=0.1 and BC: du/dr=0 at 0 and 1. I set it to solve for 5 spatial grid points while allowing the source term (0.01) at the 5th grid point only and for the rest it is zero . Even though my diffusion coefficient D=0 , why do i still get a rise of U at the 4th grid point (figure 1) . I did the same with method of lines and there were no increase in the value of U at the 4th grid point for the similar condition (figure 2). i have added the code here and also the result. Can anyone tell me where i did the mistake in pdepe ? why diffusion is occurring between the 5th and 4th grid ?
if true
function pdex
m = 1; % cyndrical coordinate
x = linspace(1,5,5); % space grid
t = linspace(0,20,1000);% time grid
opts=odeset('MaxStep', 0.4);
sol = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t,opts);
u = sol(:,:,1);
figure
% A time plot of u
plot(t,u(:,:,1));
xlabel('Distance x')
ylabel('Time t')
end
if true
% -------------description of the function
function [c,f,s] = pdex1pde(x,t,u,DuDx)
D=0; % diffusion coefficient
c = 1;
f =D*DuDx;
if (x>4)
s = 0.01;
else
s=0;
end
end
if true
if true
% -------------initial condition function u0 = pdex1ic(x) u0 = 0.1; end
end
if true
% ------- no flux boudary condition
function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)
pl = 0;
ql = 1;
pr =0;
qr = 1;
end  Kazi Haq on 21 Apr 2018
Thanks very much, Bill. I have done that before and when i increase the mesh number the solution approaches what is generated by MOL. Though this is a sample test while i wanted to use pdepe for a much larger problem where increasing the number of meshes would take significant computing resources and wanted to minimize the number of meshes. However, as you have suggested , i think a small number of meshes would not work with pdepe.

Bill Greene on 21 Apr 2018
The problem you have described is ill-posed. The source term must be prescribed over some finite length of the domain, not simply "at the 5th grid point only". Also, I think you have a misunderstanding about the role of the mesh in pdepe. The accuracy of the solution very much depends on the density of the mesh points. Especially with your discontinuity in the source term, you will need many more mesh points (e.g. 20, 50, or more) to obtain an accurate solution. You should systematically increase the number of grid points until the solution does not change significantly.