Step Function as Input Boundary Condition to MATLAB PDEPE solver
4 views (last 30 days)
Show older comments
Dear Matlab Users I have the following code to find the transient solution to the 1D heat conduction equation. The code below intends to use a rectangular heat pulse as the initial temperature but i receive an error ( Warning: Failure at t=0.000000e+00. Unable to meet integration tolerances without reducing the step size below the smallest value allowed (7.905050e-323) at time t). Can anyone tell me what is issue with the boundary condition as I set it and if there is a better way to deal with the matter.
function out = parabolic_silicon(~)
global rho cp k T_i T_o
L = 10; % [cm]
k = 1.3; % [W/cmK]
rho = 2.33; % [g/cm^3]
cp = 0.7; % [J/gK]
T_i = 295; % [K]
T_o = T_i + 5.2e-3; % [K]
t_end= 1; % [s]
m = 0;
x = linspace(0,L,500);
t = linspace(0,t_end,500);
sol = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t);
% Extract the first solution component as u.
Temperature = sol(:,:,1);
L = strcat('t = ',num2str(t(:)),' Seconds');
plot(x,sol);
% legend(L,'location','northeast');
% % --------------------------------------------------------------
% OUTPUT
out = {x,t,sol};
% % --------------------------------------------------------------
function [c,f,s] = pdex1pde(x,t,u,DuDx)
global rho cp k
c = rho*cp;
f = k*DuDx;
s = 0;
% --------------------------------------------------------------
function u0 = pdex1ic(x)
u0 = 295;
% --------------------------------------------------------------
function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)
global T_i T_o
pl = ul - (T_i + (T_o*(heaviside(t) - heaviside(t - 0.5))));
ql = 0;
pr = ur - T_i;
qr = 0;
0 Comments
Accepted Answer
Bill Greene
on 22 Mar 2018
I suggest you smooth the sharp corners in your boundary condition profile using a Heaviside function similar to the code below. Note, you can change the sharpness of the corners by varying the r variable.
function k=heavisideSmooth(t)
k1 = 0; k2 = 1; c = 0;
r = 1e3;
ecr = exp(c*r);
erx = exp(r*t);
k = (k1*ecr + k2*erx)./(ecr+erx);
end
2 Comments
Bill Greene
on 22 Mar 2018
The error message indicates the ODE solver wasn't able to take even a tiny, first step. So I assumed that was due to the sharp discontinuity in the BC at t=0.
More Answers (0)
See Also
Categories
Find more on Calculus in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!