The following code is not working
2 views (last 30 days)
Show older comments
function [x,h] = Transient()
%% Initialization
theta = 50; % degrees - angle of repose of material
n = 8; % r.p.m. - rotation speed of tube
rho = 660; % kg/m3 - bulk density of material
R = 0.075; % m - inner radius of tube
beta = 3; % degrees - inclination angle of kiln
L_tube = 2; % m,lenght of the kiln
mass_flow_rate = 0.5; % kg/min - mass flow of material
x_points=100;
t_points=100;
%% Solving initial condition
xspan_1=[0 L_tube];
h0 = 0.000001; % initial condition (bed height at the discharge end - average of particle diameter)
sol = ode45( @saeman_bed ,xspan_1 ,h0); % solver for unsteady nonliner ode
x = linspace(0,L_tube,x_points);
u0 = deval(sol,x)
%% Solving with pdepe
m = 0;
x = linspace(0,1,x_points);
t = linspace(0,2,t_points);
sol = pdepe(m,@pdex1pde,@(x0)pdex1ic(x0,x,u0),@pdex1bc,x,t);
% Extract the first solution component as u.
u = sol(:,:,1)
function [c,f,s] = pdex1pde(x,t,u,DuDx)
f_H=power((2*u/R)-power(u,2)/power(R,2),0.5);
u_T=2*pi*n*R;
c = f_H;
f = (u_T*R*cotd(theta))*power(f_H,1.5).*DuDx;
s = u_T*tand(beta)/sind(theta)*f_H.*power(1-f_H,0.5).*DuDx;
end
% --------------------------------------------------------------
function u0 = pdex1ic(x0,x,u)
u0 = interp1(x,u,x0);
% -------------------------------
end
function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)
pl = ul+0.000001 % exit
ql = 0;
pr =(3*tan(deg2rad(theta))* mass_flow_rate/rho) / (4*pi*n*power(R,3)) * power((2*ul/R) - power((ul/R),2),(-3/2)) -tan(deg2rad(beta))/cos(deg2rad(theta));
qr = -1;
end
%% function for initial condition
function [dhdx,h] = saeman_bed (z,h);
dhdx = (3*tan(deg2rad(theta))* mass_flow_rate/rho) / (4*pi*n*power(R,3)) * power((2*h/R) - power((h/R),2),(-3/2)) -tan(deg2rad(beta))/cos(deg2rad(theta));
end
end
8 Comments
Rik
on 3 Jan 2019
I see I made a copy-paste error, so indeed you should fix that line to
bcfun=@pdex1bc;
But your error seems to be a bit deeper. Can you show the output of
which pdepe -all
I'm suspecting somehing is shadowing the correct function, because your input looks like it should be valid, but it still returns this error.
Answers (1)
madhan ravi
on 4 Jan 2019
You have file named pdepe please rename it or delete it becuase it‘s shadowing the inbuilt function pdepe()
0 Comments
See Also
Categories
Find more on Geometry and Mesh 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!