周斌 , 李思维 , 陈志勇 , 张嵘
清华大学 精密仪器系, 导航工程中心, 北京 100084

收稿日期: 2016-01-19
基金项目: 总装预研基金资助项目（9140A09011514）
作者简介: 周斌(1976-), 男, 副研究员
通信作者: 陈志勇, 副研究员, E-mail:chendelta@mail.tsinghua.edu.cn

摘要：为了提高液体表面张力测量的准确性和稳定性，该文基于悬滴法中完全轮廓法计算原理，改进了二维寻优算法中变量初值的选取方法，并利用权重因子对变量误差函数进行了修正，减小了液滴轮廓偏离点对表面张力测量结果的影响。在温度为25℃，常压条件下，利用改进后的完全轮廓法对水和无水乙醇表面张力进行了测量，结果表明：改进后的完全轮廓法相比原完全轮廓法测量标准差减小了40%，充分证明改进后的完全轮廓法具有更好的测量重复性和稳定性。
关键词： 表面张力 完全轮廓法 悬滴
Full-profile fit pendent drop method for surface tension measurements
ZHOU Bin, LI Siwei, CHEN Zhiyong, ZHANG Rong
Engineering Research Center for Navigation Technology, Department of Precision Instrument, Tsinghua University, Beijing 100084, China


Abstract:The accuracy and stability of surface tension measurements are improved fitting the entire surface profile in the pendent drop method with an algorithm to determine the initial value for two-dimensional optimization. Moreover, the error function was modified by a weighting factor which reduced the influence of the droplet profile deviation point on the surface tension measurement. The full-profile fitting method was used to calculate the surface tension of water and ethanol at 25℃ and atmosphere. The final results show the standard deviation was reduced by 40% compared with the previous method, which verifies the stability and repeatability of this full-profile fitting pendent drop method.
Key words: surface tensionfull-profile fitting methodpendent drop
表面张力在化工、热力、能源等领域具有重要的意义，直接关系到材料的品质和性能。因此，对于原材料本身表面性质的认识，以及通过测量表面张力实现对材料表面性质的改善，尤为关键。
液体表面张力的测量方法有很多种，主要分为动态法和静态法^[1]。动态法主要有毛细管波法和振荡射流法；静态法主要有毛细管上升法、滴重法、最大气泡法、拉脱法^[2-3]、悬滴法^[4]、静滴法^[5]、表面波散射法^[6]、电磁悬浮法^[7]等。其中，悬滴法测量表面张力的精度较高，能够应用的温度范围较大，且测量设备简单，因此得到了广泛应用。
悬滴法在基于Young-Laplace公式的基础上，得出重力和表面张力的平衡关系。Andreas等^[9]于1938年提出了利用测量悬滴两个极限位置处的尺寸来计算液体的表面张力，即选面法。该方法主要通过测量液滴的2个特征平面来定义形状因子，并进行修正，最终求解表面张力。虽然操作简单，但计算精度较低。随着计算机和图像处理的发展，利用多点或者全部液滴轮廓点进行计算的方法，即完全轮廓法诞生并快速发展应用，这也大大提高了表面张力的测量精度。目前利用完全轮廓法测量悬滴表面张力的实际计算中，主要有一维寻优^[8]、坐标轮换法、单纯形法^[9]，以及牛顿迭代法^[10]等算法。
本文主要从完全轮廓法测量液体表面张力的基本原理入手，对该方法中涉及到的实际计算方法进行研究和改进，在二维寻优的基础上，通过改进对液滴顶点处曲率半径和形状因子初值的选取，以及利用权重因子对误差函数进行修正，从而提高表面张力的计算精度及稳定性。并利用改进后的计算方法分别对水和无水乙醇表面张力进行测量计算，来验证该方法的可行性。
1 表面张力测量方法1.1 基本理论和原理悬滴法的基本原理为：当液滴静止悬挂在毛细滴管口时，液滴的轮廓取决于液滴所受重力和表面张力的平衡。因此，通过对液滴外轮廓的采集测量，可以计算得到液体表面张力值。该方法是在Young-Laplace公式的基础上，推导得到的描述表面张力和重力平衡时的悬滴轮廓方程式^[11]。
对于任意曲面，液体内外压力差与液面曲率可表达为

$\Delta P = \gamma \left( {\frac{1}{{{R_1}}} + \frac{1}{{{R_2}}}} \right).$

(1)

式(1)即Young-Laplace方程。其中: γ表示液体的表面张力系数; R₁和R₂表示液面某一点处2个主要曲率半径值; ΔP表示液体内外压力差。
当液体只受重力作用，而没有其他外力条件时：

$\Delta P = \Delta {P_0} + \Delta \rho gz.$

(2)

其中: ΔP₀表示悬滴顶点处内外压力差; Δρ表示液体和相邻气体介质的密度差; g表示重力加速度，为恒定值; z为液滴纵坐标(如图 1所示)。对于液滴顶点处，由于对称性得：

${R_1} = {R_2}, \Delta {P_0} = 2\gamma /{R_0}.$

图 1 悬滴外形几何示意图

图选项

对于悬挂液滴，根据图 1所示关系，以及式(1)和(2)，对参数x, z, θ进行无量纲化，可以推导得出如下3个微分方程式^[12]：

$ \left\{ {\begin{array}{*{20}{c}}{\frac{{{\rm{d}}\theta }}{{{\rm{d}}S}} = 2-\beta Z-\frac{{\sin \theta }}{X}, }\\{\frac{{{\rm{d}}X}}{{{\rm{d}}S}} = \cos \theta, }\\{\frac{{{\rm{d}}Z}}{{{\rm{d}}S}} = \sin \theta .}\end{array}} \right. $

(3)

其中：β=ΔρgR₀²/γ；X=x/R₀；Z=z/R₀；S=s/R₀。坐标参数x, z, θ以及对应弧长s的描述参见图 1。由上述公式可知，当形状因子β和液滴顶点处的曲率半径R₀已知时，便可以求得表面张力系数。由于过去计算机应用和图像处理的局限性，导致对于液滴完整轮廓的提取比较困难，使得在当时通过选面法来计算液体的表面张力得到广泛应用。随着计算机和图像处理的快速发展，悬滴法逐渐由选面法这样的经验公式法过渡到了多点处理法，即完全轮廓法^[13-14]。该方法通过提取液滴外轮廓上多个点或全部点，并进行曲线拟合，从而求解表面张力。
由式(3)可知，指定一个初始值(X₀, Z₀, θ₀)和固定的步长dS，则根据Newton迭代和Runge-Kutta法，每一个形状因子β都对应一条确定的拟合曲线。令β在某个范围内变化，便可以得到一组不同的曲线拟合结果(X_i, Z_i)，如图 2所示。

图 2 β取不同值的拟合曲线

图选项

假设液滴实际轮廓点坐标为(x′_i, z′_i)，定义e_i为单点误差函数：

${e_i} = \frac{1}{2}\left[{{{\left( {{x_i}-x{'_i}} \right)}^2} + {{\left( {{z_i}-z{'_i}} \right)}^2}} \right].$

(4)

则实际轮廓与理论轮廓误差函数为

$E = \sum\limits_{i = 1}^N {{e_i}} .$

(5)

其中: x_i=R₀X_i，z_i=R₀Z_i。当式(5)作为目标函数满足某设定的收敛条件时，则得到β和R₀的解，从而求得最终的液体表面张力值。由于悬滴具有对称性，因此求解误差函数时，选择右轮廓点进行函数计算，将计算效率提高1倍。
1.2 实际轮廓提取读取液滴实际图像，利用Canny算子^[15]进行边缘检测，由于Canny算子可以有效抑制噪声，使得液滴边缘更加规整和清晰。
对数字图像从左上角开始逐行逐像素扫描，当扫描到该行第一个像素为1的点，即左轮廓点，则停止对该行的扫描并进入下一行，同时记录此像素点。由此方法，扫描结束时可得到液滴实际左轮廓。同样地，对图像从右上角逐行扫描，得到液滴实际右轮廓，从而提取出液滴完整实际轮廓。
根据滴管管口的实际尺寸(mm)和像素尺寸(pixel)对应关系，进行像素坐标和实际坐标之间的转换，最终得到液滴轮廓的实际尺寸坐标。由于数字图像坐标原点位于图像左上角，为方便后续计算，需将坐标原点变换至液滴顶点处，具体坐标变换参见图 3。

图 3 轮廓坐标变换

图选项

1.3 实验测量装置本实验测量装置如图 4所示，主要包含了功率可调的激光光源、散光片、悬挂液滴形成装置、图像采集(charge-coupled device，CCD)、以及光具减震台。CCD与液滴保持固定距离，并且聚焦在滴管出口处，使得捕获到的悬滴形状最准确。另外，增加散光片使点光源发散，以避免单点照射引起图像畸变。

图 4 实验装置

图选项

2 表面张力计算方法及改进利用完全轮廓测量液体表面张力的计算方法中，一维寻优计算简单，但测量误差较大。牛顿迭代法需要求目标函数一阶偏导、二阶偏导及Hessian矩阵，计算过程比较复杂，并且初值的选取对计算结果影响较大。坐标轮换法和单纯形法是通过判断计算结果是否满足预设的收敛条件，即目标函数是否小于某一极小值，当满足收敛条件则循环结束，得到变量值。目前坐标轮换法和单纯形法的初值选取都是依据选面法得到的。
2.1 对变量初值计算的改进利用坐标轮换法求解表面张力，需要给定合适的变量初始值，确定变量取值范围，从而减小计算量，实现快速计算。由于选面法是通过测定液滴2个极限位置处的尺寸，再利用经验公式求解β和R₀，计算精度较低。因此，对变量初值的计算进行了方法改进，以求解更加精确的变量初值，从而缩小变量取值范围。
首先根据实际液滴图像，运用最小二乘法对液滴底部实际轮廓点进行拟合，直接求解得到液滴顶点处的曲率半径R₀，为该变量的初值。并设R₀为常值，对β进行单变量寻优。令β为0.05~0.50，间隔0.001取值，利用四阶Runge-Kutta法计算得出一组理论拟合曲线。结合液滴图像实际轮廓点坐标，利用式(5)求解理论轮廓与实际轮廓间的误差函数。令误差函数最小时对应的β为该变量初值。
2.2 对目标函数的修正在提取实际液滴轮廓点过程中，可能出现个别轮廓点偏离整体轮廓线，导致目标函数值产生一定误差，造成最终变量计算结果的不准确，如图 5所示。

图 5 存在偏离点的轮廓拟合

图选项

假设中间曲线b为液滴实际轮廓线，曲线a和曲线c分别为β=0.1和β=0.2时理论拟合轮廓线。当实际轮廓线准确平滑时，利用式(5)计算误差函数，得到β=0.1时，理论拟合曲线a与实际轮廓b的误差E_a=0.050 0；β=0.2时，理论拟合曲线c与实际轮廓b的误差E_c=0.031 7，E_c < E_a。由此可得出曲线c与实际轮廓拟合结果更好，因此β=0.2为最终变量值。
当实际轮廓线有一定误差，存在个别偏离点时，利用式(5)求解误差函数值，得到E_a=0.042 3，E_c=0.043 4，E_a < E_c。因此，得到结论曲线a与实际轮廓拟合更好，β=0.1。
通过上述2种情况对比，发现实际轮廓点中个别偏离点，可能造成误差函数式(5)求解偏差，影响最终变量取值。因此，本文对误差函数式(5)通过加权进行了修正，以减小偏离点对计算结果的影响。令加权函数为^[16]

$w\left( t \right) = {\left( {1-t} \right)^2} + 0.125.$

(6)

其中：t=d_i/d_max；$ {d_i} = \sqrt {{{\left( {{x_i}-x{'_i}} \right)}^2} + {{\left( {{z_i}-z{'_i}} \right)}^2}} $；d_max为d_i最大值。(x_i, z_i)为理论轮廓坐标，(x'_i, z'_i)为实际轮廓坐标。

$ {e_i}\left( 2 \right) = \frac{1}{2}\left\{ {w\left( t \right)\left[{{{\left( {{x_i}-x{'_i}} \right)}^2} + {{\left( {{z_i}-z{'_i}} \right)}^2}} \right]} \right\}. $

(7)

利用加权后的误差公式(8)对图 5中轮廓存在偏离点的情况进行计算，得到E(2)_a=0.011 5，E(2)_c=0.011 0，E(2)_c < E(2)_a。显示曲线c与实际轮廓拟合更好，则β=0.2。该结果与实际轮廓准确平滑时得到的结论一致，证明该方法有效地减小了偏离点对计算结果的影响，使得最终变量计算更加准确。
2.3 改进后的表面张力计算方法利用节2.1中变量初值的计算方法，求得R₀和β，并根据初值条件下拟合结果产生的误差值，确定该计算中R₀的变量范围为(R₀±0.5) mm。另外，依据初值计算方法中求解得到的β值，将水和无水乙醇的β取值范围分别设为0.30~0.40和0.45~0.55。
利用坐标轮换法对R₀和β进行交替轮换单变量寻优，并先利用式(5)求解每组变量对应的误差函数值，判断是否满足预设的误差条件E < ξ。将所有误差函数满足条件对应的R₀和β存储，循环结束后，可得到一组R₀和β的对应序列。
再利用误差函数修正式(8)，对上述R₀和β序列进行新的误差计算，以减小实际轮廓中偏离点对计算结果的影响。通过计算修正后的误差函数值，得到满足修正后误差条件的R₀和β取值。最终利用公式γ=ΔρgR₀²/β求解得出液体表面张力γ。
3 实验结果利用悬滴法实验装置，对水和无水乙醇悬挂液滴进行了CCD图像采集，并读入Matlab编写的表面张力软件中，运用原完全轮廓法和改进后的完全轮廓法对其分别进行了表面张力计算，得到了如图 6所示水的表面张力测量结果和图 7中无水乙醇的表面张力测量结果。

图 6 水的表面张力计算结果

图选项

图 7 无水乙醇的表面张力计算结果

图选项

对6组水的表面张力测量数据进行分析，原完全轮廓法计算得到水的表面张力平均值为71.76 mN/m，6次实验测量值的标准差为1.36 mN/m。改进后的完全轮廓法计算得到水的表面张力平均值为71.76 mN/m，标准差为0.75 mN/m。改进后的算法与原算法相比，表面张力平均值相同，而标准差则减小了44.85%。另外，该温度条件下，水的表面张力真实值为72.00 mN/m，与该真值相比较，原完全轮廓法6次实验测量结果的标准误差为1.26 mN/m，而改进后的算法测量值与真实值之间标准误差为0.73 mN/m。可见，改进后的算法在测量准确度基本保持不变的前提下，稳定性有明显提高。
如图 7所示，对无水乙醇表面张力测量数据进行分析，原完全轮廓法计算得到无水乙醇的表面张力平均值为22.09 mN/m，标准差为1.74 mN/m；改进后的计算方法得到的表面张力平均值为22.14 mN/m，标准差为1.00 mN/m。与原算法相比，表面张力的平均值相差0.23%，而标准差减小了42.53%。与该温度下无水乙醇表面张力真实值(21.80 mN/m)相比较，原完全轮廓法得到的测量值与真实值标准误差为1.61 mN/m，改进后的算法标准误差为0.97 mN/m。再次证明了改进后的完全轮廓法与原完全轮廓法测量准确度基本一致，稳定性明显提高。
4 结论本文对完全轮廓法计算液体表面张力进行深入研究，并从变量初值选取以及误差函数修正2个方面进行了算法改进。利用水和无水乙醇表面张力测量实验进行了验证，实验结果表明：改进后的完全轮廓法相比原完全轮廓法，测量平均值基本一致，标准差减小了40%以上。充分证明改进后的计算方法具有更好的测量重复性和稳定性。

参考文献

[1]	Journal of Central South University(Science and Technology), 41(2):649-654.-->尹东霞, 马沛生, 夏淑倩. 液体表面张力测定方法的研究进展[J]. 科技通报, 2007, 23(3): 424–425.YIN Dongxia, MA Peisheng, XIA Shuqian. Progress on methods for measuring surface tension of liquids[J]. Bulletin of Science and Technology, 2007, 23(3): 424–425. (in Chinese)
[2]	Journal of Central South University(Science and Technology), 41(2):649-654.-->LI Chuanjun, CHEN Long, REN Zhongming. Application of ring method to measure surface tensions of liquids in high magnetic field[J]. Review of Scientific Instruments, 2012, 83: 1–2.
[3]	Journal of Central South University(Science and Technology), 41(2):649-654.-->王学军, 胡熙恩, 鲍卫民. 吊环法表面张力自动测量仪的研究[J]. 仪表技术与传感器, 1997, 11: 15–16.WANG Xuejun, HU Xien, BAO Weimin. Studies on automatical measurement of surface tension employing the ring method[J]. Instrument Technique and Sensor, 1997, 11: 15–16. (in Chinese)
[4]	Journal of Central South University(Science and Technology), 41(2):649-654.-->Moradian A, Mostaghimi J. High temperature surface tension measurement[J]. IEEE Transactions on Plasma Science, 2005, 33(2): 410–411. DOI:10.1109/TPS.2005.844995
[5]	Journal of Central South University(Science and Technology), 41(2):649-654.-->范建峰, 袁章福, 柯家骏. 高温熔体表面张力测量方法的进展[J]. 化学通报, 2004, 11: 804–805.FAN Jianfeng, YUAN Zhangfu, KE Jiajun. Development in measuring surface tension of high temperature molten liquid[J]. Chemistry, 2004, 11: 804–805. (in Chinese)
[6]	Journal of Central South University(Science and Technology), 41(2):649-654.-->Osada R, Hoshino T, Okada K, et al. Surface tension of room temperature ionic liquids measured by dynamic light scattering[J]. The Journal of Chemical Physics, 2009, 130: 1–5.
[7]	Journal of Central South University(Science and Technology), 41(2):649-654.-->Wang H P, Chang J, Wei B. Measurement and calculation of surface tension for undercooled liquid nickel and its alloy[J]. Journal of Applied Physics, 2009, 106: 1–3.
[8]	Journal of Central South University(Science and Technology), 41(2):649-654.--> Mohan-Nair K K. Measuring Interfacial Tension with the Pendant Drop Method[D]. Atlanta, USA:Georgia Institute of Technology, 2014.
[9]	Journal of Central South University(Science and Technology), 41(2):649-654.-->Andreas J M, Hauser E A, Tucker W B. Boundary tension by pendant drops[J]. Journal of Biomechanics, 1938, 42(8): 1001–1019.
[10]	Journal of Central South University(Science and Technology), 41(2):649-654.-->宁乔, 朱志强, 吕旭涛, 等. 图像法求液滴表面张力和接触角[J]. 空间科学学报, 2008, 28(1): 75–77.NING Qiao, ZHU Zhiqiang, Lü Xutao, et al. Determine the surface tension and contact angle of drop by image processing method[J]. Chinese Journal of Space Science, 2008, 28(1): 75–77. (in Chinese)
[11]	Journal of Central South University(Science and Technology), 41(2):649-654.-->CHANG Yaoyuan, WU Mingya, HUNG Yilin. Accurate surface tension measurement of glass melts by the pendant drop method[J]. Review of Scientific Instruments, 2011, 82: 3–4.
[12]	Journal of Central South University(Science and Technology), 41(2):649-654.-->李啸风, 陈志荣, 潘海华, 等. 结合图像数字化技术用悬滴法测动态界面张力[J]. 化工学报, 2001, 52(6): 546.LI Xiaofeng, CHEN Zhirong, PAN Haihua. Determination of dynamic interfacial tension by computer-aided pendant drop digitization[J]. Journal of Chemical Industry and Engineering, 2001, 52(6): 546. (in Chinese)
[13]	Journal of Central South University(Science and Technology), 41(2):649-654.-->Badran A A. Oscillating pendant drop:A method for the measurement of dynamic surface and interface tension[J]. Review of Scientific Instruments, 1986, 57: 259–261. DOI:10.1063/1.1138925
[14]	Journal of Central South University(Science and Technology), 41(2):649-654.-->艾芳洋, 毕胜山, 吴江涛. 全轮廓拟合悬滴法表面张力实验系统[J]. 物理实验, 2015, 35(5): 2–4.AI Fangyang, BI Shengshan, WU Jiangtao. Surface tension experiment system with full-profile fitting pendant drop method[J]. Physics Experimentation, 2015, 35(5): 2–4. (in Chinese)
[15]	Journal of Central South University(Science and Technology), 41(2):649-654.-->Canny J. A computational approach to edge detection[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1986, 8(6): 679–698.
[16]	Journal of Central South University(Science and Technology), 41(2):649-654.-->Foley T A. Scattered data interpolation and approximation with error bounds[J]. Computer Aided Geometric Design, 1986, 3(3): 163–177. DOI:10.1016/0167-8396(86)90034-8