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

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 distribution 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 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,double b, double c ,double f,double *p) /矩阵形式/b0c0 0 0 /=f0/a0b1 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)/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+)*(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 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(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=-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(serial_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); 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/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*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.);else if(serial_number1=2)/一阶迎风方法if(serial_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,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)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;iGn2+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 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; 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=-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=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);elseprintf(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; aE=0.03125*D; m=2*D*0.3125+F; n=2*D*0.3125; for(i=0;i4;i+) ai=0; ci=0; for(i=0;i5;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,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=0.03125*D+F;aE=0.03125*D;m=2*D*0.3125+F;n=2*D*0.3125;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;elsebi=-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/Gn1;F=P*u;D=T/min_x;aW=0.03125*D+F;aE=0.03125*D;m=2*D*0.3125+F;n=2*D*0.3125;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);elseprintf(No case is selected,please check out the input serial number.n);else printf(ERROR:no case or method is selected,please check out the input serial number.n);(1) With the central differencing scheme:Results: u=0.1m/s with 5 grid nodes: u=2.5m/s with 5 grid nodes: u=2.5m/s with 20 grid nodes:Figures: u=0.1m/s with 5 grid nodes:图表 1 u=2.5m/s with 5 grid nodes:图表 2 u=2.5m/s with 20 grid nodes:图表 3Discussion:A small number of grid nodes can be competent to simulate the property of fluid while the flow velocity is low. As the flow rate increases, the number of nodes needs to be a corresponding increase, otherwise the calculation results will be distorted by serious error.The central differencing scheme is not applicable to large and complex problem.(2) With the upwind differencing scheme:Results: u=0.1m/s with 5 grid nodes: u=2.5m/s with 5 grid nodes: u=2.5m/s with 20 grid nodes:Figures: u=0.1m/s with 5 grid nodes:图表 4 u=2.5m/s with 5 grid nodes:图表 5 u=2.5m/s with 20 grid nodes:图表 6Discussion:The upwind differencing scheme can be used to simulate the one-dimensional flow characteristics because of the calculation results and the theoretical analysis values are consistent. (3) With the hybrid differencing scheme:Results: u=0.1m/s with 5 grid nodes: u=2.5m/s with 5 grid nodes: u=2.5m/s with 20 grid nodes:Figures: u=0.1m/s with 5 grid nodes:图表 7 u=2.5m/s with 5 grid nodes:图表 8 u=2.5m/s with 20 grid nodes:图表 9Discussion:The hybrid differencing scheme has different calculation accuracy to deal with different problems. When flow rate increases, this method has certain advantages. (4) With the p