CN105004644A - 一种酸性天然气管线硫沉积预测方法 - Google Patents

Abstract

一种酸性天然气管线硫沉积预测方法，基于硫的溶解度模型、经典成核理论和粒子运动方程建立硫元素析出和沉积的热力学及动力学模型，并通过计算饱和距离、过饱和距离和颗粒最大迁移距离，计算酸性天然气管线的硫沉积位置、硫沉积量以及沉积硫的平均厚度，实现酸性天然气管线硫沉积预测；本发明认为饱和点到出现可见硫颗粒位置之间存在一个过饱和距离L₂，能够更精确的计算可见硫颗粒的析出位置；采用运动方程计算硫颗粒的运动轨迹，继而确定硫颗粒的最大迁移距离L₃，该方法从机理上比借用垂直井筒中硫颗粒临界流速计算模型进行硫沉积位置预测更合理，使得预测结果也更加科学、准确、可靠。

Description

一种酸性天然气管线硫沉积预测方法

技术领域

本发明涉及酸性天然气集输技术领域，特别涉及一种酸性天然气管线硫沉积预测方法。

背景技术

酸性天然气中常含有硫化氢等硫化物，在开发和集输过程中，不仅具有剧毒性和高腐蚀性，在特定条件下还会析出单质硫而引起硫元素的沉积。集输系统中的硫沉积会导致管线的过流面积减小，降低集输效率，严重时会堵塞管线和设备，造成停产。同时，沉积硫的存在也会加速集输管道和设备的腐蚀，造成管材的疲劳破坏，从而严重影响集输系统的正常运行。

国内外学者对气藏和井筒中的硫沉积问题研究较多，而集输系统中的硫沉积一直没有引起足够的重视，尚无针对集输系统硫沉积的有效预测方法。随着能源问题的日益突出，越来越多的高含硫酸性气藏被开发，而集输系统的运行压力也不断提高，更有利于硫元素的形成，这使得集输系统中的硫沉积问题越来越突出。目前，硫沉积已成为含硫酸性天然气集输过程中普遍存在的难点问题和制约含硫酸性天然气集输技术的瓶颈之一。而准确的硫沉积预测是确定酸性天然气管线硫沉积防治措施和管线中沉积硫经济清管周期的基础。

目前，已有的井筒硫沉积预测方法一般采用由硫固体颗粒垂直方向所受重力、曳力和浮力合力为零推导得到的临界流速计算模型作为硫沉积的判断准则。在地面管线中，上述三种力不在同一方向，因此该方法不能适用于地面管线的硫沉积预测。另外，传统的井筒硫沉积预测方法未考虑过饱和距离，认为只要含硫量大于其溶解度，即有固体硫颗粒析出，并采用其最大颗粒直径用于临界流速计算。而实际上，在过饱和区域，固体硫颗粒处于成核过程中，不会发生硫元素沉积；在固体硫颗粒沿管线运动过程中，硫颗粒的直径也并不是不变的，而是随着硫元素的析出逐渐增大的。传统方法并未考虑上述因素的影响，其计算精度也就不尽如人意。因此，建立更加科学可靠的酸性天然气管线硫沉积预测方法具有十分重要的意义。

发明内容

针对酸性天然气集输系统的硫沉积问题，本发明的目的在于提出一种酸性天然气管线硫沉积预测方法，基于硫的溶解度模型、经典成核理论和粒子运动方程建立硫元素析出和沉积的热力学及动力学模型，并通过计算饱和距离、过饱和距离和颗粒最大迁移距离，计算酸性天然气管线的硫沉积位置、硫沉积量以及沉积硫的平均厚度，实现酸性天然气管线硫沉积预测。

为了达到上述目的，本发明的技术方案为：

一种酸性天然气管线硫沉积预测方法，包括如下步骤：

步骤S1：建立酸性天然气热力学性质计算模型，热力学性质包括密度、动力粘度、偏差系数、焦耳-汤姆逊系数；

步骤S2：将酸性天然气管线划分为多个计算单元，利用输气管道压降、温降计算公式，建立酸性天然气管线压力温度耦合计算模型，并计算沿管线的压力温度分布；

步骤S3：建立酸性天然气中硫元素的溶解度计算模型，根据压力温度参数计算沿管线硫元素的溶解度变化；

步骤S4：计算从管道起点到硫蒸汽浓度达到其溶解度值的起始位置之间的距离，称为饱和距离；

步骤S5：根据溶解度的变化确定各区域硫元素的析出量；

步骤S6：采用经典成核理论计算从饱和点到出现可见硫颗粒位置之间的距离，称为过饱和距离；

步骤S7：通过跟踪硫颗粒的运动轨迹确定硫颗粒的最大迁移距离；

步骤S8：根据硫元素沿管线的析出量和所析出硫颗粒的最大迁移距离计算沿管线的硫沉积位置、硫沉积量以及沉积硫的平均厚度。

所述步骤S6具体为：

(1)采用经典成核理论计算需要克服吉布斯自由能的阀值：

$J={\mathrm{Ke}}^{-\frac{\mathrm{\Δ}G}{k_{B}T}}$

式中J，K，ΔG，k_B和T分别表示成核速率、预指数因子、吉布斯自由能、玻耳兹曼常数和温度，预指数因子K采用下式计算

$K={\mathrm{\θ\α}}_c{v}_s{N}^2\sqrt{\frac{2\mathrm{\σ}}{k_{B}T}}$

式中θ，α，υ，N和σ分别表示不等温因子、凝聚因子、硫的分子体积、可凝聚分子数和硫的界面张力。

可凝聚分子数采用下式计算：

N＝ρ_My_EN_A

式中ρ_M,y_E和N_A分别表示混合密度、萃取条件下硫元素的摩尔分数和阿伏加德罗常数。

硫的界面张力σ由下式计算：

$\mathrm{\σ}={p_c}^{2/3}{T}_c^{1/3}\[0.1207\×\left(1+\frac{T_{br}{\mathrm{lnp}}_c}{1-T_{br}}\right)-0.281\]{(1-T_r)}^{1.222}$

式中p_c、T_c、T_br、T_b和T_r分别表示临界压力、临界温度、对比温度、对比沸点温度和沸点温度。

天然气中存在杂质时，吉布斯自由能通过下式计算：

$\mathrm{\Δ}G=\frac{16{\mathrm{\π\σ}}^3}{3{(p_{Sat}-p_{sys})}^2}S\left(\mathrm{\θ}\right)=\frac{16{\mathrm{\π\σ}}^3}{3{(p_{Sat}-p_{sys})}^2}\frac{(2+cos\mathrm{\θ}){(1-cos\mathrm{\θ})}^2}{4}$

式中p_sat，p_sys和S(θ)分别表示饱和压力、系统压力和杂质形状因子。

(2)当成核速率超出可见成核速率的阀值时，认为首个可见的硫元素颗粒出现，从而由下式确定过饱和距离L₂，

$L_2=L_{J>J_{visible}}$

式中J_visible表示出现可见硫元素颗粒时的成核速率，建议取1个/m³d。

所述步骤S7具体为：

(1)通过求解如下运动方程计算硫颗粒的运动轨迹：

$\frac{d\stackrel{\→}{x_p}}{dt}=\stackrel{\→}{u_p}(\stackrel{\→}{x_p},t)$

式中t和分别表示硫颗粒的位置、时间和速度，硫颗粒的速度通过牛顿第二定律计算：

$m_p\left(t\right)\×\frac{d\stackrel{\→}{u_p}}{dt}=\mathrm{\Σ}\stackrel{\→}{F}$

式中m_p(t)表示瞬时硫颗粒的质量，表示总力包括曳力、重力以及布朗力，其中布朗力只有当颗粒尺寸小于2μm时才起作用。式中硫颗粒的瞬时粒径通过质量平衡计算：

$N_{i+1}\mathrm{\ρ}\frac{4{\mathrm{\πr}}_{p,i+1}^3}{3}={\mathrm{\Δm}}_i^{i+1}+N_i\mathrm{\ρ}\frac{4{\mathrm{\πr}}_{p,i}^3}{3}$

式中r_p表示颗粒半径，Δm表示计算单元中硫元素的析出量。

硫颗粒的凝聚特性采用Smoluchowski模型计算：

$\frac{N_i}{N_0}=\frac{1}{1+2{\mathrm{KN}}_0{t}_i}=\frac{1}{1+\frac{k_{B}T}{3\mathrm{\μ}}{N}_0{t}_i}$

式中N₀，K和μ分别表示初始硫颗粒数量、碰撞系数和气体粘度。

假设管壁金属表面液体无液体润湿，采用Tabakoff模型模拟粒子和管道壁之间的相互作用：

${\mathrm{\α}}_2={\mathrm{ctg}}^{-1}\left(\frac{0.95+0.00055{\mathrm{\α}}_1}{1.0-0.02108{\mathrm{\α}}_1+0.000141{\mathrm{\α}}_1^2}{\mathrm{ctg\α}}_1\right)$

$\frac{V_{p2}}{V_{p1}}=(1-0.02108{\mathrm{\α}}_1+0.000141{\mathrm{\α}}_1^2)\sqrt{\frac{1+{\mathrm{ctg}}^2{\mathrm{\α}}_2}{1+{\mathrm{ctg}}^2{\mathrm{\α}}_1}}$

式中α₁和α₂分别表示颗粒与管壁碰撞的入射角和反射角，V_p2和V_p1表示粒子碰撞前和碰撞后的速度。

(2)根据计算出的硫颗粒运动轨迹确定硫颗粒的最大迁移距离L₃。

所述步骤S8具体为：

(1)根据硫元素沿管线的析出量、所析出硫颗粒的过饱和距离和最大迁移距离计算沿管线的硫沉积位置和单位时间内硫元素的沉积量：

假设在管线入口x₁位置处，析出硫元素y₁g/s，其过饱和距离为L₂₁，最大迁移距离为L₃₁，则在x₁+L₂₁+L₃₁处，有y₁g/s硫元素沉积；

同理，在管线入口x₂位置处，析出硫元素y₂g/s，其过饱和距离为L₂₂，最大迁移距离为L₃₂，则在x₂+L₂₂+L₃₂处，有y₂g/s硫元素沉积；

同理，在管线入口x₃位置处，析出硫元素y₃g/s，其过饱和距离为L₂₃，最大迁移距离为L₃₃，则在x₃+L₂₃+L₃₃处，有y₃g/s硫元素沉积；

……；

同理，在管线入口x_n位置处，析出硫元素y_n g/s，其过饱和距离为L_2n，最大迁移距离为L_3n，则在x_n+L_2n+L_3n处，有y_n g/s硫元素沉积；

……；

管线某位置所有计算出的硫元素沉积量之和，即为此处单位时间内硫元素的沉积量，管线中硫沉积量不为0的区域即为硫沉积位置。

(2)计算酸性天然气管线中硫元素的总沉积量和沉积硫的平均厚度：

管线中某一时间段内硫元素的总沉积量：

$M=\underset{0}{\overset{L}{\∫}}(m\×t)dx$

式中t为时间，x为离管线入口的距离，L为管线总长，M为t时间内的管线中硫元素的总沉积量。

管线中沉积硫的平均厚度：

$\mathrm{\δ}=\frac{M}{\mathrm{\π}D\·L_s}$

式中δ为沉积硫平均厚度，M为管线中硫元素的总沉积量，D为管线内径，L_s为硫沉积区域长度，M为t时间内的管线中硫元素的总沉积量。

本发明的有益效果是，本发明认为饱和点到出现可见硫颗粒位置之间存在一个过饱和距离L₂，能够更精确的计算可见硫颗粒的析出位置；采用运动方程计算硫颗粒的运动轨迹，继而确定硫颗粒的最大迁移距离L₃，该方法从机理上比借用垂直井筒中硫颗粒临界流速计算模型进行硫沉积位置预测更合理，使得预测结果也更加科学、准确、可靠。

附图说明

图1为本发明构建预测模型的流程图。

图2为本发明中管线压力温度分布耦合计算流程图。

图3为本发明中硫固体颗粒的运动轨迹示意图。

具体实施方式

下面结合附图对本发明作进一步详细描述。

一种酸性天然气管线硫沉积预测方法，包括以下步骤：

参照图1，步骤S1：建立酸性天然气密度、动力粘度、偏差系数、焦耳-汤姆逊系数等热力学性质计算模型。

(1)酸性天然气密度由下式计算：

$\mathrm{\ρ}=\frac{{\mathrm{\Σy}}_i{M}_i}{{\mathrm{\Σy}}_i{V}_i}$

式中ρ，y_i，M_i和V_i分别表示混合气体的密度以及i组分的摩尔分数、相对分子质量和摩尔容积。

(2)酸性天然气的动力粘度计算：

①Dempsey方法

$\begin{array}{c}\ln \left(\frac{{\mathrm{\μT}}_{pr}}{{\mathrm{\μ}}_1}\right)=A_0+A_1{p}_{pr}+A_2{p}_{pr}^2+A_3{p}_{pr}^3+T_{pr}(A_4+A_5{p}_{pr}+A_6{p}_{pr}^2+A_7{p}_{pr}^3)\\ +T_{pr}^2(A_8+A_9{p}_{pr}+A_{10}{p}_{pr}^2+A_{11}{p}_{pr}^3)+T_{pr}^3(A_{12}+A_{13}{p}_{pr}+A_{14}{p}_{pr}^2+A_{15}{p}_{pr}^3)\end{array}$

其中

μ₁＝(1.709×10^-5-2.062×10^-6γ_g)(1.8T+32)+8.188×10^-3-6.15×10^-3lg(γ_g)

式中p_pr和T_pr分别表示天然气的拟对比压力和拟对比温度，μ₁表示标准大气压和给定温度下单组分气体的动力粘度，γ_g表示天然气的相对密度，相关参数见表1。

表1 Dempsey法参数取值表

②Lee关系式

μ＝10^-4K exp(Xρ_g ^γ)

其中

$K=\frac{(9.4+0.02M_g){\left(1.8T\right)}^{1.5}}{209+19M_g+1.8T}$

$X=3.5+\frac{986}{1.8T}+0.01M_g$

Y＝2.4-0.2X

式中μ表示在天然气的粘度，T表示气体绝对温度，M_g表示气体分子量；ρ_g表示天然气的密度，其计算式为：

${\mathrm{\ρ}}_g=3.4844\frac{{\mathrm{\γ}}_{g}p}{ZT}$

受H₂S、CO₂、N₂等非烃气体存在的影响，酸性天然气的粘度比常规气体粘度高，需要对其粘度进行非烃校正：

①Stanging校正

${{\mathrm{\μ}}_1}^{\′}={\mathrm{\μ}}_1+{\mathrm{\μ}}_{{\mathrm{co}}_2}+{\mathrm{\μ}}_{h_{2}s}+{\mathrm{\μ}}_{N_2}$

其中

${\mathrm{\μ}}_{h_{2}s}=y_{h_{2}s}\·(8.49\×{10}^{-3}\lg ({\mathrm{\γ}}_g)+3.73\×{10}^{-3})$

${\mathrm{\μ}}_{{\mathrm{co}}_2}=y_{{\mathrm{co}}_2}\·(9.08\×{10}^{-3}\lg ({\mathrm{\γ}}_g)+6.24\×{10}^{-3})$

${\mathrm{\μ}}_{N_2}=y_{N_2}\·(8.48\×{10}^{-3}\lg ({\mathrm{\γ}}_g)+9.59\×{10}^{-3})$

式中分别表示H₂S粘度校正值、CO₂粘度校正值和N₂粘度校正值，γ_g表示天然气的相对密度，分别表示H₂S、CO₂、N₂的摩尔含量。

②杨继盛校正

μ＝10^-4K'exp(Xρ_g ^γ)

其中

$K^{\′}=\frac{(9.4+0.02M_g){\left(1.8T\right)}^{1.5}}{209+19M_g+1.8T}+K_{H_{2}S}+K_{{\mathrm{CO}}_2}+K_{N_2}$

式中和分别表示天然气中H₂S、CO₂和N₂存在时引起的附加粘度校正系数。

对于0.6＜γ_g＜1.0的天然气：

$K_{H_{2}S}=Y_{H_{2}S}(0.000057{\mathrm{\γ}}_g-0.000017)\×{10}^4$

$K_{{\mathrm{CO}}_2}=Y_{{\mathrm{CO}}_2}(0.000050{\mathrm{\γ}}_g+0.000017)\×{10}^4$

$K_{N_2}=Y_{N_2}(0.000050{\mathrm{\γ}}_g+0.000047)\×{10}^4$

对于1.0＜γ_g＜1.5的天然气：

$K_{H_{2}S}=Y_{H_{2}S}(0.000029{\mathrm{\γ}}_g+0.0000107)\×{10}^4$

$K_{{\mathrm{CO}}_2}=Y_{{\mathrm{CO}}_2}(0.000024{\mathrm{\γ}}_g+0.000043)\×{10}^4$

$K_{N_2}=Y_{N_2}(0.000023{\mathrm{\γ}}_g+0.000074)\×{10}^4$

式中和分别表示天然气中H₂S、CO₂和N₂体积百分数。

(3)酸性天然气偏差系数计算：

①DPR(Dranchuk-Purvis-Robinson)方法

$Z=1+(A_1+\frac{A_2}{T_{pr}}+\frac{A_3}{T_{pr}^3}){\mathrm{\ρ}}_r+(A_4+\frac{A_5}{T_{pr}}){\mathrm{\ρ}}_r^2+\left(\frac{A_5{A}_6}{T_{pr}}\right){\mathrm{\ρ}}_r^5+\frac{A_7}{T_{pr}^3}{\mathrm{\ρ}}_r^2(1+A_8{\mathrm{\ρ}}_r^2)\exp (-A_8{\mathrm{\ρ}}_r^2)$

其中

ρ_r＝0.27p_pr/(ZT_pr)

式中p_pr，T_pr和ρ_r分别表示天然气的拟对比压力、拟对比温度和密度，相关参数见表2。

表2 DPR法参数取值表

②DAK(Dranchuk-Abu-Kassem)方法

$\begin{array}{c}Z=1+(A_1+\frac{A_2}{T_{pr}}+\frac{A_3}{{T_{pr}}^3}\frac{A_4}{{T_{pr}}^4}+\frac{A_5}{{T_{pr}}^5}){\mathrm{\ρ}}_r+(A_6+\frac{A_7}{T_{pr}}+\frac{A_8}{{T_{pr}}^2}){{\mathrm{\ρ}}_r}^2-A_9(\frac{A_7}{T_{pr}}+\frac{A_8}{{T_{pr}}^2}){{\mathrm{\ρ}}_r}^5\\ +A_{10}\left(1+A_{11}{{\mathrm{\ρ}}_r}^2\right)\left(\frac{{{\mathrm{\ρ}}_r}^2}{{T_{pr}}^3}\right)EXP\left(-A_{11}{{\mathrm{\ρ}}_r}^2\right)\end{array}$

其中

${\mathrm{\ρ}}_r=\frac{0.27p_{pr}}{{\mathrm{ZT}}_{pr}}$

式中各参数取值见表3。

表3 DAK方法参数取值表

③HY(Hall-Yarborough)方法

$Z=\frac{1+y+y^2-y^3}{{(1-y)}^3}-(14.76t-9.76t^2+4.58t^3)y+(90.7t-242.2t^2+42.4t^3)y^{(1.18+2.82t)}$

$Z=\frac{0.06125p_{pr}t}{y}\exp \[-1.2{(1-t)}^2\]$

其中

$t=\frac{1}{T_{pr}}=\frac{T_{pc}}{T}$

式中y表示对比密度，在已知p、T时，联立两式求解y，再将y代入式中就可求出偏差系数Z值。

酸性天然气中含有H₂S、CO₂、N₂等非烃气体，需要对其偏差系数进行非烃校正：

①Wicher-Azia校正

a.临界参数校正：

T′_ci＝T_ci-ε

p′_ci＝p_ciT_ci/T_ci

其中

ε＝15(M-M²)+4.167(N^0.5-N²)

式中T_ci，p_ci，T_ci′和p_ci′分别表示i组分的临界温度、临界压力、校正临界温度和校正临界压力，M表示酸性天然气中H₂S与CO₂的摩尔分数之和，N表示酸性天然气中H₂S的摩尔分数。

b.温度修正：

T′＝T+1.94(p/2760-2.1×10^-8p²)

②CKB(Car-Kobayshi-Burrows)校正

$T_{pc}^{\′}=T_{pc}-44.4y_{{\mathrm{co}}_2}+72.2y_{H_{2}S}-138.9y_{N_2}$

$p_{pc}^{\′}=p_{pc}+3.034y_{{\mathrm{CO}}_2}+4.137y_{H_{2}S}-1.172y_{N_2}$

式中和分分别表示天然气中H₂S、CO₂和N₂的摩尔分数。

(4)焦耳-汤姆逊系数采用PR方程计算：

$C_J=\frac{R}{C_p}\frac{\left(2r_A-r_{B}T-2r_{B}BT\right)Z-\left(2r_{A}B+r_{B}AT\right)}{\[3Z^2-2\left(1-B\right)Z+\left(A-2B-3B^2\right)\]T}$

其中

A＝r_Ap/R²/T²

B＝r_Bp/R/T

$r_A=0.457235{\mathrm{\α}}_1{R}^2{T}_{pci}^2/p$

r_B＝0.077796RT_pci/p_pci

${\mathrm{\α}}_i={\[1+m_i\left(1-T_{pri}^{0.5}\right)\]}^2$

$m_i=0.3746+1.5423{\mathrm{\ω}}_i-0.2699{\mathrm{\ω}}_i^2$

上述公式中T_pri，p_pri，P_pci和ω_i分别表示组分i的临界温度、对比温度、临界压力和偏心因子。

参照图2，步骤S2：将酸性天然气管线划分为多个计算单元，利用输气管道压降计算公式(或伯努利方程)和温降计算公式，建立酸性天然气管线压力温度耦合计算模型，并计算沿管线的压力温度分布。

(1)管段压降采用下列方法计算：

①输气管道压降计算公式

$p_{in}^2-p_{out}^2=\frac{M^{2}ZRT}{A^2}(\mathrm{\λ}\frac{L}{D}-2ln\frac{p_{in}}{p_{out}})$

式中p_in，p_out，M，D，λ，Z，R，T，L和A分别表示输气管线计算管段的入口和出口压力、天然气质量流量、管道内径、管道水力摩阻系数、天然气压缩系数、天然气气体常数、天然气平均温度、计算管段长度和管道流通面积。式中水力摩阻系数可采用Coledrook-White公式计算：

$\frac{1}{\sqrt{\mathrm{\λ}}}=-2\lg (\frac{k}{3.7D}+\frac{2.51}{\mathrm{Re}\sqrt{\mathrm{\λ}}})$

式中λ，k，D和Re分别为管道水力摩阻系数、内壁粗糙度、内径和雷诺数。

②伯努利方程

$\frac{p_1}{\mathrm{\ρ}g}+\frac{V_1^2}{2g}+Z_1=\frac{p_2}{\mathrm{\ρ}g}+\frac{V_2^2}{2g}+Z_2+(1+\underset{i=1}{\overset{N}{\Σ}}{K}_i)f\frac{L}{D}\frac{V^2}{2g}$

式中p，ρ，g，V，Z，K_i，f，D和L分别表示压力、密度、比重、速度、高度、局部阻力系数、水力摩阻系数、直径和长度，下标1和2分别表示计算管段的入口和出口。式中水力摩阻系数可采用下式计算：

$f=\frac{0.25}{{\[{\log }_{10}(\frac{k}{3.7D}+\frac{5.74}{{\mathrm{Re}}^{0.9}})\]}^2}$

(2)管段温降采用下式计算：

$T_x=T_0+(T_Q-T_0)e^{-ax}-D_i\frac{p_Q-p_Z}{aL}(1-e^{ax})$

其中

$a=\frac{K\mathrm{\π}D}{{\mathrm{MC}}_p}$

式中p_Q和p_Z分别表示输气管线计算管段的起点和终点压力，T_Q表示计算管段中天然气的起点温度，T₀表示土壤或环境温度，x表示计算段任意一点至起点的距离，T_x表示计算管段距离起点x处的温度，K表示总传热系数，D表示管路外径，M表示天然气的质量流量，C_p表示天然气的定压比热，L－表示计算管段长度，D_i表示焦耳汤姆逊系数。

步骤S3：建立酸性天然气中硫元素的溶解度计算模型，根据压力温度参数计算沿管线硫元素的溶解度变化。

固体硫在酸性气体中的溶解度可采用以下两种方法计算：

①Chrastil方法

C_r＝ρ^kexp(a/T+b)

式中C_r，ρ，T分别表示硫的溶解度、气流密度和温度，k，a，b为溶解度系数，取值见表4。

表4 硫在混合酸气中的溶解度系数

②缔合模型

$y_s=\frac{p_s^0\exp \[V^s(p-p_s^0)/\left(RT\right)\]}{{\mathrm{p\φ}}_s^g}$

式中y_s，V^s，p，R和T分别表示固体硫单质的摩尔溶解度、饱和蒸气压、固态摩尔体积、气相逸度系数、气体常数和天然气的热力学温度。

步骤S4：计算从管道起点到硫蒸汽浓度达到其溶解度值的起始位置之间的距离(称为饱和距离)L₁。

步骤S5：根据溶解度的变化确定各区域硫元素的析出量。对于某个计算单元，若入口位置天然气中的含硫量大于出口位置天然气的溶解度，则其差值即为该单元中析出硫元素的量。反之，则认为该计算单元中没有硫元素析出。

步骤S6：采用经典成核理论计算从饱和点到出现可见硫颗粒位置之间的距离(称为过饱和距离)L₂。

(1)采用经典成核理论计算需要克服吉布斯自由能的阀值：

$J={\mathrm{Ke}}^{-\frac{\mathrm{\Δ}G}{k_{B}T}}$

式中J，K，ΔG，k_B和T分别表示成核速率、预指数因子、吉布斯自由能、玻耳兹曼常数和温度，预指数因子K采用下式计算

$K={\mathrm{\θ\α}}_c{v}_s{N}^2\sqrt{\frac{2\mathrm{\σ}}{k_{B}T}}$

式中θ，α，υ，N和σ分别表示不等温因子、凝聚因子、硫元素的分子体积、可凝聚分子数和硫元素的界面张力。

可凝聚分子数采用下式计算：

N＝ρ_My_EN_A

式中ρ_M,y_E和N_A分别表示混合密度、萃取条件下硫元素的摩尔分数和阿伏加德罗常数。

硫元素的界面张力σ由下式计算：

$\mathrm{\σ}={p_c}^{2/3}{T_c}^{1/3}\[0.1207\×\left(1+\frac{T_{br}{\mathrm{lnp}}_c}{1-T_{br}}\right)-0.281\]{(1-T_r)}^{1.222}$

式中p_c、T_c、T_br、T_b和T_r分别表示临界压力、临界温度、对比温度、对比沸点温度和沸点温度。

天然气中存在杂质时，吉布斯自由能通过下式计算：

$\mathrm{\Δ}G=\frac{16{\mathrm{\π\σ}}^3}{3{\left(p_{Sat}-p_{sys}\right)}^2}S\left(\mathrm{\θ}\right)=\frac{16{\mathrm{\π\σ}}^3}{3{\left(p_{Sat}-p_{sys}\right)}^2}\frac{\left(2+\cos \mathrm{\θ}\right){\left(1-\cos \mathrm{\θ}\right)}^2}{4}$

式中p_sat，p_sys和S(θ)分别表示饱和压力、系统压力和杂质形状因子。

(2)当成核速率超出可见成核速率的阀值时，认为首个可见的硫元素颗粒出现，从而由下式确定过饱和距离L₂：

$L_2=L_{J>J_{visible}}$

式中J_visible表示出现可见硫元素颗粒时的成核速率，建议取1个/m³d。

参照图3，步骤S7：通过跟踪硫颗粒的运动轨迹确定硫颗粒的最大迁移距离L₃。

(1)通过求解如下运动方程计算硫颗粒的运动轨迹：

$\frac{d\stackrel{\→}{x_p}}{dt}=\stackrel{\→}{u_p}(\stackrel{\→}{x_p},t)$

式中和分别表示硫颗粒的位置、时间和速度。硫颗粒的速度通过牛顿第二定律计算：

$m_p\left(t\right)\×\frac{d\stackrel{\→}{u_p}}{dt}=\mathrm{\Σ}\stackrel{\→}{F}$

式中m_p(t)表示瞬时硫颗粒的质量，表示总力包括曳力、重力以及布朗力，其中布朗力只有当颗粒尺寸小于2μm时才起作用。式中硫颗粒的瞬时粒径通过质量平衡计算：

$N_{i+1}\mathrm{\ρ}\frac{4{\mathrm{\πr}}_{p,i+1}^3}{3}={\mathrm{\Δm}}_i^{i+1}+N_i\mathrm{\ρ}\frac{4{\mathrm{\πr}}_{p,i}^3}{3}$

式中r_p表示颗粒半径，Δm表示计算单元中硫元素的析出量。

硫颗粒的凝聚特性采用Smoluchowski模型计算：

$\frac{N_i}{N_0}=\frac{1}{1+2{\mathrm{KN}}_0{t}_i}=\frac{1}{1+\frac{k_{B}T}{3\mathrm{\μ}}{N}_0{t}_i}$

式中N₀，K和μ分别表示初始硫颗粒数量、碰撞系数和气体粘度。

假设管壁金属表面液体无液体润湿，采用Tabakoff模型模拟粒子和管道壁之间的相互作用：

${\mathrm{\α}}_2={\mathrm{ctg}}^{-1}\left(\frac{0.95+0.00055{\mathrm{\α}}_1}{1.0-0.02108{\mathrm{\α}}_1+0.000141{\mathrm{\α}}_1^2}{\mathrm{ctg\α}}_1\right)$

$\frac{V_{p2}}{V_{p1}}=(1-0.02108{\mathrm{\α}}_1+0.000141{\mathrm{\α}}_1^2)\sqrt{\frac{1+{\mathrm{ctg}}^2{\mathrm{\α}}_2}{1+{\mathrm{ctg}}^2{\mathrm{\α}}_1}}$

式中α₁和α₂分别表示颗粒与管壁碰撞的入射角和反射角，V_p2和V_p1分别表示粒子碰撞前和碰撞后的速度。

(2)根据计算出的硫颗粒运动轨迹确定硫颗粒的最大迁移距离L₃。

步骤S8：根据硫元素沿管线的析出量和所析出硫颗粒的最大迁移距离计算沿管线的硫沉积位置、硫沉积量以及沉积硫的平均厚度。

(1)根据硫元素沿管线的析出量、所析出硫颗粒的过饱和距离和最大迁移距离计算沿管线的硫沉积位置和单位时间内硫元素的沉积量：

例如，假设在管线入口x₁位置处，析出硫元素y₁g/s，其过饱和距离为L₂₁，最大迁移距离为L₃₁，则在x₁+L₂₁+L₃₁处，有y₁g/s硫元素沉积；

同理，在管线入口x₂位置处，析出硫元素y₂g/s，其过饱和距离为L₂₂，最大迁移距离为L₃₂，则在x₂+L₂₂+L₃₂处，有y₂g/s硫元素沉积；

同理，在管线入口x₃位置处，析出硫元素y₃g/s，其过饱和距离为L₂₃，最大迁移距离为L₃₃，则在x₃+L₂₃+L₃₃处，有y₃g/s硫元素沉积；

……；

同理，在管线入口x_n位置处，析出硫元素y_n g/s，其过饱和距离为L_2n，最大迁移距离为L_3n，则在x_n+L_2n+L_3n处，有y_n g/s硫元素沉积；

……；

管线某位置所有计算出的硫元素沉积量之和，即为此处单位时间内硫元素的沉积量，管线中硫沉积量不为0的区域即为硫沉积位置。

(2)计算酸性天然气管线中硫元素的总沉积量和沉积硫的平均厚度：

管线中某一时间段内硫元素的总沉积量：

$M=\underset{0}{\overset{L}{\∫}}(m\×t)dx$

式中t为时间，x为离管线入口的距离，L为管线总长，M为t时间内的管线中硫元素的总沉积量。

管线中沉积硫的平均厚度：

$\mathrm{\δ}=\frac{M}{\mathrm{\π}D\·L_s}$

式中δ为沉积硫平均厚度，M为管线中硫元素的总沉积量，D为管线内径，L_s为硫沉积区域长度，M为t时间内的管线中硫元素的总沉积量。

上述实施方式是示例性的，不能理解为对本发明的限制，本领域的技术人员在不脱离本发明基本原理的情况下可以对上述实施方式进行变化、修改、替换和变形。本发明的保护范围由所附权利要求极其等同限定。

Claims (4)

1.一种酸性天然气管线硫沉积预测方法，其特征在于，包括如下步骤：

步骤S1：建立酸性天然气热力学性质计算模型，热力学性质包括密度、动力粘度、偏差系数、焦耳-汤姆逊系数；

步骤S2：将酸性天然气管线划分为多个计算单元，利用输气管道压降、温降计算公式，建立酸性天然气管线压力温度耦合计算模型，并计算沿管线的压力温度分布；

步骤S3：建立酸性天然气中硫元素的溶解度计算模型，根据压力温度参数计算沿管线硫元素的溶解度变化；

步骤S4：计算从管道起点到硫蒸汽浓度达到其溶解度值的起始位置之间的距离，称为饱和距离；

步骤S5：根据溶解度的变化确定各区域硫元素的析出量；

步骤S6：采用经典成核理论计算从饱和点到出现可见硫颗粒位置之间的距离，称为过饱和距离；

步骤S7：通过跟踪硫颗粒的运动轨迹确定硫颗粒的最大迁移距离；

步骤S8：根据硫元素沿管线的析出量和所析出硫颗粒的最大迁移距离计算沿管线的硫沉积位置、硫沉积量以及沉积硫的平均厚度。

2.根据权利要求1所述的一种酸性天然气管线硫沉积预测方法，其特征在于：

所述步骤S6具体为：

(1)采用经典成核理论计算需要克服吉布斯自由能的阀值：

$J={\mathrm{Ke}}^{-\frac{\mathrm{\Δ}G}{k_{B}T}}$

式中J，K，ΔG，k_B和T分别表示成核速率、预指数因子、吉布斯自由能、玻耳兹曼常数和温度，预指数因子K采用下式计算

$K={\mathrm{\θ\α}}_c{v}_s{N}^2\sqrt{\frac{2\mathrm{\σ}}{k_{B}T}}$

式中θ，α，υ，N和σ分别表示不等温因子、凝聚因子、硫的分子体积、可凝聚分子数和硫的界面张力。

可凝聚分子数采用下式计算：

N＝ρ_My_EN_A

式中ρ_M,y_E和N_A分别表示混合密度、萃取条件下硫元素的摩尔分数和阿伏加德罗常数；

硫的界面张力σ由下式计算：

$\mathrm{\σ}={p_c}^{2/3}{T}_c^{1/3}\[0.1207\×(1+\frac{T_{br}{\mathrm{lnp}}_c}{1-T_{br}})-0.281\]{(1-T_r)}^{1.222}$

式中p_c、T_c、T_br、T_b和T_r分别表示临界压力、临界温度、对比温度、对比沸点温度和沸点温度；

天然气中存在杂质时，吉布斯自由能通过下式计算：

$\mathrm{\Δ}G=\frac{16{\mathrm{\π\σ}}^3}{3{(p_{Sat}-p_{sys})}^2}S\left(\mathrm{\θ}\right)=\frac{16{\mathrm{\π\σ}}^3}{3{(p_{Sat}-p_{sys})}^2}\frac{(2+\cos \mathrm{\θ}){(1-\cos \mathrm{\θ})}^2}{4}$

式中p_sat，p_sys和S(θ)分别表示饱和压力、系统压力和杂质形状因子；

(2)当成核速率超出可见成核速率的阀值时，认为首个可见的硫元素颗粒出现，从而由下式确定过饱和距离L₂，

$L_2=L_{J>J_{visible}}$

式中J_visible表示出现可见硫元素颗粒时的成核速率，建议取1个/m³d。

3.根据权利要求1所述的一种酸性天然气管线硫沉积预测方法，其特征在于：

所述步骤S7具体为：

(1)通过求解如下运动方程计算硫颗粒的运动轨迹：

$\frac{d\stackrel{\→}{x_p}}{dt}=\stackrel{\→}{u_p}(\stackrel{\→}{x_p},t)$

式中t和分别表示硫颗粒的位置、时间和速度，硫颗粒的速度通过牛顿第二定律计算：

$m_p\left(t\right)\×\frac{d\stackrel{\→}{u_p}}{dt}=\mathrm{\Σ}\stackrel{\→}{F}$

式中m_p(t)表示瞬时硫颗粒的质量，表示总力包括曳力、重力以及布朗力，其中布朗力只有当颗粒尺寸小于2μm时才起作用。式中硫颗粒的瞬时粒径通过质量平衡计算：

$N_{i+1}\mathrm{\ρ}\frac{4{\mathrm{\πr}}_{p,i+1}^3}{3}={\mathrm{\Δm}}_i^{i+1}+N_i\mathrm{\ρ}\frac{4{\mathrm{\πr}}_{p,i}^3}{3}$

式中r_p表示颗粒半径，Δm表示计算单元中硫元素的析出量。

硫颗粒的凝聚特性采用Smoluchowski模型计算：

$\frac{N_i}{N_0}=\frac{1}{1+2{\mathrm{KN}}_0{t}_i}=\frac{1}{1+\frac{k_{B}T}{3\mathrm{\μ}}{N}_0{t}_i}$

式中N₀，K和μ分别表示初始硫颗粒数量、碰撞系数和气体粘度；

假设管壁金属表面液体无液体润湿，采用Tabakoff模型模拟粒子和管道壁之间的相互作用：

${\mathrm{\α}}_2={\mathrm{ctg}}^{-1}\left(\frac{0.95+0.00055{\mathrm{\α}}_1}{1.0-0.02108{\mathrm{\α}}_1+0.000141{\mathrm{\α}}_1^2}{\mathrm{ctg\α}}_1\right)$

$\frac{V_{p2}}{V_{p1}}=(1-0.02108{\mathrm{\α}}_1+0.000141{\mathrm{\α}}_1^2)\sqrt{\frac{1+{\mathrm{ctg}}^2{\mathrm{\α}}_2}{1+{\mathrm{ctg}}^2{\mathrm{\α}}_1}}$

式中α₁和α₂分别表示颗粒与管壁碰撞的入射角和反射角，V_p2和V_p1表示粒子碰撞前和碰撞后的速度；

(2)根据计算出的硫颗粒运动轨迹确定硫颗粒的最大迁移距离L₃。

4.根据权利要求1所述的一种酸性天然气管线硫沉积预测方法，其特征在于：

所述步骤S8具体为：

(1)根据硫元素沿管线的析出量、所析出硫颗粒的过饱和距离和最大迁移距离计算沿管线的硫沉积位置和单位时间内硫元素的沉积量：

假设在管线入口x₁位置处，析出硫元素y₁g/s，其过饱和距离为L₂₁，最大迁移距离为L₃₁，则在x₁+L₂₁+L₃₁处，有y₁g/s硫元素沉积；

同理，在管线入口x₂位置处，析出硫元素y₂g/s，其过饱和距离为L₂₂，最大迁移距离为L₃₂，则在x₂+L₂₂+L₃₂处，有y₂g/s硫元素沉积；

同理，在管线入口x₃位置处，析出硫元素y₃g/s，其过饱和距离为L₂₃，最大迁移距离为L₃₃，则在x₃+L₂₃+L₃₃处，有y₃g/s硫元素沉积；

……；

同理，在管线入口x_n位置处，析出硫元素y_n g/s，其过饱和距离为L_2n，最大迁移距离为L_3n，则在x_n+L_2n+L_3n处，有y_n g/s硫元素沉积；

……；

管线某位置所有计算出的硫元素沉积量之和，即为此处单位时间内硫元素的沉积量，管线中硫沉积量不为0的区域即为硫沉积位置；

(2)计算酸性天然气管线中硫元素的总沉积量和沉积硫的平均厚度：

管线中某一时间段内硫元素的总沉积量：

$M=\underset{0}{\overset{L}{\∫}}(m\×t)dx$

式中t为时间，x为离管线入口的距离，L为管线总长，M为t时间内的管线中硫元素的总沉积量；

管线中沉积硫的平均厚度：

$\mathrm{\δ}=\frac{M}{\mathrm{\π}D\·L_s}$

式中δ为沉积硫平均厚度，M为管线中硫元素的总沉积量，D为管线内径，L_s为硫沉积区域长度，M为t时间内的管线中硫元素的总沉积量。