粘性流体力学—平面驻点流动(西门茨流动)

1、平面驻点流动（西门茨流动）平面驻点流动（西门茨流动）指导教师：指导教师：报报 告告 人人 ： 0 01 1问题概述问题概述0202方程推导与求解方程推导与求解目录目录0404工程实际应用工程实际应用0303利用利用Matlab求数值解求数值解ayVaxU平面驻点流动，在平面有势流动中平面驻点流动，在平面有势流动中x、y两个方向的流速分别为：两个方向的流速分别为：第一节第一节 问题概述问题概述（1-1）（1-2）于是由伯努利方程，可得压强分布为于是由伯努利方程，可得压强分布为式中，式中，p0为驻点为驻点O的压强，的压强， . 以上以上3个式子的流速与压强均满足势流方程，并且是不可压缩粘个式子的流

2、速与压强均满足势流方程，并且是不可压缩粘性流动运动方程的精确解。但是不可压缩粘性流动的运动方程中多一性流动运动方程的精确解。但是不可压缩粘性流动的运动方程中多一粘性项粘性项 。对于势流。对于势流 ， 为流速的势函数。则为流速的势函数。则u2u0)()(222u 则则N-S方程中的粘性项对于势流而言恒等于零。但方程中的粘性项对于势流而言恒等于零。但势流解却不能势流解却不能满足满足“无滑移无滑移”这个粘性流动的边界条件这个粘性流动的边界条件.（1-3）)(21220yFxapp2)(yyF 对于平面驻点流动，为了满足粘性流动的无滑移的条件，西门茨对于平面驻点流动，为了满足粘性流动的无滑移的条件，西

3、门茨给出其精确解给出其精确解. 现假定：现假定： 式中式中p0表示驻点表示驻点O（x=0，y=0）处的压强；）处的压强；p为任意点（为任意点（x，y）处的压强；处的压强；a为常数。为常数。 可证明，该假定将自动满足连续性方程。可证明，该假定将自动满足连续性方程。第二节第二节 方程推导和求解方程推导和求解 由平面运动的由平面运动的N-S方程可以确定方程可以确定f与与F两个函数，将上式带入恒定两个函数，将上式带入恒定的的N-S方程方程)1)122222222yvxvypyvvxvuyuxuxpyuvxuu可得：可得：fFaf ffaf ff 22221（2-1）（2-2） 边界条件边界条件: 由以

4、上三个边界条件可得由以上三个边界条件可得: 首先解首先解 ，作变量置换，令，作变量置换，令 3233222)()(AyAyffAyAyffAyAyff)()(,Ayfy可得：可得：带入式可得：带入式可得： 32222)(AaAfaf ff 22（2-3）（2-4）（2-5） 如果式如果式2-5中中aAaaA/3222AaA式中，式中，a为由势流解得到的常数，为由势流解得到的常数，为流体的运动粘度，为已知量。为流体的运动粘度，为已知量。012 ，方程将大大简化，则需要满足，方程将大大简化，则需要满足由此方程可改写为由此方程可改写为)()(,ayfya则方程可简化为则方程可简化为边界条件为边界条件

5、为0000（2-6）（2-7）（2-8）)()(,Ayfy方程方程012 仍然是非线性的，难以求得解析解。仍然是非线性的，难以求得解析解。可采用数值解法，得下表：可采用数值解法，得下表：0 00.40.40.80.81.21.21.61.62 22.42.42.82.83.23.23.63.60 00.08810.08810.31240.31240.6220.6220.97980.97981.3621.3621.75531.75532.1532.1532.55232.55232.9542.9540 00.41450.41450.68590.68590.84670.84670.93230.932

6、30.97320.97320.99050.99050.9970.9970.99920.99920.9990.9991.23261.23260.84630.84630.52510.52510.29380.29380.14740.14740.06580.06580.02650.02650.0090.0090.00280.00280 0yaUudd22dd 则则)(UaaUyAaUu)(Uu 在在=2.4左右，左右，? = 0.99，即此时，即此时粘性流动的流速已接近势流流速，只差粘性流动的流速已接近势流流速，只差百分之一。可以此点距固体壁面的距离百分之一。可以此点距固体壁面的距离作为边界层的厚度作

7、为边界层的厚度，则，则aa4 . 2弗勒塞林求解的平面驻点流动和轴对称驻点流动解：弗勒塞林求解的平面驻点流动和轴对称驻点流动解：（2-9）（2-10）解式解式 ，可得压强，可得压强 ，对其积分可得，对其积分可得fFaf f 221p)2(1)(22ffayF联立式联立式2-2、2-7可得：可得：aaaayp )( ,在边界层内在边界层内 都只是都只是1的数量级，的数量级，因而沿壁面法线的压强梯度因而沿壁面法线的压强梯度aayp当当 很小时，压强梯度也很小。很小时，压强梯度也很小。yaapp0此外此外表明流动过程中压力逐渐增大至表明流动过程中压力逐渐增大至0p（2-11）（2-12）012 ,1

8、0000，ya平面驻点流动方程：平面驻点流动方程：边界条件：边界条件： 其中其中:在在0,4的范围内的范围内的数值解。的数值解。 求解求解Uudd第三节第三节 利用利用MATLAB求数值解求数值解1 1、问题描述：、问题描述：第一步：将方程化为一阶常微分方程组。第一步：将方程化为一阶常微分方程组。1)2()3() 1 ()3()3()2()2() 1 () 1 (2fffdffdffdff边界条件：边界条件： 1)2(0)2(00) 1 (0fff，2 2、求解过程、求解过程第二步：建立第二步：建立ode.m和和lbc.m两个两个M文件。文件。 程序说明：程序说明： ode.m文件描述一阶常微

9、分方程组，由于方程是三阶微分方程，所以需要三个一阶微分方程来描述。 lbc.m文件描述边界条件，f0表示初始值，finf表示末端值，本问题告诉我们的是f(1)和f(2)的初始值与f(2)的末端值，其中finf(2)-1表示finf(2)-1=0，其它以此类推。ode.m文件程序如下： function dfdx=ode(x, f) dfdx=f(2);f(3);-f(1)*f(3)+f(2)2-1; lbc.m 文件程序如下： function res=lbc(f0,finf) res=f0(1);f0(2);finf(2)-1; 第三步：求解方程。第三步：求解方程。 程序说明：程序说明： (

10、1)solinit是被指定为x和f域的范围。x是初始网格点，f表示在节点solinit.x(i)处f(x(i)的初始值猜测解solinit.f(:,1)，一般用bvpinit实现。 (2)bvp4c是MATLAB 7.0软件求解一阶常微分方程组的库函数，调用格式如程序所述(bvp4c的调用格式有三种，程序中的只是其中一种)。但由于sol不能直接输出数值解，所以要用bvp4c的配置函数deval。在MATLAB 7.0工作窗口输入程序：infinity=4;solinit=bvpinit(0:0.4:infinity,0 0 0);sol=bvp4c(ode,lbc,solinit);x=0:0

11、.4:infinityf=deval(sol,x)plot(x,f(1,:),ob,x,f(2,:),rp,x,f(3,:),b*) /*绘图命令*/xlabel(轴it eta);ylabel(轴it phi) legend(平面驻点流动phi曲线,平面驻点流动dphi/deta曲线,平面驻点流动d2phi/deta2曲线)title(平面驻点流动的数值解)3、平面驻点的数值解和相应的曲线、平面驻点的数值解和相应的曲线 从结果中可以看出从结果中可以看出 在在 之后就趋近于之后就趋近于1，在，在 时时 =0.99，说明流动的流，说明流动的流速已经达到来流速度的速已经达到来流速度的99%。8 . 24 . 2第四节第四节 工程实际应用工程实际应用 盆地尺度的地下水流动过程中，多个水流系统交汇可以形成滞流盆地尺度的地下水流动过程中，多个水流系统交汇可以形成滞流区，是油气聚集、沉积矿产形成的重要部位区，是油气聚集、沉积矿产形成的重要部位 。 驻点可以分为盆地内部驻点和盆地底界驻点两大类，盆地内部驻驻点可以分为盆地内部驻点和盆地底界驻点两大类，盆地内部驻点（点（SP 1、SP 2）位于逆向局部水流系统的下方，是四个水流系统同）位于逆向局部水流系统的下方，是四个水流系统同时发生汇聚和发散