cfd 空间5阶weno、时间3阶tvd r-k格式.cpp

CFD各种差分格式,ENO/Lax-Friedriches/Beam-Worming/Lax实例

#include<iostream.h>
#include<stdio.h>
#include<math.h>
#define PI 3.14159
void main()
{
	int j,n;
	int G=30;
	double H=0.01;
    double T=0.005;
    double a0,a1,a2;
    double a01,a11,a21;
	double a02,a12,a22;
	double w0,w1,w2;
	double w01,w11,w21;
	double w02,w12,w22;
//	double u[1000][1000],u1[1000][1000],u2[1000][1000];
//	double L[1000][1000],L1[1000][1000],L2[1000][1000];
//	double H0[1000][1000],H1[1000][1000],H2[1000][1000];	
	double  **u=new double* [1000];
	double  **u1=new double* [1000];
	double  **u2=new double* [1000];
	double  **L=new double* [1000];
	double  **L1=new double* [1000];
	double  **L2=new double* [1000];
	double  **H0=new double* [1000];
	double  **H1=new double* [1000];
	double  **H2=new double* [1000];
    for(int i=0;i<1000;i++) 
	{
        u[i]=new  double[1000];   
		u1[i]=new  double[1000]; 
		u2[i]=new  double[1000]; 
		L[i]=new  double[1000]; 
		L1[i]=new  double[1000]; 
		L2[i]=new  double[1000]; 
		H0[i]=new  double[1000]; 
		H1[i]=new  double[1000]; 
		H2[i]=new  double[1000]; 
	}

//  for(j=0;j<501;j++)         u[0][j]=1.0;     //初始条件离散;
//	for(j=501;j<1000;j++)      u[0][j]=0.0;

    for(j=0;j<301;j++)         u[0][j]=0.0;
    for(j=301;j<601;j++)
	{
		u[0][j]=sin(PI*(0.01*j-3));
	}
	for(j=601;j<801;j++)
        u[0][j]=0.75*0.01*j-4.5;   
    for(j=801;j<1000;j++)      u[0][j]=-1.0;

//   for(i=0;i<1000;i++)   
//      delete   []u[i];              
//  delete   []u;   

	for(n=0;n<60;n++)
	{
		for(j=9*n+2;j<998-6*n;j++)
		{
			a0=0.3/(10e-6+13.0/12.0*(u[n][j]-2*u[n][j+1]+u[n][j+2])*(u[n][j]-2*u[n][j+1]+u[n][j+2])+0.25*(3*u[n][j]-4*u[n][j+1]+u[n][j+2])*(3*u[n][j]-4*u[n][j+1]+u[n][j+2]))
				*(10e-6+13.0/12.0*(u[n][j]-2*u[n][j+1]+u[n][j+2])*(u[n][j]-2*u[n][j+1]+u[n][j+2])+0.25*(3*u[n][j]-4*u[n][j+1]+u[n][j+2])*(3*u[n][j]-4*u[n][j+1]+u[n][j+2]));
			a1=0.6/(10e-6+13.0/12.0*(u[n][j-1]-2*u[n][j]+u[n][j+1])*(u[n][j]-2*u[n][j+1]+u[n][j+2])+0.25*(u[n][j-1]-u[n][j+1])*(u[n][j-1]-u[n][j+1]))
				*(10e-6+13.0/12.0*(u[n][j-1]-2*u[n][j]+u[n][j+1])*(u[n][j]-2*u[n][j+1]+u[n][j+2])+0.25*(u[n][j-1]-u[n][j+1])*(u[n][j-1]-u[n][j+1]));  
			a2=0.1/(10e-6+13.0/12.0*(u[n][j-2]-2*u[n][j-1]+u[n][j])*(u[n][j-2]-2*u[n][j-1]+u[n][j])+0.25*(u[n][j-2]-4*u[n][j-1]+3*u[n][j])*(u[n][j-2]-4*u[n][j-1]+3*u[n][j]))
				*(10e-6+13.0/12.0*(u[n][j-2]-2*u[n][j-1]+u[n][j])*(u[n][j-2]-2*u[n][j-1]+u[n][j])+0.25*(u[n][j-2]-4*u[n][j-1]+3*u[n][j])*(u[n][j-2]-4*u[n][j-1]+3*u[n][j]));
			w0=a0/(a0+a1+a2);
			w1=a1/(a0+a1+a2);
			w2=a2/(a0+a1+a2);
			H0[n][j-2]=w0*(1.0/3.0*u[n][j]+5.0/6.0*u[n][j+1]-1.0/6.0*u[n][j+2])+w1*(-1.0/6.0*u[n][j-1]+5.0/6.0*u[n][j]+1.0/3.0*u[n][j+1])
				+w2*(1.0/3.0*u[n][j-2]-6.0/7.0*u[n][j-1]+6.0/11.0*u[n][j]);
		}
		for(j=9*n+3;j<998-6*n;j++)
		{
			L[n][j-3]=-1.0/H*(H0[n][j-2]-H0[n][j-3]);
			u1[n][j]=u[n][j]+T*L[n][j-3];
		}
  
		
		for(j=9*n+5;j<996-6*n;j++)
		{
			a01=0.3/(10e-6+13.0/12.0*(u1[n][j]-2*u1[n][j+1]+u1[n][j+2])*(u1[n][j]-2*u1[n][j+1]+u1[n][j+2])+0.25*(3*u1[n][j]-4*u1[n][j+1]+u1[n][j+2])*(3*u1[n][j]-4*u1[n][j+1]+u1[n][j+2]))
				*(10e-6+13.0/12.0*(u1[n][j]-2*u1[n][j+1]+u1[n][j+2])*(u1[n][j]-2*u1[n][j+1]+u1[n][j+2])+0.25*(3*u1[n][j]-4*u1[n][j+1]+u1[n][j+2])*(3*u1[n][j]-4*u1[n][j+1]+u1[n][j+2]));
			a11=0.6/(10e-6+13.0/12.0*(u1[n][j-1]-2*u1[n][j]+u1[n][j+1])*(u1[n][j]-2*u1[n][j+1]+u1[n][j+2])+0.25*(u1[n][j-1]-u1[n][j+1])*(u1[n][j-1]-u1[n][j+1]))
				*(10e-6+13.0/12.0*(u1[n][j-1]-2*u1[n][j]+u1[n][j+1])*(u1[n][j]-2*u1[n][j+1]+u1[n][j+2])+0.25*(u1[n][j-1]-u1[n][j+1])*(u1[n][j-1]-u1[n][j+1]));  
			a21=0.1/(10e-6+13.0/12.0*(u1[n][j-2]-2*u1[n][j-1]+u1[n][j])*(u1[n][j-2]-2*u1[n][j-1]+u1[n][j])+0.25*(u1[n][j-2]-4*u1[n][j-1]+3*u1[n][j])*(u1[n][j-2]-4*u1[n][j-1]+3*u1[n][j]))
				*(10e-6+13.0/12.0*(u1[n][j-2]-2*u1[n][j-1]+u1[n][j])*(u1[n][j-2]-2*u1[n][j-1]+u1[n][j])+0.25*(u1[n][j-2]-4*u1[n][j-1]+3*u1[n][j])*(u1[n][j-2]-4*u1[n][j-1]+3*u1[n][j]));
			w01=a01/(a01+a11+a21);
			w11=a11/(a01+a11+a21);
			w21=a21/(a01+a11+a21);
			H1[n][j-2]=w01*(1.0/3.0*u1[n][j]+5.0/6.0*u1[n][j+1]-1.0/6.0*u1[n][j+2])+w11*(-1.0/6.0*u1[n][j-1]+5.0/6.0*u1[n][j]+1.0/3.0*u1[n][j+1])
				+w21*(1.0/3.0*u1[n][j-2]-6.0/7.0*u1[n][j-1]+6.0/11.0*u1[n][j]);
		}
		for(j=9*n+6;j<996-6*n;j++)
		{
			L1[n][j-3]=-1.0/H*(H1[n][j-2]-H1[n][j-3]);
			u2[n][j]=0.75*u[n][j]+0.25*u1[n][j]+0.25*T*L1[n][j-3];
		}
    
		
		for(j=9*n+8;j<994-6*n;j++)
		{
			a02=0.3/(10e-6+13.0/12.0*(u2[n][j]-2*u2[n][j+1]+u2[n][j+2])*(u2[n][j]-2*u2[n][j+1]+u2[n][j+2])+0.25*(3*u2[n][j]-4*u2[n][j+1]+u2[n][j+2])*(3*u2[n][j]-4*u2[n][j+1]+u2[n][j+2]))
	    		*(10e-6+13.0/12.0*(u2[n][j]-2*u2[n][j+1]+u2[n][j+2])*(u2[n][j]-2*u2[n][j+1]+u2[n][j+2])+0.25*(3*u2[n][j]-4*u2[n][j+1]+u2[n][j+2])*(3*u2[n][j]-4*u2[n][j+1]+u2[n][j+2]));
			a12=0.6/(10e-6+13.0/12.0*(u2[n][j-1]-2*u2[n][j]+u2[n][j+1])*(u2[n][j]-2*u2[n][j+1]+u2[n][j+2])+0.25*(u2[n][j-1]-u2[n][j+1])*(u2[n][j-1]-u2[n][j+1]))
				*(10e-6+13.0/12.0*(u2[n][j-1]-2*u2[n][j]+u2[n][j+1])*(u2[n][j]-2*u2[n][j+1]+u2[n][j+2])+0.25*(u2[n][j-1]-u2[n][j+1])*(u2[n][j-1]-u2[n][j+1]));  
			a22=0.1/(10e-6+13.0/12.0*(u2[n][j-2]-2*u2[n][j-1]+u2[n][j])*(u2[n][j-2]-2*u2[n][j-1]+u2[n][j])+0.25*(u2[n][j-2]-4*u2[n][j-1]+3*u2[n][j])*(u2[n][j-2]-4*u2[n][j-1]+3*u2[n][j]))
				*(10e-6+13.0/12.0*(u2[n][j-2]-2*u2[n][j-1]+u2[n][j])*(u2[n][j-2]-2*u2[n][j-1]+u2[n][j])+0.25*(u2[n][j-2]-4*u2[n][j-1]+3*u2[n][j])*(u2[n][j-2]-4*u2[n][j-1]+3*u2[n][j]));
			w02=a02/(a02+a12+a22);
			w12=a12/(a02+a12+a22);
			w22=a22/(a02+a12+a22);
			H2[n][j-2]=w02*(1.0/3.0*u2[n][j]+5.0/6.0*u2[n][j+1]-1.0/6.0*u2[n][j+2])+w12*(-1.0/6.0*u2[n][j-1]+5.0/6.0*u2[n][j]+1.0/3.0*u2[n][j+1])
				+w22*(1.0/3.0*u2[n][j-2]-6.0/7.0*u2[n][j-1]+6.0/11.0*u2[n][j]);
		}
		for(j=9*n+9;j<994-6*n;j++)
		{
			L2[n][j-3]=-1.0/H*(H2[n][j-2]-H2[n][j-3]);
			u[n+1][j]=1.0/3.0*u[n][j]+2.0/3.0*u2[n][j]+2.0/3.0*T*L2[n][j-3];
		}
	}

	FILE *fp;
	fp=fopen("d:\\WENO1.plt","w");
	
	for(j=9*G;j<1000-6*G;j++)
	{   
	    fprintf(fp,"%f%c",j*H,' ');
	    fprintf(fp,"%f%c",u[G][j],'\n');
	}
/*	for(j=9*G;j<1000-6*G;j++)
	{   

		fprintf(fp,"%f%c",j*0.01+0.15,' ');
		fprintf(fp,"%f%c",u[0][j],'\n');
	}*/
	fclose(fp);

	for(i=0;i<1000;i++) 
	{
        delete []u[i];   
		delete []u1[i]; 
		delete []u2[i]; 
		delete []L[i]; 
		delete []L1[i]; 
		delete []L2[i]; 
		delete []H0[i]; 
		delete []H1[i]; 
		delete []H2[i]; 
	}
	delete []u;   
	delete []u1; 
	delete []u2; 
	delete []L; 
	delete []L1; 
	delete []L2; 
	delete []H0; 
	delete []H1; 
	delete []H2; 
}