MATLAB: Is the boundary condition for the heat transfer problem not satisfied using the PDEPE function in MATLAB

boundaryconditionsdifferentialheatMATLABneumannotpartialpdepesatisfiedslopetransferviolated

The left boundary of the heat transfer problem that I am trying to solve is not correct. I imposed a Neuman boundary condition on the left bound with "du/dx=0". The solution should look flat with a slope of zero on the left boundary but the solution plot shows a positive slope.

Best Answer

We have verified that there is a limitation in MATLAB 7.0 (R14)+ in the way that the PDEPE function handles boundary conditions with certain values.
PDEPE imposes Boundary Conditions in terms of flux "f()":
    p() + q()*f() = 0
At the left bound, "x=0". If the flux term is in a form such that it is zero for any specified boundary condition, the boundary condition will be ignored.
For example, "f() = (A+B*x)*du/dx" where "A" and "B" are constants. If "A = 0", flux is "(0+B*x)*du/dx".
You specify the boundary condition on the left bound to be "du/dx =0".
On the left bound, "B*(0)*du/dx = 0". This boundary condition is always satisfied regardless the value of "du/dx" and MATLAB will disregard the previously specified boundary condition.
Hence this particular combination:
     f() = x*du/dx
     0 <= x 
     Boundary Condition at 0,  du/dx = 0
gives an incorrect boundary result.
To work around this issue, shift the value of "x" from 0 to a very small number called delta that is close to 0. Define another point next to delta as follow:
x = [delta, 1.01*delta, linspace(2*delta,1000,100)];
The boundary condition on the left side should also be changed from:
pl = 0; 
ql = 1;
pr = 0; 
qr = 1;
to
pl = 0; 
ql = 1/(A+Bx*xl);
pr = 0; 
qr = 1;
Related Question