重庆大学高等流体作业流体作业2

-1-研究生课程考核试卷科目：高等流体力学教师：何川姓名：苗闪闪学号：20111002060专业：动力工程及工程热物理类别：学术考生成绩：卷面成绩平时成绩课程综合成绩阅卷评语：阅卷教师（签名）重庆大学研究生院制-2-四、设的气体以的速度以零攻角的速度定常绕流长度为的大平板，试用数值解讨论边界层内的流动规律。解：设以匀速U运动的流体沿切线方向绕过一静止平板，由于流体的粘性作用，板面上的流体流速被降为零，板壁附近极薄区域的流体沿板的法向存在很大的速度梯度，该区域被称为流动边界层。，为大雷诺数绕流。沿板面取x坐标，板的发现方向取y坐标舍边界层的厚度为δ=δ(x)，因厚度很薄，既有。V∞Yx0根据问题的特性，忽略质量力，边界层内的速度矢量,压强，流体力学基本方程具体化为：C.E.(4.1)M.E.(4.2)(4.3)下面，应数量级比较的方法分析方程中各项的数量级关系。-3-首先，选特征参数U、L，奖方程组中个变物无量纲化，并在无量纲化过程中注意使用各无量纲变量的取值范围为(0，1)或取1的数量级（记作(~1)）。即令：注：大雷诺数流动时，压力采用流动造成的动压作为参照具有相当的数量级关系。将各无量纲变量代入式（4.1），并比较其各项的数量级关系，有：将各无量纲变量代入式（4.2），并比较其各项的数量级关系，有：1将各无量纲变量代入式（4.3），并比较其各项的数量级关系，有：1根据物理依据和工程事实，在连续型方程中，认为A=B；在x方向的大量方程中，去掉数量级较低的第一项，保留第二项，即边界层中的粘性作用，认为该项于压力梯度项具有相同的数量级，则有；分析y方向的动量方程，对-4-流项与粘性扩散项具有相同的数量级，而压力梯度项却具有高得多的数量级。得到简化的近似方程组：C.E.(4.4)M.E.(4.5)(4.6)由式（4.6）可以认为，在任一过断流面上，边界层内各点的压力与其外边界上的流势压力pe相等，即，而边界层外势流区满足这里，即有边界层方程即可表示为C.E.(4.4)M.E.(4.7)对于平板绕流，Ve=U=c，则平板然刘的边界层发成可以简化为C.E.(4.4)M.E.(4.8)其定界边界条件为：y=0:u=0,v=0(4.9a)y→∞:u=U(4.9b)-5-下面采用无量纲相似性解法。i引入流函数，将未知函数的数量由两个化为一个令是连续性方程自动满足，动量方程变形为(4.10)下面用但参数群论发证明上述方程所对应的边界层流存在相似性解：令，带入上述方程：即与元方程对比，当时，与原方程有相同的形式，因此，该边界层流动存在相似性解。ii作相似性变换引入无量纲相似性变换则-6-有式（4.10）可变换为常微分方程(4.11)边界条件相应变为：f(0)=0，f’(0)=0，f’(∞)=1(4.12)iii将初边值问题变换为初值问题令f’’(0)=A，这里，A为非零常数。引入新变量则有：则式（4.11）变为(4.13)-7-相应的边界条件变为：f”(0)=A，F”(0)=1定界条件相应变换为：F(0)=0，F’(0)=0，F”(0)=1(4.14)，iv龙格库塔解方程（1）标准四阶龙格-库塔法解微分方程组的原理对于下面的一阶常微分方程组初值问题：,,,,,,,,,dxfxyztdtdygxyztdtdzhxyztdt，其中，当0tt时000000()()()xtxytyztz取步长t，逐步做递进式运算112341123411234122612261226nnnnnnxxKKKKyyLLLLzzmmmm其中123412341234,,,,,,,,,,,KKKKLLLLmmmm分别为：111,,,,,,,,,nnnnnnnnnnnnKtfxyztLtgxyztmthxyzt111211121112,,,2222,,,2222,,,2222nnnnnnnnnnnnKLmtKtfxyztKLmtLtgxyztKLmtmthxyzt-8-222322232223,,,2222,,,2222,,,2222nnnnnnnnnnnnKLmtKtfxyztKLmtLtgxyztKLmtmthxyztV433343334333,,,,,,,,,nnnnnnnnnnnnKtfxKyLzmttLtgxKyLzmttmthxKyLzmtt（2）对三阶非线性常微分方程求解根据龙格-库塔法求解常微分方程的原理，要先将本问题中的原三阶非线性常微分方程化为一阶常微分方程组的形式然后再求解。对由方程0FFF，我们可设11222'''''''FFFFFFFFF(1)我们再设()()'()''xtFytFztF则（1）式可变为dxydtdyzdtdzxzdt，再由初值条件00F，00F，01F可推得t=0时，(0)0x，(0)0y，(0)1z于是对于该问题的一阶常微分方程组dxydtdyzdtdzxzdt，当0tt时，000,,xxyyzz取步长为h则有：-9-1213243*1()21()2()khykhylkhylkhyl，1213243*1()21()2()lhzlhzmlhzmlhzm，1211322433*((*))11(())())2211(()())22(()())mhxzmhxkzmmhxkzmmhxkzm则原方程化简为1123411234112340000226(22)6(22)60,0,0,1nnnnnnhxxKKKKhyyllllhzzmmmmtxyz设置一定的迭代停滞条件，这样逐步迭代就可以算出F，从而得到A。同理得到A后就可以对000,00,0ffffffA进行化简得到：1123411234112343200001226(22)6(22)6'0,0,0,,nnnnnnxxKKKKhyyllllhzzmmmmtxyzAAF进一步用龙格-库塔法迭代求解。（3）编写matlab程序h=0.005;%步长h为0.005F=zeros(3,1000);F(3,1)=1;K=zeros(3,4);forn=1:1000,K(1,1)=h*F(2,n);K(2,1)=h*F(3,n);K(3,1)=h*(-F(1,n)*F(3,n));K(1,2)=h*(F(2,n)+K(2,1)/2);K(2,2)=h*(F(3,n)+K(3,1)/2);-10-K(3,2)=h*(-(F(1,n)+K(1,1)/2)*(F(3,n)+K(3,1)/2));K(1,3)=h*(F(2,n)+K(2,2)/2);K(2,3)=h*(F(3,n)+K(3,2)/2);K(3,3)=h*(-(F(1,n)+K(1,2)/2)*(F(3,n)+K(3,2)/2));K(1,4)=h*(F(2,n)+K(2,3));K(2,4)=h*(F(3,n)+K(3,3));K(3,4)=h*(-(F(1,n)+K(1,3))*(F(3,n)+K(3,3)));F(:,n+1)=F(:,n)+(K(:,1)+2*K(:,2)+2*K(:,3)+K(:,4))/6;end;plot(0:h:1000*h,F(1,:),'r-',0:h:1000*h,F(2,:),'g-',0:h:1000*h,F(3,:),'m-');%画三条曲线在一张图上title('F及其一阶二阶导数的变化曲线');%表头命名为“F及其一阶二阶导数的变化曲线xlabel('变量T');%坐标轴x代表变量Ttext(3.3,3.5,'F的值');%命名Ftext(3.2,1.99,'F的一阶导数');%命名F’text(3.2,0.3,'F的二阶导数');%命名F’’gridon;%开启网格线%求A值A=power(F(2,1000),-3/2)%=0.4696%（2）求f及其一阶导数和二阶导数t=0.005;%步长t为0.005f=zeros(3,1000);f(3,1)=A;K=zeros(3,4);forn=1:1000,K(1,1)=h*f(2,n);K(2,1)=h*f(3,n);K(3,1)=h*(-f(1,n)*f(3,n));K(1,2)=h*(f(2,n)+K(2,1)/2);-11-K(2,2)=h*(f(3,n)+K(3,1)/2);K(3,2)=h*(-(f(1,n)+K(1,1)/2)*(f(3,n)+K(3,1)/2));K(1,3)=h*(f(2,n)+K(2,2)/2);K(2,3)=h*(f(3,n)+K(3,2)/2);K(3,3)=h*(-(f(1,n)+K(1,2)/2)*(f(3,n)+K(3,2)/2));K(1,4)=h*(f(2,n)+K(2,3));K(2,4)=h*(f(3,n)+K(3,3));K(3,4)=h*(-(f(1,n)+K(1,3))*(f(3,n)+K(3,3)));f(:,n+1)=f(:,n)+(K(:,1)+2*K(:,2)+2*K(:,3)+K(:,4))/6;end;plot(0:h:1000*h,f(1,:),'r-',0:h:1000*h,f(2,:),'g-',0:h:1000*h,f(3,:),'m-');%画三条曲线在一张图上title('f及其一阶二阶导数的变化曲线');%表头命名为“f及其一阶二阶导数的变化曲线xlabel('变量T');%坐标轴x代表变量Ttext(3.2,2.5,'f的值');%命名ftext(3.2,1.2,'f的一阶导数');%命名f’text(3.2,0.1,'f的二阶导数');%命名f’’gridon;%开启网格线%（3）求解边界层厚度x=0:0.01:1;h=4.855*sqrt(2*1.1/100000/5.7*x);plot(x,h)gridon%（4）速度变化x=0.3;%分别改为0.85、0.3、0.02T=0:t:1000*t;y=T.*sqrt(2*1.1/100000*x/5.7);u=5.7*f(2,:);-12-subplot(2,1,1);%表示将两个速度分布绘制在一个窗口中plot(y,u);title('x方向速度u的变化曲线');xlabel('y-距板面的距离');ylabel('u-x方向的速度');gridon;v=sqrt(1.1*5.7/100000/2/x)*(T.*f(2,:)-f(1,:));subplot(2,1,2);plot(y,v);title('y方向速度v的变化曲线');xlabel('y-距板面的距离');ylabel('v-y方向的速度');gridon;h=4.855*sqrt(2*1.1/100000/5.7*x)%（5）排挤厚度和动量厚度x=0:0.01:1;h1=sqrt(2*1.1/100000/5.7*x).*(sum(1-f(:,2)))*0.05;%计算排挤厚度subplot(2,1,1);plot(x,h1);title('排挤厚度随x的变化');xlabel('x');ylabel('排挤厚度h1');gridon