引言

昆明市人工煤气的复杂组成介质导致了管道多相流动和腐蚀问题。由于人工煤气中含有一定量的水、CO$_2$等腐蚀性介质，以及焦油、灰尘、萘和氨等杂质，在长期运行过程中，管道内存在较为严重水沉积、固体颗粒颗粒沉积和腐蚀等问题。因此，如何分析多相流动条件下管道的腐蚀特征，评价腐蚀速率，了解安全运行能力，是制定天然气替代方案中首先需要解决的基础性问题。

对管道流体混合机理的系统研究开始于Taylor，其研究主要以传质学、传热学和流体力学相关理论作为基础，在混气机理的研究中，主要研究内容集中在扩散系数与混气量计算模型上。Taylor、Taylor-CW、GRI扩散系数计算方法在国内外虽然得到了应用，但其适用范围有限，Taylor混合系数在低雷诺数情况下对管线置换的混合系数能够较好地预测，但高雷诺数下的误差较大；Taylor-CW、GRI扩散系数计算方法对高雷诺数的工况能较好地预测，但却不适用于低雷诺数工况^[1-3]。混气量计算模型主要为一维和二维模型，并且大多是将水力计算与对流扩散计算分开，存在一定的局限性^[4-6]。Ekambara等针对管道层流流动扩散建立了数学模型，并采用数值方法对混气规律进行分析，但其求解方法仍然分为两步，首先进行水力计算，获得流速等参数后再对混气规律进行分析，其计算结果与实验值吻合的较好^[7-10]。基于Ekambara等的成果，以FLUENT为模拟平台，进一步采用低雷诺数$k$$-$$\varepsilon$模型对输气管道气体置换过程进行二维数值模拟^[11-13]。

因目前混气影响规律分析大多是基于单条短管，未对长距离管线和管网进行分析。对于实际管线，管道还可能存在局部变形、管内存在杂质、积液以及管道存在上坡和下坡等情况，需要进一步对这些混气影响因素进行分析，建立多相流动腐蚀速率预测模型计算所有燃气管道的腐蚀情况，进而为管道剩余寿命及安全能力评价提供依据。本文基于双流体模型建立多相流输送管道数值模拟模型，分析了整个管网的运行情况，得到了管网内的压力、流量、持液率等参数；以此为基础，利用Pipephase建立燃气管网仿真模型，输入双流体模型求解的边界条件，进一步分析持液率，识别判断具备腐蚀条件的管段，根据持液率、CO$_2$分压、温度、pH值等影响重要度，并基于全动态多相流模拟(简称OLGA)多相流软件的Dewarrd腐蚀速率预测模型，对各管段的腐蚀速率，得到各积水管段的腐蚀等级。针对目前采用的腐蚀预测模型不适用于人工煤气管道的实际情况，基于腐蚀热力学原理建立了人工煤气管网二氧化碳腐蚀速率预测模型，为昆明市人工煤气管网的腐蚀速率分析提供条件。

1 双流体模型

人工煤气成分较为复杂，因含水在流动过程中状态主要为气、液两相，低洼管段容易积水，人工煤气在流动过程中会出现液滴现象，通过双流体模型求解得到管道内的压力、温度、气液相流量和流速等关键性数值^[14-17]。

1.1 质量守恒方程

针对管道中的气相，根据质量守恒定理，可列出如下方程

$ \dfrac{\partial }{{\partial t}}\left( {{V_{\rm{g}}}{\rho _{\rm{g}}}} \right) = - \dfrac{1}{A}\left[{\dfrac{\partial }{{\partial x}}\left( {A{V_{\rm{g}}}{\rho _{\rm{g}}}{v_{\rm{g}}}} \right)} \right] + {\psi _{\rm{g}}} + {G_{\rm{g}}} $

(1)

针对管道中的液相，方程如下

$ \dfrac{\partial }{{\partial t}}\left( {\mathop V\nolimits_{\rm{L}} \mathop \rho \nolimits_{\rm{L}} } \right) = - \dfrac{1}{A}\left[{\dfrac{\partial }{{\partial x}}\left( {A\mathop V\nolimits_{\rm{L}} \mathop \rho \nolimits_{\rm{L}} \mathop v\nolimits_{\rm{L}} } \right)} \right] - \mathop \psi \nolimits_{\rm{g}} \dfrac{{\mathop V\nolimits_{\rm{L}} }}{{\mathop V\nolimits_{\rm{L}} + \mathop V\nolimits_{\rm{D}} }} - \mathop \psi \nolimits_{\rm{e}} + \mathop \psi \nolimits_{\rm{d}} + \mathop G\nolimits_{\rm{L}} $

(2)

考虑其中的液相夹带，方程如下

$ \dfrac{\partial }{{\partial t}}\left( {\mathop V\nolimits_{\rm{D}} \mathop \rho \nolimits_{\rm{L}} } \right) = - \dfrac{1}{A}\left[{\dfrac{\partial }{{\partial x}}\left( {A\mathop V\nolimits_{\rm{D}} \mathop \rho \nolimits_{\rm{L}} \mathop v\nolimits_{\rm{D}} } \right)} \right] - \mathop \psi \nolimits_{\rm{g}} \dfrac{{\mathop V\nolimits_{\rm{L}} }}{{\mathop V\nolimits_{\rm{L}} + \mathop V\nolimits_{\rm{D}} }} - \mathop \psi \nolimits_{\rm{e}} + \mathop \psi \nolimits_{\rm{d}} + \mathop G\nolimits_{\rm{D}} $

(3)

式中：

$t$—时间，s；

$V$—体积分率，%；

$\rho$—密度，kg/m$^3$；

$A$—流过管线的横断面积，m$^2$；

$x$—管道长度，m；

$v$—流速，m/s；

$\psi$—气液相间的质量传递速率，液相蒸发为正，kg/s；

$G$—质量源，kg；

$\psi_{\rm{e}}$、$\psi_{\rm{d}}$—液滴夹带及沉积速率，kg/s；

$V_{\rm{g}}$、$V_{\rm{L}}$、$V_{\rm{D}}$—气相、液膜、液滴体积分率，且$V_{\rm{g}}$+$V_{\rm{L}}$+$V_{\rm{D}}$=1。

下标：

g—气相；

L—液膜；

D—液滴。

注：气相指流动的人工煤气；液膜指附着在管道内壁的溶解各种成分的水膜；液滴指流动过程中形成的水滴；沉积速率指液滴单位时间内沉积的质量；体积分率指单位体积管段中气相、液膜、液体所占百分比，总和为1。

1.2 动量守恒方程

考虑各相本身的动量变化以及各相间动量传递引起的动量变化，根据动量守恒原理可列出气相和液滴相的联合动量守恒方程如下

$ \dfrac{\partial }{{\partial t}}\left( {{V_{\rm{g}}}{\rho _{\rm{g}}}{v_{\rm{g}}} + {V_{\rm{D}}}{\rho _{\rm{L}}}{v_{\rm{D}}}} \right) = - \dfrac{1}{A}\dfrac{\partial }{{\partial x}}\left( {A{V_{\rm{g}}}{\rho _{\rm{g}}}v_{\rm{g}}^2 + \\A{V_{\rm{D}}}{\rho _{\rm{L}}}v_{\rm{D}}^2} \right) - \left( {{V_{\rm{g}}} + {V_{\rm{D}}}} \right)\dfrac{{\partial p}}{{\partial x}} - \dfrac{1}{2}{\lambda _{\rm{g}}}{\rho _{\rm{g}}}\left| {{v_{\rm{g}}}} \right|{v_{\rm{g}}}\dfrac{{{S_{\rm{g}}}}}{{4A}} - \dfrac{1}{2}{\lambda _{\rm{i}}}{\rho _{\rm{g}}}\left| {{v_{\rm{r}}}}\right|{v_{\rm{r}}}\dfrac{{{S_{\rm{i}}}}}{{4A}} +\left( {{V_{\rm{g}}}{\rho _{\rm{g}}} \\+{V_{\rm{D}}}{\rho _{\rm{L}}}} \right){\rm{g}}\cos \alpha + {\psi _{\rm{g}}}\dfrac{{{V_{\rm{L}}}}}{{{V_{\rm{L}}} + {V_{\rm{D}}}}}{v_{\alpha}} + {\psi _{\rm{i}}}{v_{\rm{i}}} - {\psi _{\rm{d}}}{v_{\rm{D}}} $

(4)

同理，可列出液相动量守恒方程为

$ \dfrac{\partial }{{\partial t}}\left( {{V_{\rm{L}}}{\rho _{\rm{L}}}{v_{\rm{L}}}} \right) = - \dfrac{1}{A}\dfrac{\partial }{{\partial x}}\left( {A{V_{\rm{L}}}{\rho _{\rm{L}}}v_{\rm{L}}^2} \right) - {V_{\rm{L}}}\dfrac{{\partial p}}{{\partial x}} -\\ \dfrac{1}{2}{\lambda _{\rm{L}}}{\rho _{\rm{L}}}\left| {{v_{\rm{L}}}} \right|{v_{\rm{L}}}\dfrac{{{S_{\rm{L}}}}}{{4A}} - \dfrac{1}{2}{\lambda _{\rm{i}}}{\rho _{\rm{g}}}\left| {{v_{\rm{r}}}} \right|{v_{\rm{r}}}\dfrac{{{S_{\rm{i}}}}}{{4A}} + \left( {{V_{\rm{g}}}{\rho _{\rm{g}}} + {V_{\rm{D}}}{\rho _{\rm{L}}}} \right){\rm{g}}\cos \alpha -\\ {\psi _{\rm{g}}}\dfrac{{{V_{\rm{L}}}}}{{{V_{\rm{L}}} + {V_{\rm{D}}}}}{v_{\rm{ { α} }}} - {\psi _{\rm{e}}}{v_{\rm{i}}} + {\psi _{\rm{d}}}{v_{\rm{D}}} - {V_{\rm{L}}}d\left( {{\rho _{\rm{L}}} - {\rho _{\rm{g}}}} \right){\rm{g}}\dfrac{{\partial {V_{\rm{L}}}}}{{\partial x}}\sin \alpha $

(5)

式中：

$p$—压力，Pa；

$\lambda$—摩阻系数，无因次；

$v_{\rm{r}}$—气液两相的相对速度，m/s；

g—重力加速度，m/s$^2$；

$\alpha$—管轴线与垂线夹角，(°)；

$S_{\rm{g}}$、$S_{\rm{L}}$、$S_{\rm{i}}$—气相、液膜、气液相间界面的湿周，m；

$v_{\alpha}$—各相间的质量传递速度，kg/s；对于液膜蒸发，$\psi_{\rm{g}}>0$，$v_{\alpha}=v_{\rm{L}}$；对于液滴蒸发，$\psi_{\rm{g}}>0$，$v_{\alpha}=v_{\rm{D}}$；对于气体冷凝，$\psi_{\rm{g}}<0$，$v_{\alpha}=v_{\rm{g}}$；

$d$—管道直径，mm；

下标：

i—气液相间的界面。

1.3 能量守恒方程

考虑气相、液相、夹带液滴相间的热量传递，以及管道内流体与外界环境之间的热量传递，可列出上述三相的混合能量方程，如下

$ \dfrac{\partial }{{\partial t}}\left[{{m_{\rm{g}}}\left( {{E_{\rm{g}}} + \dfrac{1}{2}V_{\rm{g}}^2 + {\rm{g}}h} \right) + {m_{\rm{L}}}\left( {{E_{\rm{L}}} + \dfrac{1}{2}V_{\rm{L}}^2 + {\rm{g}}h} \right) + {m_{\rm{D}}}\left( {{E_{\rm{D}}} + \dfrac{1}{2}V_{\rm{D}}^2 +\\ {\rm{g}}h} \right)} \right] = - \dfrac{\partial }{{\partial x}}\left[{{m_{\rm{g}}}{v_{\rm{g}}}\left( {{H_{\rm{g}}} + \dfrac{1}{2}V_{\rm{g}}^2 + {\rm{g}}h} \right) + \\ {m_{\rm{L}}}{v_{\rm{L}}}\left( {{H_{\rm{L}}} + \dfrac{1}{2}V_{\rm{L}}^2 + {\rm{g}}h} \right) + {m_{\rm{D}}}{v_{\rm{D}}}\left( {{H_{\rm{D}}} + \dfrac{1}{2}V_{\rm{D}}^2 + {\rm{g}}h} \right)} \right] $

(6)

式中：

${m_{\rm{f}}} = {V_{\rm{f}}}{\rho _{\rm{f}}}$，其中f=g，L，D；

$E$—单位质量的内能，J；

$h$—高程，m；

$H$—质量源的焓，J。

在动量守恒、质量守恒和能量守恒方程的基础上，一个完整的双流体计算模型还应包括热物性参数计算模型、气液两相流型判别模型，以及描述相间质量、动量和热量传递关系的闭合关系式。通过求解上述双流体模型，可得到昆明市人工煤气管网示意图各管段的压力、温度和气液相流速等数值(表 1)。

2 昆明市人工煤气管网持液率分析

分析持液率的目的是判断管道积水段和具备管道腐蚀的环境条件以及识别腐蚀区域。利用Pipephase建立昆明市人工煤气燃气管网仿真模型。设置管网入口节点1为压力和温度边界条件，其余节流点为流量边界条件；人工煤气的组分如表 2所示，同时，根据历史检测数据，设置人工煤气中水的体积分数为0.9%。根据中国天然气网(http://-www.wea-ther.com.cn/-cityintro/101290101.-shtml)提供的数据，昆明市的年均气温为16.4 ℃，因此，设定管道入口温度和环境温度均为14.9 ℃。根据这些参数，计算得到的稳态条件下昆明市人工煤气管网内的持液率分布如图 1所示。

由图 1可知，人工煤气管网的各个管段内都存在不同程度的积液，其中持液率小于0.100的管道为42条，持液率0.100$\sim$0.200的管道为39条，其余管道的持液率为0.200$\sim$0.559。持液率的变化主要与管道的管径、输量、气液相流速、管道温度相关。

昆明市人工煤气管网持液率大于0.400的管段共有24段，如图 2所示，持液率最高的管道为L78，持液率值为0.559。现场开挖过程中，也证实了积液段最易腐蚀。中压L78在2014年6月管段曾发生过腐蚀泄漏事故。这说明持液率预测结果能够从一定程度上反映现场的腐蚀情况，最易导致管道积液和发生内腐蚀的管段如图 3所示。

以上述平均气温求解的持液率数据为基础，进一步结合昆明市的气候特点，设定春秋季平均气温为13.2 ℃、夏季平均气温为19.6 ℃，冬季为9.3 ℃，分析不同温度下管道内的持液率，结果如图 4所示。

图 4表明，管道内的持液率是与环境温度密切相关的，在冬季环境温度为9.3 ℃时，管道内的持液率和积液量要远远大于夏季环境温度为19.6 ℃时的值。这也说明，在冬季管道内要注意及时排出积液，降低管道内的积液量，从而提高管道的输送效率和缓解腐蚀情况。

3 腐蚀速率预测 3.1 腐蚀速率预测模型

使用OLGA中的De Waard95模型来进行腐蚀速率预测，需考虑Fe$^{2+}$浓度、温度、pH值、CO$_2$分压、与流速相关的传质过程以及腐蚀产物膜以及与流速无关的腐蚀动力学过程等因素，在双流体研究的基础之上，可以进一步结合腐蚀预测模型计算管道的腐蚀速率。昆明市人工煤气管网的运行温度介于0$\sim$30 ℃，De Waard95模型是壳牌公司提出的一种半经验模型，在温度低于80 ℃时具有较高的精度。进行多相流腐蚀计算时选择De Waard95模型预测管道的腐蚀速率^[18-24]。

De Waard 95模型的表达式如下

$ \dfrac{1}{{ V_{{\rm{corr}}} }} = \dfrac{1}{{ V_{\rm{r}} }} + \dfrac{1}{{ V_{\rm{m}} }} $

(7)

$ \lg V_{\rm{r}} = 4.93 - \dfrac{{1119}}{{T+ 273}} + 0.58\lg p_{{\rm{C}}{{\rm{O}}_2}} -\\\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;0.34\left( {p {\rm{H}}_{{\rm{act}}} - p {\rm{H}}_{{\rm{C}}{{\rm{O}}_2}} } \right) $

(8)

$ V_{\rm{m}} = 2.45\dfrac{{U^{0.8} }}{{ d^{0.2} }} p_{{\rm{C}}{{\rm{O}}_2}} $

(9)

$ p {\rm{H}}_{{\rm{C}}{{\rm{O}}_2}} = 3.82 + 0.00384T - 0.5\lg p_{{\rm{C}}{{\rm{O}}_2}} $

(10)

式中：

$V_{{\rm{corr}}}$—腐蚀速率，mm/a；

$V_{\rm{r}}$—反应速率，mm/a；

$V_{\rm{m}}$—传质速率，mm/a；

$T$—介质温度，℃；

$p_{{\rm{C}}{{\rm{O}}_2}}$—CO$_2$分压，Pa；

$p{\rm{H}}_{{\rm{act}}}$—实际pH值；

$p{\rm{H}}_{{\rm{C}}{{\rm{O}}_2}}$—CO$_2$饱和溶液的pH值；

$U$—介质的液相流动速率，m/s。

钢制管道的腐蚀速率主要与反应腐蚀速率和传质腐蚀速率两方面有关，腐蚀速率随温度和CO$_2$分压上升而变大；传质速率随流速和CO$_2$分压的上升而变大，但随管径的增大而减小。腐蚀速率的主要影响因素是CO$_2$分压、介质流速以及温度等。

3.2 腐蚀速率预测结果分析

考虑到L78管段内积液量较高，曾发生过腐蚀导致的泄漏事故。因此，通过在双流体模型基础上再利用De Waard 95的腐蚀模型计算腐蚀速率。L78管段的腐蚀情况如图 5所示。

图 5表明，计算得到的平均腐蚀速率为0.054 3 mm/a，这与室内动态腐蚀实验结果动态腐蚀速率0.277 4 mm/a和管道服役条件下(现场壁厚测试)的腐蚀速率0.190 9 mm/a存在明显差异，OLGA多相流下的腐蚀预测结果仅仅是实际值的28.4%，因此需要根据实验测试结果对该腐蚀预测模型与方法进行改进，建立适于人工煤气管网的腐蚀预测模型，从而为管网的运行能力分析提供依据与条件。

3.3 人工煤气管网腐蚀速率预测模型建立

人工煤气管道腐蚀过程中的电化学和化学反应与CO$_2$含量有关，化学反应速率取决于溶液pH值。反应释放的Fe$^{2+}$、CO$_3^{2-}$和HCO$_3^-$可沉淀至管道表面而形成FeCO$_3$的腐蚀产物。

$ {\rm{Fe}}_{{\rm{(aq)}}}^{{\rm{2 + }}}{\rm{ + CO}}_{{\rm{3(aq)}}}^{{\rm{2 - }}} \to {\rm{FeC}}{{\rm{O}}_{{\rm{3(g)}}}} $

(11)

在含有CO$_2$的反应环境中，管道的腐蚀速率取决于扩散控制或激活控制的反应^[25-30]。

化学反应速率常数随温度变化关系满足阿伦尼乌斯公式(Arrhenius Equation)

$ k = {a^{ - {E_{\rm{a}}}/{\rm{R}}{T}}} $

(12)

式中：

$k$—反应速率常数，无因次；

$a$—指前因子，无因次；

$E_{\rm{a}}$—表观活化能，J/mol；

R—摩尔气体常量，R=8.314 J/(mol$\cdot$K)。

为了更准确地预测腐蚀速率，根据人工煤气管道的输送工况，考虑不同的CO$_2$分压、温度条件进行动态反应釜实验，以确定CO$_2$分压对内腐蚀速率的影响规律，具体实验数据见图 6。

在同样温度条件下，随着CO$_2$分压的增加，溶解于水中的CO$_2$增加，0.050 MPa CO$_2$分压下的管材腐蚀速率明显大于0.005 MPa，在20 ℃时，腐蚀速率分别是0.850 0 mm/a与0.180 0 mm/a；在同样的CO$_2$分压条件下，腐蚀速率随着温度的增加而显著增加，在0.005 MPa的CO$_2$分压时，管材的腐蚀速率从20 ℃的0.180 0 mm/a增加到50 ℃的0.990 0 \mbox{mm/a}。根据不同分压、温度下的腐蚀速率实验测试结果数据绘制ln $k\sim1/T$直线图，可以得到阿伦尼乌斯公式在不同CO$_2$分压下的腐蚀公式参数，即由斜率可得表观活化能$E_{\rm{a}}$，由截距确定指前因子$a$的大小(在反应温度变化不大时$E_{\rm{a}}$可视为常数)，将微分式作定积分可得不同温度下的反应速率常数与其相对应的温度之间的关系。

$ V_{{\rm{corr}}} = a {{\rm{e}}^{ - Q/{\rm{R}}T}} p_{{\rm{C}}{{\rm{O}}_2}}^{k} $

(13)

式中：$Q$—反应热量，J/mol。

根据实验结果测试的腐蚀速率数据，确定相应$a$、$Q$参数，其中$a$=250 640，$Q$=33 600 J/mol，$k$=0.11，建立${p_{{\rm{C}}{{\rm{O}}_2}}} \leqslant 0.005 {\rm{MPa}}$条件下的腐蚀速率模型为

$ V_{{\rm{corr}}} = 250640 \times {{\rm{e}}^{ - 33600/{\rm{R}}T}} \times p_{{\rm{C}}{{\rm{O}}_2}}^{0.11} $

(14)

在$0.005 {\rm{MPa}} < {p_{{\rm{C}}{{\rm{O}}_2}}} \leqslant 0.050 {\rm{MPa}}$条件下的腐蚀速率预测模型，其中$a$=917 250，$Q$=33 550 J/mol，$k$=0.11。其腐蚀预测公式为

$ V_{{\rm{corr}}}= 917250 {{\rm{e}}^{ - 33550/{\rm{R}}T}} p_{{\rm{C}}{{\rm{O}}_2}}^{0.11} $

(15)

通过De Waard95模型求解与实验测试结果对比(见表 3)，建立的腐蚀速率预测模型误差平均值为8.9%，最小误差为1.8%，最大误差为18.8%，明显优于OLGA中多相流腐蚀预测结果，可用于人工煤气管道的腐蚀速率预测，为管网的安全运行能力评价提供准确的腐蚀速率数据。

昆明市人工煤气管网各管段存在明显的含水现象，通过管网的多相流持液率计算其最大的持液率为0.559，根据昆明市人工煤气管网各管段持液率情况，利用腐蚀预测模型，进行昆明市人工煤气管网的腐蚀速率计算，根据腐蚀等级划分依据，从而可以确定昆明市人工煤气管网的腐蚀等级(图 7)，其中“高”等级的管段共14段，主要分布于管网前端压力较高与末端持液率高的管段，利用确定的腐蚀速率值便可实现目标管网的安全能力评价。

目前昆明人工煤气管道有相当数量从未进行过通球内检测，而且由于服役时间比较长，管道内部情况比较复杂(可能存在结垢、积蜡或管道变形等不利于实施通球检测的因素)，无法实施内检测。对这种既没有历史内检测数据、又无法实施内检测或实施内检测风险很高的管道，只有采取上述的的安全评估方法，既可避免内检测作业又能间接掌握管道上腐蚀缺陷发展程度。在基本掌握具有腐蚀缺陷的管道安全状况，可以制定具有针对性的管道管理策略和措施。

4 结论

(1) 基于双流体模型进行人工煤气管网的多相流动分析，管网中持液率小于0.100的管道为42条，持液率0.100$\sim$0.200的管道为39条，持液率大于0.400的管段共有24段，持液率最高的管道为L78，其持液率值为0.559。

(2) 不同环境温度下管道的持液率存在变化，这主要与管道对的管径、输量、气液相流速、管道温度相关。

(3) 在人工煤气管网多相流分析模型的基础上，分析了管网的持液率分布规律，采用OLGA多相流腐蚀预测模型分析了最高持液率管段的腐蚀速率。

(4) 针对已有腐蚀模型预测误差问题，利用描述化学反应速率常数随温度变化关系的阿伦尼乌斯公式(Arrhenius Equation)，结合不同二氧化碳分压条件下的动态反应釜测试实验，建立了人工煤气管网腐蚀速率预测模型，对比实验测试结果，腐蚀速率预测模型误差平均值为8.9%，明显优于OLGA多相流腐蚀预测结果。

(5) 基于腐蚀热力学原理建立人工煤气管网二氧化碳腐蚀速率预测模型，直接为昆明市人工煤气管网的腐蚀速率分析提供条件，进而为管道剩余安全能力评价提供依据，判断各管道抢维修优先次序。