CN105004466A - 一种高精度非接触气动摩擦阻力测量方法及测量装置 - Google Patents

Abstract

本发明公开了一种基于粒子图像测速(PIV)原理的高精度非接触气动摩擦阻力测量方法及测量装置，用于在空气动力学风洞试验中通过非接触光学方法，快速测量模型表面气动摩擦阻力的一维分布，测量仪包括：双脉冲激光器和片光系统、接收和记录粒子散射光的CCD相机、同步控制器与计算机。在风洞试验中使用该测量仪测量得到模型近壁面的粒子图像对序列，其后采用单像素系综互相关算法对其进行图像处理计算，能够得到空间分辨率较高的近壁区平均流向速度场，根据牛顿内摩擦定理，由壁面处的速度梯度计算出精度较高的摩擦力及其一维分布。本发明具有精度高，对设备要求较低，测量时间快，对温度、湿度等外界条件不敏感等优点。

Description

一种高精度非接触气动摩擦阻力测量方法及测量装置

技术领域

本发明涉及激光测速技术领域中的粒子图像测速(Particle Image Velocimetry，简称PIV)技术，具体来说，是一种高精度非接触式气动摩擦阻力测量方法及测量装置。

背景技术

对于航空航天飞行器而言，由气流相对运动引起的壁面摩擦阻力是影响飞行器空气动力学特性的重要力学指标，是在地面风洞试验中需要准确测量的物理量。测量壁面摩擦阻力的方法包括：牛顿内摩擦定律方法，边界层动量积分方法，Preston管法，Stanton管法，热比拟方法，油膜干涉法和液晶涂层法等。这些方法各有利弊，简述如下：直接测量法基于牛顿内摩擦定律，通过测量壁面处平均流向速度的法向梯度来计算壁面摩擦阻力，虽然容易实现，但受测速仪器分辨率的限制，测量精度较低；边界层动量积分方法通过对速度型(平均流向速度的法向分布)积分可以有效消除部分测量误差，从而提高精度，但是测量工作量较大，至少需要测量两个流向站点的边界层速度型；Preston管和Stanton管实际上是Pitot管的一种发展和延伸，通过测量壁面外某点的平均流速来估计壁面处的速度梯度，虽然操作容易，但是测速管道的引入会对流场产生一定的干扰，测量精度不高；基于MEMS技术的微型传感器尺寸小、惯性小，一定程度上提高了测量的时间和空间分辨率，但传感器加工困难，成本较高；油膜干涉法利用薄油膜在剪切应力作用下厚度发生变化的特性实现对剪切应力的测量，但标定复杂，且对光路系统有很高要求；液晶涂层法虽然操作方便，但测试精度很低。此外，上述方法都引入了测试探头或涂层与流场进行接触测量，不可避免的对测试流场产生干扰，从而影响测试精度。近年发展起来的PIV技术以非接触光学测量的方法测量靠近壁面的流场速度，能够解决探头/涂层干扰问题。在PIV测速的基础上应用牛顿内摩擦定律，即可获得壁面摩擦阻力。由于壁面摩擦阻力计算依赖于壁面处流向速度的法向梯度，因此提高摩擦阻力测量精度的关键是得到空间分辨率较高的速度场。

传统的PIV图像处理算法因为查询窗口的存在而面临空间分辨率不足的难题，具体方式为：

应用传统PIV互相关算法计算跨帧图像对中示踪粒子的运动速度时，通常采用窗口互相关方法，查询窗口大小一般选择32像素×32像素。假设查询窗口大小为2M×2N，则相关函数为：

$R_{{\mathrm{II}}^{\′}}(\mathrm{\Δ}r,\mathrm{\Δ}s)=\underset{i=-M}{\overset{M}{\Σ}}\underset{j=-N}{\overset{N}{\Σ}}I(i,j)I^{\′}(i+\mathrm{\Δ}r,j+\mathrm{\Δ}s)---\left(1\right)$

式中，I为t时刻图像灰度，I'为t+Δt时刻图像灰度。查询窗口是一个2M×2N大小的矩形，(i，j)表示在矩形中位置坐标为(i,j)处的像素，假设坐标原点设在矩形的中心处；(Δr,Δs)表示t+△t时刻图像上与(i,j)坐标点相距(Δr,Δs)处的像素。

通过在相关区(根据来流速度，实验的视场范围以及微距镜头与测试平面的距离等条件确定)寻找使得相关函数取得最大值时的(x,y)，得到示踪粒子组在Δt时间内的位移Δs，根据速度等于位移除时间，即便可得到查询窗口区域内速度的平均值，此时的空间分辨率为32像素。

但由于在近壁区速度梯度急剧变化，查询窗口的存在会导致速度估值偏高，从而造成壁面处速度梯度计算不准。为了将PIV测速原理应用到非接触气动摩擦阻力测量上，必须解决PIV图像处理算法空间分辨率不高的问题。

发明内容

本发明提供一种高精度非接触气动摩擦阻力测量方法及装置，用于在空气动力学风洞试验中通过非接触光学方法，快速测量模型表面气动摩擦阻力的一维分布。

本发明一种高精度非接触气动摩擦阻力测量方法，通过下述步骤实现：

步骤1：在流场中播撒示踪粒子。

步骤2：照明流场。

步骤3：待流场稳定。；

步骤4：采集示踪粒子跨帧图像对。

对于典型风洞试验，采样时间需要超过2分30秒，采集的粒子跨帧图像对需要超过2000对；每对示踪粒子跨帧图像对为t与t+Δt时刻记录下的粒子图像，Δt为跨帧时间。

步骤5：对步骤4得到示踪粒子跨帧图像对进行处理，采用单像素系综互相关算法得到模型表面近壁区平均流向速度场及速度梯度。

A、对粒子跨帧图像对进行分组；

将步骤4中采集的全部示踪粒子跨帧图像对分开成两个组，第一组图像为每对示踪粒子跨帧图像对中在t时刻记录下的示踪粒子图像，记为set1。第二组图像为每对示踪粒子跨帧图像对中在t+Δt时刻记录下的示踪粒子图像，记为set2。

B、set1与set2中的图像匹配；

以灰度强度显示set1中每幅图像中的每一个示踪粒子，并采用单像素系综互相关方法，在set2中找到与set1中图像对应的图像中的示踪粒子进行匹配，具体方式为：

令示踪粒子跨帧图像对中，set1中的示踪粒子跨帧图像中像素点坐标为(i，j)，set2中对应的像素点邻域内相距坐标(Δr,Δs)的像素点坐标为(i+Δr，j+Δs)。

则N_f个示踪粒子跨帧图像对的系综相关函数用公式表示为：

$R_{\mathrm{\Δ}r,\mathrm{\Δ}s}(i,j)=\frac{1}{N_f-1}\×\frac{\underset{n=1}{\overset{N_f}{\Σ}}\left(I_{i,j}^n-\stackrel{\‾}{I_{i,j}}\right)\left(I_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}^n-\stackrel{\‾}{I_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}}\right)}{{\mathrm{\σ}}_{i,j}{\mathrm{\σ}}_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}}$

式中，和分别为第n个示踪粒子跨帧图像对中位于set1中图像与set2中图像像素点的灰度值；和是N_f个示踪粒子跨帧图像对中位于set1中图像与set2中图像像素点的灰度强度系综平均值，σ表示灰度强度的标准差：

${\mathrm{\σ}}_{i,j}=\sqrt{\frac{1}{N_f-1}\underset{n=1}{\overset{N_f}{\Σ}}(I_{i,j}^n-\stackrel{\‾}{I_{i,j}})}.$

得到N_f个示踪粒子跨帧图像对的系综相关函数之后，利用互相关算法找到系综相关函数R_Δr,Δs(i,j)的峰值，进而得到示踪粒子在set1和set2间的位移，最终得到示踪粒子的运动速度，且利用中心差分法得到某一法向高度y处的速度梯度。

步骤6：根据牛顿内摩擦定律计算出飞行器模型表面的摩擦阻力。

针对上述高精度非接触气动摩擦阻力测量方法的测量装置，包括风洞、激光器、片光系统、CCD相机、同步控制器与计算机；

其中，风洞用来产生并且控制气流；激光器的发射端安装有片光系统，片光系统用来将激光器发射的激光束扩散成片光；片光与流场的速度方向平行且垂直于飞行器模型需要进行摩擦阻力测量的表面；CCD相机的光轴与片光平面相垂直；通过CCD相机采集一定时间内沿飞行器模型需要进行摩擦阻力测量的表面外法线方向一定高度内的示踪粒子的跨帧图像对，转变成数字信号并存入控制计算机；同步控制器具有三个接口，分别通过连接线与激光器、CCD相机以及计算机相连接；通过同步控制器接收计算机发送的数字信号，同时用来触发CCD相机和激光器工作，使得CCD相机可在双曝光模式下采集示踪粒子的跨帧图像对。

1、本发明高精度非接触气动摩擦阻力测量方法及测量装置，在图像采集过程中采用长焦微距镜头能够以较高的分辨率记录粒子图像，利用单像素系综相关算法进行图像处理可以得到每一个像素上的时间平均速度，因此具有很高的空间分辨率，为准确获得壁面摩擦阻力奠定了基础；

2、本发明高精度非接触气动摩擦阻力测量方法及测量装置，能够获得更高精度的壁面摩擦阻力一维分布，并且对测量设备要求较低，测量时间快，对温度、湿度等外界条件不敏感。

附图说明

图1为本发明粒子图像测速方法流程图。

图2为单像素精度算法示意图。

图3为根据牛顿内摩擦定律计算壁面摩擦阻力的原理示意图。

图4为本发明在单像素精度算法实施例中计算得出的零压力梯度下典型二维湍流边界层的流向平均速度分布U⁺(y⁺)和流向速度脉动强度分布u_rms ⁺(y⁺)。

图5为本发明在单像素精度算法实施例中于低速风洞试验中测量得到的典型粒子图像。

图6a为本发明在单像素精度算法实施例中采用单像素系综相关算法和传统窗口相关算法计算的速度绝对误差的壁面法向分布。

图6b为本发明在单像素精度算法实施例中采用单像素系综相关算法计算的平均速度随机误差σ(e_U)的壁面法向分布。

图7为本发明在单像素精度算法实施例中采用单像素系综相关算法和传统窗口相关算法计算的壁面摩擦系数c_f相对误差的对比。

图8为本发明高精度非接触气动摩擦阻力测量仪结构示意图。

图中：

1-风洞     2-激光器      3-片光系统

4-CCD相机  5-同步控制器  6-计算机

7-示踪粒子

具体实施方式

下面结合附图对本发明作进一步详细说明。

本发明高精度非接触气动摩擦阻力测量方法，如图1所示，通过下述步骤实现：

步骤1：在流场中播撒示踪粒子；

在PIV实验中，示踪粒子大小及浓度的选择一定程度上影响着测量精度，故在进行实验前需要根据实验条件及要求选取适当大小的示踪粒子及浓度，通常情况下气体环境中粒径在1-5μm时综合性能较好；而由本发明通过验证实验已经得知，利用单像素系综相关算法的最优粒子直径为3～4像素，具体选择还要根据具体实验的视场范围、微距镜头与测试平面的距离、激光器输出光强等参数进一步确定。

步骤2：采用激光器照明流场。

步骤3：开启风洞待流场稳定。

步骤4：通过CCD相机，采集示踪粒子跨帧图像对；

对于典型风洞试验，CCD相机采样时间需要超过2分30秒，采集的粒子跨帧图像对需要超过2000对；每对示踪粒子跨帧图像对为CCD相机曝光t与t+Δt时刻记录下的粒子图像，Δt为跨帧时间。上述CCD相机曝光的间隔时间和跨帧时间的选择要根据流场的流速以及视野范围来选择，选择的依据就是使得示踪粒子在跨帧时间内移动的像素要大于粒子在图像中所占据的像素个数，小于10个像素大小左右。

步骤5：将由步骤4得到示踪粒子跨帧图像对，导入计算机进行数据处理，采用单像素系综互相关算法进行计算，得到模型表面近壁区平均流向速度场及速度梯度；

A、对粒子跨帧图像对进行分组；

将由步骤4中由CCD相机采集的全部示踪粒子跨帧图像对分开成两个组，如图3所示，第一组图像为每对示踪粒子跨帧图像对中CCD相机在曝光t时刻记录下的示踪粒子图像，记为set1；而第二组图像为每对示踪粒子跨帧图像对中CCD相机在曝光t+Δt时刻记录下的示踪粒子图像，记为set2。

B、set1与set2中的图像匹配；

以灰度强度显示set1中每幅图像中的每一个示踪粒子，并采用单像素系综互相关方法，在set2中找到与set1中图像对应的图像中的示踪粒子进行匹配，具体方式为：

令示踪粒子跨帧图像对中，set1中的示踪粒子跨帧图像中像素点坐标为(i，j)，set2中对应的像素点邻域内相距坐标(Δr,Δs)的像素点坐标为(i+Δr，j+Δs)，如图2所示；

则全部N_f个示踪粒子跨帧图像对的系综相关函数用公式表示为：

$R_{\mathrm{\Δ}r,\mathrm{\Δ}s}(i,j)=\frac{1}{N_f-1}\×\frac{\underset{n=1}{\overset{N_f}{\Σ}}\left(I_{i,j}^n-\stackrel{\‾}{I_{i,j}}\right)\left(I_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}^n-\stackrel{\‾}{I_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}}\right)}{{\mathrm{\σ}}_{i,j}{\mathrm{\σ}}_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}}---\left(2\right)$

式中，和分别为第n个示踪粒子跨帧图像对中位于set1中图像与set2中图像像素点的灰度值；和是N_f个示踪粒子跨帧图像对中位于set1中图像与set2中图像像素点的灰度强度系综平均值，σ表示灰度强度的标准差：

${\mathrm{\σ}}_{i,j}=\sqrt{\frac{1}{N_f-1}\underset{n=1}{\overset{N_f}{\Σ}}(I_{i,j}^n-\stackrel{\‾}{I_{i,j}})}.---\left(3\right)$

得到N_f个示踪粒子跨帧图像对的系综相关函数之后，利用互相关算法找到系综相关函数R_Δr,Δs(i,j)的峰值，即为set1中位于像素点(i,j)的示踪粒子在set2中对应示踪粒子间的时间平均空间位置。为了提高精度，采用二维高斯回归方法将系综相关函数R_Δr,Δs(i,j)的峰值的识别精度提高到亚像素量级。由此可以得到示踪粒子在set1和set2间的位移，即示踪粒子在t时刻与t+Δt时刻间的位移。最后，由于跨帧时间Δt已知，由此可得到示踪粒子的运动速度；且利用中心差分法可以得到某一法向高度y处的速度梯度。

步骤6：根据牛顿内摩擦定律计算出飞行器模型表面的摩擦阻力。

摩擦应力是单位面积上的摩擦阻力，二者之间关系是摩擦应力乘以面积等于摩擦阻力。

如图3所示，摩擦应力计算公式为：

${\mathrm{\τ}}_w=\mathrm{\μ}\frac{\∂u}{\∂y}{|}_{y=0}---\left(4\right)$

其中，μ为流体的动力粘性系数，为飞行器模型中需要进行摩擦阻力测量的表面近壁区平均流向速度的法向梯度，表示单位法向高度上速度的变化，如图3所示；u为流场速度，y是飞行器模型中需要进行摩擦阻力测量的表面外法线坐标。

上述方法求得的摩擦应力是单位面积上的摩擦阻力，二者之间关系是摩擦应力乘以面积等于摩擦阻力。

本发明中采用的单像素系综相关算法区别于一般PIV互相关算法的地方是：牺牲时间分辨率来换取更高的空间分辨率，故而单像素系综相关算法只能在每一个像素上获得时间平均速度。

实施例：

以零压力梯度下典型二维湍流边界层为例，通过对比单像素系综相关算法与传统相关算法在计算边界层的速度误差，来说明单像素系综相关算法在计算近壁区速度场时具有更高的精度和空间分辨率。

本例模拟了低速风洞1实验中自由来流速度U_∞＝10m/s的二维湍流边界层。如图4所示,给出了其流向平均速度分布U⁺(y⁺)和流向速度脉动强度分布u_rms ⁺(y⁺)的理论值。如图5所示为给出低速风洞试验中测量得到的典型粒子图像。

利用单像素系综相关算法和传统PIV互相关算法处理图像以得到速度场，图6a所示为两种方法计算的速度绝对误差的法向分布，是流向平均速度在不同法向高度上的平均。传统PIV互相关算法查询窗口尺寸为32×32像素，重叠率为50％，速度估计误差在在0.04像素附近。而单像素系综相关算法计算得到的从壁面处的0.005像素增加到y⁺＝5处的0.025像素，可见后者对近壁区速度型的估计精度更高。单像素系综相关算法的不确定度可用平均速度的测量不确定度(随机误差)来表征，如图6b所示，对比速度脉动型u_rms ⁺(y⁺)可以发现，测量不确定度与速度脉动水平呈正相关，由于在真实的湍流流动中，近壁区的湍流强度相对较小，所以测量不确定性的影响并不明显。

图7为两种方法计算出的壁面摩擦阻力系数相对误差的对比，可以看出，传统PIV互相关算法得到的|e_cf|约为0.1％，且与跨帧粒子图像对的样本数量无关。单像素系综相关算法对样本数量的依赖性较强，样本空间数量时速度估算的不确定度较大，但当跨帧图像对样本超过2000对后，单像素系综相关算法的摩擦阻力测量误差|e_cf|即已经小于传统PIV互相关算法；而当样本数量的达到10000对时，|e_cf|减小到约为0.01％，即摩擦阻力测量精度提高一个量级。由此证实，本发明所使用的测量仪和单像素系综相关算法，能够显著提高壁面摩擦力的测量精度。

以上实施例说明，本发明测量仪和单像素系综相关算法，在测量近壁面平均速度上具有很高的精度，从而能够显著提高壁面摩擦力的测量精度。相比于传统PIV方法，摩擦阻力的测量精度提高一个量级。

针对上述方法本发明还提出一种高精度非接触气动摩擦阻力测量装置，包括风洞1、激光器2、片光系统3、CCD相机4、同步控制器5与计算机6，如图8所示。

其中，风洞1用来产生并且控制气流，用来模拟飞行器或实体周围气体的流动情况，并可量度气流对实体的作用效果以及观察物理现象，它是进行空气动力实验最常用、最有效的工具之一。激光器2采用双脉冲Nd:YAG激光器2，激光器2的发射端安装有片光系统3，片光系统3用来将激光器2发射的激光束扩散成片光。由于需测量流向平均速度沿飞行器模型表面外法线方向分布，因此片光应与流场的速度方向平行且垂直于飞行器模型需要进行摩擦阻力测量的表面。CCD相机4是一种半导体器件，能够把光学影像转化为数字信号，并存入计算机6。由于边界层厚度尺寸较小，粒子浓度相对较低，为了拍摄到分辨率较高的图像，CCD相机4选择Nikon AF-S VR 105mm f/2.8G IF-ED长焦微距镜头，使CCD相机4的光轴与片光平面相垂直；通过CCD相机4采集一定时间内沿飞行器模型需要进行摩擦阻力测量的表面外法线方向一定高度内的示踪粒子7的跨帧图像对，转变成数字信号并存入控制计算机6。同步控制器5具有三个接口，分别通过连接线与激光器2、CCD相机4以及计算机6相连接。通过同步控制器5接收计算机6发送的数字信号，同时用来触发CCD相机4和激光器2工作，使得CCD相机4可在双曝光模式下采集示踪粒子7的跨帧图像对。

Claims (5)

1.一种高精度非接触气动摩擦阻力测量方法，其特征在于：通过下述步骤实现：

步骤1：在流场中播撒示踪粒子；

步骤2：照明流场；

步骤3：待流场稳定；

步骤4：采集示踪粒子跨帧图像对；

对于典型风洞试验，采样时间需要超过2分30秒，采集的粒子跨帧图像对需要超过2000对；每对示踪粒子跨帧图像对为t与t+Δt时刻记录下的粒子图像，Δt为跨帧时间；

步骤5：对步骤4得到示踪粒子跨帧图像对进行处理，采用单像素系综互相关算法得到模型表面近壁区平均流向速度场及速度梯度；

A、对粒子跨帧图像对进行分组；

将步骤4中采集的全部示踪粒子跨帧图像对分开成两个组，第一组图像为每对示踪粒子跨帧图像对中在t时刻记录下的示踪粒子图像，记为set1；第二组图像为每对示踪粒子跨帧图像对中在t+Δt时刻记录下的示踪粒子图像，记为set2；

B、set1与set2中的图像匹配；

以灰度强度显示set1中每幅图像中的每一个示踪粒子，并采用单像素系综互相关方法，在set2中找到与set1中图像对应的图像中的示踪粒子进行匹配，具体方式为：

令示踪粒子跨帧图像对中，set1中的示踪粒子跨帧图像中像素点坐标为(i，j)，set2中对应的像素点邻域内相距坐标(Δr,Δs)的像素点坐标为(i+Δr，j+Δs)；

则N_f个示踪粒子跨帧图像对的系综相关函数用公式表示为：

$R_{\mathrm{\Δ}r,\mathrm{\Δ}s}(i,j)=\frac{1}{N_f-1}\×\frac{\underset{n=1}{\overset{N_f}{\Σ}}\left(I_{i,j}^n-\stackrel{\‾}{I_{i,j}}\right)\left(I_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}^n-\stackrel{\‾}{I_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}}\right)}{{\mathrm{\σ}}_{i,j}{\mathrm{\σ}}_{i+\mathrm{\Δ}r,j+\mathrm{\Δ}s}}$

式中，和分别为第n个示踪粒子跨帧图像对中位于set1中图像与set2中图像像素点的灰度值；和是N_f个示踪粒子跨帧图像对中位于set1中图像与set2中图像像素点的灰度强度系综平均值，σ表示灰度强度的标准差：

${\mathrm{\σ}}_{i,j}=\sqrt{\frac{1}{N_f-1}\underset{n=1}{\overset{N_f}{\Σ}}(I_{i,j}^n-\stackrel{\‾}{I_{i,j}})}.$

得到N_f个示踪粒子跨帧图像对的系综相关函数之后，利用互相关算法找到系综相关函数R_Δr,Δs(i,j)的峰值，进而得到示踪粒子在set1和set2间的位移，最终得到示踪粒子的运动速度，且利用中心差分法得到某一法向高度y处的速度梯度；

步骤6：根据牛顿内摩擦定律计算出飞行器模型表面的摩擦阻力。

2.如权利要求1所述一种高精度非接触气动摩擦阻力测量方法，其特征在于：所述步骤1中示踪粒子直径为3～4像素。

3.如权利要求1所述一种高精度非接触气动摩擦阻力测量方法，其特征在于：所述步骤4中示踪粒子在跨帧时间内移动的像素大小小于10个像素。

4.如权利要求1所述一种高精度非接触气动摩擦阻力测量方法，其特征在于：步骤5的B中，采用二维高斯回归方法将系综相关函数R_Δr,Δs(i,j)的峰值的识别精度提高到亚像素量级。

5.针对权利要求1所述高精度非接触气动摩擦阻力测量方法的测量装置，其特征在于：包括风洞、激光器、片光系统、CCD相机、同步控制器与计算机；

其中，风洞用来产生并且控制气流；激光器的发射端安装有片光系统，片光系统用来将激光器发射的激光束扩散成片光；片光与流场的速度方向平行且垂直于飞行器模型需要进行摩擦阻力测量的表面；CCD相机的光轴与片光平面相垂直；通过CCD相机采集一定时间内沿飞行器模型需要进行摩擦阻力测量的表面外法线方向一定高度内的示踪粒子的跨帧图像对，转变成数字信号并存入控制计算机；同步控制器具有三个接口，分别通过连接线与激光器、CCD相机以及计算机相连接；通过同步控制器接收计算机发送的数字信号，同时用来触发CCD相机和激光器工作，使得CCD相机可在双曝光模式下采集示踪粒子的跨帧图像对。