- **Main CFD Forum**
(*http://www.cfd-online.com/Forums/main/*)

- - **start-up flow poiseuille w/finite difference**
(*http://www.cfd-online.com/Forums/main/80968-start-up-flow-poiseuille-w-finite-difference.html*)

start-up flow poiseuille w/finite differencehello,
I'd appreciate if someone with matlab knowledge could help me. I'm trying to solve the above mentioned pde in matlab. I have the following 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 velocityN 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))) Same with the initial condition say a constant 10 atm/cm is applied at t=0 thanks |

All times are GMT -4. The time now is 02:33. |