CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   Main CFD Forum (http://www.cfd-online.com/Forums/main/)
-   -   SIMPLE algorithm for solving blasius flow (http://www.cfd-online.com/Forums/main/95238-simple-algorithm-solving-blasius-flow.html)

 tanmoy December 10, 2011 14:46

SIMPLE algorithm for solving blasius flow

i've tried solving 2 dimensional N-S equations for blasius flow bt wd sum absurd results. i 2k d equations from d book on cfd by j.d anderson. i m facing difficulty with how to give the boundary conditions,nebody plz help me. d prob is dat,d code runs well bt d reuslts r far away from reality :(
m using SIMPLE algorithm with staggered grid................. plzzzzzzzzzzzzz help me :confused:............. thnx :)

 mazhar December 18, 2011 01:22

There are few typo mistakes in Anderson book. If you can share your code perhaps I can help you.
My Email ID is iqbalmazhar107@yahoo.com

 tanmoy December 22, 2011 05:44

hey,dis is my code. i've written it in C. plz help me.............

#include <stdio.h>
#include <conio.h>
#include<math.h>
void main()
{
/* declaration of variables*/

float U,rho,meu,L,Re,udt,uddt,vdt,vddt,a1,a2,a3,a4,b1,b2 ,b3,b4,P;
float dt,dx,A,B,dy,u[11][22],v[12][23],p[11][21],pdash[11][21],a,b,c,d;
int pr,pc,ur,uc,vr,vc,k; //counters

FILE *f=fopen("output1.txt","w");
FILE *g=fopen("output2.txt","w");
FILE *s=fopen("output3.txt","w");
clrscr();
printf ("\nEnter the free stream velocity: ");
scanf ("%f",&U);
printf ("\nEnter the length of the plate: ");
scanf ("%f",&L);
printf ("\nEnter the pressure at the inlet:");
scanf ("%f",&P);
printf ("\nEnter the value of density of the fluid: ");
scanf ("%f",&rho);
printf ("\nEnter the value of coefficient of viscosity: ");
scanf ("%f",&meu);
Re=(rho*U*L)/meu;
printf ("\nThe reynold's number is: %f",Re);
//initializing the matrices
for (ur=0;ur<=10;ur++)
{
for (uc=0;uc<=21;uc++)
{
u[ur][uc]=0;
}
}
for (vr=0;vr<=11;vr++)
{
for (vc=0;vc<=22;vc++)
{
v[vr][vc]=0;
}
}
for (pr=0;pr<=10;pr++)
{
for (pc=0;pc<=20;pc++)
{
p[pr][pc]=0;
pdash[pr][pc]=0;
}
}
//giving the velocity spike
v[2][4]=(U)/2;
dx=L/20;
dy=0.001378;
dt=.001;
for (k=0;k<=250;k++)
{
for (uc=0;uc<=21;uc++)
{
u[10][uc]=U;
u[0][uc]=0;
}
for (vc=0;vc<=22;vc++)
{
v[0][vc]=-v[1][vc];
v[11][vc]=-v[10][vc];
}
for (ur=0;ur<=10;ur++)
{
u[ur][0]=2*U-u[ur][1];
}

for (vr=0;vr<=11;vr++)
{
v[vr][1]=0;
v[vr][0]=0;
}
for (pr=0;pr<=10;pr++)
{
pdash[pr][0]=0;
pdash[pr][20]=0;
}
for (pc=0;pc<=20;pc++)
{
pdash[0][pc]=0;
pdash[10][pc]=0;
}
//calculating the u-velocity at all the internal nodes
for (ur=1,pr=1,vr=1;ur<=9;ur++,pr++,vr++)
{
for (uc=1,pc=1,vc=1;uc<=20;uc++,pc++,vc++)
{
vdt=(v[vr+1][vc]+v[vr+1][vc+1])/2;
vddt=(v[vr][vc]+v[vr][vc+1])/2;
a1=(rho*(u[ur][uc+1]*u[ur][uc+1]-u[ur][uc-1]*u[ur][uc-1]))/(2*dx);
a2=(rho*(u[ur+1][uc]*vdt-u[ur-1][uc]*vddt))/(2*dy);
a3=(u[ur][uc+1]-2*u[ur][uc]+u[ur][uc-1])/(dx*dx);
a4=(u[ur+1][uc]-2*u[ur][uc]+u[ur-1][uc])/(dy*dy);
A=-(a1+a2)+meu*(a3+a4);
u[ur][uc]=u[ur][uc]+((A*dt)/rho)-(dt/(rho*dx))*(p[pr][pc]-p[pr][pc-1]);
}
}
//calculating the v-velocity at all the internal nodes
for (vr=2,ur=1,pr=1;vr<=11;vr++,ur++,pr++)
{
for (vc=2,uc=1,pc=1;vc<=21;vc++,uc++,pc++)
{
udt=(u[ur-1][uc+1]+u[ur][uc+1])/2;
uddt=(u[ur-1][uc]+u[ur][uc])/2;
b1=(rho*(v[vr-1][vc+1]*udt-v[vr-1][vc-1]*uddt))/(2*dx);
b2=(rho*(v[vr][vc]*v[vr][vc]-v[vr-2][vc]*v[vr-2][vc]))/(2*dy);
b3=(v[vr-1][vc+1]-2*v[vr-1][vc]+v[vr-1][vc-1])/(dx*dx);
b4=(v[vr][vc]-2*v[vr-1][vc]+v[vr-2][vc])/(dy*dy);
B=-(b1+b2)+meu*(b3+b4);
v[vr-1][vc]=v[vr-1][vc]+((B*dt)/rho)-(dt/(rho*dx))*(p[pr][pc]-p[pr-1][pc]);
}
}
for (ur=0;ur<=10;ur++)
{
u[ur][21]=u[ur][20];
}
for (vr=0;vr<=22;vr++)
{
v[vr][22]=v[vr][21];
}
for (pr=1,ur=1,vr=2;pr<=9;pr++,ur++,vr++)
{
for (pc=1,uc=1,vc=2;pc<=19;pc++,uc++,vc++)
{
a=2*((dt)/(dx*dx)+(dt)/(dy*dy));
b=-(dt/(dx*dx));
c=-(dt/(dy*dy));
d=(1/dx)*(rho*(u[ur][uc+1]-u[ur][uc]))+(1/dy)*(rho*(v[vr][vc]-v[vr-1][vc]));
pdash[pr][pc]=-(1/a)*(b*(pdash[pr][pc+1]+pdash[pr][pc-1])+c*(pdash[pr+1][pc]+pdash[pr-1][pc])+d);
}
}
for (pr=0;pr<=10;pr++)
{
for (pc=0;pc<=20;pc++)
{
p[pr][pc]=p[pr][pc]+0.1*pdash[pr][pc];
}

}
}
for (ur=0;ur<=10;ur++)
{
for (uc=0;uc<=21;uc++)
{
printf("\nu[%d][%d]=%f",ur,uc,u[ur][uc]);
fprintf(f,"\nu[%d][%d]=%f",ur,uc,u[ur][uc]);
}
}

getch();
}

 mazhar January 1, 2012 12:07

Sorry for the delay, but I am not familiar with C++. I have written a program for Couttee flow, as discussed by Anderson, in MATLAB. If you are interested I can post it here.

 vishnukmd7 January 3, 2012 11:53

If you are trying to solve the flat plate problem numerically, its a waste of time to try and solve the NS. Solve the simplified and final equation, which is the blasius equation for a flat plate. Its one dimentional( on eta) though non-linear, can be solved easily by runge kutta( any order ) shooting method or traditional (non-linear) difference method.

I would recommend Applied Numerical Analysis by Brian Bradie. The treatment of theory of convergence and order is good and the style of explanation excellent.

There is a slight problem in implementing shooting method, since the initial condition of f(0) is unknown.But a predictor corrector would do the trick.

 tanmoy January 9, 2012 12:19

dear mazhar,plz do it na...... wl be greatful to u for ur help.ty :)

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