跨江输油管道成品油泄漏扩散影响分析

杨大慎，陈 强，谢 成，陈晓华，姜红涛，赵彦琳

(1.国家石油天然气管网集团有限公司华南分公司，广东 广州 510620；2.中国石油大学(北京)机械与储运工程学院，北京 102249)

近年来，随着社会经济的发展，对能源的需求不断增长，管道运输在成品油远距离运输中具有不可替代的作用。运输管道的油品泄漏会对企业经济和生态环境带来巨大的损失和破坏[1]。2009—2019年海底输油管道事故统计信息显示，我国共发生100余起泄漏事件，仅2001年发生的东海油田管道泄漏事故造成的直接经济损失就超过2 000万元。因此，油气运输管道的安全运行越来越受到人们的重视，而开展泄漏油品的轨迹和扩散范围的研究可为事故的紧急处理提供及时有效的参考。

针对油品泄漏扩散特征，国内外学者进行了一些相关研究，例如BEMPORD[2]考虑水流的影响，对分层流工况下圆孔的浮射流轨迹进行了研究，建立了稳定分层流动液体中水下浮力圆形射流的数值模拟模型，用于分析稳定或中性分层情况下不同流动条件时潜水浮力射流的行为。JOHVANSEN[3]提出了一种基于拉格朗日概念的海底井喷模型，研究了分层水柱中上升气泡在水中的溶解，得到模拟结果与野外试验观测结果较吻合的结论。UDITHA等[4]研究了气泡尺寸分布在深水天然气或石油运输中的作用机制，通过数值模拟和实验对比证明了气泡破裂和聚结在泄漏模型中的重要性。KANDEL等[5]通过实验合成水溶性阴离子分散剂，研究了不同分散剂对抑制成品油污染扩散的影响。高龙骧等[6]改进了经典溢油模型，建立了内河溢油动态模型，实现了对溢油过程的仿真与预测。刘冰等[7]针对不同的河流水体和河流环境溢流的特点，结合基本控制方程与水动力数学模型，改进Fay模型、简化Navy漂移模型，建立了新的河流溢油环境影响评价模型。祝慧娜等[8]结合水环境的脆弱性和对人体的危害，建立了河流水环污染风险模糊综合评价模型。

综上可以看出，输油管道成品油泄漏研究大多数是从实验研究和数值模拟两个方面开展，而且主要针对海水管道、井下泄漏等油品污染的扩散情况进行研究，较少涉及河流跨越段输油管道泄漏后成品油在河流流域扩散速度与扩散范围的研究。因此，笔者基于北盘江流域成品油管道跨越段，在不同工况下对不同泄漏位置、不同泄漏时间和不同泄漏方式下成品油在江中的扩散速度与扩散范围进行数值模拟，获得不同影响因素作用下油品扩散的轨迹、速度和污染范围，以期为管道运输安全、工程抢险救灾、水环境保护提供技术参考。

1 研究原理及方法

1.1 研究方法概述

研究区选取某一河段流域，通过Getdata软件获得河岸坐标点，再导入矢量作图软件SolidWorks中构建几何模型，并设置油品污染物的泄漏位置；通过CFD前处理软件对构建的几何结构模型进行网格划分，生成结构性网格；根据实际问题进行不同程度的归一理想化，调整各物性参数、边界条件，建立污染扩散模型，对计算结果进行对比分析，整理出不同条件下该流域的污染扩散规律。具体研究思路如图1所示。

图1 研究思路

1.2 网格划分

采用ICEM 软件来生成网格。考虑到结构化网格的优点，如容易实现区域的边界拟合，适于流体和表面应力集中等方面的计算，网格生成速度快、质量好，数据结构简单等，故将计算域模型画成结构化网格，并在计算过程中不断加密网格，找出网格无关化解。

1.3 离散方法

采用ANSYS R19.2进行数值模拟，选择有限体积法求解，采用QUICK离散格式，压力和速度耦合求解采用SIMPLE(semi-implicit method for pressure linked equations)算法[9]。基于SIMPLE算法，以压力为基本求解变量，通过密度延迟修正方法将SIMPLE算法扩展到可压缩流动的计算[10]。使用非交错网格使程序更加简洁通用，运用动量插值法消除由非交错网格所产生的非物理性压力振荡现象[11]。使用标准κ-ε湍流模型对粘性流场进行模拟，对于固体壁面附近的区域采用壁面函数法处理[12]。SIMPLE算法是一种压力修正法，采用预测-修正计算方法，加快了每个时间步长的收敛速度，效率较高。动量空间离散采用二阶迎风格式，体积分数的空间离散应用Geo-Reconstruct格式，湍动能和比耗散率的空间离散均采用一阶迎风格式[13]。

1.4 湍流模型控制方程

按照假设，流体为不可压缩的牛顿流体，其基本控制方程采用连续性方程与不可压缩流体的雷诺N-S方程:

▽·u=0

(1)

(2)

式中：u为流体速度；ρ为流体密度；p为流体压力；v为流体运动黏度；fi为作用在流体上的外力。

1.5 VOF多相流模型

VOF多相流模型基于流体体积函数进行计算[14]，即利用许多单个网格中一种流体体积与网格体积的比值f来进行界面追踪[15]。定义流体体积分数fq为网格内第q相流体所占体积与该网格总体积之比[16]。若fq=1，则表示网格内均为第q相流体，不含其他流体；若fq=0，则表示网格内不存在第q相流体；若0

(3)

∑fq=1

(4)

在二维计算域中，假设fo、fw分别为成品油和水在一个网格中所占的体积分数，则fo、fw的计算公式为：

fo=Vo/V

(5)

fw=Vw/V

(6)

其中，V、Vo、Vw分别为网格体积、网格中成品油的体积和网格中水的体积。

结合式(3)～式(6)，可推导出计算域内成品油和水两相混合流动过程的连续性方程：

(7)

(8)

此外，在VOF方法中，多相流的物性参数φ由计算区域内各相的物性参数和体积分数共同决定，即：

φ=∑φqfq

(9)

1.6 模型假设

在数值模拟中，建立模型需要针对现实情况做出相应的简化和假设，忽略次要因素，考虑主要影响，以提高计算效率。

(1)所建立的模型只包含两种流体，油品和水，且二者在流动过程中完全不相溶，二者均不可压缩，密度、粘度等物理性质均保持不变；泄漏油品以油膜的形式随河流表面向下游扩散迁移。

(2)在管道发生完全断裂时，上下游的阀门关闭后，管道中的残余油品经0.5 h就完全流入河流中，在山坡上无残留油品，河流两侧泄漏口泄漏时完全沿着山坡流入河道，油品泄漏时风速保持恒定。

(3)在管道发生小孔泄漏时，从上下游阀门关闭到封堵成功，管道内压力和小孔泄漏流量基本保持不变。

(4)油品和水不发生任何化学反应，油水混合物挥发量不计。

2 北盘江计算模型

2.1 北盘江概述

北盘江为典型的山区性河流，洪枯水季节分明。每年6月至次年5月为枯水期。洪水期径流主要由降雨形成，水量约占全年的80%，水位暴涨；枯水期径流主要由地下水补给，流量约占全年的20%，水位相对平稳。北盘江多年平均流量为309 m3/s，多年平均径流量为98.18亿m3。由于北盘江跨越段所处地区，山区、河谷交错，成品油管道一旦发生泄漏，水流弥散作用将造成泄漏油品迅速向下游迁移，短时间内即对下游造成大面积污染，给当地的生态、渔业、旅游、农业等造成严重影响。因此，对北盘江跨越段成品油管道的泄漏扩散特征进行研究是非常必要的。

研究区选取放马坪阀室至小水井阀室的北盘江跨越管段，采用单跨悬索跨越方案，管道通过支架固定到桥上，管道总长203 m，共设置支架21个。第一个支架距离桥端3 m，支架之间距离大约为10 m。现场调研发现，距离桥端33 m处的支架因雨水进入较多，腐蚀相对严重，故将其设置为典型泄漏计算点。

2.2 流域计算模型的建立

通过CAD建模，选取有限的区域段建立相应的计算域。以北盘江跨越段为起点，以直线距离20.2 km的位置为终点。东西方向(x方向)距离跨度约为9.4 km，南北方向(y方向)距离跨度约为16.0 km。对北盘江跨越段下游河两岸进行数据分析处理，进一步细化所需河道信息，顺序提取河岸两侧轮廓，获得河道坐标数据。将河道坐标数据点导入矢量建模软件中，并将北盘江跨越段在河流中的起点设为坐标原点(0,0)。通过绘图软件将地理位置坐标进行离散化，获得北盘江的下游河道计算模型的坐标参数。北盘江河岸两侧轮廓及河道坐标数据点分布如图2所示。

图2 北盘江河岸两侧轮廓及河道坐标数据点分布

以图2的边界数据为基础，进行计算域的建模。将河道两岸的散点坐标导入icem软件，得到河流两岸边界的轨迹线，在起始位置和出口位置形成闭合区域，再对封闭曲线进行内部填充，从而得到二维的计算域模型图，计算域y方向长度为12 014 m，x方向长度为9 634 m。入口边界条件设为velocity-inlet；出口边界条件设为pressure-outlet；两侧河岸的边界条件设为wall，无滑移边界条件；成品油进入河道区域的边界长度设为10 m，边界条件设为mass-flow-inlet。

在模拟计算中，北盘江跨越段的河道信息及相关流域的水文资料均被考虑进模型中，由于北盘江是季节性河流，河水的流速和流量均受季节影响，故笔者将北盘江分为枯水期、平水期和丰水期3种，以探讨不同时期泄漏溢油的污染状况和扩散范围。北盘江跨越段不同时期水流流量和流速的平均值如表1所示。

表1 北盘江跨越段不同时期水流流量与流速的平均值

同时，将泄漏工况按照破裂严重性分为两大类：管道发生完全断裂和局部小孔泄漏。管道发生完全断裂可能是因为山区意外地质灾害、人为破坏或者年久失修管道锈蚀，不同位置发生的完全断裂所造成的溢油扩散范围可能不同，故将泄漏点分为管道南侧、管道中间和管道北侧3处，以探讨不同泄漏点位置对溢油扩散情况的影响。当管道破损不大，仅发生小孔泄漏时，成品油泄漏量相对于管道完全断裂来说比较小，因此仅考虑成品油在跨越段管道中间泄漏的情况，成品油从小孔泄漏后直接落入北盘江中。

2.3 泄漏量计算

北盘江跨越段的管道长度L为203 m，管道直径D为406.4 mm，壁厚为12.7 mm，设计压力为10 MPa。管道材质为X60钢，制管方式为LSAW，输送介质为97#、93#汽油和0#柴油等成品油。跨越段桥梁距离水面的垂直距离约为18 m。泄漏检测系统可以在成品油泄漏发生20 s内做出正确的泄漏报警，阀门由全开至全关需20 s。假设该段成品油最大输量Q年为1 000万t/年，成品油密度ρ为850 kg/m3。

2.3.1 管道完全断裂的泄漏量

当管道发生完全断裂事故时，从事故开始到阀门关闭这段时间溢油量Q1的估算公式为：

(20+20)≈14.92 m3

(10)

式中：T为一年的周期长度；t为泄漏时间。

紧急关闭阀门后，管道残余量一般为上下游阀室之间所剩余的成品油量。估算时按照泄漏的最大情况来分析，所以留在管道内的成品油泄漏量Q2为：

(11)

则最大成品油泄漏量Q为：

Q=Q1+Q2=41.24 m3

(12)

2.3.2 管道小孔泄漏的泄漏量

当管道发生小孔泄漏，设定小孔的孔径d为10%D(按照整理的有关数据，管道发生泄漏的最大小孔泄漏尺寸为10%D)。小孔的通流长度即管道的壁面厚度l为12.7 mm，则l/d=12.7/40.6=0.31<0.5，符合流体力学中关于薄壁小孔的定义要求。

通过流体力学中的伯努利方程，可得小孔泄漏的连续性方程为：

(13)

式中：ρ为溢油的密度；v1、v2分别为泄漏小孔内侧截面和外侧截面处溢油的速度；p1、p2分别为小孔内侧截面和外侧截面处的压力，p1取管道压力，即10 MPa，p2取环境压力，即101 325 Pa；ζ为动能修正系数；Cv为小孔流速系数。

则管道小孔泄漏的成品油流量q为：

(15)

式中：A为泄漏小孔外侧截面的面积；Cc为小孔收缩系数，即最小通流截面的面积与泄漏小孔外侧截面面积的比值；Cd为小孔流量系数，由实验测得。完全收缩时，流体在泄漏孔处的雷诺数较大，Cd的取值范围为0.61～0.62；不完全收缩时，Cd的取值范围为0.7～0.8。本次计算为完全收缩，故Cd的取值为0.62。经计算，小孔泄漏的成品油流量q为0.11 m3/s。

当管道突发小孔泄漏事故时，会经历泄漏警报和阀门关闭过程。虽然阀门关闭后泄漏的成品油流量相对比较小，但是泄漏持续时间较长，根据经验推算紧急堵漏的时间为2 h。则从泄漏发生至阀门关闭再到紧急封堵泄漏这段时间的溢油量Q′为：

Q′=q×(20+20+7 200)=796.4 m3

(16)

3 模拟结果分析

参考北盘江跨域段管道的数据与水文资料确定模拟计算中所需参数，根据所提供的参数条件设置计算条件，如北盘江江水的密度为998.2 kg/m3，动力黏度为1.003 mPa·s；成品油密度设为850.0 kg/m3，动力黏度为0.08 Pa·s。

3.1 管道完全断裂扩散结果分析

工况①～工况③分别为北盘江枯水期管道南侧、管道中间、管道北侧发生完全断裂时成品油的泄漏工况，枯水期3个工况不同时刻泄漏成品油扩散的数值模拟结果分别如图3～图5所示。计算结果显示，泄漏0.5 h后，成品油完全流进入北盘江，泄漏成品油刚开始集中在河流中央，成块状聚集。泄漏发生2 h后，油品逐渐被冲散，块状主体主要集中在跨越段下游4.1～4.4 km；泄漏4 h后，成品油扩散范围更大，油水混合物还附着在岸边，主要集中在跨越段下游7.3～8.2 km。

图3 工况①不同时刻的成品油扩散情况

图4 工况②不同时刻的成品油扩散情况

图5 工况③不同时刻的成品油扩散情况

工况④～工况⑥分别为北盘江平水期管道南侧、管道中间、管道北侧发生完全断裂时成品油的泄漏工况，平水期3个工况不同时刻泄漏成品油扩散的数值模拟结果分别如图6～图8所示。在平水期3个工况，同样是泄漏0.5 h后成品油全部流入北盘江，泄漏成品油最初集中在跨越段管道附近河道内，呈大范围块状聚集；泄漏2 h后，成品油主要集中在跨越段下游9.6～11.8 km处的河道中；泄漏4 h后，成品油大多已经成为油水混合细带状，进一步污染到下游14.6～15.7 km处。

图6 工况④不同时刻的成品油扩散情况

图7 工况⑤不同时刻的成品油扩散情况

图8 工况⑥不同时刻的成品油扩散情况

工况⑦～工况⑨分别为北盘江丰水期管道南侧、管道中间、管道北侧发生完全断裂时成品油的泄漏工况，丰水期3个工况不同时刻泄漏成品油扩散的数值模拟结果如图9～图11所示。由于丰水期水流速度较快，泄露4 h后工况⑦～工况⑨的油品扩散范围已超出计算域，因此选取了t=3 h的扩散情况进行对比分析。丰水期3个工况泄漏的成品油在北盘江中扩散得较为迅速，按照泄漏0.5 h后成品油全部流入北盘江，此时成品油富集于管道附近大部分河岸。泄漏2 h后，油膜集中在跨越段下游12.2～13.8 km的河道范围内；泄漏3 h后，成品油扩散成小碎块状，且有部分滞留于河流弯头处，在下游15.6～17.8 km处有较大块油膜存在。

图9 工况⑦不同时刻的成品油扩散情况

图10 工况⑧不同时刻的成品油扩散情况

图11 工况⑨不同时刻的成品油扩散情况

3.2 管道小孔泄漏扩散结果分析

图12 工况⑩不同时刻的成品油扩散情况

图13 工况不同时刻的成品油扩散情况

4 结论

笔者基于有限元体积法和κ-ω模型，采用追踪多相流界面的VOF方法，建立了一种跨越段管道泄漏污染扩散模型，得到以下结论：

(1)模型计算结果可得出不同工况下成品油扩散范围与泄漏时间的关系，丰水期管道完全断裂泄漏扩散范围最广，平水期次之，枯水期最小。

(2)管道泄漏点位置在泄漏前期对成品油扩散影响较大，但随着扩散时间的增加，在泄漏后期成品油扩散范围主要由水流速度和河岸的地形决定。在成品油污染扩散前期，泄漏的成品油大部分聚集在一起，尚未被水流冲散，适合收集处理。因此，在管道泄漏点附近，宜设置围油栏等应急处理装置进行油品拦截回收。

(3)河道中的弯道容易产生漩涡，漩涡附近会造成油品滞留，也会使经过的大块油膜被冲散，部分油品会附着在此处的河岸上。应急抢险时可在弯道处安装拦油装置，有利于节约时间和回收成本。

(4)小孔泄漏造成的成品油扩散泄漏流量小、泄漏持续时间长，不会产生大块的油膜。丰水期水流速度快，小孔泄漏的油品基本被水流冲散，污染扩散快且不利于收集。但无论是枯水期还是丰水期，泄漏前期的应急处理都很关键，对降低泄漏造成的流域污染损失起到重要作用。

(5)VOF方法追踪互不相容且物性参数相差较大的两种或多种流体的多相流流动界面时效果较好，且可根据离散方法改进计算精度。