一维对流扩散问题求解.docx

1、一维对流扩散问题求解Subjects：A property is transported by means of convection and diffusion through the one-dimensional domain. The governing equation is; boundary conditions are at x=0 and at x=L. Using five equally spaced cells and the central differencing scheme for convection and diffusion calculate the d

2、istribution of as a function of x for(i)Case 1: u=0.1 m/s, (ii) Case2: Case 1: u=2.5 m/s, and compare the results with the analytical solution.(iii)Case 3: recalculate the solution for u=2.5m/s with 20 grid nodes and compare the results with the analytical solution, the following date apply: length

3、L=1.0 m,.Solution:Programming language:File name:二维稳态导热模拟.cppCompile software: Microsoft Visual C+ 6.0Code:#include #include#includeusing namespace std;#define Gn1 5 /网格节点数.#define Gn2 20 /网格节点数.#define P 1 /流体密度，单位：kg/m/m/m.#define ProA 1 /A端传递特性.#define ProB 0 /B端传递特性./FDD/追赶法过程void TDMA1(double a

4、,double b, double c ,double f,double *p) /矩阵形式 /b0 c0 0 0 /=f0 /a0 b1 c1 0 /=f1 /. . . . . /0 0 . b(N-1) c(N-2) /=f(N-1) double dGn1-1,uGn1,llGn1-1, yGn1,XGn1; int i; for(i=0; i=Gn1-2;i+) di = ci; u0 = b0; for(i=1; i=Gn1-1;i+) lli-1=ai-1/ui-1; ui=bi-lli-1*ci-1; y0=f0; for(i=1; i=0;i-) Xi=(yi-ci*Xi+1

5、)/ui; for (i=0; i=Gn1-1;i+) *(p+i)=Xi; void TDMA2(double a,double b, double c ,double f,double *p) double dGn2-1,uGn2,llGn2-1, yGn2,XGn2; int i; for(i=0;i=Gn2-2;i+) di=ci; u0 = b0; for(i=1;i=Gn2-1;i+) lli-1=ai-1/ui-1; ui=bi-lli-1*ci-1; y0 = f0; for(i=1;i=0;i-) Xi=(yi-ci*Xi+1)/ui; for (i=0;i=Gn2-1;i+

6、) *(p+i)=Xi; void main() double u,L=1,T=0.1,min_x,F,D,aW,aE,m,n; double aGn2,bGn2,cGn2,fGn2,XGn2; int serial_number1, serial_number2,i; printf(There are four methods to solve these cases:n); printf(1.The central differencing scheme.n); printf(2.The upwind differencing scheme.n); printf(3.The hybrid

7、differencing scheme.n); printf(4.The power-law scheme.n); printf(Please input the serial number of the method which you want to use:); scanf(%d,&serial_number1); printf(There are three cases in this subject:n); printf(1. u=0.1m/s with 5 grid nodes.n); printf(2. u=2.5m/s with 5 grid nodes.n); printf(

8、3. u=2.5m/s with 20 grid nodes.n); printf(Please input the serial number of the case which you want to solve:); scanf(%d,&serial_number2); if (serial_number1=1)/中心差分法 if (serial_number2=1)/第一题开始 u=0.1; min_x=L/Gn1; F=P*u; D=T/min_x; aW=D+F/2; aE=D-F/2; m=2*D+F; n=2*D-F; for(i=0;iGn1-1;i+) ai=-aW; ci

9、=-aE; for(i=0;iGn1;i+) if(i=0) bi=-ci+m; fi=m*ProA; else if(i=Gn1-1) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; void TDMA1(double a,double b, double c ,double f,double *p); TDMA1(a,b,c,f,X); printf(Result：X); printf(%d,Gn1); printf(=n); for (i=0;i=Gn1-1;i+) printf( ); printf(%lfn,Xi); else if(se

10、rial_number2=2)/第二题开始 u=2.5; min_x=L/Gn1; F=P*u; D=T/min_x; aW=D+F/2; aE=D-F/2; m=2*D+F; n=2*D-F; for(i=0;iGn1-1;i+) ai=-aW, ci=-aE; for(i=0;iGn1;i+) if(i=0) bi=-ci+m; fi=m*ProA; else if(i=Gn1-1) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; TDMA1(a,b,c,f,X); printf(Result：X); printf(%d,Gn1); print

11、f(=n); for (i=0;i=Gn1-1;i+) printf( ); printf(%lfn,Xi); else if(serial_number2=3)/第三题开始 u=2.5; min_x=L/Gn2; F=P*u; D=T/min_x; aW=D+F/2; aE=D-F/2; m=2*D+F; n=2*D-F; for(i=0;iGn2-1;i+) ai=0; ci=0; for(i=0;iGn2;i+) bi=0; fi=0; for(i=0;iGn2-1;i+) ai=-aW; ci=-aE; for(i=0;iGn2;i+) if(i=0) bi=-ci+m; fi=m*P

12、roA; else if(i=Gn2-1) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; TDMA2(a,b,c,f,X); printf(Result：X); printf(%d,Gn2); printf(=n); for (i=0;i=Gn2-1;i+) printf( ); printf(%lfn,Xi); else printf(NO case is selected,please check out the input serial number.); else if(serial_number1=2)/一阶迎风方法 if(serial

13、_number2=1)/第一题开始 u=0.1; min_x=L/Gn1; F=P*u; D=T/min_x; aW=D+F; aE=D; m=2*D+F; n=2*D; for(i=0;iGn1-1;i+) ai=0; ci=0; for(i=0;iGn1;i+) bi=0; fi=0; for(i=0;iGn1-1;i+) ai=-aW; ci=-aE; for(i=0;iGn1;i+) if(i=0) bi=-ci+m; fi=m*ProA; else if(i=Gn1-1) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; TDMA1(a,b

14、,c,f,X); printf(Result：X); printf(%d,Gn1); printf(=n); for (i=0;i=Gn1-1;i+) printf( ); printf(%lfn,Xi); else if(serial_number2=2)/第二题开始 u=2.5; min_x=L/Gn1; F=P*u; D=T/min_x; aW=D+F; aE=D; m=2*D+F; n=2*D; for(i=0;iGn1-1;i+) ai=-aW; ci=-aE; for(i=0;iGn1;i+) if(i=0) bi=-ci+m; fi=m*ProA; else if(i=Gn1-1

15、) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; TDMA1(a,b,c,f,X); printf(Result：X); printf(%d,Gn1); printf(=n); for (i=0;i=Gn1-1;i+) printf( ); printf(%lfn,Xi); else if(serial_number2=3)/第三题开始 u=2.5; min_x=L/Gn2; F=P*u; D=T/min_x; aW=D+F; aE=D; m=2*D+F; n=2*D; for(i=0;iGn2;i+) ai=0; ci=0; for(i=0;i

16、Gn2+1;i+) bi=0; fi=0; for(i=0;iGn2-1;i+) ai=-aW; ci=-aE; for(i=0;iGn2;i+) if(i=0) bi=-ci+m; fi=m*ProA; else if(i=Gn2-1) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; TDMA2(a,b,c,f,X); printf(Result：X); printf(%d,Gn2); printf(=n); for (i=0;i=Gn2-1;i+) printf( ); printf(%lfn,Xi); else printf(No case

17、is selected,please check out the input serial number.n); else if(serial_number1=3)/混合方法 if(serial_number2=1)/第一题开始 u=0.1; min_x=L/Gn1; F=P*u; D=T/min_x; aW=F; aE=0; m=2*D+F; n=2*D; for(i=0;iGn1-1;i+) ai=-aW; ci=-aE; for(i=0;iGn1;i+) if(i=0) bi=-ci+m; fi=m*ProA; else if(i=Gn1-1) bi=-ai-1+n; fi=n*ProB

18、; else bi=-ci-ai-1; fi=0; TDMA1(a,b,c,f,X); printf(Result：X); printf(%d,Gn1); printf(=n); for (i=0;i=Gn1-1;i+) printf( ); printf(%lfn,Xi); else if(serial_number2=2)/第二题开始 u=2.5; min_x=L/Gn1; F=P*u; D=T/min_x; aW=F; aE=0; m=2*D+F; n=2*D; for(i=0;iGn1-1;i+) ai=-aW; ci=-aE; for(i=0;iGn1;i+) if(i=0) bi=

19、-ci+m; fi=m*ProA; else if(i=Gn1-1) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; TDMA1(a,b,c,f,X); printf(Result：X); printf(%d,Gn1); printf(=n); for (i=0;i=Gn1-1;i+) printf( ); printf(%lfn,Xi); else if (serial_number2=3)/第三题开始 u=2.5; min_x=L/Gn2; F=P*u; D=T/min_x; aW=F; aE=0; m=2*D+F; n=2*D; for(i=

20、0;iGn2-1;i+) ai=-aW; ci=-aE; for(i=0;iGn2;i+) if(i=0) bi=-ci+m; fi=m*ProA; else if(i=Gn2-1) bi=-ai-1+n; fi=n*ProB; else bi=-ci-ai-1; fi=0; TDMA2(a,b,c,f,X); printf(Result：X); printf(%d,Gn2); printf(=n); for (i=0;i=Gn2-1;i+) printf( ); printf(%lfn,Xi); else printf(No case is selected,please check out the input serial number.n); else if(serial_number1=4)/指数方法 if(serial_number2=1)/第一题开始 u=0.1; min_x=L/Gn1; F=P*u; D=T/min_x; aW=0.03125*D+F;