# mesh size & urf SIMPLE FVM MATLAB

 Register Blogs Members List Search Today's Posts Mark Forums Read

 August 14, 2012, 08:52 mesh size & urf SIMPLE FVM MATLAB #1 Member   Michael Moor Join Date: May 2012 Location: Ireland Posts: 30 Rep Power: 13 Hello everyone, I have developed a laminar 2D solver for Poiseulle flow using FVM and SIMPLE, and the results match nicely to the Hagen-Poiseulle velocity profile, but if i increase the mesh size, the solution blows up... My results are coming out nicely, but i cannot seem to be able to increase the mesh size beyond 36*36... I have to decrease the pressure under-relaxation to 0.01 to get it to iterate, but then it just oscillates at a continuity residual of around e-2... do you know the reason for this and why it happens? or how to fix it? I am trying to maybe implement a method such that: if MassResidual <= 10^-2 alphap = 0.05 disp('alphap is now 0.1') end if MassResidual <= 10^-3 alphap = 0.1 disp('alphap is now 0.3') end in an effort to try and march the process along? is this a correct approach? I am also solving the equations using Jacobi: kk=1; ip=0; while kk<=200; pcheck = pdash; %for i=2:nXp-1; %for j=2:nYp-1; %pdash(i,j)=(bdash(i,j)+saE(i,j)*pdash(i+1,j)+saW( i,j)*pdash(i-1,j)+saN(i,j)*pdash(i,j+1)+saS(i,j)*pdash(i,j-1))/saP(i,j); %end %end kk=kk+1 ; mx = max(max(max(abs(pcheck-pdash)))); ip=ip+1; if (kk>50&&mx<10^-6), break, end end and for residuals: MassResidual=0.0; for i=2:nXp-1; for j=2:nYp-1; MassResidual=max([MassResidual abs(bdash(i,j))]); end end MassResidual=abs(MassResidual); I'd appreciate any assistance! Best Regards, Michael 0565113@gmail.com

 August 14, 2012, 10:12 #2 Member   Michael Moor Join Date: May 2012 Location: Ireland Posts: 30 Rep Power: 13 Apologies, I did not add all of the relevant information... This is the momentum equation that i solve for, including under-relaxation: kk=1; while kk<=100; %Solve iteratively using the Jacobi Method to yield u, which will now %be the most up to date value, i.e. u* for i=3:nX-1; for j=2:nYp-1; u(i,j)=((biJ(i,j)+uaE(i,j)*u(i+1,j)+uaW(i,j)*u(i-1,j)+uaN(i,j)*u(i,j+1)+uaS(i,j)*u(i,j-1) + (1-alphau)*uaP(i,j)*u(i,j)/alphau))*alphau/uaP(i,j); end end kk=kk+1; % Counter end and i correct the pressure with under-relaxation: %% Correct the Pressure, including under-relaxation pnew(2:nXp,2:nYp-1) = pstar(2:nXp,2:nYp-1) + alphap*pdash(2:nXp,2:nYp-1);