# Boundary condition for 2d lid driven cavity using ghost cells

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

 April 7, 2009, 03:12 Boundary condition for 2d lid driven cavity using ghost cells #1 Senior Member   TWB Join Date: Mar 2009 Posts: 119 Rep Power: 9 Hi, I'm trying to use ghost cells to implement boundary condition for the 2d lid driven cavity problem. I'm using staggered grid and I wonder if my implementation is correct because the ans diverged after 20+ time steps. My BC are: u(1,:)=-u(-1,:) where u/v(-1,:) is the ghost cell for left wall v(1,:)=-v(-1,:) p(1,:)=-p(-1,:) u(size_x,:)=-u(size_x+1,:) where u/v(size_x+1,:) is the ghost cell for right wall v(size_x,:)=-v(size_x+1,:) p(size_x,:)=-p(size_x+1,:) u(:,1)=-u(:,-1) where u/v(:,-1) is the ghost cell for bottom wall v(:,1)=-v(:,-1) p(:,1)=-p(:,-1) u(:,size_y)=-u(:,size_y+1)+2. where u/v(:,size_y+1) is the ghost cell for top wall v(:,size_y)=-v(:,size_y+1) p(:,size_y)=-p(:,size_y+1) There is a "+2" becos at the top wall, u=1. At other walls, dp/dn=0, For left/right wall du/dx=0,dv/dx=0 For bottom wall du/dy=0,dv/dy=0 Is my implementation correcto? tks!

 April 9, 2009, 03:03 Pressure condition #2 New Member   Join Date: Mar 2009 Posts: 5 Rep Power: 9 Hallo, velocities look correct, but try no negative value of pressure anywhere e.g. p(1,=p(-1, instead of p(1,=-p(-1,. For upper wall(the moving one) I used simple velocity prescription e.g. u=1, v=0. Good luck. L.

 April 12, 2009, 01:18 #3 Senior Member   TWB Join Date: Mar 2009 Posts: 119 Rep Power: 9 Oh ya, the pressure was a typo error. it's just dp/dn=0. For the upper wall, due to the staggered grid implementation,the wall is actually positioned in between the u(:,size_y) and u(:,size_y+1), hence the implementation is u(y)+u(y+1)/2 = 1.0, where 1.0 is the wall velocity. Also, I'm using the fractional method. Hence, there's the momentum eqn and poisson eqn. I also need to update the velocity at some point in time. Does it mean that I must update the ghost pts constantly?

 April 12, 2009, 14:17 #4 Member   private Join Date: Mar 2009 Posts: 74 Rep Power: 9 Don't you mean [u(y+) + u(y-)]/2 = u(boundary) (1.0 for your problem)?

 April 13, 2009, 10:51 #5 Senior Member   TWB Join Date: Mar 2009 Posts: 119 Rep Power: 9 Ya tt's right otd! Anyway, I managed to get it working now, after correcting a small careless mistake. Tks!

 August 22, 2011, 17:16 #6 New Member   Vignesh Suresh Join Date: Jul 2011 Posts: 10 Rep Power: 7 Hey quarkz....can you please tell me what the mistake was since I am having the same problem. My residuals keep increasing after the 15th iteration. Thanks

 August 22, 2011, 18:04 #7 Senior Member   TWB Join Date: Mar 2009 Posts: 119 Rep Power: 9 Hi illuminati5288, I am sorry I can't remember since it's 2 yrs ago. I am sure you will find your mistake if you looked carefully. Good luck!

 August 23, 2011, 15:46 #8 New Member   Vignesh Suresh Join Date: Jul 2011 Posts: 10 Rep Power: 7 Hey quarkz can you remember if it was due to some mistake in the boundary conditions or did have something to do with the signs? I am just unable to figure out why my code diverges. I have followed the steps as illustrated in some of the papers and in the books. If you dont mind can you please check my code and tell me where I am going wrong. Thanks. #include #include"fstream" #include"conio.h" #include"stdlib.h" #include"math.h" #include"iomanip" using namespace std; double Re=100.0; double nu; int grid_size=32; int counter,i,j; double dx,dy,dt,old,residual_L2,residual,g,pressure_resid ual,u_residual,v_residual,velocity_residual,x,y; double ux,vy,Fxx,Gyy,Hx_xy,Hy_xy; double Lx=1.0; double Ly=1.0; double U=1.0; double w=1.0; double factor=pow(10.0,-9.0); double factor2=pow(10.0,-4.0); double factor3=pow(10.0,-4.0); double diffusion_time,advection_time; double time=0.0; double u[1000][1000],v[1000][1000],p[1000][1000],F[1000][1000],G[1000][1000],Hx[1000][1000],Hy[1000][1000]; double ff,gf,hx,hy; double old_p[1000][1000],old_u[1000][1000],old_v[1000][1000]; double p_ref,q,qh,phi; double nodal_u[1000][1000],nodal_v[1000][1000],nodal_p[1000][1000]; double stream[1000][1000]; int flag=0; void input(); void ic_u(); void ic_v(); void ghost_u(); void ghost_v(); void F_QUICK(); void G_QUICK(); void Hx_QUICK(); void Hy_QUICK(); void time_step(); void update_uv(); void pressure_ic(); void pressure_bc(); void interior_pressure_poisson(); void pressure_correction(); void initial(); void calculation(); void nodal_velocities(); void stream_fucntion(); void output(); ofstream fout ("OUTPUT.dat"); int main() { initial(); calculation(); output(); } void initial() { input(); time_step(); pressure_ic(); ic_u(); ghost_u(); ic_v(); ghost_v(); F_QUICK(); G_QUICK(); Hx_QUICK(); Hy_QUICK(); pressure_bc; } void calculation() { system("CLS"); cout<<"SOLVING PRESSURE POISSON EQUATION"<v_residual) velocity_residual=u_residual; else velocity_residual=v_residual; if(velocity_residual<=factor3) flag=1; } cout<<"."; }while(flag!=1); cout<=0;i--) { q=(u[i][j]+u[i][j+1])/2.0; if(q<0) { phi=(1.0/8.0)*((3.0*u[i][j])+(6.0*u[i][j+1])-(u[i][j+2])); ff=q*phi; } else if(q>0) { phi=(1.0/8.0)*((3.0*u[i][j+1])+(6.0*u[i][j])-(u[i][j-1])); ff=q*phi; } else { ff=0.0; } qh=(u[i][j+1]-u[i][j])/dx; F[i][j]=((ff)-(nu*qh)); } } } void G_QUICK() { for(j=0;j=0;i--) { q=(v[i][j]+v[i+1][j])/2.0; if(q<0) { phi=(1.0/8.0)*((3.0*v[i+1][j])+(6.0*v[i][j])-(v[i-1][j])); gf=q*phi; } else if(q>0) { phi=(1.0/8.0)*((3.0*v[i][j])+(6.0*v[i+1][j])-(v[i+2][j])); gf=q*phi; } else { gf=0.0; } qh=(v[i][j]-v[i+1][j])/dy; G[i][j]=(gf)-(nu*qh); } } } void Hx_QUICK() { for(j=1;j=1;i--) { q=(v[i][j]+v[i][j-1])/2.0; if(q<0) { phi=(1.0/8.0)*((3.0*u[i][j])+(6.0*u[i-1][j])-(u[i-2][j])); hx=q*phi; } else if(q>0) { phi=(1.0/8.0)*((3.0*u[i-1][j])+(6.0*u[i][j])-(u[i+1][j])); hx=q*phi; } else { hx=0.0; } qh=(v[i][j]-v[i][j-1])/dx; Hx[i][j]=(hx)-(nu*qh); } } } void Hy_QUICK() { for(j=1;j=1;i--) { q=(u[i][j]+u[i-1][j])/2.0; if(q<0) { phi=(1.0/8.0)*((3.0*v[i][j-1])+(6.0*v[i][j])-(v[i][j+1])); hy=q*phi; } else if(q>0) { phi=(1.0/8.0)*((3.0*v[i][j])+(6.0*v[i][j-1])-(v[i][j-2])); hy=q*phi; } else { hy=0.0; } qh=(u[i-1][j]-u[i][j])/dy; Hy[i][j]=(hy)-(nu*qh); } } } void interior_pressure_poisson() { for(counter=1;counter<=1000000;counter++) { residual_L2=0.0; for(j=1;j=1;i--) { old=p[i][j]; ux=(u[i][j+1]-u[i][j]); vy=(v[i][j]-v[i+1][j]); Fxx=(F[i][j-1]-(2.0*F[i][j])+F[i][j+1]); Hx_xy=(((Hx[i][j+1]-Hx[i][j]))-((Hx[i+1][j+1]-Hx[i+1][j]))); Hy_xy=(((Hy[i][j+1]-Hy[i][j]))-((Hy[i+1][j+1]-Hy[i+1][j]))); Gyy=(G[i-1][j]-(2.0*G[i][j])+G[i+1][j]); g=((1.0/dt)*((ux*dx)+(vy*dy)))-(Fxx+Hx_xy+Hy_xy+Gyy); residual=w*(((1.0/5.0)*(p[i][j-1]+p[i][j+1]+p[i-1][j]+p[i+1][j]))+((1.0/20.0)*(p[i-1][j-1]+p[i+1][j-1]+p[i-1][j+1]+p[i+1][j+1]))-(old)-((6.0/20.0)*g*dx*dy)); p[i][j]=(old+residual); residual_L2+=pow(residual,2.0); } } residual_L2=pow((1.0/(grid_size*grid_size))*residual_L2,0.5); if(residual_L2<=factor) { break; } } } void pressure_correction() { p_ref=p[grid_size][1]; for(i=1;i=1;i--) { u[i][j]=(u[i][j]-((dt/dx)*((F[i][j]-F[i][j-1])+(p[i][j]-p[i][j-1])+(Hx[i][j]-Hx[i+1][j])))); } } for(j=1;j=2;i--) { v[i][j]=(v[i][j]-((dt/dy)*((G[i-1][j]-G[i][j])+(p[i-1][j]-p[i][j])+(Hy[i][j+1]-Hy[i][j])))); } } ghost_u(); ghost_v(); F_QUICK(); G_QUICK(); Hx_QUICK(); Hy_QUICK(); } void nodal_velocities() { for(i=1;i=1;i--) { for(j=1;j=1;i--) { for(j=2;j

 August 26, 2011, 14:35 #9 New Member   Vignesh Suresh Join Date: Jul 2011 Posts: 10 Rep Power: 7 hey quarkz, I have kind of narrowed down the cause for divergence. When a moving lid is horizontal then the divergence occurs in the x-velocity but no divergence in the y-velocity. The vice-versa occurs when the moving lid is vertical. Any ideas why this is happening. Thanks

 January 20, 2013, 06:54 Reynolds Effects #10 New Member   Aryan Join Date: Jan 2013 Posts: 8 Rep Power: 5 Hi Vignesh It might be occured because of the high Re regime. I suggest check your program for the Re = 1 Good Luck

 Thread Tools Display Modes Linear Mode

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts BB code is On Smilies are On [IMG] code is On HTML code is OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post murali CFX 5 August 3, 2012 08:56 sudha Main CFD Forum 6 November 21, 2010 06:11 Adrian Main CFD Forum 0 October 24, 2008 18:41 student Main CFD Forum 1 July 20, 2008 13:50 Karl FLUENT 7 May 11, 2002 22:12

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