液力偶合器漩涡流场PIV试验测量影响因素分析*

柴博森 项玥 马文星† 王志豪 魏亚宵 李振男

(1.吉林大学 机械科学与工程学院， 吉林 长春 130022； 2.吉林省产品质量监督检验院， 吉林 长春 130103)

液力偶合器流场是多种流动状态和物理效应并存的非定常多尺度漩涡流场[1- 2].精确的流场试验测量结果是研究液力偶合器流体动力学特性的重要载体，对于其性能改进和结构优化设计极为重要.随着流场测量理论与方法的快速发展，光学非接触式流动测量技术逐渐成为液力元件流场测量的主流技术，其中以粒子图像测速技术(PIV)和激光多普勒测速技术(LDV)应用较为广泛.德累斯顿工业大学Christen教授和福伊特公司Kernchen[3]采用LDV技术对部分充液率下液力偶合器的流场结构分布进行了研究，并与流场数值模拟结果进行了对比验证.Keunchul等[4]采用LDV技术对液力变矩器的主流区域流场进行了较高精度的试验测量.相对于LDV技术，PIV技术在全流场区域同步性测量方面具有明显优势.Kunisaki等[5]将PIV技术与电子计算机断层扫描(CT)技术相结合，测量了液力变矩器导轮流场结构.Kubota等[6]应用PIV技术测量了液力变矩器闭锁离合器区域流场结构分布.笔者所在课题组[7- 8]基于PIV技术测量了液力偶合器流场，并分析了不同工况下的流场结构演化规律.王玉岭[9]应用LDV技术测量了液力变矩器的流场分布，并对流场特征进行了降阶处理研究.从相关文献及技术调研来看，国内外在液力偶合器流场测量研究方面以定性研究居多，关于流场试验量化测量与误差分析的研究却鲜见报道.

在液力偶合器流场试验测量中，试验模型设计、示踪粒子选取、原始图像采集、流场数据后处理等均会不同程度地引入误差.这些误差的出现与叠加将严重地影响PIV试验测量精度，导致流场测量结果失真[10].文中通过理论建模与流场试验结果分析研究PIV流场试验测量中的误差来源，详细分析示踪粒子选择、图像采集和流场计算3个方面引起误差的原因，并提出降减和消除误差的途径与方法，以期为液力偶合器非定常漩涡流场PIV试验的精确测量提供参考.

1 试验系统

试验系统主要由机械部分、激光片光系统和图像采集部分构成，如图1(a)所示.机械部分包括林普机电有限公司生产的YS7124 型号三相异步变频调速电机、长春通用机械厂生产的2CY型号齿轮泵负载装置及CLZ型号联轴器.激光片光系统由西安远讯光电科技有限公司制造，可输出1～2 mm厚的激光片光.图像采集部分由元奥仪器生产的FR-340数字相机及配套采集软件组成.采用有机玻璃材料制造循环圆有效直径为230 mm的透明型液力偶合器样机，泵轮外壁面和涡轮外壁面之间的距离为90 mm，容积为2 710 mL，泵轮叶片数为13，涡轮叶片数11，且叶片厚度均为4 mm，如图1(a)所示.

图1 PIV试验系统及液力偶合器样机Fig.1 PIV system and prototype of hydrodynamic coupling

2 示踪粒子

2.1 跟随性力学模型

假设投入流场中的示踪粒子为球形刚性粒子，粒子跟随流体运动，在Stokes阻力FD、压力梯度力FP、虚拟质量力Fv、Basset力FB、外部势力F的作用下，粒子满足Lagrange运动方程

(1)

式中：dp为示踪粒子直径，m；vp为粒子运动速度，m/s；ρp为粒子密度，kg/m3；μ为流体黏性系数，kg/(m·s)；vf为流体质点速度，m/s.

2.2 跟随时间

假定ρp等于流体密度ρf，当ρp足够小时，Stokes黏性阻力起主导作用，忽略其他受力，只考虑粒子受Stokes黏性阻力作用，由式(1)得

(2)

t=0时，vp0=0，解方程(2)得

vp=vf(1-e-t/)

(3)

t0=-

(4)

常温(25℃)下取传动介质水的密度ρf=1.00 g/cm3，取水动力黏性系数μ=0.839×10-3Pa·s，采用PSP示踪粒子，取ρf=1.03 g/cm3.为了定量说明粒径对跟随性的影响，通过冷场发射电子显微镜实测3种粒子的物理形状和相对直径大小，并采集3种粒子的流场图像，如图2所示.

图2 粒径相对大小及流场图像Fig.2 Relative size of particle diameter and flow field images

1号粒子直径大小不一致，较大粒子等效直径的统计值为14.198 μm，较小粒子的等效直径为1.240 μm，流场图像中不同直径大小的粒子图像灰度差异性较大，对于基于图像灰度特征不变性的流场计算将引入误差；2号粒子等效直径为11.523 μm；3号粒子等效直径为1.437 μm.根据式(4)计算3种粒子的跟随时间t0，结果如表1所示.1号粒子中不同直径大小的粒子跟随时间相差数十倍，粒子跟随性无法判断；相对于2号粒子的跟随时间41.70 μs，3号粒子的跟随时间仅为0.65 μs，流动跟随性最好.

表1 示踪粒子跟随时间Table 1 Follow time of tracer particles

2.3 粒子浓度

粒子浓度由源密度NS和像密度Ni界定[12],

(5)

(6)

式中：C为单位体积内粒子含量数；Δz0为激光片光厚度，mm；M0为图像放大率，像素/mm；dt为粒子图像直径，像素；AI为查询域大小，像素2.

(7)

根据PIV测速原理，为了采集最好的图像效果，便于后续处理，要求NS≪1且Ni≫15[13].由式(5)-(7)可确定投入液力偶合器流场中球形粒子的最佳质量.试验结果表明，容积为2 710 mL的液力偶合器内加入1.0 g 3号粒子，能够满足PIV试验的粒子浓度要求.在制动工况下，液力偶合器内部循环流量最大，涡轮流场图像具有明显的漩涡流动特征，为了便于后续分析与结果讨论，采集涡轮径向切面连续2帧流场图像，如图3所示.

图3 制动工况下涡轮的连续2帧流场图像

Fig.3 Two consecutive flow field images of turbine in brake condition

3 图像采集

3.1 离轴角误差

假设CCD相机光轴与激光片光垂直，粒子图像经过透镜中心，由主光路投射到成像平面上.目标平面(x-y平面)上粒子初始位置在1处，初始坐标为(x0,y0,0)，根据透镜成像原理，1处粒子在像平面上点X1处成像，其坐标为(-M1x0,-M1y0,0)，图像放大倍数M1=di/d0，如图4所示.激光光片厚度Δz0=1～2 mm，物距d0=260 mm，由于Δz0相对于d0很小，因此可忽略不计，在后续分析中将只考虑初始坐标在x-y平面上的粒子运动情况.

经过Δt时间后，粒子由位置1运动到位置2，如图4所示.

图4 粒子位移引起的CCD像平面测量误差

Fig.4 Measurement error on CCD image plane caused by particle displacement

此时粒子位移与CCD相机光轴不垂直，且粒子运动脱离了目标平面，粒子最终坐标为 (x0+Δx,y0+Δy,Δz).位置2处粒子将在像平面上点X2处成像，坐标为(-M2(x+Δx),-M2(y+Δy),0)，图像放大倍数M2=di/(d0-Δz).根据透镜成像原理，位置2处粒子和目标平面内位置3处粒子在像平面内成同一像，由成像关系可得

(8)

其中，

ΔX=X2-X1=

示踪粒子从点1到点2的真实位移和目标平面内从点1到点3的虚拟位移之间将产生一个相对误差ε，它与离轴角θ有关，如图4所示.

(9)

式中，θx-z和θy-z分别是离轴角θ在x-z和y-z平面上的投影角，().如果Δz增大，将引起离轴角θ增大，此时图像透视误差将显著地影响流场测量精度.为了降低透视误差，需要根据待测流场特性合理策划2D-PIV视场布置方案，保证离轴角θ取最小数值，从而有效地降低相对误差ε.

3.2 光路折射误差

在激光片光照射下，示踪粒子反射出的光线依次通过水介质、有机玻璃和空气后，最终被CCD相机采集.理想情况下，粒子成像点P、成像中心和真实粒子点3点共线，此时成像点P对应的真实粒子点应是P2点，如图5所示.

图5 光路折射误差Fig.5 Refraction error of optical path

光线传播介质的折射率不同，将导致粒子成像线路发生改变，与成像点P真正对应的粒子点应该是P1，如图5所示.此时将产生粒子点的位置误差ΔP，根据几何关系得

ΔP=tg(tanθ1-tanθ2)+ZL(tanθ1-tanθ3)

(10)

式中：ΔP为流体介质折射产生的位置误差；ZL为激光切面到偶合器内壁距离，mm；tg为偶合器有机玻璃壁面厚度，mm；θ1为光线入射角，(°)；θ2为光线通过有机玻璃的折射角，(°)；θ3为光线通过水的折射角，(°).

根据射线光学理论[14- 16]可得距离ZL为

ZL=αz0+βtg

(11)

其中，

na=1.000、nw=1.333、ng=1.520分别表示空气、水介质和有机玻璃的折射率.

根据上述关系式可知，误差值ΔP与ZL、tg和θ1这3个因素有关.ZL越远，则ΔP越大；厚度tg越厚，则ΔP越大；粒子离光轴距离越远，则ΔP也越大.为了降低光路折射误差，将液力偶合器样机的图像采集区域设计为大型平面，如图1(b)所示.在保证样机强度的前提下应尽量减小玻璃壁厚tg；同时需要限制CCD相机的空间位置，保证相机光轴与激光片光垂直.

4 PIV计算

4.1 算法选择及计算结果

算法的合理选择是获得正确流场计算结果的关键.以图3流场图像为例，选择传统互相关算法和递归互相关算法分别计算.传统互相关算法流程如图6所示，流速提取步骤如下：

步骤1 将第1帧粒子图像划分成若干大小尺寸一致的查询区域(IW)，查询区域大小为n×n像素，以灰度函数f(x,y)表示，区域中心点坐标为(x1，y1).

步骤2 在第2帧图像中选取一个m×m像素大小的搜索区域(SW)，其中m>n；在SW区域范围内，以逐个像素点为中心选取n×n像素的相关区域，相关区域大小和第1步选取的查询区域大小一致，第2帧图像中相关区域用灰度函数g(x,y)表示，区域中心点坐标为(x2,y2).

步骤3 按一定的搜索路径计算第2帧图像中搜索区域内包含的若干相关区域与第1帧中查询区域的互相关系数Ci(f,g)，其中(i=0，1，2，…，m2-1)；在Ci(f,g)的集合内选择最大互相关系数Cmax，此时相关域图像g(x,y)和查询域图像f(x,y)最匹配；通过锁定互相关计算的最高峰值位置，提取f(x,y)和g(x,y)之间的位移，结合连续2帧图像间的时间间隔Δt，获得查询区域内的平均速度.

步骤4 移动第1帧图像中查询区域的位置，重复上述步骤直至覆盖完第1帧图像中的每一个查询区域，即可获得全流场速度矢量分布.

图6 传统互相关算法Fig.6 Outline of conventional cross-correlation algorithm

基于传统互相关算法对图3连续两帧粒子图像进行计算，查询区域大小分别取为64×64像素、32×32像素和16×16像素，计算结果如图7所示.

正值表示流速沿顺时针方向，负值表示沿逆时针方向

Fig.7 Calculation results of traditional PIV cross-correlation algorithm

递归互相关算法流程如图8所示，流速提取步骤如下.

步骤1 按传统互相关算法设置64×64像素的初始查询域，计算查询域f1(x,y)和相关域g1(x,y)的互相关系数C1(x,y)，通过锁定互相关系数最大峰值和时间间隔Δt，获得第1阶段流场速度.

步骤2 参照第1阶段流场速度，根据查询区域中心计算点周围8个不同的速度矢量，改变搜索区域形状，如图8中扭曲变形实线框所示，将第1阶段查询区域缩放为4个32×32像素的新查询区域，计算新查询区域f2(x,y)和相关区域g2(x,y)的互相关系数C2(x,y)，获得第2阶段流场速度.

步骤3 重复上述操作步骤，继续缩放查询区域大小，以16×16像素作为第3阶段计算的查询区域大小，获得第3阶段流场速度.

图8 递归互相关PIV算法流程Fig.8 Outline of recursive PIV cross-correlation algorithm

基于递归互相关算法对图3中的连续2帧图像进行3次迭代计算，计算结果如图9所示.

正值表示流速沿顺时针方向，负值表示流速沿逆时针方向

Fig.9 Calculation results of recursive PIV cross-correlation algorithm

4.2 结果讨论与分析

传统互相关算法中,查询区域大小对计算结果的精度和可靠性有重要影响.查询区域越小越能体现流场的局部特征，但是查询区域取得过小，区域粒子数量减少，将降低图像匹配计算的准确性，导致流速误矢量增加；查询区域过大，不仅会降低图像空间分辨率，而且还会增加运算量，导致计算效率降低.如图7(a)左图所示，选择64×64像素进行查询区域计算后获得的流速场信息较少，流场结构无法体现制动工况下液力偶合器内的整体漩涡流动规律，流速范围为0.10～0.50 m/s；涡量场结构能够体现出主流区域的漩涡流动和局部涡旋结构，涡量范围为-200～300 s-1，如图7(a)右图所示.选择32×32像素进行查询区域计算的结果如图7(b)左图所示，流速范围为0.05～0.35 m/s；涡量场中布满较多小尺度涡旋结构，但是主流区域上未能体现出整体漩涡流动特征，涡量范围为-200～200 s-1，如图7(b)右图所示.选择16×16 像素进行查询区域计算的结果如图7(c)左图所示，流场结构信息较为丰富，流速范围为0.02～0.18 m/s，流速梯度较大，多尺度涡旋结构明显，高流速区域以小尺度涡旋结构分布于整体流场结构中；涡量场中涡量范围为-200～150 s-1，2个低涡量数值区域靠近轮毂偏向下侧叶片，如图7(c)右图所示.

基于递归互相关算法，以64×64像素作为查询区域起始值，计算结果如图9所示，流速范围为0.10～0.70 m/s.在主流区域上存在较为明显的整体漩涡流动，见图9(a)区域1，在叶片与轮毂的交汇区域存在小尺度涡旋结构，见图9(a)区域2；涡量场中涡量数值范围为-200～350 s-1，在主流区域上整体漩涡流动沿顺时针方向，涡量数值为100～350 s-1，见图9(b)区域3，在轮毂和叶片交汇区域存在与主流区域循环流动相反方向的角涡结构，涡量数值为-200～-50 s-1，见图9(b)区域4，该角隅区域上出现的小尺度涡旋将对主流区域的整体漩涡流动造成延迟阻碍，引起液力能量损耗.

通过流场结构对比分析可看出，对于制动工况下的液力偶合器漩涡流场计算，采用递归互相关算法时，通过自适应改变查询区域的大小能够得到更加准确的流场结构，尤其是对于主流区域上漩涡流动特征和角隅区域上涡旋结构特征的捕捉能力更强；而采用固定查询区域的传统互相关算法无法适应漩涡流场结构的多变性，不能准确计算流速梯度较大的漩涡流场结构.

5 结论

文中利用电子显微镜实测示踪粒子物理形状及粒径大小，结合流场图像采集效果，认为采用等效直径为1.5 μm的PSP示踪粒子应用于水介质液力偶合器流场试验更为合适.推导了适合于液力偶合器PIV流场试验的粒子浓度检验公式，根据该公式确定容积为2 710 mL的液力偶合器流场中加入1.0 g直径为1.5 μm的PSP示踪粒子能够满足粒子浓度要求.建立了图像采集过程中的离轴角误差和光路折射误差的数学模型，提出将液力偶合器采集区域设计为大型平面并控制激光片光厚度及空间位置的采集策略，可有效地降低相关误差.针对制动工况下液力偶合器流场图像，采用传统互相关算法，使用固定查询区域的分析方法，无法捕捉到多尺度漩涡流场结构；而采用递归互相关算法，使用自适应查询区域分析方法，能够提取多尺度下漩涡流场结构特征，尤其是在主流区域整体漩涡流动特征和角隅区域局部涡旋结构特征提取方面具有明显优势.

参考文献：

[1] 初长祥，马文星.工程机械液压与液力传动系统(液力卷) [M].北京：化学工业出版社，2015.

[2] 刘应诚.液力偶合器设计制造与使用维修 [M].北京：化学工业出版社，2016.

[3] CHRISTEN M,KERNCHEN R.Fluid velocity in constant fill turbo couplings:Measurements using laser Doppler velocimetry [J].Antriebstechnik,2001,40(1):71- 74.

[4] KEUNCHUL C L,SEOUNG-CHOOL Y ,HAROLD J S.Quantification of primary flows of a torque converter using laser doppler velocimetry [C]∥Fisita World Automotive Congress.Korea:F2000A102,Seoul,2000:1- 5.

[5] KUNISAKI Y,KOBAYASHI T,SAGA T,et al.A study on internal flow field of automotive torque converter three-dimensional flow analysis around a stator cascade of automotive torque converter by using PIV and CT techniques [J].JSAE,2001,22(4):559- 564.

[6] KUBOTA T,KOBAYASHI T,SAGA T,et al.Application study of PIV measurement of flow field around lock-up clutch of automotive torque converter [J].JSAE,2003,24(4):426- 430.

[7] 柴博森，刘春宝．基于粒子图像测速技术的液力偶合器漩涡流动特性研究 [J].农业工程学报，2013，29(23)：86- 92．

CHAI Bo-sen,LIU Chun-bao.Study on vortex flow characteristics of hydrodynamic coupling based on particle image velocimetry technique [J].Transactions of the Chinese Society of Agricultural Engineering,2013,29(23):86- 92.

[8] 柴博森，项玥，马文星，等．制动工况下液力偶合器流场湍流模型分析与验证 [J].农业工程学报，2016，32(3)：34- 40．

CHAI Bo-sen,XIANG Yue,MA Wen-xing,et al.Analysis and experimental verification of turbulence models in flow simulation for hydrodynamic coupling under braking condition [J].Transactions of the Chinese Society of Agricultural Engineering,2016,32(3):34- 40.

[9] 王玉岭．液力变矩器流场分布特征及降阶处理研究 [D].北京：北京理工大学，2015．

[10] 李丹勋，曲兆松，禹明忠，等．粒子示踪测速技术原理与应用 [M]．北京：科学出版社，2012．

[11] CHRISTENSEN K.Review of particle image velocimetry:a practical guide,second edition [J].Aiaa Journal,2015,46(11):2974- 2975．

[12] ADRIAN R J,WESTERWEEL J．Particle image velocimetry [M].Cambridge:Cambridge University Press,2011．

[13] WESTERWEEL J,ELSINGA G E,ADRIAN R J．Particle image velocimetry for complex and turbulent flows [J].Annual Review of Fluid Mechanics,2013,45(1):409- 436.

[14] 崔宏滨，李永平，康学亮．光学 [M]．2版.北京：科学出版社，2015．

[15] CHARONKO J J,VLACHOS P P.Estimation of uncertainty bounds for individual particle image velocimetry measurements from cross-correlation peak ratio [J].Measurement Science and Technology,2013,24(6):1- 16.

[16] BERNHARD Wieneke.PIV uncertainty quantification from correlation statistics [J].Measurement Science and Technology,2015,26(7):1- 10.