重庆大学高等流体作业流体计算

1、1.齐次方程的变量代换将上述定解问题中的高阶常微分方程表示为一阶常微分方程组：即有：因此可令：相应边界条件为：。将上述常微分方程的边值问题变换为初值问题令：则有：代入原方程有：由于A为非零常数，则有：。初始条件：即转化的初值问题为： 这样原问题转化为：为求得值，须计算%四阶龙格库塔法求解微分方程matlab%设f(0)=1 %分别用x，y，z表示（f=x f=y f=z）方程f+0.5ff=0可以分别记为：% f=x=y; % f=y=z; % f=y=z=-0.5xz;x=0;y=0;z=1;h=0.1;n=100;t=0;U=22 szj=t,x,y,z; for i=1:n-1 % i=

2、1:n定义龙格库塔程序 k1=y; l1=z; m1=-x*z; k2=y+l1*h/2; l2=z+m1*h/2; m2=-(x+k1*h/2)*(z+m1*h/2); k3=y+l2*h/2; l3=z+m2*h/2; m3=-(x+k2*h/2)*(z+m2*h/2); k4=y+l3*h; l4=z+m3*h; m4=-(x+k3*h)*(z+m3*h); x=x+h*(k1+2*k2+2*k3+k4)/6; y=y+h*(l1+2*l2+2*l3+l4)/6; z=z+h*(m1+2*m2+2*m3+m4)/6; t=t+h; szj=szj;t,x,y,z; %f=szj(:,2)

3、 f=szj(:,3) f=szj(:,4) endplot(szj(:,1),szj(:,2), k-+)%fhold onplot(szj(:,1),szj(:,3), b-*)%fhold onplot(szj(:,1),szj(:,4), r-o)%fxlabel(自变量n);legend(f因变量,f一阶导函数,f二阶导函数);title(四阶龙格库塔法求解方程数值解);2.打靶修正法：#include#includemain()double k1,k2,k3,k4,l1,l2,l3,l4,m1,m2,m3,m4;double a,b=0,u,v,m,n,y0,y1,x=0,y=0,

4、z=1,h=0.1,h1=0; /定义变量和步长，赋初值/do y0=y; k1=y; l1=z; m1=-0.5*x*z; k2=y+l1*h/2; l2=z+m1*h/2; m2=-0.5*(x+k1*h/2)*(z+m1*h/2); k3=y+l2*h/2; l3=z+m2*h/2; m3=-0.5*(x+k2*h/2)*(z+m2*h/2); k4=y+l3*h; l4=z+m3*h; m4=-0.5*(x+k3*h)*(z+m3*h); x=x+h*(k1+2*k2+2*k3+k4)/6;/打靶法求解初边值问题（四阶龙格库塔法）/ y=y+h*(l1+2*l2+2*l3+l4)/6;

5、 z=z+h*(m1+2*m2+2*m3+m4)/6; printf(%fn, y);/f=y/ while(fabs(y0-y)/y)0.); a=pow(y,-1.5); printf(%fn, a); /f(0)=A 求A= / x=0; y=0; z=a; m=0.*22/2/0.02; n=2*0.*0.02/22; /对x= 0.02处/ do y0=y; u=22*y;/水平方向的速度u/ v= (h1*y-x)*pow (m, 0.5); /垂直方向的速度v/ y1=h1*pow(n,0.5);/边界层厚度/ printf (%f %f %f %fn, h1, y1, u, v

6、); l1=z; k1=y; m1=-0.5*x*z; l2=z+m1*h/2; k2=y+l1*h/2; m2=-0.5*(x+k1*h/2)*(z+m1*h/2); l3=z+m2*h/2; k3=y+l2*h/2; m3=-0.5*(x+k2*h/2)*(z+m2*h/2); l4=z+m3*h; k4=y+l3*h; m4=-0.5*(x+k3*h)*(z+m3*h); x=x+h*(k1+2*k2+2*k3+k4)/6; z=z+h*(m1+2*m2+2*m3+m4)/6; y=y+h*(l1+2*l2+2*l3+l4)/6; /四阶龙格库塔法求x=0.01m处的u, v/ h1=h

7、1+h; while (fabs (y-y0)/y)0.); x=0;do n=2*0.*x/22; /求边界层厚度/ h=h1*pow (n, 0.5); printf (%f %fn, x, h); x=x+0.01; while(x=1.0); 3流场分析平板绕流边界层无量纲速度分布%（2）流场分析 figure%x方向速度a=0.10 0.20 0.30 0.40 0.50 0.60 0.70 0.80 0.90 1.00 1.10 1.20 1.30 1.40 1.50 1.60 1.70 1.80 1.90 2.00 2.10 2.20 2.30 2.40 2.50 2.60 2.

8、70 2.80 2.90 3.00 3.10 3.20 3.30 3.40 3.50 3.60 3.70 3.80 3.90 4.00 4.10 4.20 4.30 4.40 4.50 4.60 4.70 4.80 4.90 5.00 5.10 5.20 5.30 5.40 5.50 5.60 5.70 5.80 5.90 6.00 6.10 6.20 6.30 6.40 6.50 6.60 6.70 6.80 6.90 7.00 7.10 7.20 7.30 7.40 7.50 7.60 7.70 7.80 7.90 8.00 8.10 8.20 8.30 8.40 8.50 8.60 8.

9、70 8.80 ;b=0.73 1.46 2.19 2.92 3.65 4.38 5.10 5.82 6.54 7.26 7.96 8.66 9.36 10.04 10.71 11.37 12.01 12.64 13.26 13.85 14.43 14.99 15.52 16.04 16.53 16.99 17.44 17.85 18.25 18.61 18.96 19.27 19.57 19.84 20.09 20.31 20.52 20.70 20.87 21.02 21.15 21.27 21.38 21.47 21.55 21.62 21.68 21.73 21.78 21.81 21

10、.85 21.87 21.90 21.92 21.93 21.94 21.96 21.96 21.97 21.98 21.98 21.99 21.99 21.99 21.99 21.99 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00 22.00; plot(a,b, *b-)%fxlabel(自变量n);ylabel(速度u(m/s);title(x方向的速度);%流场分析，y方向速度fig

11、ure %y方向速度vc=0.0.0.0.0.001020.0.0.0.0.0.0.005770.0.0.0.0.0.0.0.0.0.0.0190.0.021690.023020.0.0.0.0.0.0.031430.0.0.0.0.0.0.0.03790.0.0.0.039790.0.0.0.0.0.0.0.0.0.041870.0.042020.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.042320.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.00150.0.0.0.0.003890.0.0.0.006460.0.0.0.0.0

12、.010970.0.0.013290.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.024260.0.0.0.0.0.024390.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.002580.0.0.003960.0.0.0.0.0.0.0.0.00910.00970.0.010880.0.0.0.0.0.0.0.0.0.0.0.0.016680.016950.0.0.0.0.0.0.0.0.0.01850.018570.0

13、18630.0.0.0.0.0.0.0.0.0.0.01890.0.0.0.0.018920.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.0.;plot(a,c(1,:), k-+)hold onplot(a,c(2,:),b-*)hold onplot(a,c(3,:),r-o)xlabel(自变量n);ylabel(v(m/s);title(速度v);legend(x=0.2,x=0.6,x=1);figure %流场分析，边界层厚度l=0.010.020.030.040.050.060.070.080.090.10.110.120.130.140.150.160.

14、170.180.190.20.210.220.230.240.250.260.270.280.290.30.310.320.330.340.350.360.370.380.390.40.410.420.430.440.450.460.470.480.490.50.510.520.530.540.550.560.570.580.590.60.610.620.630.640.650.660.670.680.690.70.710.720.730.740.750.760.770.780.790.80.810.820.830.840.850.860.870.880.890.90.910.920.930.940.950.960.970.980.99;r=0.000890.0.0.001780.001990.002180.0.0.002670.0.0.0.0.003330.0.003560.003670.0.0.003980.0.0.0.004360.004450.0.0.0.0.0.0.0.0.005190.0.005340.0.0.0.