Eng-Tips is the largest engineering community on the Internet

Intelligent Work Forums for Engineering Professionals

  • Congratulations waross on being selected by the Tek-Tips community for having the most helpful posts in the forums last week. Way to Go!

Start-up Poiseuille Flow

Status
Not open for further replies.

longbow75

Chemical
Aug 22, 2006
7
0
0
US
Hi,
I'm trying to solve start-up flow in a hollow cylinder with matlab and finite difference method. I'm still new in matlab and I have this code:

clear;clc;
h = 1;
Dr = 0.01; Dt = 0.01;
nmax = 4000; epsilon = 0.0001;
r = [0:Dr:h];
t = [0:Dt:h];
beta=0.1;
n1 = fix(h/Dr)+1;
m1 = fix(h/Dt)+1;
% Initialize solution matrices
velocity = zeros(n1,m1)
velocityN =velocity
% Iterative process for the solution
for k = 1:nmax
% boundary conditions:
velocityN:),n1) = zeros(n1,1);

% Calculate psiN in interior points
for i = 2:n1-1
for j = 1:m1
velocityN(i,j+1)= (beta*(1-(1/(2*i)))*velocity(i-1,j)+(1-2*beta)*velocity(i,j)+beta*(1+1/(2*i))*velocity(i+1,j));
end;
end;
% Check convergence
nfail = 0;
for i = 1:n1
for j = 1:m1
if abs(velocityN(i,j)-velocity(i,j))>epsilon
nfail = nfail+1;
end;
end;
end;

if nfail > 0
velocity = velocityN;
else
return;
end;
end;
fprintf('No convergence after %i iterations.',k);
[X,Y] = meshgrid(r,t);
contour(X,Y,velocity',10)

however, for r=0 I have dv/dr=0 as the second boundary condition. How do I incorporate this into the program. dV/dr should be:
velocityN(0,j+1)=((1-4*beta)*velocity(0,j)+4*beta(velocity(1,j)))

 
Status
Not open for further replies.
Back
Top