船舱油液泄漏入水模拟方法

顾乾乾，庄学强，林添金，王永坚，高亚丽

(1.集美大学轮机工程学院，福建 厦门 361021；2.福建省船舶与海洋工程重点实验室，福建 厦门 361021)

0 引言

据国际油轮船东污染组织(ITOPF)统计，1970—2020年全球范围内发生7 t及以上船舶溢油事故1845 起，总溢油量约586万t，低于7 t的溢油事故更是不计其数[1]。海上溢油事故频发不仅使自然环境、生态系统受到严重损害，经济蒙受损失，而且还严重危害人体健康[2-3]。为了对船舶溢油进行防控，世界上很多国家在防污法律法规的制定，溢油防污技术的开发，以及溢油模拟预测等方面做了大量的工作。其中，溢油模拟预测能够预测油类污染物在海上的行为过程，为油污清理以及生态环境损害评估提供科学支撑信息[4-5]。溢油模拟预测包括溢油油液泄漏入水过程模拟、油膜漂移扩散模拟、溢油风化过程模拟、溢油对生态损害模拟等[6]。油膜漂移扩散、溢油风化等方面的模拟研究较多，但是在溢油油液泄漏入水过程模拟的研究较少，进展也十分缓慢。究其原因主要是该过程的准确模拟涉及到复杂船体结构建模和事故发生时的复杂海况条件耦合运算，模拟起来繁琐困难。目前业界普遍做法是对整个泄漏过程情景参数做不同程度的简化，但这势必影响泄漏过程模拟结果的准确性。船舶溢油事故中油液泄漏入水过程模拟结果是后续的溢油漂移、扩散、风化行为归宿模拟计算和溢油环境损害预测计算的重要输入源项，源项的准确性将会直接影响后续模拟预测的准确性，进而影响应急措施的科学性和损害索赔的客观性。因此，溢油油液泄漏入水过程模拟研究对整个溢油模拟以至于溢油应急技术的研究都有十分重要的意义。

1 油液泄漏入水行为模拟方法

1.1 基于溢油量估算的泄漏过程极简化模拟

溢油量是衡量船舶溢油事故严重程度最重要的参数之一，人们经常把溢油泄漏行为描述简化处理为溢油量估算，这是一种只考虑泄漏结果，不考虑过程的极简化模拟方法。目前关于泄漏量估算方法有HFO(heavy fuel oil)变换量估算法、油膜监测估算法、光学估计法等等。

HFO变换量估算法是一种基础的估算方法，依据船舶日常油量记录数据，HFO总库存量变化与HFO设备运转消耗量的减少来推算出溢油量，其计算公式为[7]：

G=Z-C-R。

(1)

式中:G为船舶破损后估计的溢油量；Z为破损前船舶总储存油量；C为船舶记录总油量时刻至船舶破损时这期间船舶的总耗油量；R为船舶破损后油舱剩余油量。

溢油油膜监测估算法是在船舶污染事故发生后，通过现场观测、卫星遥感等多种监测技术获取油膜面积[8-9]，根据油膜颜色与油膜厚度的对应关系获得油膜厚度，其海面溢油量计算公式[10]为:

(2)

式中:G为溢油估计量；Si为第i种颜色的油膜面积；Hi为第i种颜色油膜的厚度；ρ为油品密度；n为油膜颜色分区数。

光学估计法是利用粒子图像测速技术，通过对2个连续的视频进行观察，分析视频帧之间的距离，得出溢油的瞬时速度或是平均速度，再与破损面积、油品密度和溢油时间相乘便得到溢油估计量[11-12]，其计算公式为：

G=S×V×ρ×t。

(3)

式中:G为溢油估计量；S为破损处横截面积；V为溢油流速；ρ为油品密度,t为溢油泄漏时间。

4) 自磨刃强化处理工艺对降低功耗、提高刀片使用寿命和节约钢材有重要的意义，不远的将来一定会得到推广。

HFO变换量估算法、油膜监测估算法仅仅关注了泄漏结果，只得到一个静态的泄漏量数据，简单的将溢油过程简化为连续定常泄漏或者瞬时泄漏过程，都没有考虑泄漏过程中液位、溢油速度等参数的动态变化对溢油泄漏的影响。

1.2 基于伯努利方程的小孔射流模拟方法

孔洞射流模式是随液面高度和泄漏时间变化的泄漏源强模型，更加接近于船舶溢油实际泄漏行为[13-14]。船舶发生破损，穿孔处与外部大气或是海水相通，造成液体流出，破孔位置在水线上和水线下造成的不同泄漏情况，如图1、图2所示。破孔位置在水线上的泄漏行为简单，主要影响因素是液体重力[15]。图1为油船侧壁破孔在水线之上横截面示意图。其中：AB表示泄漏口流通面积；At为油舱的上表面积；h为从破裂中心线到油舱油的自由表面的垂直距离；Pg为舱内气体空间压力；用时间标度td来表征油舱破裂后溢油的持续时间。

(4)

破孔位置在水线下的泄漏行为复杂，如图2所示。涉及到的影响因素较多，采用孔洞射流模型，该模型由伯努力方程(Bernoulli equation)推导而得，油液泄漏率为：

Q=CDABρL[(2(Pg-Pout)/ρL)+2gh]1/2。

(5)

式中：CD为泄漏口的流量系数，m2；Pg为油舱上方的压力，Pa；Pout为破损口舷外压力(根据破损口位置，可能是大气压或者是水压)，Pa；ρL为油液密度，kg/m3。

式(5)是一个理想化的方程，仅仅考虑了泄漏过程中随着油液下降造成的压力变化，这与泄漏发生时的复杂海况条件不相符合。船体运动或是海浪波动的影响会打破油舱内外压力平衡，从而发生进一步泄漏：1)初始阶段，持续到通过断裂处建立静水压力平衡；2)由于船体运动改变静水压力平衡而产生泄漏[16-17]。

(6)

(7)

第二阶段破孔处的瞬时速度为：

(8)

式中：Ll为舱内液位；LB为从船底至船体破孔处距离；Vg为舱内气体空间体积；ρl为舱内液体的密度；ρw为水的密度；uB为破孔处液体的流速；CD为流量系数；g为引力常数；a是正弦振荡的振幅；T为运动周期。

当水线下出现泄漏时，将溢油分为两个阶段，第一阶段由于油品和海水密度以及液位不同，破损孔油品侧压高于海水侧压力，存在压力差造成的泄漏，这个阶段持续到油舱内容物与周围环境之间通过断裂处建立流体静水平衡；第二阶段为船舶运动或是波浪造成的压力波动而导致的泄漏，该阶段有水的流入以及水和油液的流出，只要船舶在运动，出流和入流就会一直持续，当破孔两侧被水覆盖后出流和入流都将是水[18-21]。船舶破损发生泄漏，泄漏孔的形状多是不规则的，且射流流量系数基本上是依靠经验模型得来，可见，将普适性系数应用于特定泄漏情况，计算结果与实际泄漏将产生误差，精确度较低。

1.3 基于CFD的精确模拟方法

计算流体力学模型(CFD 模型)是基于纳维-斯托克斯方程组(N-S 方程)的三维流体守恒定律模型[22]，能克服射流模型中流量系数难于确定、复杂船体结构无法考虑，以及事故、风、潮、流情景无法叠加等一系列问题，可以提供油流从油舱通过破损口泄漏入水详细流动行为的描述[23]。其基本控制方程包括连续性方程、动量方程、能量方程等[5]。

1)连续性方程(即质量守恒方程)为：

(∂ρ)/(∂t)+(∂(ρu))/(∂x)+(∂(ρv))/(∂y)+(∂(ρw))/(∂z)=0。

(9)

其中:u、v、w分别为x、y、z方向上的速度分量，m/s；ρ为流体密度，kg/m3。

2)动量守恒方程为：

(10)

(11)

(12)

其中：μ为流体的动力粘度，N·s/m2；P为作用在流体微元体上的压力，Pa；Su、Sv、Sw为动量守恒方程的广义源项。

3)能量守恒方程为：

(13)

其中：ST为流体的内热源及由于粘性作用流体机械能转换为热能的部分；Cp为比热容；T为温度；k为流体的传热系数。

CFD模型比其他方法能更加直观、准确地模拟出每时刻溢油的动态，能叠加风、浪、流、船舶动态，船体结构等影响因素，更真实地模拟溢油泄漏情景，是该领域研究的发展方向。

2 实验与计算方法对比验证

为了对比验证上述计算方法的准确性与适用性，本文引入双壳油舱模型，开展多组单壳和双壳油舱模型的泄漏模拟对比实验。实验装置包括三维立体油舱破孔水槽、计时器、高清照相机和微量电子秤。出于对环境污染以及安全的考虑，实验选用自来水代替油品，并用黑色染色剂对水体进行染色，方便观察记录。实验模拟温度设定为室温，水面油舱破孔泄漏持续时间以及不同时刻的泄漏状态采用XL—013A专业秒表计时器和高清照相机记录。

2.1 单壳结构泄漏验证

极简化模拟和小孔射流计算方法是数值计算，通过时间的变化计算泄漏速度、泄漏量的数值变化，输出的是数值结果，无法给出泄漏流动的动态过程。在单壳中，开始阶段液体均在压力的作用下从容器内快速喷出，形成一条弯曲的喷射水柱。随着液面不断降低，压力逐渐减小，水柱喷出距离随之减小，最后逐渐沿壁面向下流动，持续到液体液位低于破损孔，泄漏停止。液体流动状态随时间变化关系如图3所示，左侧为实验验证，右侧为CFD仿真模拟。这说明可模拟仿真结果和实验观测结果在时间点和液体溢出状态上基本一致，且相较于其他只能输出计算数据的计算方法，CFD仿真计算模拟能够直接输出泄漏动态过程，具有视觉上的直观性。

进一步对比分析实验验证数据和三种模拟计算结果，得到单壳油舱水上泄漏情况的极简化模拟、小孔射流计算模拟、CFD仿真模拟和实验观测数据曲线对比图，如图4所示。从图4中可以看出，小孔射流、CFD仿真模拟和实验观测结果流速趋势基本一致，极简化模拟是依据实验泄漏量以及泄漏持续时间得出的泄漏平均速度，只能作为参考，与实际泄漏情况相差较大。

实验泄漏过程：开始阶段，由于破损瞬间液体内部存在黏性以及液体和壁面存在摩擦，液体流动受到影响；泄漏持续阶段，在惯性的作用下液体流动平稳，惯性占据主导作用；结束阶段，惯性作用减弱，液体内部粘性以及液体和壁面的摩擦对液体流动影响变大，泄漏速度快速减小，液体沿壁面向下流动。小孔射流模拟计算结果与实验观测结果在开始阶段(0～5 s)和结束阶段(15～35 s)相差较大，仅在泄漏持续阶段(5～15 s)结果较接近，这是由于小孔射流模拟计算主要依据破损孔处液体高度，忽略了液体黏性和摩擦的影响作用，将泄漏计算整个过程视为稳定变化过程，泄漏系数是固定值，所以小孔射流模拟计算在液体粘性和摩擦影响较大的开始和结束阶段结果偏差较大，惯性作用影响较大的泄漏持续阶段结果偏差较小，造成了计算结果局部阶段偏差小，整体泄漏过程偏差大，计算精度不高。CFD仿真模拟结果与实验观测结果在泄漏发生的整个过程数据基本一致，数据贴合度较高。

2.2 双壳结构泄漏验证

相较于单壳油舱模型，双壳油舱模型外壳破孔比内壳破孔大，破孔中心处于同一水平线上。双壳结构对于泄漏速度没有影响，与单壳结构泄漏速度一致，但是对泄漏流动状态以及泄漏总量有影响。双壳油舱液体泄漏状态如图5所示。从图5可以看出，双层壳油舱破损舱内液体在压力差的作用下瞬间以柱状喷出，随着压力差的减小，喷出距离逐渐缩短，缩短到一定程度后，一部分液体流入双层壳，一部分液体喷出壳外，喷出距离进一步减小，液体将不再喷出外壳，全部流向双壳内。油舱外壳破孔处内壁下沿Water体积分数图(采集处Water的体积占采集处液体体积的比重)如图6所示。从图6可以看出，在泄漏发生第12秒时，液体开始流入双层壳底部，第15秒时，液体全部流入双层壳内。图7为双壳油舱总泄漏量对比图。从图7可以看出，极简化模拟为均匀泄漏，泄漏积累量均匀增加，由于无法计算流入双层壳空间内的液体量，最终泄漏量高于实验结果；小孔射流计算方法计算泄漏积累量曲线前12 s趋势与实验结果基本一致，后续泄漏由于液体流入双壳空间，小孔射流计算方法无法得出双壳结构截留液体的时刻以及双壳结构截留液体的总量，从而导致泄漏总量计算偏大，偏离实际结果；CFD模拟计算方法泄漏积累量在0～12 s时以一定比例逐步增加，第12秒由于液体开始流入双壳空间内，曲线出现大幅转折，曲线变缓，在第15秒喷出液体全部流入双壳空间内，泄漏量不再增加。实验结果表明，CFD模拟计算与实验观测泄漏量曲线趋势和数据贴合度基本一致，能准确模拟出双壳结构对整个泄漏过程以及泄漏参数的影响变化，具有简洁、高效和精确度高的特点。

综上所述，CFD模拟计算方法能直观显示液体流动状态，使用更加方便、简洁；CFD模拟计算结果更加贴合实验数据，精确度高于其他计算方法；CFD模拟计算方法能够计算复杂油舱结构，甚至叠加复杂海况下的泄漏情况，这也是其他计算方法所不具备的。因此CFD计算模拟方法具有更强的实用性，以及更精确、直观的计算结果。

3 结论

本文分析了几种常用的船舶溢油泄漏计算方法，并进一步引入双壳船舶模型，开展实验分析对比验证。极简化模拟方法适用于溢油量估算，计算精度低，与实际相差较大；小孔射流计算方法能够计算简单结构下溢油泄漏瞬时速度、泄漏体积量等参数变化，但是无法计算双壳结构对溢油泄漏的影响作用以及参数变化；CFD模拟计算方法能够直观显示出泄漏动态过程，计算曲线与实际泄漏曲线基本一致，计算精度高，并且能够计算复杂船体结构的泄漏情况。因此，CFD模拟计算方法能耦合复杂船体结构，更准确地模拟油舱油液泄漏入水行为，是该方面研究的发展方向。