水系统气泡形成过程的分子动力学模拟

velocity-verlet算法求得气泡系统的平衡状态python案例文章标题：探究气泡系统平衡状态：基于velocity-verlet算法的Python案例在科学研究领域中，气泡系统的平衡状态一直是一个备受关注的话题。

通过基于velocity-verlet算法的Python案例，我们可以更加深入地理解气泡系统平衡状态的相关理论和方法。

本文将从简单到复杂，由浅入深地讨论气泡系统平衡状态的研究过程，共享对这一主题的全面理解和个人观点。

一、气泡系统的基本原理气泡系统是由气泡和液相组成的多相体系，其平衡状态可以通过气泡的分布、形状和动力学行为来描述。

在研究气泡系统的平衡状态时，我们需要考虑气泡间的相互作用、液相的粘性和密度等因素。

通过velocity-verlet算法，我们可以模拟气泡系统在不同条件下的平衡状态，并进一步探讨其动态特性。

二、velocity-verlet算法的原理和应用velocity-verlet算法是一种常用的分子动力学模拟算法，通过同时更新位置和速度的方式来模拟粒子系统的运动。

该算法具有较好的数值稳定性和计算效率，在研究气泡系统平衡状态时得到了广泛的应用。

通过velocity-verlet算法，我们可以模拟气泡系统的动力学行为，探究气泡与液相之间的相互作用和平衡态的演化过程。

三、基于Python的气泡系统平衡状态模拟在进行气泡系统平衡状态的模拟时，我们可以借助Python编程语言中的相关库和工具来实现。

通过编写相应的程序和算法，我们可以对气泡系统的平衡状态进行模拟和分析。

在实际操作中，我们可以通过定义系统的能量函数、计算粒子间的相互作用力以及实现velocity-verlet算法来模拟气泡系统在不同条件下的平衡态，并对其进行深入的研究和分析。

四、案例分析与结论以一个具体的案例为例，我们将基于velocity-verlet算法的Python 模拟结果进行分析，并得出相应的结论。

通过对气泡系统的平衡状态进行模拟和分析，我们可以更加深入地理解气泡间的相互作用、液相的动力学行为以及平衡态的演化过程。

基于水动力学计算的气泡上升规律研究缴健;何生荣;李玉婷;罗玉峰【摘要】The relationship between rising time of air bubbles with different diameter and river depth is the basis for discharge measurement by rising air bubble technique.Motion of rising air bubbles is studied in this paper.The motion is classified into accelerat%不同直径气泡上升时间与水深之间的关系是气泡浮子法测流的基础。

针对静水中气泡上升的加速与匀速过程进行了计算与分析,将气泡上升过程划分为加速与匀速段,基于流体力学原理,求得气泡匀速上升速度与直径之间的关系、加速时间及其高度占总上升过程的比例,以及不同直径气泡冒出时间与水深之间的关系。

结果表明,气泡直径在[0.5,2.0]mm区间时最终速度的流态处于过渡区,在（2.0,4.0]mm区间时最终速度的流态处于紊流区,且最终速度随直径增大而增大;加速上升时间及其高度占总过程的比例极小,故可对气泡上升全过程近似做匀速处【期刊名称】《水利信息化》【年(卷),期】2011(000)003【总页数】4页(P41-44)【关键词】浮子测流法;气泡;上升速度;上升规律;水深;水动力学【作者】缴健;何生荣;李玉婷;罗玉峰【作者单位】河海大学水利水电学院,江苏南京210098;水利部南京水利水文自动化研究所,江苏南京210012;河海大学水利水电学院,江苏南京210098;河海大学水利水电学院,江苏南京210098【正文语种】中文【中图分类】TV131.20 引言河道流速是水文测验中的重要观测项目，目前一般采用流速仪法进行观测。

该方法需要缆道或测船，费工费时且难以实现自动化和连续测流。

实验二：分子动力学模拟-水分子扩散系数一、前言分子动力学模拟的基本思想是将物质看成是原子和分子组成的粒子系统（many-body systems ），设置初始位能模型，通过分析粒子的受力状况，计算粒子的牛顿运动方程，得到粒子的空间运动轨迹，可以求得复杂体系的热力学参数以及结构和动力学性质。

分子动力学模拟的理论是统计力学中的各态历经假说(Ergodic Hypothesis)，即保守力学系统从任意初态开始运动，只要时间足够长，它将经过相空间能量曲面上的一切微观运动状态，系统力学量的系综平均等效力学量的时间平均，因此可以通过计算系综的经典运动方程来得到力学量的性质。

比如，由N 个粒子组成的系综的势能计算函数为：int U U U VDW += (1-1)VDW U 表示粒子内和粒子之间的Van der Waals 相互作用；int U 表示粒子的内部势能（键角弯曲能，键伸缩能、键扭转能等）；根据经典力学方程，系统中第i 个粒子的受力大小为：U k z j y i x U F i i i i i ⎪⎪⎭⎫ ⎝⎛∂∂+∂∂+∂∂-=-∇= (1-2) 那么第i 个粒子的加速度可以通过牛顿第二定律得到：()()ii i m t F t a = (1-3) 由于体系有初始位能，每个粒子有初始位置和速度，那么加速度对时间进行积分，速度对时间积分就可以获得各个任意时刻粒子的速度和位置： i i i a v dt d r dtd ==22 (1-4) t a v v i i i +=0 (1-5)20021t a t v r r i i i i ++= (1-6) i r 和v 分别是系统中粒子t 时刻的位置和速度，0i r 和0i v 分别是系统中粒子初始时刻的位置和速度。

依据各态历经假说，可获得任意物理量Q 的系综平均，因此得到体系的相关性质：()()[]dt t r Q t t Q Q t t ⎰∞→==01lim (1-7) 分子动力学模拟能够计算体系的能量，粒子间的相互作用，角动量，角度以及二面角分布，剪切粘度，结构参数，压力参数，热力学参数，弹性性质，动力学性质等。

液体界面现象的分子动力学模拟在研究生期间，我的课题跟计算模拟有关，故本文利用分子动力学计算机模拟的方法,研究了液体的界面现象。

主要分三部分来说明，这三部分分别是分子动力学模拟方法，分子动力学模拟过程和具体的结果与讨论。

1.分子动力学模拟方法要研究气液界面性质，在模拟体系中必须要有气相、液相以及气液界面的存在，于是我们采用如图1所示的模拟盒，让液相和气相处于模拟盒的中部，固相处于模拟盒的两侧，这样在整个模拟体系中将有多种界面。

对气液界面现象进行分子动力学模拟研究时，以氩和甲烷为例，其分子间的相互作用采用L-J12-6 模型，分子i 和j 相互作用的势函数表示为：采用直角坐标系0xyz，分子在初始时刻为饱和液相，分子以边长为L的正方体FCC结构分布于模拟盒中部，其质心在原点（0，0，0）处。

模拟过程中，系统采用正则系综，用Verlet 算法对运动微分方程求解，并不断调整体系质心使之处于原点位置。

其中的运动方程和运动速度分别表示为：式中：Uij为分子i和j之间的势能，rij为它们之间的距离，fi为第i个分子所受到的作用力，mi为分子i的质量，vi为其速度，t为时间，ri为分子i的位移。

由于我们采用的是正则（NVT）系综，需要进行温度修正。

模拟盒在x、y方向的长度Lx=Ly=L，为了保证足够的气相空间，取在z方向的长度Lz=3L。

本模拟在三个方向上都采用周期性边界条件。

在计算机分子动力学模拟过程中，将所用到的各物理量进行无量纲化处理，以确保模拟体系各物理量之间的数量级较为接近，减小模拟过程中出现的误差。

2．分子动力学模拟过程在进行分子动力学模拟之初，每个分子的具体位置及其速度的大小和方向需要赋予初值。

在该研究中，采用直角坐标系0xyz，分子在初始时刻以边长为L的正方体FCC结构分布于模拟盒中部，其质心在原点（0，0，0）处；接着随机地给每一个分子速度赋予初值（本文采用Boltzmann分布对速度赋值，以节约机时），这样就可以跟踪每一个分子的运动状态。

《计算材料学》实验讲义实验二：分子动力学模拟-水分子扩散系数一、前言分子动力学模拟的大体思想是将物质看成是原子和分子组成的粒子系统（many-body systems ），设置初始位能模型，通过度析粒子的受力状况，计算粒子的牛顿运动方程，取得粒子的空间运动轨迹，能够求得复杂体系的热力学参数和结构和动力学性质。

分子动力学模拟的理论是统计力学中的各态历经假说(Ergodic Hypothesis)，即保守力学系统从任意初态开始运动，只要时刻足够长，它将通过相空间能量曲面上的一切微观运动状态，系统力学量的系综平均等效劳学量的时刻平均，因此能够通过计算系综的经典运动方程来取得力学量的性质。

比如，由N 个粒子组成的系综的势能计算函数为：int U U U VDW += (1-1)VDW U 表示粒子内和粒子之间的Van der Waals 彼此作用；int U 表示粒子的内部势能（键角弯曲能，键伸缩能、键扭转能等）；依照经典力学方程，系统中第i 个粒子的受力大小为：Uk z j y i x U F i i i i i ⎪⎪⎭⎫⎝⎛∂∂+∂∂+∂∂-=-∇= (1-2) 那么第i 个粒子的加速度能够通过牛顿第二定律取得： ()()i i i m t F t a= (1-3)由于体系有初始位能，每一个粒子有初始位置和速度，那么加速度对时刻进行积分，速度对时刻积分就能够够取得各个任意时刻粒子的速度和位置：ii i a v dt d r dtd ==22 (1-4) ta v v i i i+=0 (1-5)20021t a t v r r i i i i++= (1-6)i r 和v 别离是系统中粒子t 时刻的位置和速度，0i r 和0i v 别离是系统中粒子初始时刻的位置和速度。

依据各态历经假说，可取得任意物理量Q 的系综平均，因此取得体系的相关性质：()()[]dt t r Q t t Q Q tt ⎰∞→==01lim (1-7)分子动力学模拟能够计算体系的能量，粒子间的彼此作用，角动量，角度和二面角散布，剪切粘度，结构参数，压力参数，热力学参数，弹性性质，动力学性质等。

水被加热产生气泡的过程仿真摘要：一、引言二、水加热产生气泡的物理原理1.水分子结构2.加热过程中的水分子运动3.气泡形成的原因三、气泡产生的过程仿真1.仿真的方法与技术2.仿真结果分析3.结果与实际现象的对比四、仿真在科学研究中的应用1.对气泡生长过程的理解2.对其他相变过程的启示3.仿真技术的发展前景五、结论正文：一、引言水被加热产生气泡这一现象在日常生活中十分常见，但关于这一过程的物理原理以及气泡产生的具体过程却鲜为人知。

本文将通过对这一过程的仿真，来揭示其中的奥秘。

二、水加热产生气泡的物理原理1.水分子结构：水分子由一个氧原子和两个氢原子组成，呈V 形结构。

在水分子间存在着较强的氢键作用力，使水分子聚集在一起形成水滴。

2.加热过程中的水分子运动：当水被加热时，水分子的运动速度会加快，氢键作用力逐渐减弱。

随着温度的升高，水分子的热运动越来越剧烈，直至水分子逐渐转化为水蒸气。

3.气泡形成的原因：在水加热的过程中，由于水分子间的作用力减弱，使得水分子中的氢原子和氧原子逐渐分离。

当水分子中的氢氧分离达到一定程度时，气泡就会形成。

三、气泡产生的过程仿真1.仿真的方法与技术：本文采用分子动力学方法进行仿真，模拟水分子在加热过程中的运动状态。

通过计算机对水分子之间的相互作用力进行计算，并随着温度的变化调整水分子的运动速度。

2.仿真结果分析：经过一段时间的仿真计算，我们观察到了气泡的形成过程。

从仿真结果中可以看出，气泡的形成是一个逐渐发展的过程，从水分子逐渐分离，到形成稳定的气泡。

3.结果与实际现象的对比：通过观察实际水加热产生气泡的现象，我们可以发现仿真结果与实际情况非常接近。

这说明我们的仿真方法是有效的，可以用来研究气泡产生的过程。

四、仿真在科学研究中的应用1.对气泡生长过程的理解：通过仿真，我们可以更加直观地了解气泡产生的过程，有助于我们深入理解气泡生长机制。

2.对其他相变过程的启示：仿真技术不仅可以应用于气泡产生过程的研究，还可以用于其他相变过程的研究，例如冰的融化、油的蒸发等。

纳米尺度下气泡核化生长的分子动力学研究佚名【摘要】采用分子动力学方法模拟纳米尺度下液体在固体壁面上发生核化沸腾的过程,主要研究壁面浸润性对气泡初始核化过程和气泡生长速率的影响以及固-液界面效应在液体核化沸腾的能量传递过程中所起到的作用.研究结果发现:壁面浸润性越强,气泡在固壁处越容易核化.该结果与经典核化理论中“疏水壁面易于产生气泡”的现象产生了明显的区别.其根本原因是在纳米尺度下,固-液界面热阻效应不能被忽略.一方面,在相同的壁温下,通过增强固-液相互作用,可以显著降低界面热阻,使得热量传递效率提高,导致靠近壁面处的流体温度升高,气泡核化等待时间缩短,有利于液体沸腾核化.另一方面,气泡的生长速率随着壁面浸润性的增强而明显升高.当气泡体积生长到一定程度时,会在壁面处形成气膜,从而导致壁面传热性能恶化.因此,通过壁面的热流密度呈现出先增大后减小的规律.【期刊名称】《物理学报》【年(卷),期】2018(067)023【总页数】11页(P159-169)【关键词】纳米尺度;气泡核化;浸润性;界面效应【正文语种】中文1 引言随着微纳米技术在不同领域中的广泛应用,微纳设备中因高热流密度的产生而导致各种微纳器件的热管理问题成为亟待解决的核心难题.工质核化相变原理因具有高效换热特性,而作为一种理想方案被应用于微纳制造领域[1−3].利用相变原理进行换热将导致气泡产生及两相流的形成.因此,对微纳尺度下流体的相变诱导气泡形核生长及脱离等现象的研究成为传热传质领域一项具有挑战性的研究课题.深入研究其物理机制并优化微系统内部的热质传递过程,是节约资源和提高能源利用效率的重要手段.众所周知,纳米尺度下固-液界面效应对于微纳系统传热传质及工质相变现象具有显著影响.已有研究[4−6]发现,在微米通道的单相流动中,壁面浸润性影响固-液界面结构及流体原子与固体原子的黏附情况,进而影响界面能量与动量的传递.同时,界面效应也会对纳米尺度相变过程产生影响.因此,有必要对微纳设备中涉及的异质核化的机理及影响因素进行深入的分析与探讨.壁面浸润性对于异质核化具有重要的影响,根据经典核化理论,表面沸腾中液体在光滑固体加热表面上形成气泡所需的活化能与过热度及表面接触角有关.当壁面接触角越大时,气泡核化所需的活化能越少,发生沸腾的极限过热度越低.此结论已被大量实验数据及理论分析证实.Bourdon等[7]通过实验研究固体表面浸润性对池沸腾起始核化点的影响规律,发现起始核化过热度随表面接触角的增大呈单调递减趋势.Jo 等[8]认为在低热流密度区域,虽然疏水壁面更加有利于核化的发生,但是临界热流密度(CHF)却很低.此外,异质浸润性壁面与均质浸润性壁面相比,更有利于核化沸腾.Quan等[9]对固壁上异质核化沸腾进行理论分析,发现当近壁区液体存在温度梯度时,在异质核化的初始核化阶段,其临界核化半径随着接触角的增大而减小,可以促进气泡核化.Xu和Qian[10]采用范德瓦耳斯理论分析单气泡异质成核过程,得到与之类似的结论:当壁面为均质浸润性时,增大接触角或者过热度将促使气泡在加热壁面处膨胀扩展,增大气泡脱离直径,易于产生膜态沸腾.在常规尺度范围,传热传质以连续介质力学为理论依据,而微观尺度则以纳米级微观粒子理论作为依据,界面特性影响尤为突出,使得微纳相变过程明显区别于常规尺度.目前,已有学者对微纳尺度的相变传热问题展开研究.由于尺寸限制,在纳米尺度下开展实验研究存在很大的困难,而分子动力学方法直接面向原子层面,是研究原子尺度物理现象的有效手段.在均质核化方面,Kinjo和Matsumoto[11]采用分子动力学方法研究负压情况下的气泡空化过程,发现核化速率比经典核化理论预测值大8个数量级.Kimura和Maruyama[12]研究平板上的异质核化规律,模拟结果与经典核化理论预测的一致.Kinjo等[13]采用分子动力学方法研究受限微通道内流体的核化过程,发现壁面浸润性导致流体出现3种不同的核化过程,即弱吸附力壁面处的气膜形成过程、中等吸附力通道内的异质核化过程以及强吸附力通道内的均质核化过程.Mao和Zhang[14]模拟气泡均质成核过程,结果显示气泡生长规律与用RaleighPlesset方程预测的结果不一致.研究发现采用经典核化理论预测的极限过热度、临界核化半径、核化速率及气泡生长速率等数据与纳米尺度下获取到的实验结果相比,在某些情况下会相差很大的数量级[15].目前学术界对于该问题产生的深层原因仍然没有统一的定论.在异质核化方面,Nagayama等[16]认为气泡核化位置受壁面浸润性影响,在亲水性微通道内发生均质核化,气泡在主流液体区域产生;在疏水性微通道内发生异质核化,气泡在固体壁面处产生;此外,对于超疏水壁面,则不会在壁面处产生气泡,而是在固-液之间形成一层气膜.Bai和Li[17]研究发现对于浸润性较强的壁面,近壁区中能量较高的分子向主流液体区的运动受到限制,能量传递效率较低,从而导致界面接触热阻增大,固体界面温差也大幅增加,因此不利于核化发生.Novak等[18]通过模拟流体氩在固体壁面处的异质核化过程,发现气泡在浸润性弱的壁面处发生异质核化的时间有所减少,其原因是随着固-液相互作用增强,流体在近壁区呈类固体的形式排列,从而导致气泡在类固体排列流体区域上产生.Carey和Wemhof f[19]修正了Redlich-Kwong流体状态方程,用于预测受固壁影响的近壁区流体物性.近壁区流体由于受固壁的作用,导致该区域内的流体压力增高,流体发生核化所需的临界过热度显著提高.因此,距离固体壁面有一定距离的液体反而先达到核化条件开始沸腾,进而发生均质核化.然而也有部分文献研究得出相反的结果,Hens等[20]认为在非浸润壁面上气泡难以形成,亲水性壁面为气泡核化提供有利场所并促进了气膜的形成.Yamamoto和Matsumoto[21]推断固-液相互作用影响界面能量传递,从而调控沸腾核化行为.当壁面浸润性越强时,流体吸收的能量越多,因而更加易于核化的发生.综上所述,微纳尺度界面效应影响突出,使得微纳相变具有显著区别于常规尺度的行为特点.尽管对于纳米尺度下核化相变的研究已经取得了一定成果,但是目前的研究手段多数都是在初始时刻给定系统一个过饱和状态,最终达到一个稳定型态,整个研究过程类似等温系统中的空化反应,对于由壁面传热产生热纳米气泡的核化生长过程,仍然缺乏完整深入的研究,未能全面揭示固-液界面效应在沸腾相变过程中所起的作用.此外,不同学者对壁面浸润性与气泡核化之间的内在联系经常会持有相反的意见.因此,关于壁面浸润性对液体核化沸腾的影响还有待进一步验证分析.本文利用分子动力学方法模拟了纳米尺度下液体在不同浸润性壁面发生异质核化沸腾的完整过程,并着重分析探讨了固体壁面浸润性对气泡核化生长的作用机制.分析两者之间的关系,有助于加强对纳米尺度下气泡核化机理的理解,同时能为实际应用提供可靠的理论支撑.本模拟采用开源分子动力学模拟软件LAMMPS实现,原子位型实现可视化采用VMD软件.2 物理模型及模拟细节图1为模拟系统分别在初始、核化以及终了状态的分子模型xz平面分布图,y方向与xz平面垂直.本文模拟的是二维气泡的核化过程,模拟体系尺寸为Lx× Ly×Lz=172.9σ×5.76σ×461.2σ(σ为流体氩原子之间的尺寸参数),x和z方向的尺度显著大于y方向.系统沿着x,y方向均采用周期性边界条件,z方向为固定边界条件.固体壁面原子按照面心立方(FCC)晶格排列,晶格常数为1.15σ,其(100)晶面与流体原子接触,共有29100个固体原子.壁面厚度d=4.61σ,纳米凹槽的深度h=17.29σ,宽度w=18.45σ.固体壁面最外两层原子固定不动,作为边界壁来维持系统的稳定.对其余固体原子施加弹簧力作用,使其在初始位置附近振动,本文中采用的弹簧系数为3249.1ε·σ−2(ε为流体氩原子之间的能量参数).沿着z方向,流体氩原子分为两部分被置于固壁间,在上下侧液膜之间有一段真空区域.靠近下壁面流体液膜的厚度为115.29σ,上侧液膜厚度为11.52σ.初始氩原子按照FCC晶格排列方式布置,晶格常数为1.72σ,共有95035个流体原子.图1 模拟系统在(a)初始、(b)核化和(c)终了状态的分子模型x-z平面图Fig.1 .Molecular distribution in x-z plan of system atthe(a)beginning,(b)nucleation,and(c)f i nal period.流体氩原子之间的相互作用采用Lennard-Jones(L-J)势能模型,表达式为式中r为原子间的距离;流体氩原子之间的尺寸参数σ=0.3405 nm,能量参数ε=1.67×10−21J;原子质量m=6.69×10−23g.固体原子之间的相互作用也采用L-J 势能模型,σs=0.2475 nm,εs=8.35× 10−20J.固-液之间的势能作用通过对L-J势能模型进行修正[22]:式中由Lorentz-Berthelot守则得到能量参数εsl=尺寸参数σsl=(σs+ σl)/2, 下标s表示固体,l表示流体.调节能量系数α与尺寸系数β的值可以得到不同的浸润性.本文算例参数设置分别为α=0.14,β=0.6,0.7,0.8,0.9,1.0及α=0.2,β=1.0,在模拟过程中势能截断半径为3.5σ.超过此距离的原子,其相互作用忽略不计.采用Velocity-Verlet算法求解运动方程,时间步长取为0.0023τ,其中为特征时间,m为质量.在模拟过程中首先对整个系统采用正则系综(NVT),维持其温度为恒定值kB为玻尔兹曼常数,运行100万步,使系统达到平衡后撤掉系统整体温度热浴,仅对壁面施加温度控制,将下壁面温度升高到上壁面温度降低至流体原子采用微正则系综(NVE).根据固-液壁面浸润性,算例运行时间在40—600万步范围内,观察流体的核化轨迹.本文为了直观展示流体被壁面加热发生核化的过程,忽略了控制初始温度阶段的时间步统计,直接从加热阶段开始记录时间.在本文中,为了观察气泡在核化过程中的密度以及温度变化,将流体核化区域沿着z 方向划分为n层,第k个切片层(1 6 k 6 n)在第JStart到第JEnd内时间段内的无量纲粒子密度为式中Nk为第k层的液体原子数目;z为各液体分层的厚度;A为液体计算区域在xy 平面的面积,A=Lx×Ly.第k个切片层(16 k 6 n)在第JStart到第JEnd内时间段内的温度为式中α代表x,y和z三个方向;为粒子沿着α方向的热运动速度.根据原子携带的能量和其具有的速度可获得通过流体的热流密度,热流密度的微观表达采用(5)式计算:式中V为系统的容积;N为流体原子数目;ei为每个原子的能量,包括动能和势能的总和;Fij为第i个原子受到来自与第j个原子之间的相互作用力;vi为第i个原子的速度;rij为第i个原子与第j个原子之间的距离.3 结果与讨论通过调节固-液势能参数来改变壁面浸润性,图2表示两个粒子间的无量纲势能与相互作用力随无量纲距离的变化关系,Φ为势能,F为相互作用力.图2(a)中势阱深度反映固体原子对流体原子的束缚强度,势阱越深流体原子越不容易挣脱固体原子的束缚.当吸收的热量不足以克服势能壁垒的限制时,流体原子将在其平衡位置附近振动,呈现出类固体排列方式.图2(b)中无量纲力为正值时,固-液原子之间相互排斥;无量纲力为负值时,固-液原子之间相互吸引.由图2可见,当α一定时,β越小,固-液之间的作用力越弱,壁面浸润性更趋向于疏水性;当β一定,α越大,固-液之间的作用力越强,壁面浸润性越趋向于亲水性.上述结论与文献[22]的模拟结果一致.本文研究对比不同浸润性壁面处气泡的核化生长过程,可以用纯物质T-ρ相图分布中的状态点来观察液体氩的核化状态.从气泡核化过程可知,气泡核化位置沿着z方向主要分布在0—90σ之间.因而,将此空间范围设置为成核区域,统计该区域内流体的平均数密度和温度.每1000步输出一次模拟结果.由图3得知,成核区域的原子数密度随着流体温度的升高逐渐减小,从初始气液饱和态逐渐过渡到稳定过热态,表明本文的成核条件符合核化动力学规律.此外,需要统计远离液体核化区域的蒸汽相空间内的温度与压力,以便确定沸腾发生的状态.监测沿着z方向250σ—300σ气相区域内的温度与压力,发现气相温度一直维持在左右,压力为其相应的饱和压力0.007654εσ−3.因此,本文着重研究过饱和沸腾.图2 两个粒子间的无量纲(a)势能与(b)相互作用力随无量纲距离的变化关系Fig.2 .Dimensionless(a)potential energy and(b)interaction force depends on the distance between two particles.图3 核化区域热力学状态点及饱和线与包络线相分布图Fig.3 .Thermodynamic state point and phase diagram with coexistence curve and spinodal curve in nucleation region.3.1 气泡核化过程图4 不同壁面浸润性下异质核化过程(a)α=0.14,β=0.6;(b)α=0.14,β=0.7;(c)α=0.14,β=0.8;(d)α=0.14,β=0.9;(e)α=0.1 4,β=1.0;(f)α=0.20,β=1.0Fig.4 .Heterogeneous nucleati on process with dif f erent wallwettability:(a)α=0.14,β=0.6;(b)α=0.14,β=0.7;(c)α=0.14,β=0.8;(d)α=0.14,β=0.9;(e)α=0.14,β=1.0;(f)α=0.20,β=1.0.图4(a)—(f)表示气泡在不同浸润性壁面处发生异质核化的过程.由图可见,给壁面加热一段时间后,气泡胚核首先在凹槽内部形成,其体积随着壁面传热量的增多而增大.当气泡足够大露出穴面后,开始向两个方向生长,一方面向上膨胀,气泡高度增大;另一方面气泡向侧面铺展扩张,三相接触线逐渐向外延伸.由于沿着x,y方向为周期性边界条件,当气泡生长到一定程度会引起边界聚合,在壁面处形成气膜,推动液膜向上运动.由图4可知,随着α与β的增大,固-液相互作用增强,达到气泡初始核化所需过热度的时间缩短,即核化等待时间减少.当固-液势能参数为α=0.14,β=0.6时,固-液相互作用最弱,气泡发生核化的等待时间约为7360τ;固-液势能参数取α=0.2,β=1.0时,壁面浸润性最强,核化等待时间仅276τ.根据图3不同浸润性壁面处核化区域过热状态的演化规律,体系进入包络线区域的快慢程度有所差异,即核化等待时间不同.当势能参数为α=0.14,β=0.6时,状态曲线在亚稳态过热区紧贴饱和线缓慢上升后才进入稳态过热区;当势能参数为α=0.14,β=1.0时,则迅速进入稳态过热区域,气泡发生核化更快.图5给出不同研究结果中气泡约化核化时间随壁面浸润性的变化关系,CA表示接触角.气泡核化时间包括核化等待时间与气泡生长时间.由于不同尺度下,气泡核化时间的数量级存在很大的差别,本文为方便比较不同尺度下气泡核化时间随壁面浸润性的变化趋势,均采用约化比值.Phan等[23]与Gong和Cheng[24]的研究结果分别以接触角为43.8o与28.1o时的气泡核化时间为参考值,对气泡核化时间进行约化.本文的模拟结果则将势能参数为α=0.14,β=0.6的气泡核化时间作为参考值,约化核化时间取为1.Phan等[23]采用实验方法研究水在100 mm×5 mm纳米涂层壁面上的沸腾核化,发现随着壁面浸润性的增强,气泡脱离直径增大,核化时间增长,而气泡脱离频率却降低.Gong和Cheng[24]使用格子玻尔兹曼方法探究介观尺度下固体表面浸润性对液体饱和沸腾传热特性的影响机制.综上可知,在常规尺度范围,壁面接触角越大,即壁面浸润性越弱,气泡核化时间越短,越容易核化.在纳米尺度下,α与β值越大,壁面浸润性越强,核化等待时间与气泡生长时间均缩短,核化周期减小,反而更有利于核化发生.这是由于尺度效应的存在,导致核化规律在不同尺度下呈现出截然相反的变化趋势.图5 气泡约化核化时间与壁面浸润性的关系Fig.5 .Relationship between reduced nucleation time and wall wettability.图6 流体沿着z方向的密度分布图 (a)势能参数为α=0.14,β=1.0;(b)势能参数为α=0.14,β=0.6Fig.6 .Density prof i le of argon along zcoordinate:(a)Potential energy parameters is α =0.14,β =1.0;(b)potential energy parameters is α =0.14,β =0.6.在统计流体密度分布时,将整个系统内流体区域沿着z方向均匀划分层,每层的高度是2.3σ.图6表示不同浸润条件下流体沿z方向的密度分布随加热时间的变化趋势.图中红色虚线表示固-液界面位置.不同浸润条件下,其流体数密度分布规律类似.在初始时刻,液体区域的密度ρσ3维持在0.78附近.随着加热时间增加,气泡胚核在近壁区核化生长,核化区域的流体数密度逐渐减小.直到气泡消失,在壁面处形成气膜且液体层脱离壁面后,靠近高温壁面处的流体数密度接近气相密度.对比图6(a)与图6(b)流体数密度变化规律,可以看出壁面浸润性越强,核化区域靠近密度减小越快,反映出气泡核化生长速率随壁面浸润性的增强而增大.为了对比不同浸润条件下的核化状态,选取同一个时刻的流体数密度分布.图7表示时间为690τ时流体沿着z方向的数密度分布.结合上述核化过程原子位型图,近壁核化区域的流体有3种状态,即过热液体态、气泡核化态、过热气体态.势能参数为α=0.14,β=0.6时核化区域流体氩的数密度最大.保持α=0.14不变,增大β值,当β=0.7和0.8时,紧贴壁面处流体数密度有轻微的减小,几乎维持不变;当β=0.9时,近壁处流体数密度呈大幅度下降趋势,由图4可知,此时气泡胚核露出穴面;当固-液势能参数为α=0.2,β=1.0时,流体的密度已经接近气相的密度值,液膜脱离壁面,在壁面附近形成气膜.综上所述,在浸润性强的壁面处液体核化时间短,气泡核化速率快,气泡胚核体积更大.因此,浸润性强的表面更有利于气泡成核.图7 时间为690τ时氩流体沿着z方向的数密度分布Fig.7 .Number density profi le of argon along z coordinate at 690τ.3.2 气泡生长规律通过计算核化气泡的等效半径来定量分析壁面浸润性对气泡核化生长的作用.由于本文模拟的是二维气泡,因此气泡体积相当于在xz平面上的面积,为了计算气核的等效半径,将液氩核化区域在xz平面内划分成二维小网格,每个网格的尺寸为2σ×2σ,每100步统计一次网格内的流体原子数密度.将液体数密度与蒸汽数密度之和的1/2作为判断网格是气态还是液态的阈值[25].当网格内的无量纲原子数密度小于0.39时,将此网格标记为气态;当网格内的无量纲原子数密度大于0.39时,则认为此网格为液体.然后将核化区域内所有气态网格叠加起来得到气泡的二维面积S∗,计算气泡的等效半径[26]:图8表示气泡等效半径随加热时间的变化关系.由图8可见,不同浸润性壁面处的气泡生长规律类似,在初始加热阶段,液体没有发生核化,气泡的等效半径为0;当流体吸收的热量聚集到一定程度时,近热壁区发生核化,气泡等效半径开始增大.从图8可以看出,当气泡核化进入稳定生长阶段,气泡等效半径几乎遵循线性增长规律,此过程中气泡等效半径的生长速率随着壁面浸润性的增强而显著提升,势能参数为α=0.20,β=1.0时的气泡等效半径约等于势能参数为α=0.14,β=0.6的3.27倍.值得注意的是,本文将气膜在壁面处形成的时间记录为气泡停止生长的时间,即气泡等效半径变化曲线的统计结束时间.由图8可知,图中每条曲线终止时刻的气泡等效半径值有所差异,其随着壁面浸润性的增强反而减小.结合图4(e)—(f)分析其原因:一方面是壁面浸润性越强,气泡胚核不止在凹槽产生,还会在平壁面处形成,进而多个气核迅速聚合形成气膜;另一方面,流体原子在浸润性强的壁面形成气泡后,蒸汽与固壁间的液体膜更容易蒸发相变,三相接触线沿着x方向迁移,气泡铺展直径增大.当壁面浸润性弱时,气泡反而更倾向于向上生长.因此,浸润性越强的壁面在液膜脱离壁面时气泡的体积与等效半径越小.图8 气泡等效半径随加热时间的变化关系Fig.8 .Relationship between bubble equivalent radius and heating time.3.3 核化传热机理为研究壁面浸润性对气泡核化生长的作用机理,需要从固-液能量传递方面来进行分析.图9表示靠近壁面处的流体沿着z方向的温度分布随时间的变化关系.热量通过固-液界面扩散到流体区域,近壁区的液温先升高,在液体区域形成温度梯度,直至流体内部的温度呈现线性分布.当固-液势能参数为α=0.14,β=1.0时,近壁区流体温度几乎在920τ就达到线性分布,固-液界面温度约为,固-液温度阶跃约为. 而固-液势能参数为α=0.14,β=0.6时,其需要约2760τ才达到线性分布,此时,固-液界面处的流体温度值约为,固-液温度阶跃约为. 据此推断,界面热阻随固-液势能作用的减弱而增大,导致通过固-液界面传递热量的效率降低,进而使得流体温度升高达到线性分布的时间越来越长,对气泡的核化生长时间产生了显著影响.需要说明的是,图9中在−20σ—0范围内流体温度随加热的进行出现减小现象,这是由于此范围为凹槽中的液体原子核化变为气体,粒子数大量减少,造成统计误差.根据核化动力学,当活化分子团能量积累到临界活化能后,才能形成稳定的气泡并开始生长.由图9可知,浸润性强的壁面传热效率高,流体温度快速升高,明显高于浸润性弱的情形,使得核化区域液体能量迅速累积达到临界活化能,满足核化条件.而且气泡在生长过程中周围过热液体温度越高,越容易蒸发相变.该结论也可从核化区域流体总势能变化规律得到验证.图10表示核化区域流体总势能随加热时间的变化规律.值得注意的是,图中势能值为负值.初始加热阶段,气泡尚未形成稳定的胚核,液体原子仍紧贴壁面.当势能参数为α=0.2,β=1.0时,固-液作用最强,核化区域的流体总势能值最小.与之相反,当势能参数为α=0.14,β=0.7时,其总势能值最大.随着加热进行,前者的势能值显著增大,反而高于后者,促进了气泡的核化生长.此外,壁面势能参数α与β越大,能垒越深,固-液相互作用越强,液体原子越难挣脱固体壁面束缚,导致气泡胚核的初始核化位置随浸润性变化而有所差异.如图4所示,气泡胚核在纳米凹槽内初始形成阶段,壁面浸润性越弱,越趋于在壁面处形成;而浸润性越强,则气核在偏离固体壁面的液体中央形成.在气泡生长后期,贴近壁面处的液膜受热蒸发.对于浸润性强的壁面,其热阻小,近壁面处液膜温度高,液体更容易相变成气体.综上所述,尽管能垒越深,紧贴壁面的液体原子越难以挣脱束缚,但是近壁区流体吸收能量增多,明显利于形成气泡胚核,只是紧贴壁面的流体原子受到固体壁面的束缚程度不同,导致初始形核位置有所差异.图9 核化区域氩流体沿着z方向的温度分布 (a)势能参数为α=0.14,β=1.0;(b)势能参数为α=0.14,β=0.6Fig.9 .Temperature prof i le of argon along zcoordinate in nuclear boiling region:(a)Potential energy parameters is α =0.14,β=1.0;(b)poten tial energy parameters is α =0.14,β =0.6.为进一步探究核化过程中热量传递的规律,图11给出气泡核化区域流体的温度随时间的变化规律.由图可知,壁面浸润性越强,通过界面传递的热量越多,流体温度升高越快.在气泡核化等待阶段,底板传递的热量直接加热近壁区液体,便于能量聚集使其发生核化.当液体发生核化,气泡生长到一定程度,部分底板与气体直接接触,气固直接传热效率较低,气泡生长核化的热量主要来源于周围高温液体.如图11温度曲线所示,当核化区域流体温度上升到一定程度开始下降,这是由于底板传给液体的热量不足以弥补气泡边界液体蒸发相变及气泡体积膨胀所需能量,核化区域温度反而降低.当势能参数为α=0.2,β=1.0时,核化区域温度几乎没有下降过程,在1023τ时刻温度又以不同速率迅速升高;而α=0.14,β=1.0,0.9,0.8,0.7时则分别在。

水里气泡的原理水里气泡的形成原理涉及到液体的物理性质以及气体溶解与析出的过程。

当液体中存在气体时，由于气体的分子热运动与压强差异，会导致气体分子从液体中溶解和析出，形成气泡。

首先，液体的表面张力与气液界面张力是水里气泡形成的基本原理。

液体分子间存在相互作用力，使其分子呈现尽量小接触面积的趋势。

这种现象被称为表面张力。

在液体中存在气泡时，气液界面周围的液体分子会受到周围液体分子的牵引，导致气泡表面产生一定凸起，即气液界面处有向内凹陷的趋势，使得液体表面张力对气泡产生的作用力方向指向气泡内部。

这个作用力使得气泡具有收缩的特性，使气泡内部的气体被约束在气泡内。

其次，气体溶解与析出的过程是水里气泡形成的重要原理。

气体在液体中的溶解与析出是由于气体溶解度受温度、压强和溶液成分等因素的影响而发生。

当液体中存在饱和溶解气体时，溶解气体分子与液体分子之间的相互作用力达到一定平衡，液体内部气体分子的动能趋于平均值。

然而，当液体受到机械或热力作用时，如气泡形成、振动、加热等，会导致液体内部的气体分子的动能增大，溶解气体分子从液体中析出形成气泡。

这些气泡在液体中的增大与形成与液体分子的热运动、压强差异等因素密切相关。

此外，气泡的形成原理还与传质过程有关。

传质过程是指气体分子在液体中的扩散与对流的过程。

当液体中存在气体溶解时，溶液内部的气体分子会通过扩散方式从高浓度区向低浓度区传递。

当外部因素改变液体内的气体浓度分布时，扩散过程会被加强或减弱，从而影响气体的溶解度和气泡的生成。

此外，如果液体中存在气泡的运动或扩散，会引起液体的对流，对流过程会加速气泡的形成与传输。

总结起来，水里气泡的形成原理主要涉及液体表面张力的作用以及气体溶解与析出的过程。

液体表面张力对气泡的收缩与稳定具有重要作用，而气体的溶解度和传质过程也影响着气泡的形成与扩散。

这些因素的相互作用使得水里气泡形成成为可能。

气泡动力学特性的研究与应用随着科技的发展，气泡动力学逐渐成为了研究和应用的重要方向。

气泡既是一种普遍存在于自然界中的物质，又是一种可用于工业生产和科学实验研究的重要手段。

气泡的动力学特性研究既有理论意义，也有实际应用价值，本文将就此探讨。

一、气泡动力学特性的基本概念气泡是一种空气或其它气体包裹在水（或其它液体）中的球形或半球形体。

气泡通常是由于振荡、撞击、渗漏等原因形成。

在自然界中，气泡广泛应用于海洋、人体生理、大气、地表水、燃烧和环保等领域。

此外，气泡也是科学实验和工业生产中常用的物质。

气泡动力学特性研究的目的是解析气泡所受到的运动和外力作用的物理特性，如气泡在液体中的流动、振荡、破裂、生长等过程。

气泡在液体中的运动主要受到重力、表面张力、动量和浮力等力的作用。

气泡大小和形状、液体性质、气泡运动速度等因素都对气泡运动和特性产生影响。

依据不同研究对象和方法，气泡动力学特性研究可以分为理论分析、实验和数值模拟三种不同形式。

二、气泡动力学特性的研究方法（一）理论分析气泡动力学特性的理论分析主要通过数学物理方程模型建立，通过求解方程得到特定气泡的运动和特性。

气泡运动与物理特性耦合的物理方程组主要包括Navier-Stokes方程、质量守恒方程、气泡表面张力方程、以及边界条件等方程式。

通过对方程解析求解，可以得到气泡育形、壁压、速度、流场等运动参数和字符参数。

理论分析的优势在于可以给出简洁而通用的模型，能够预测和探索气泡特定运动特性，还可以为实验和数值模拟提供参数参考。

不过，理论分析方法的不足之处在于常常需要解答很多数学问题来获得分析和预测结果，这需要特定的数学技术，难以解决实际工程和生产中的某些问题。

（二）实验气泡动力学特性的实验研究可以通过光学实验、水力学模型实验、压力实验等方式进行。

常见的实验设备包括气泡发生器、气泡观测装置、高速相机、光学显微镜等。

实验能够定量获取气泡的运动速度、形态、壁压、生长和破裂循环等动态信息，具有无可替代的优势。

有限水域水中爆炸气泡脉动的数值模拟随着人们对海洋资源的开发和利用的不断加深，越来越多的工业污染和人为干扰正在导致海洋环境的恶化。

其中，有限水域是容易受到这些干扰的海区之一。

在这些区域中，水体的物理和化学性质可能会发生重大变化，从而影响到生态系统的平衡和健康。

因此，对于这些海洋区域的研究和监测显得尤为重要。

其中，水中爆炸气泡脉动的数值模拟是一种有效的手段，可以帮助我们更好地了解这些海域的物理特征和环境变化。

在水中爆炸气泡脉动的数值模拟中，主要就是对水波的物理特征进行模拟。

水波一般包括两种类型：远场波和近场波。

远场波是指成簇的波，它们频繁出现在爆炸点附近，伴随着高强度的涡旋和涡街。

近场波则是指趋向于中心和远离中心的波动，其大小和形状可以通过数学模拟计算获得。

在数值模拟中，我们可以使用 Navier-Stokes 方程，即流体力学基本方程，来描述水波的物理特征。

这些方程可以帮助我们计算出水波中每个点的运动和压力等参数，并可以通过计算机模拟得到这些参数随时间、空间和深度变化的趋势。

通常，我们使用计算流体力学（CFD）来进行这种数值模拟。

在这个过程中，我们需要为模拟设置合适的初始参数和运动参数，如初始流速、初始压力、边界条件等。

这些参数的设置直接影响着数值计算的准确性和结果的分析。

通过数值模拟，我们可以深入了解水中爆炸气泡脉动的机制和原理，从而指导海洋环境的监测和管理。

我们可以据此预测爆炸气泡造成的水波和浪涌等影响、预判海洋环境的变化趋势、以及制定相应的环保和资源保护政策。

总之，水中爆炸气泡脉动的数值模拟是一种有效的工具，可用于研究海洋环境的物理和其它特征，从而为我们更好地了解这些海域、保护海洋资源和维护生态平衡提供有益的参考和指导。

在水中爆炸气泡脉动的数值模拟中，液体和气体的物理参数是进行模拟的关键。

下面列出一些相关数据，并进行分析：1. 液体的密度：水的密度约为1000 kg/m3。

当液体密度越大时，水波的速度就会越快，水的表面张力也会增强。

分子动力学模拟计算水分子扩散系数分子动力学模拟是一种计算物质在原子水平上运动行为的方法，通过模拟大量的原子或分子通过经典力场相互作用的过程来研究物质的宏观性质。

水分子是地球上最重要的分子之一，其在自然界中的运动行为与生命活动密切相关。

水分子的扩散系数是描述水分子扩散速率的重要参数，它在环境污染、生化反应等领域中有着重要应用价值。

本文将介绍使用分子动力学模拟计算水分子扩散系数的方法。

在分子动力学模拟中，水分子通常采用经典的力场模型进行建模。

该模型使用Coulomb势函数描述电荷之间的相互作用，使用Lennard-Jones 势函数描述范德华力的作用。

当然，也有一些更高级的模型如远程碰撞模型等可以进行更精确的描述。

通过定义原子或分子间的力场，可以计算分子的受力情况。

在模拟中，通常需要设置一定的初始条件和模拟参数。

初始条件包括系统的初始位置、速度和受力情况等。

模拟参数包括模拟的时间步长、模拟的总时间等。

通常可以根据实际情况进行参数的选择。

为了减少模拟的误差，通常需要对系统进行一段时间的热化来使系统达到平衡状态，然后再进行实际的模拟计算。

在模拟过程中，通过计算分子之间的相互作用力，可以得到每个分子的运动轨迹。

通过分析分子的平均位置和速度变化等参数，可以得到水分子的扩散系数。

水分子的扩散系数可以通过计算分子的平均自由程和平均运动时间来获得。

平均自由程是分子在单位时间内所能行进的平均距离，平均运动时间是分子通过平均速度所需要的时间。

从这个角度来看，水分子的扩散系数与分子的速度和行进距离密切相关。

当然，实际的分子动力学模拟是非常复杂和耗时的。

模拟中需要考虑大量的分子和大量的相互作用，这就需要大量的计算和存储能力。

此外，还需要有效的算法和技术来加快计算的速度和改进模拟的精度。

因此，分子动力学模拟通常需要在高性能计算机中进行，以确保计算结果的准确性和可靠性。

总之，分子动力学模拟是一种计算水分子扩散系数的重要方法。

水中放电气泡脉动过程的数值模拟研究水中放电气泡脉动过程是水中放电分析中重要的一个模型，具有实际工程应用价值和科学研究价值。

在最近几年来，由于计算机技术的发展，使得计算机模拟成为了水中放电气泡脉动过程的重要研究手段。

以数值模拟的方法提供了一种可行的方法来研究水中放电气泡脉动过程。

首先，通过计算机的模拟，研究人员可以清楚的了解水中放电气泡的脉动过程，计算出不同的参数下水中放电气泡的形状和脉动特性，如气泡的大小、个数、周期、变化率等，从而更精确地模拟水中放电气泡脉动过程。

其次，水中放电气泡脉动过程的数值模拟可以提供有用的结果，便于研究人员了解水中放电气泡脉动过程的特性，从而为实际工程应用提供重要信息，如污水处理、机器冷却系统等。

再者，数值模拟法可以用来研究水中放电气泡的脉动的影响因素，如放电条件、温度、环境流动等，可以更清晰地了解水中放电气泡的脉动特性，从而实现对水中放电气泡脉动过程的精确模拟。

最后，通过数值模拟，研究人员可以更深入地了解水中放电气泡脉动过程，可以更准确的预测水中放电气泡的动态特性，为科学研究和实际应用提供充分的参考。

本研究采用了PISO算法来数值模拟水中放电气泡脉动过程，PISO算法是一种基于边界的迭代算法，可以有效地求解水力学方程。

首先，构建了体积三角剖分的模型，以确定水中放电气泡脉动过程数值模拟的计算范围；然后，利用质量守恒方程和能量守恒方程，建立压强场和温度场；接着，通过压强场和温度场对模型空间和时间变量进行分别求解；最后，在得到最终求解结果后，与实验结果进行比对，得到结论。

本文研究的水中放电气泡脉动过程的数值模拟，可以提供有用的实际应用信息，如污水处理和机械冷却系统。

此外，通过模拟，研究人员还可以对水中放电气泡的脉动的影响因素等特性进行系统的研究，为科学研究和实际应用提供有价值的信息。

总之，本文通过使用PISO方法，研究了水中放电气泡脉动过程的数值模拟。

研究结果表明，数值模拟方法可以准确模拟水中放电气泡脉动过程，为实际工程的设计和科学研究提供有用的信息。

气泡在水中上升运动的数值模拟朱仁庆;李晏丞;倪永燕;侯玲【摘要】基于流体体积函数(VOF)模型,借助Fluent软件,数值模拟了气泡在水中上升运动.考虑不同初始位置以及气泡大小对气泡在水中运动的影响,监测气泡在不同时刻的变形,分析了速度随时间的变化,并考察了气泡在不同密度比和粘度比的酒精流场和乙醚流场中运动.结果表明:直径大的气泡在上升过程中速度变化较大,上下表面速度差较大,大气泡较不稳定.气泡运动中,底部射流区域的速度先达到最大,然后降低,降低到一定程度会反弹.外部流体与气泡粘度比、密度比、表面张力系数对气泡运动有较大影响.【期刊名称】《江苏科技大学学报（自然科学版）》【年(卷),期】2010(024)005【总页数】7页(P417-422,451)【关键词】气泡;数值模拟;上升速度;流体体积函数法【作者】朱仁庆;李晏丞;倪永燕;侯玲【作者单位】江苏科技大学,船舶与海洋工程学院,江苏,镇江,212003;江苏科技大学,船舶与海洋工程学院,江苏,镇江,212003;江苏科技大学,船舶与海洋工程学院,江苏,镇江,212003;江苏科技大学,船舶与海洋工程学院,江苏,镇江,212003【正文语种】中文【中图分类】U661.1水中浮泡运动常见于船舶与海洋工程实际中,如:螺旋桨空化,水下爆炸引起的气泡,波浪破碎发生卷吸而引起的空泡等.气泡在流体中运动是强非线性的,运动时界面变形较大,因此气泡运动数值模拟越来越受国内外学者的关注,而气泡运动界面追踪是研究重点.目前已发展多种界面追踪技术并应用于气泡运动数值模拟,并且取得了一定的成果.界面模拟方法有:边界积分法[1-5]，VOF法[6-8]，Level Set法[9-10]，Lattice-Boltzmann法[11]，Front Tracking法[12-13].文献[14] 采用了Front Tracking法对粘性流体中气泡进行数值模拟,并分析了气泡上升运动速度随时间的变化规律.文献[15]采用边界积分法分析了二维气泡在无粘流体中上升运动.文献[16] 采用Lattice Boltzmann法对单个气泡运动,以及2个气泡和3个气泡运动进行了数值模拟,获得气泡运动的速度等值线图和速度随时间变化曲线图,取得一定的成果. 本文基于VOF技术中的PLIC界面重构方法,采用速度和压力耦合方法求解运动方程,对单个气泡在水中的运动进行了数值模拟,追踪了液界面变化,同时分析了不同气泡直径和气泡的初始高度对气泡上升时运动速度的影响.综合考虑了气泡在不同外流场中运动,分析了由密度比、粘性比及表面张力系数对气泡上升运动的影响.1 数值模型1.1 控制方程1) 考虑表面张力的动量方程(1)式中，v为速度矢量;ρ为流体密度;μ为粘性系数;p为压强;F为表面张力源项.2) 不可压缩流体连续性方程(2)3) 采用VOF法追踪界面的相函数输运方程(3)式中，aq为第q项体积分数.对于两相流方程(1)中ρ和μ由体积分数决定ρ=ρ1aq+(1-aq)ρ2(4)μ=μ1aq+(1-aq)μ2(5)式中，ρ1,ρ2，μ1,μ2分别为2种不同流体的密度和粘度.1.2 表面张力计算本文所用的表面张力模型是由文献[17]提出的连续表面力模型.采用CSF模型计算表面张力时,首先要计算界面的曲率和界面法向.定义aq为第q相体积分数,借助于体积分数分布,可得界面法向矢量n(6)表面曲率其中单位法向矢量(7)若一个单元只有两相,故(8)2 几何模型与计算条件为了消除固壁对气泡运动产生的影响[14],本文选取计算区域大于10D(D为直径),为0.1m×0.2m,通过Gambit软件划分网格,网格间距为5×10-4m,计算边界均为无滑移边界条件，计算几何模型见图1.气泡初始时刻在水中保持静止,初始压强和速度均为0,其形状为圆形(二维).气泡密度为 1.22kg/m3；粘度系数为1.789×10-5N·s/m2.水的密度为9.982×102kg/m3;粘度系数为1×10-3N·s/m2；表面张力系数为0.0728N·s/m2.图1 计算几何模型Fig.1 Computational geometry model描述气泡特性常用的无量纲参数主要有Morton数、Reynolds数、Weber数、密度比ρf/ρb和粘度比μf/μb,下标f和b分别代表外部流场和气泡.本文考虑的气泡运动场为低雷诺数的流场,其密度比为814.5，粘度比为55.9.3 结果分析与讨论3.1 单个气泡动力学特性本文模拟了直径D=10mm气泡在水中上升运动,观察气泡在运动过程中的变形.并对气泡的运动速度和压强变化进行监测.在表面张力作用下,保持了气泡内部压强和外部流体压强的平衡,保证了气泡稳定.同时由于表面张力作用在气泡表面,气泡的内部压强要大于外部流场压强.初始时刻气泡上下表面存在一个压力差,其下表面所受的压力梯度较大，在上下表面的压力差作用下气泡向上运动.在压力差与气泡表面发展出的涡片共同诱导出一个从下方推向气泡的射流.初期的射流并不能穿透气泡上表面,只是促使气泡底部向上凹陷.射流不断向气泡顶部发展,当射流长度达到一定程度,仍不能穿透气泡表面,射流开始向气泡横向发展,并形成马蹄状气泡[18].单个气泡在静止流场上升过程中,气泡的外形变化如图3～7所示,数值模拟结果与文献[19]实验结果一致(图2).图2 水中气泡上升运动(实验结果)Fig.2 The rising of bubble in the water (experimental results)气泡在水中运动,上表面的速度随时间逐渐增大,增大到一定程度后速度保持微小增幅,继续上升,直至与自由表面接触发生破碎(图3).气泡在t=0.01s时刻的速度等值线图,气泡仍保持圆形,此时气泡在界面附近处的速度U，V(单位：m/s)最大(图4).经过0.05s,气泡射流作用下下表面发生凹陷,形成月牙状(图4a)).气泡在底部y方向的速度V较大,在气泡凹陷形成的一对脚处,水平速度U比较大.在t=0.1s时(图5),底部射流发展为横向,抹平了气泡对脚,形成扁平帽子形状.此时气泡的各方向速度已经平稳.当气泡上升到自由表面处,由于考虑了表面张力作用,气泡顶部被自由表面的表面张力束缚,导致气泡上升受阻,气泡在浮力作用下继续上升,速度变小,在压力和自由表面张力共同作用下,气泡在水平方向发生拉伸，直至在t=0.26s时,气泡突破自由表面的束缚,发生破碎.而气泡下表面仍保持惯性继续上升,同时由于气泡破碎产生较大的压强梯度,导致自由表面上升(图6～7).图3 t=0.01s时相函数分布和U，V速度等值线Fig.3 Phase function distribution and contour of U and V velocity at t=0.01s图4 t=0.05s时相函数分布和U，V速度等值线Fig.4 Phase function distribution and contour of U and V velocity at t=0.05s图5 t=0.1s时相函数分布和U，V速度等值线Fig.5 Phase function distribution and contour of U and V velocity at t=0.1s图6 t=0.25s时相函数分布和U，V速度等值线Fig.6 Phase function distribution and contour of U and V velocity at t=0.25s图7 t=0.26s时相函数分布和U，V速度等值线Fig.7 Phase function distribution and contour of U and V velocity at t=0.26s3.2 气泡大小和初始位置对气泡运动的影响本文就直径为6，8mm气泡分别在相同的初始位置(指距自由水面高度,初始自由水面高为0.8m),考虑气泡上升运动过程中的速度随时间变化，针对气泡上下表面的速度进行分析.直径较大的气泡在水中运动时较难保持形状稳定,变形较大,而且上升速度和小气泡上升速度相比较大.直径小的气泡在水中容易保持其稳定形态,其发生变形时间比大气泡晚些.气泡与自由表面接触时,直径较大的气泡产生射流较强,导致自由液面抬升要高于小气泡.图8为气泡直径为6mm,在不同初始位置气泡运动速度随时间的变化曲线.图8a)，b)初始位置分别为0.03，0.05m.初始时刻气泡底部在射流作用下速度(Vbot)在很短时刻内达到一个峰值,在运动过程中气泡下表面速度逐渐减小,此时上表面速度逐渐增加.在t=0.05s时,上下表面速度近似相平衡,此时气泡上下表面速度保持动态平衡,射流发展为气泡横向,此时气泡的形状近似稳定.图8 D=6mm气泡在不同初始位置时的速度变化曲线Fig.8 Bubble velocity versus time when initial position is 0.03 and 0.05m(D=6mm)气泡顶部的速度(Vtop)在初始时刻也有较大的增幅,在t=0.05s以后增幅减小,上下表面速度近似相等.保持一定的振幅，气泡接近自由表面时,由于自由表面在表面张力的作用下对气泡上升运动起到阻碍作用,在t=0.26s之后气泡上下表面速度都发生降低,直至气泡破裂.在气泡破裂时,上表面速度在压力梯度作用下突然增大,随后速度降低(图8b)).图9为直径8mm,初始位置分别为0.03，0.05m时气泡速度随时间变化曲线.在初始时刻，直径较大的气泡底部产生射流速度要比直径小的气泡大,而且气泡上下表面的速度随时间变化，上下振荡的幅度比直径为6mm气泡振荡幅度要大.气泡下表面产生射流导致气泡下表面速度发生周期性变化,呈衰减趋势(图9a)).由分析可知,不同初始位置对相同初始直径的气泡运动速度影响不是很大.直径大的气泡在初期产生的射流强度要大于小气泡产生的射流强度.小气泡在水中运动比大气泡要稳定.大气泡的上下速度振荡较大,容易产生较大变形,所以大气泡在水中运动易破裂.图9 D=8mm,初始位置为0.03和0.05 m时的速度变化曲线Fig.9 Bubble velocity versus time when initial position is 0.03 and 0.05m (D=8mm)3.3 外流场发生变化时对气泡运动的影响气泡在流体中运动时受到表面张力、粘性力、浮力、重力和压力梯度力等作用.为了考察各个力对气泡运动规律的影响,分别考虑了气泡在外流场为酒精和乙醚时的上升运动规律.水、酒精和乙醚参数见表1.表1 外流场的流体参数Table 1 Parameters of the ambient liquid流体密度/(kg·m-3)粘度/(N·s·m-2)密度比(ρf/ρb)粘度比(μf/μb)表面张力系数/(N·m-1)水998.20.001814.555.90.0728酒精7900.0012644.967.10.023乙醚8040.00395656.3220.80.0165图10为直径8mm气泡在酒精中上升运动时,初始射流导致速度达到一个峰值,随后速度逐渐衰减,从峰值到最小值周期为0.025s.气泡在水中上升时,底部射流导致达到峰值的速度衰减到最小值周期为0.05s.分析可知,密度比减小,气泡速度衰减的周期减小.由图10a)，b)可知,密度比相差不大情况下,气泡在流场中上下表面速度衰减趋势相同,在粘度比较大的乙醚流体中,气泡的上表面达到一定速度后保持恒定速度上升.表面张力系数较小时,气泡初期产生的射流速度较大,同时气泡运动靠近自由液面时,由于表面张力系数作用,对气泡的运动影响减小，速度趋势趋于平缓.4 结论1) 采用VOF法获得了单个气泡在水中运动的时刻历程,追踪气泡运动时界面变化,较清晰反应了气泡界面运动的规律,分析了气泡上升运动对自由液面影响.2) 通过分析单个气泡在自由液面水中上升运动时的速度场,得到气泡运动速度分布图,气泡界面处的底部速度和气泡在射流凹陷处速度最大.图10 D=8mm,初始位置为0.03 m，外流场分别为酒精和乙醚时气泡速度变化曲线Fig.10 Bubble velocity versus time when initial position is 0.03 m,D=8mm, external flow field are alcohol and ether3) 通过比较直径不同和初始高度不同的气泡在水中的运动规律,直径大的气泡运动时较易产生大的变形,初始高度越大的气泡产生的射流速度越大.4) 不同外部流场的粘度比、密度比、表面张力系数对气泡运动有较大影响,密度比对气泡底部射流有影响,密度比越大影响就越明显.粘度比对气泡上升过程保持稳定有影响,粘度系数较大,气泡的运动速度越趋近于一个恒定值.表面张力系数对气泡产生射流速度有影响,表面张力系数越大,对射流影响越大;同时气泡靠近自由液面时,表面张力对气泡上升运动有阻碍作用.参考文献(References)[1] Lorstad D, Francois M, Shyy W, et al. Assessment of volume of fluid and immersed boundary methods for droplet calculations[J]. International Journal for Number, Methods in Fluids, 2004,46(2):109-125.[2] Tryggvason G, Bunner B, Esmaeeli A, et al. A front-tracking method for the computations of multiphase flow[J]. Journal of Computational Physics, 2001,169:708-759.[3] 宗智,何亮,张恩国.水中结构物附近三维爆炸气泡的数值模拟[J].水动力学研究与进展： A辑,2007,22 (5):592-602.Zong Zhi, He Liang, Zhang Enguo. Numerical simulation of a three-dimensional underwater explosion bubble near a structure[J]. Journal of Hydrodynamics:Ser A, 2007, 22 (5) : 592-602. (in Chinese)[4] 张阿漫,姚熊亮.近自由面水下爆炸气泡的运动规律研究[J].物理学报,2008,57(1):339-352.Zhang Aman, Yao Xiongliang. The law of the underwater explosion bubble motion near free surface[J]. Acta Physica Sinica, 2008,57(1):339-352. (in Chinese)[5] 张阿漫,姚熊亮.单个三维气泡的动力学特性研究[J].应用力学学报,2008,25(1):107-111.Zhang Aman, Yao Xiongliang. Dynamics for single three-dimensional bubble[J]. Chinese Journal of Applied Mechanics, 2008,25(1):107-111. (in Chinese)[6] Lorstad D. Numerical modeling of deforming bubble transport related to cavitating hydraulic turbines[D].Sweden: Department of Heat and Power Engineering, Lund University, 2003.[7] Puckett E G, Almgren A S, Bell J B, et al. A high-order projection method for tracking fluid interfaces in variable density incompressible flows[J]. Journal of Computional Physics, 1997,130:269-282.[8] 端木玉,朱仁庆.流体体积方程的求解方法[J].江苏科技大学学报：自然科学版，2007,21(2):10-15.Duan Muyu, Zhu Renqing. Method for solving fluid volwm equation[J]. Journal of Jiangsu University of Science and Technology: Natural Science Edition, 2007,21(2):10-15. (in Chinese)[9] Osher S, Fedkiw R P. Level set methods: an overview and some recent results[J]. Journal of Computional Physics, 2001,169:463-502.[10] Sussman M, Smereka P, Osher S. A level set approach for computing solutions to incompressible two-phase flow[J]. Journal of Computional Physics, 1994,114:146-159.[11] Watanabe T, Ebihara K. Numerical simulation of coalescence and breakup of rising droplets[J]. Computational Fluids, 2003,32:823-834. [12] 陈斌， Kawamura T, Kodama Y. 静止水中单个上升气泡的直接数值模拟[J].工程热物理学报, 2005,26(6):980-982.Chen Bin, Kawamura T, Kodama Y. Direct numerical simulations of a singlebubble rising in still water[J]. Journal of Engineering Thermophysics, 2005, 26(6)：980-982. (in Chinese)[13] 陈斌.倾斜壁面附近上升气泡的直接数值模拟[J]. 工程热物理学报,2007,28(26)：965-967.Chen Bin. Direct numerical simulation of a single bubble rising along an inclination wall[J]. Journal of Engineering Thermophysics, 2007,28(26)：965-967. (in Chinese)[14] Hua Jinsong, Lou Jing. Numerical simulation of bubble rising in viscous liquid [J]. Journal of Computational Physics， 2007,222:765-769. [15] Robinson P B, Boulton-Stone J M, Blake J R. Application of boundary integral method to the interaction of rising two-dimensional, deformable gas bubbles [J]. Journal of Engineering Mathematics, 1995,29:393-412. [16] Gupta A, Kumar R. Lattice Boltzmann simulation to study multiple bubble dynamics[J]. International Journal of Heat and Mass Transfer, 2008,51:5192-5203.[17] Brackbill J U, Kothe D B, Zemach C. A continuum method for modeling surface tension[J]. Journal of Computational Physics, 1992,100:335-354. [18] 张淑君,吴锤结.气泡之间相互作用的数值模拟[J].水动力学研究与进展, 2008,23(6):683-684.Zhang Shujun, Wu Chuijie. Numerical simulation of the interactions between two three-dimensional deformable bubbles[J]. Chinese Journal of Hydrodynamics, 2008, 23(6 ):683-684. (in Chinese)[19] Walters J K, Davidson J F. The initial motion of a gas bubble formed in an inviscid liquid, Part 1: The two dimensional bubble [J]. Journal of FluidMechanics, 1962,12:408-416.。