基于图像处理的颗粒流动层模型研究

洪西洋，杨 晖，李 然，郑 刚

(1.上海理工大学 光电信息与计算机工程学院，上海 2000932；2. 上海理工大学 医疗器械与食品学院，上海 2000932)



基于图像处理的颗粒流动层模型研究

洪西洋1，杨 晖1，李 然1，郑 刚2

(1.上海理工大学 光电信息与计算机工程学院，上海 2000932；2. 上海理工大学 医疗器械与食品学院，上海 2000932)

如何方便且准确地测量滚筒中颗粒的运动速率分布是研究的难点。文中将给出一种基于数字图像处理算法的实验方法来解决此难题。通过找到最合适的极短曝光时间，可捕获在此曝光时间内颗粒运动轨迹的图像，然后利用数字图像处理算法计算出感兴趣区域的颗粒运动路程，即可计算出颗粒的运动速率分布。通过精确的实验数据分析可看出，本实验方法与理论仅有<8%的误差，有力证实了通过图像法测量颗粒速率分布可展开有效的研究工作。

图像处理；连续流模型；颗粒运动速率

HONG Xiyang1, YANG Hui1, LI Ran1, ZHENG Gang2

(1. School of Optical-Electrical and Computer Engineering,University of Shanghai for Science and Technalogy,Shanghai 2000932, China; 2.School of Medical Devices and Food Sciences,University of Shanghai for Science and Technalogy, Shanghai 2000932, China)

颗粒流越来越多的受到关注，目前的研究工作主要受技术需求驱动，对物理的认知偏少。颗粒材料的表面流动研究对工业生产和自然现象至关重要[1]。例如工业上旋转窖、搅拌机和喂养等系统中材料的运输、处理和存储等；沙丘、熔岩流、雪崩和河流沉积物等自然现象[1]。本文介绍连续流分析的基础框架，并给出滚筒中表面颗粒流的简化描述数学模型。通过图像处理方法进行实验数据的测量，实验结果表明，转动滚筒中一般颗粒流模型的可靠性，基于图像法的实验数据测量误差率<10%，适合滚筒中颗粒流的一般实验研究。

1 滚筒中表面流一般模型

图1(a)为一般连续流模型图。假设静态层上面的流动层颗粒流是单向的，并对两个颗粒层交界面的曲率影响可忽略。垂线平均连续流公式和x方向动量平衡公式被简化用于以下的假设；由于扩张的流动颗粒相对运动缓慢的颗粒数量较小，故认为流动层的密度近似常数。流动层速率呈线性vx=2u(y/δ)[1-4]，其中，u(x,t)是流动层的垂线平均速率；δ是流动层厚度。忽略颗粒缓慢的塑形形变，交界面的剪切应力为

图1 两种连续流模型示意图

(1)

其中，d为颗粒的直径；tanβs是有效动摩擦系数；βs为静态休止角；β为交界面切线与水平面的夹角；ρ为流动层密度；f(ρ)为密度相关常系数；g为重力加速度；δ为流动层厚度。根据实验性证据取f(ρ)=cδ/d(c≈1.5)。颗粒流动各物理量之间的关系方程为[4]

(2)

(3)

其中，Γ是从流动层流入静态层的流量[2-5]，假设静摩擦力完全作用于堆层交界处，由摩尔-库伦准则有[6-7]

Txy|y=0=-ρgδcosβtanβm

(4)

其中，tanβm是有效静摩擦系数[8]；βm为动态最大角；β为交界面切线与水平面的夹角；ρ为流动层密度；g为重力加速度；δ为流动层厚度。

图1(b)为转动滚筒中的连续流模型图。最简单情况，滚筒被填充50%。此时认为表面几乎是平整的，则交界面处的流量为Γ=ωx。带入式(2)得

(5)

其中，ω为转动滚筒的转动角速度；L在此处为转动滚筒半径；u为流动层颗粒垂线平均速率。

将式(5)带入式(3)简化得

(6)

由式(5)和式(6)联立消除流动层厚度δ得

(7)

(8)

其中，A=sin(β-βs)/cosβs。

对于滚筒这一特殊的对称模型；本文更关心滚筒中表面流的模型，由于剪切率随β的变化微小，为方便实验，取中点x=0的剪切率，此时cosβm≈cosβ由式(8)可得简化模型

(9)

其中，弗劳德数Fr=ω2L/g；尺寸比率s=d/L；A=sin(βm-βs)/cosβs。

2 实验方法

图2 实验装置示意图

2.1 实验装置

本实验中二维滚筒半径为L=0.15 m，厚度0.01 m；球形玻璃珠颗粒半径d=0.009 m；滚筒中填充接近1/2的颗粒；采用海顿步进电机，通过改变方波的频率控制滚筒转速，使颗粒的运动达到相对稳定的状态。

用休止角测量仪测出颗粒的静态休止角βs=33.1°。控制电机使颗粒运动达到稳定连续流状态，此时滚筒的转速为ω=0.86 rad/s。

图3 崩塌面和铅垂线拟合

为测得准确的βm，在昏暗的室内实验环境下调节照明系统，使拍摄的图像方便角度测量；使图像的一边与铅垂线平行。通过检测崩塌面最上面的亮点拟合出如图3所示的拟合线，根据拟合数据得出角度的正切值，即可计算出动态角度为βm=35.2°。

图4 图像预处理

保持测量βm时滚筒的转速不变，在黑暗的实验环境下，调节照明系统以达到合适的光照场景，采集文中所需的实验图像；通过调节相机的曝光时间采集到相对理想的图像，如图4(a)所示。此时相机的曝光时间为Δt=30 ms。

对于图4(a)的实验图像，文中是无法直接得到所需数据的。必须对图4(a)进行数字图像处理。

通过对图4(a)的分析，需要对图像进行一些预处理。首项对图像进行灰度化处理，将图像转化为灰度图方便后期的图像处理工作。然后根据直方图进行灰度变换，使得图像信息能更加明显，凸显信息部分；以便能更加直观的处理图像中包含的信息。经灰度变换后，得到图4(b)。

2.2 图像处理核心算法

从获得的实验图像可看出，其中包含了大量的背景噪声，这对获得图像信息是不利的，甚至掩盖了想要得到的图像信息。所以必须清除图像的背景噪声。

顶帽变换：用原图像减去上述经过开运算后的图像可去除面积较为规则的背景图像，保留尺寸较小的目标图像；常用于解决图像受不均匀光照导致的背景亮度不统一而无法直接做阈值变换的情况。对于本实验图像顶帽变换效果显著；顶帽变换算法定义如下

imtophat(A,B)=A-A·B

(10)

经过顶帽变换后的图像，如图4(c)所示；经过顶帽变换后的图像虽去除了大部分的背景噪声，可图像的对比度也有较大的削弱。对上图进行二值化处理后的结果图4(d)所示。

最小外接矩形算法：对于本实验，想要获得的数据是图中比较规则线段的长度，对图4(d)进行分析，最终选择利用最小外接矩形算法来求解线段的长度。本算法通过凸包算法求出线段点集的凸包，然后用枚举法找到凸包的最小外接矩形即线段的最小外接矩形。通过算法找到线段的最小外接矩形后，计算其对角线长度即线段长度。得出的初步实验结果如图5所示。

图5 实验结果图

3 数据和讨论

根据图像法测出的结果得出中点(x=0)垂线附近的颗粒运动速率，如图6(a)所示。

图6 颗粒流速率分布

由图6(a)的数据计算得出平均速率u′=0.25 m/s。将实验参数βs=33.1°、βm=35.2°、ω=0.86 rad/s、L=0.15 m、d=0.009 m带入深度平均速度理论式(9)

得到理论深度平均速率u=0.23 m/s

比较深度平均速率的理论值和实际值可知，误差为η=8%。

根据实验数据， 方向上速率的变化趋势如图6(b)所示，可知水平方向速率先递增后递减，与理论上速率的变化趋势基本吻合[5]。

4 结束语

由本文的理论与实验数据的有效验证可知，滚筒中颗粒流简化模型适用于滚筒中颗粒流的理论和实验研究，以及对滚筒中颗粒流的深入研究。通过实验的验证，滚筒的颗粒流简化模型理论与实际相吻合。这也在一定程度上说明了，滚筒中一般颗粒流模型的可靠性，同时也证实了实验中的一般现象。基于图像法的实验数据测量和验证相对方便简单，对实验器材和实验环境没有过度依赖，其误差率<8%，适合滚筒中颗粒流的一般实验研究。

[1] Yang H, Li R, Kong P,et al. Avalanche dynamics of granular materials under the slumping regime in a rotating drum as revealed by speckle visibility spectroscopy[J].Physical Review E,2015, 9(4):42-56.

[2] Orpe A V,Khakhar D V.Segregation of granular materials in rotating cylinders[J]. Physica A:Statistical Mechanics and its Applications,2003,318(1):129-136.

[3] Boutreux T,Raphael E,Gennes P G. Surface flows of granular materials: a modified picture for thick avalanches[J]. Physical Review E Statistical Physics Plasmas Fluids & Related Interdisciplinary Topics,1998,58(4):4692-4700.

[4] Hill K M,Gioia G,Tota V V. Structure and kinematics in dense free-surface flow[J].Physical Review Letters,2003,91(6):8535-8547.

[5] Hutter K,Koch T,Savage S.The dynamics of avalanches of granular materials from initiation to runout, part II:experiments[J].Acta Mechanical,1995,109(1): 127-143.

[6] Jop P, Forterre Y,Pouliquen O. Crucial role of side walls for granular surface flows:consequences for the rheology[J]. Journal of Fluid Mechanics,2005,541(6): 167-192.

[7] Khakhar D V, McCarthy J J, Shinbrot T,et al. Transverse flow and mixing of granular materials in a rotating cylinder[J]. Physics of Fluids,1997,9(9):31-43.

[8] Jain N, Ottino J M,Lueptow R M. An experimental study of the flowing granular layer in a rotating tumbler[J].Physics of Fluids,2002,14(2): 572-582.

[9] Khakhar D V, Orpe A V, Andresen P,et al. Surface flows of granular materials: model and experiments in heap formation[J].Journal of Fluid Mechanics,2001,441(1): 255-264.

[10] Orpe A V,Khakhar D V. Solid-fluid transition in a granular shear flow[J]. Physical Review Letters,2004,93(6) : 68-79.

[11] MiDi G D R. On dense granular flows[J]. European Physical Journal E,2004,14(4):341-365.

Particle Flow Model Based on Image Processing

This paper proposes an experimental method based on digital image processing algorithm for the accurate and convenient measurement of the velocity distribution of particles of the rotating cylinder in the roller. The most suitable short exposure time is determined first for the particle trajectory images to be captured within this exposure time. The movement distance of the region of interest is then calculated by the algorithm of the digital image processing to get the distribution of particle velocity. The accurate analysis of the experiment data show that the method and theory of this experiment has a mere error of less than 8%, indicating that measuring the particle distribution by image processing can expand the research work effectively.

image processing; continuum flow; velocity of particle

2015- 12- 29

国家自然科学基金资助项目(11572201)

洪西洋(1992-)，男，硕士研究生。研究方向：测试系统和图像处理。杨晖(1981-)，男，副教授。研究方向：颗粒技术。

10.16180/j.cnki.issn1007-7820.2016.10.026

TP 391.41

A

1007-7820(2016)10-089-04