MATLAB: Hi, is this the right code for the pde (∂u/∂t= 2 + (∂^2 u)/(∂x^2 )) with boundary conditions of 0.5 and 3.5 and initial condition of x^2 – 0.5 where x ranges from -1 to 2

Question

function[c,f,s] = pde_description(x,t,u,dudx)
% pde_description outputs the descrition of the partial differential
% equation to be solved in terms of the variables c,f and s given the 
% input variables x,t,u and dudx
%


% Version no:1 Date of Creation: 03/05/2018 09:29  Author: Naomi Akinyede


c = 1;
f = 2*x + dudx;
s = 0;
end 
function u0 = pde_initial_cond(x)
% pde_initial_cond outputs the initial conditions of pde_description given
% the in the input variable x
%
% Version no:1 Date of Creation: 03/05/2018 09:29  Author: Naomi Akinyede
u0 = x^2 - 0.5;
end
function [pl,ql,pr,qr] = pde_boundary_cond(xl,ul,xr,ur,t)
% pde_boundary cond outputs the boundary conditions of pde_descrition given
% the xl,xr,ul and ur corresponding to the left and right boundary
% conditions and t corresponding to some time vector
%
% Version no:1 Date of Creation: 03/05/2018 09:29  Author: Naomi Akinyede
pl = ul - 0.5;
ql = 0;
pr = ur - 3.5;
qr = 0;
end
xx = linspace(-1,2,10);
tt = linspace(0,0.015,15);
sol = pdepe(0,@pde_description,@pde_initial_cond,@pde_boundary_cond,xx,tt);

Answer 1

It's more usual to define

c = 1;
f = dudx;
s = 2;

instead of

c = 1;
f = 2*x + dudx;
s = 0;

but the code looks correct to me.

For t-> Inf, you should get

u(x) = -x^2+2*x+3.5

as stationary solution for your problem. This is an indication whether your code is correct or not.

Best wishes

Torsten.