北航计算流体力学程序

#include

#include

# define N 10

void main()

{

float

a0[N][2*N+1],b0[2*N+1][ 3*N+1],a1[N][ 2*N+1],b1[2*N+1][ 3*N+1],a[N][ 2*N+ 1],b[2*N+1][ 3*N+1];

float c[2*N],d[2*N];

int i,j;

float m=1,n=1;

for(i=0;i<=N-1;i++)

for(j=0;j<=2*N-1;j++)

a0[i][j]=float(i)/N;

for(i=0;i<=N-1;i++)

a0[i][2*N]=0;

for(i=1;i<=2*N;i++)

for(j=0;j<=3*N;j++)

b0[i][j]=float((i+N))/N;

for(j=0;j<=2*N-1;j++)

b0[0][j]=1;

for(j=2*N;j<=3*N;j++)

b0[0][j]=0;

b0[2*N][ 3*N]=3;

while(m>0.00001&&n>0.00001)

{

for(i=1;i<=N-2;i++)

for(j=1;j<=2*N-1;j++)

a1[i][j]=(a0[i-1][j]+a0[i+1][j]+a0[i][j-1]+a0[i][j+1])/4;

for(j=1;j<=2*N-1;j++)

a1[N-1][j]=(b0[0][j]+a0[N-2][j]+a0[N-1][j-1]+a0[N-1][j+1])/4; for(i=1;i<=N-1;i++)

a1[i][0]=float(i)/N;

for(i=1;i<=N-1;i++)

a1[i][2*N]=0;

for(j=0;j<=2*N;j++)

a1[0][j]=0;

for(i=1;i<=2*N;i++)

for(j=1;j<=3*N-1;j++)

b1[i][j]=(b0[i-1][j]+b0[i+1][j]+b0[i][j-1]+b0[i][j+1])/4;

for(j=1;j<=2*N-1;j++)

b1[0][j]=(b0[0][j-1]+b0[0][j+1]+a0[N-1][j]+b0[1][j])/4;

for(i=0;i<=2*N;i++)

b1[i][0]=float((i+N))/N;

for(j=2*N;j<=3*N;j++)

b1[0][j]=0;

for(i=1;i<=2*N-1;i++)

b1[i][3*N]=b1[i][3*N-1];

for(j=0;j<=3*N;j++)

b1[2*N][j]=3;

for(i=0;i<=N-1;i++)

for(j=0;j<=2*N;j++)

a[i][j]=fabs(a1[i][j]-a0[i][j]);

for(i=0;i<=2*N;i++)

for(j=0;j<=3*N;j++)

b[i][j]=fabs(b1[i][j]-b0[i][j]);

m=a[0][0];

for(i=0;i<=N-1;i++)

for(j=0;j<=2*N;j++)

if(a[i][j]>m)

m=a[i][j];

n=b[0][0];

for(i=0;i<=2*N;i++)

for(j=0;j<=3*N;j++)

if(b[i][j]>n)

n=b[i][j];

for(i=0;i<=N-1;i++)

for(j=0;j<=2*N;j++)

a0[i][j]=a1[i][j];

for(i=0;i<=2*N;i++)

for(j=0;j<=3*N;j++)

b0[i][j]=b1[i][j];

}

for(i=1;i<=2*N-1;i++)

c[i]=N*(b0[i+1][3*N]-b0[i-1][3*N])/2;

c[0]=(b0[1][3*N]-b0[0][3*N])*N;

c[2*N]=(b0[2*N][3*N]-b0[2*N-1][3*N])*N; for(i=0;i<=N-1;i++)

{for(j=0;j<=2*N;j++)

printf("%f ",a0[i][j]); printf("\n");}

for(i=0;i<=2*N;i++)

{ for(j=0;j<=3*N;j++)

printf("%f ",b0[i][j]); printf("\n");}

for(i=0;i<=2*N;i++)

{d[i]=1+float(i)/float(N);

printf("%f,%f",c[i],d[i]);

printf("\n");}

}