0 引言

大牛地气田位于鄂尔多斯盆地伊陕斜坡东北段，属于典型的低渗致密砂岩气藏，气层有效厚度薄，剩余未动用储量品位差，必须采用水平井多级压裂才能实现经济开发。如何同时考虑气液两相流和压裂施工及气井产水对渗透率的影响，获得多因素影响下的压裂水平井产能，现有文献还未见相关的报道。目前关于压裂水平井产能的计算，国内外学者^[1-6]已经做了大量的研究工作，路爽等^[7]基于偏心直井产能公式，利用汇源反应和叠加原理建立了压裂水平井产能公式，但没有考虑气水两相渗流的影响；袁琳等^[8]、李旭成等^[9]和孙恩慧等^[10]虽然建立了气水同产压裂水平井产能计算方法，但是没有考虑渗透率各向异性及水平段分布方向对产能的影响，且没有考虑压裂施工和气井产水对基质有效渗透率的影响，而用来作产能分析的商业化试井软件，所基于的理论基础也都是单相流动，没有考虑气液两相渗流的影响^[11-13]，同时致密气藏试井压力恢复需要的关井时间一般较长，而对于大牛地气田这种低压产水气井，长时间关井存在井筒积液造成气井水淹的风险。

基于保角变换理论和气水两相渗流理论，同时将基质有效渗透率作为变量来考虑压裂施工和气井产水对储层基质渗透率的影响，利用势函数叠加原理建立压裂水平井气液两相产能方程，分析泄气边界、水气比、渗透率各向异性程度、水平井井筒与渗透率主值方向夹角对气井产能和生产压差的影响，以期为气井产能预测和生产控制提供依据。

1 模型的建立及求解 1.1 模型的假设

建立的模型有以下假设条件：①压裂水平井位于顶底封闭的无限大的气藏中心；②假设垂直裂缝完全穿过产层，流体先从基质流入裂缝，再沿裂缝进入水平井筒；③不考虑流体在井筒内部流动的压力损失，裂缝垂直于井筒且沿井筒方向等间距分布，相互平行；④考虑应力敏感、压裂施工和地层产水对渗透率的影响；⑤考虑渗透率各向异性对气井产能的影响。

1.2 渗透率各向异性的处理

压裂水平井保角变换^[14-15]如图 1所示，为了计算方便，首先须将各向异性储层W平面转化为各向同性储层W1平面^[16-17]，在W平面内考虑渗透率主值方向为沿y轴方向，x轴方向渗透率为K_x，y轴方向渗透率为K_y，压裂水平井长度为L，井筒与渗透率主值方向的夹角为θ₁，裂缝半长为Lf，裂缝宽度为Wf，则在W平面内气体流动的渗流微分方程为

$ {K_x}\frac{{{\partial ^2}\psi }}{{\partial {x^2}}} + {K_x}\frac{{{\partial ^2}\psi }}{{\partial {y^2}}} = {\mu _{\rm{g}}}\varphi {C_{\rm{t}}}\frac{{\partial \psi }}{{\partial t}} $

(1)

式中：K_x为x轴方向渗透率，mD；K_y为y轴方向渗透率，mD；μ_g为天然气黏度，mPa·s；φ为储层平均孔隙度；C_t为综合压缩系数，MPa^-1。

令$X = x\sqrt {\frac{K}{{{K_x}}}} , Y = y\sqrt {\frac{K}{{{K_x}}}} , K = \sqrt {{K_x}{K_y}} $，带入式（1）可得

$ \frac{{{\partial ^2}\psi }}{{\partial {x^2}}} + \frac{{{\partial ^2}\psi }}{{\partial {y^2}}} = \frac{1}{\eta }\frac{{\partial \psi }}{{\partial t}} $

(2)

其中：

$ \eta = \frac{{\sqrt {{K_x}{K_y}} }}{{{\mu _{\rm{g}}}\varphi {C_{\rm{t}}}}} $

(3)

这样就把W平面上各向异性的气藏转换成为了W1平面上各向同性的气藏，转换后的压裂水平井长度L₁为

$ {L_1} = 2\sqrt {{{\left( {\sqrt {\frac{K}{{{K_x}}}} \frac{L}{2}\sin {\theta _1}} \right)}^2} + {{\left( {\sqrt {\frac{K}{{{K_x}}}} \frac{L}{2}\cos {\theta _1}} \right)}^2}} $

(4)

由于人工裂缝始终垂直于井筒，所以转换后的裂缝宽度Wf₁为

$ W{f_1} = Wf\frac{{{L_1}}}{L} $

(5)

转换后的裂缝半长Xf₁为

$ X{f_1} = \sqrt {{{\left( {\sqrt {\frac{K}{{{K_x}}}} \frac{{Xf}}{2}\cos {\theta _1}} \right)}^2} + {{\left( {\sqrt {\frac{K}{{{K_x}}}} \frac{{Xf}}{2}\sin {\theta _1}} \right)}^2}} $

(6)

转换后井筒与y轴的夹角θ₂为

$ {\theta _2} = \arctan \left( {\frac{{\sqrt {\frac{K}{{{K_x}}}} \frac{L}{2}\sin {\theta _1}}}{{\sqrt {\frac{K}{{{K_y}}}} \frac{L}{2}\cos {\theta _1}}}} \right) $

(7)

引入保角变换系数

$ z = {W_1}{{\rm{e}}^{{\rm{i}}{\theta _2}}} $

(8)

令$ z = \zeta + {\rm{ i}}\eta , {W_1} = X + {\rm{ i}}Y$，代入式（1）可得W1平面与Z平面的对应关系为

$ \left\{ \begin{array}{l} \zeta = X\cos {\theta _2} - Y\sin {\theta _2}\\ \eta = X\sin {\theta _2} + Y\cos {\theta _2} \end{array} \right. $

(9)

各向异性储层通过2次变换后变为图 1所示的Z平面上水平井井筒在η轴上的渗透率各向同性储层，W平面上任意一点与Z平面上对应点的关系为

$ \left\{ \begin{array}{l} \zeta = x\sqrt {\frac{K}{{{K_x}}}} \cos {\theta _2} - y\sqrt {\frac{K}{{{K_y}}}} \sin {\theta _2}\\ \eta = x\sqrt {\frac{K}{{{K_x}}}} \sin {\theta _2} + y\sqrt {\frac{K}{{{K_y}}}} \cos {\theta _2} \end{array} \right. $

(10)

1.3 气藏到裂缝的流动

对各向同性的均质气藏，水平井位于顶底封闭、水平方向无限大的气藏中心生产时，考虑气液两相流和应力敏感的影响，定义气水两相广义拟压力为^[18-19]

$ \varphi \left( p \right) = \int_0^p {\left( {\frac{{{k_{{\rm{rg}}}}{\rho _{\rm{g}}}}}{{{\mu _{\rm{g}}}}} + \frac{{{k_{{\rm{rw}}}}{\rho _{\rm{w}}}}}{{{\mu _{\rm{w}}}}}} \right){{\rm{e}}^{ - \alpha \left( {{p_i} - p} \right)}}{\rm{d}}p} $

(11)

式中：ρ_g为气相密度，g/cm³；ρ_w为水相密度，g/cm³；α为应力敏感系数，MPa^-1；μ_w为地层水黏度，mPa·s；p_i为原始地层压力，MPa；

考虑缝间干扰的影响，则N条裂缝同时生产时，Z平面上地层内任意一点的气水两相广义拟压力分布为^[8]

$ \begin{array}{l} \psi \left( p \right) = \frac{{1.842\left( {{\rho _{{\rm{gsc}}}} + {R_{{\rm{wg}}}}{\rho _{{\rm{wsc}}}}} \right)}}{{Kh}}\sum\limits_{i = - {N_0}}^{{N_0}} {{q_{{\rm{gscfj}}}}\arccos h\frac{1}{{\sqrt 2 }}} \left\{ {1 + \frac{{{\xi ^2}}}{{Xf_1^2}} + {{\left( {\frac{{\eta - id}}{{Xf}}} \right)}^2} + } \right.\\ \;\;\;\;\;\;\;\;\;\;\;{\left. {\sqrt {{{\left[ {1 + \frac{{{\xi ^2}}}{{Xf_1^2}} + \frac{{{{\left( {\eta - id} \right)}^2}}}{{Xf_1^2}}} \right]}^2} - 4\frac{{{\xi ^2}}}{{Xf_1^2}}} } \right\}^{0.5}} \end{array} $

(12)

$ d = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} {L_1}/N,\\ {L_1}/\left( {2N} \right), \end{array}&\begin{array}{l} N\;为奇数\\ N\;为偶数 \end{array} \end{array}} \right. $

(13)

$ j = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} i + {N_0} + 1,\\ \frac{{i + {N_0}}}{2} + 1, \end{array}&\begin{array}{l} N\;为奇数\\ N\;为偶数 \end{array} \end{array}} \right. $

(14)

$ {N_0} = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} \left( {N - 1} \right)/2,\\ N - 1, \end{array}&\begin{array}{l} N\;为奇数\\ N\;为偶数 \end{array} \end{array}} \right. $

(15)

式中：q_gscfj为标准状况下第j条裂缝的气相流量，m³/d；d为裂缝间距或半间距，m；N为裂缝条数。

裂缝条数为奇数时，i从-N₀开始以1为步长递增；裂缝条数为偶数时，i从-N₀开始以2为步长递增，Z平面内第j条裂缝所在位置坐标为（0，md），将坐标带入式（12）可得第j条裂缝处的气水两相广义拟压力为

$ \begin{array}{l} \psi \left( {0,md} \right) = \frac{{1.842\left( {{\rho _{{\rm{gsc}}}} + {R_{{\rm{wg}}}}{\rho _{{\rm{wsc}}}}} \right)}}{{Kh}} \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\sum\limits_{i = - {N_0}}^{{N_0}} {\arccos \sqrt {1 + \frac{{\left( {md - id} \right)}}{{Xf_1^2}}} } + {C_1} \end{array} $

(16)

$ m = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} - {N_0} + j - 1,\\ - {N_0} + 2\left( {j - 1} \right), \end{array}&\begin{array}{l} N\;为奇数\\ N\;为偶数 \end{array} \end{array}} \right. $

(17)

令W平面内沿水平井井筒方向的泄气边界为r_e，则由式（7）和式（10）可得对应Z平面上的泄气边界R_e为

$ \begin{array}{l} {R_{\rm{e}}} = {r_{\rm{e}}}\sin {\theta _1}\sqrt {\frac{K}{{{K_x}}}} \sin \left[ {\arctan \left( {\frac{{\sqrt {\frac{K}{{{K_x}}}} \frac{L}{2}\sin {\theta _1}}}{{\sqrt {\frac{K}{{{K_y}}}} \frac{L}{2}\cos {\theta _1}}}} \right)} \right] + \\ \;\;\;\;\;\;\;{r_{\rm{e}}}\cos {\theta _1}\sqrt {\frac{K}{{{K_y}}}} \cos \left[ {\arctan \left( {\frac{{\sqrt {\frac{K}{{{K_x}}}} \frac{L}{2}\sin {\theta _1}}}{{\sqrt {\frac{K}{{{K_y}}}} \frac{L}{2}\cos {\theta _1}}}} \right)} \right] \end{array} $

(18)

Z平面内在泄气边界（0，R_e）处，气水两相广义拟压力为

$ \begin{array}{l} \psi \left( {0,{R_{\rm{e}}}} \right) = \frac{{1.842\left( {{\rho _{{\rm{gsc}}}} + {R_{{\rm{wg}}}}{\rho _{{\rm{wsc}}}}} \right)}}{{Kh}} \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\sum\limits_{i = - {N_0}}^{{N_0}} {\arccos \sqrt {1 + \frac{{\left( {{R_{\rm{e}}} - id} \right)}}{{Xf_1^2}}} } + {C_1} \end{array} $

(19)

式中：R_e为Z平面上的泄气半径，m；ρ_gsc为标准状况下气相密度，g/cm³；ρ_wsc为标准状况下水相密度，g/cm³；R_wg为生产水气比；Xf₁位Z平面上的裂缝半长，m。

利用化简公式$ {\rm{arccos}}\sqrt {1{\rm{ }} + {x^2}} {\rm{ }} = {\rm{ ln}}(x + \sqrt {1{\rm{ }} + {x^2}} )$，式（19）与式（18）相减可得Z平面内任意一条裂缝边界处压力满足：

$ \begin{array}{l} {\psi _{\rm{e}}}{\psi _{{\rm{fj}}}} = \frac{{1.842\left( {{\rho _{{\rm{gsc}}}} + {R_{{\rm{wg}}}}{\rho _{{\rm{wsc}}}}} \right)}}{{Kh}}\sum\limits_{i = - {N_0}}^{{N_0}} {{q_{{\rm{gscfj}}}}} \cdot \\ \;\;\;\;\;\;\;\;\;\;\ln \frac{{\left( {\frac{{{R_{\rm{e}}}}}{{X{f_1}}} - \frac{{id}}{{X{f_1}}}} \right) + \sqrt {1 + \frac{{{{\left( {{R_{\rm{e}}} - id} \right)}^2}}}{{Xf_1^2}}} }}{{\left| {\frac{{md}}{{X{f_1}}} - \frac{{id}}{{X{f_1}}}} \right| + \sqrt {1 + \frac{{{{\left( {md - id} \right)}^2}}}{{Xf_1^2}}} }} \end{array} $

(20)

1.4 裂缝到井筒的流动

对于裂缝内部的气液两相流，由于裂缝半长远大于井筒半径，流体从裂缝边缘向井筒内的流动可以等效为边界压力为P_fj及井底压力为P_wfj的平面径向流（图 2）。

等效渗流半径R为

$ R = \sqrt {2hX{f_1}} $

(21)

定义裂缝内部气水两相拟压力为

$ m\left( p \right) = \int_{{p_i}}^p {\left( {\frac{{{k_{{\rm{rg}}}}{\rho _{\rm{g}}}}}{{{\mu _{\rm{g}}}}} + \frac{{{k_{{\rm{rw}}}}{\rho _{\rm{w}}}}}{{{\mu _{\rm{w}}}}}} \right){\rm{d}}p} $

(22)

考虑裂缝内部高速非达西渗流时的气水两相运动方程分别为^[8]

$ \frac{{{k_{{\rm{rg}}}}{\rho _{\rm{g}}}}}{{{\mu _{\rm{g}}}}}{\rm{d}}p = \frac{{1.842{q_{{\rm{gscf}}}}{\rho _{{\rm{wsc}}}}}}{{{k_{\rm{f}}}w}}\frac{{{\rm{d}}r}}{r} + \frac{{2.59 \times {{10}^{ - 4}}}}{{k_{{\rm{rg}}}^{0.5}{\mu _{\rm{g}}}}}\frac{{\rho _{{\rm{gsc}}}^2q_{{\rm{gscf}}}^2}}{{k_{\rm{f}}^{1.5}{w^2}}}\frac{{{\rm{d}}r}}{{{r^2}}} $

(23)

$ \frac{{{k_{{\rm{rw}}}}{\rho _{\rm{w}}}}}{{{\mu _{\rm{w}}}}}{\rm{d}}p = \frac{{1.842{q_{{\rm{wsc}}}}{\rho _{{\rm{wsc}}}}}}{{{k_{\rm{f}}}w}}\frac{{{\rm{d}}r}}{r} $

(24)

在对应区间上积分，得到气水两相渗流条件下裂缝-井筒渗流模型为

$ \begin{array}{l} {m_{{\rm{fj}}}} - {m_{{\rm{wfj}}}} = \frac{{1.842\left( {{q_{{\rm{gscf}}}}{\rho _{{\rm{wsc}}}}{R_{{\rm{wg}}}}} \right)}}{{{k_{\rm{f}}}w}}\ln \frac{R}{{{r_{\rm{w}}}}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\frac{{2.59 \times {{10}^{ - 4}}\rho _{{\rm{gsc}}}^2q_{{\rm{gscfj}}}^2}}{{k_{\rm{f}}^{1.5}w_g^2}}\int_{{r_{\rm{w}}}}^R {\frac{{f\left( p \right)}}{{{r^2}}}{\rm{d}}r} \end{array} $

(25)

$ f\left( p \right) = \frac{1}{{k_{{\rm{rg}}}^{0.5}{\mu _{\rm{g}}}}} $

(26)

1.5 气藏-裂缝-井筒的耦合

联立式（20）和式（26）可得，考虑渗透率各向异性和水平井井筒与渗透率主值方向的夹角时，致密砂岩气藏压裂水平井任意一条裂缝处的产能方程为

$ \begin{array}{l} {\psi _{\rm{e}}} - {\psi _{{\rm{wfj}}}} = \frac{{1.842\left( {{\rho _{{\rm{gsc}}}} + {R_{{\rm{wg}}}}{\rho _{{\rm{wsc}}}}} \right)}}{{Kh}}\sum\limits_{i = - {N_0}}^{{N_0}} {{q_{{\rm{gscfj}}}}} \cdot \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\ln \frac{{\left( {\frac{{{R_{\rm{e}}}}}{{X{f_1}}} - \frac{{id}}{{X{f_1}}}} \right) + \sqrt {1 + \frac{{{{\left( {{R_{\rm{e}}} - id} \right)}^2}}}{{Xf_1^2}}} }}{{\left| {\frac{{md}}{{X{f_1}}} - \frac{{id}}{{X{f_1}}}} \right| + \sqrt {1 + \frac{{{{\left( {md - id} \right)}^2}}}{{Xf_1^2}}} }} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\frac{{1.842\left( {{q_{{\rm{gscfj}}}} + {\rho _{{\rm{wsc}}}}{R_{{\rm{wg}}}}} \right)}}{{{k_{\rm{f}}}w}}\ln \frac{R}{{{r_{\rm{w}}}}} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\frac{{2.59 \times {{10}^{ - 4}}\rho _{{\rm{gsc}}}^2q_{{\rm{gscfj}}}^2}}{{k_{\rm{f}}^{1.5}w_{\rm{g}}^2}}\int_{{r_{\rm{w}}}}^R {\frac{{f\left( p \right)}}{{{r^2}}}{\rm{d}}r} \end{array} $

(27)

2 压裂施工及地层产水的影响

致密气藏生产过程中普遍存在产水现象^[20-22]，一方面是压裂施工过程中压裂液侵入地层引起的，另一方面是气藏本身在生产过程中存在可动水。压裂液或可动水在近井地带的聚集，以及压裂施工本身产生的次级裂缝，必然会造成压裂改造区储层有效渗透率发生改变，因此在利用产能公式求解气井产能时，须要利用井底流压、产气量测试数据，设置合适的有效渗透率范围，对气液两相流产能方程设置迭代求解精度进行牛顿迭代求解，通过插值拟合法得到受压裂施工和气井产水影响下的储层有效渗透率，然后再通过牛顿迭代法求解气井产能，具体求解思路如图 3所示。

3 实例计算及产能影响因素分析

大牛地气田某气井实际储层参数和气井参数如下：实钻水平段长度为1 200 m，水平井井筒半径为0.1 m，气层有效厚度为12 m，气藏温度为363.5 K，天然气相对密度为0.7，试验测得应力敏感指数为0.02 MPa^-1。实际压裂裂缝12条，裂缝半长为150 m，裂缝宽度为5 mm，裂缝渗透率为4 000 mD，压裂返排14天后进行流压测试，折算地层中部流压为14.73 MPa，地层静压力为25.98 MPa，平均产气量为8.39万m³/d。

气井生产水气比取邻井压裂返排后正常生产期间的水气比1:1 m³/万m³，受压裂施工和地层产水的影响，设置储层基质有效渗透率为0.001~ 2.000 mD，渗透率拟合步长为0.005 mD，通过牛顿迭代求解和插值拟合法得到储层有效渗透率为0.34 mD，气井无阻流量为12.67万m³/d（图 4）。根据IPR指示曲线法进行配产，IPR曲线上偏离直线段的点即为合理配产点（图 4），气井初期合理产气量为3.2万m³/d，气井生产曲线如图 5所示，初期配产4万m³/d时气井压降速率较快，生产40天后配产调整为3.2万m³/d左右，压降速率明显变缓，并实现稳产2年左右，说明按照本文方法计算得到的无阻流量可靠性高。

3.1 泄气边界

由于致密气藏压力传递较慢，随着生产时间的增加，气井的泄气边界逐渐增加，同时受砂体边界的约束，泄气边界不可能是无限大的，因此定义无因次泄气边界为沿水平段方向泄气边界到水平段中点的距离与压裂水平井水平段长度之比，分别设置无因次泄气边界为0.46，0.50，0.55，0.60，0.65，0.70，0.75，研究泄气边界和气井产能对生产压差的影响，计算结果如图 6所示，可以看出：①生产初期无因次泄气边界小于0.55时，气井产能主要受压裂改造区控制，随着气井产能的增加，生产压差几乎呈线性增加；②当无因次泄气边界大于0.55时，受储层供气能力不足的影响，随着气井产能的增加，生产压差呈下凹型快速增加。因此，在气井生产过程中，为避免不必要的压力损失，生产初期在保证不出砂的前提下可适当提高配产生产，而在生产的中后期由于压力损失随配产增加呈下凹型快速增长，宜适当降低配产生产。

3.2 水气比

分别设置水气比为1:1 m³/万m³，2:1 m³/万m³，3:1 m³/万m³，4:1 m³/万m³，5:1 m³/万m³，计算不同产气量条件下的气井生产压差，得到图 7所示图版，可以看出：①相同水气比条件下，气井产能越大所需生产压差越大；②相同气井产能条件下，水气比越大，气井所需生产压差越大。这是由于气井产水量的增加，必然会造成储层中两相渗流阻力的增加，从而导致附加的压力损失增加，因此气井生产过程中在保证携液能力的同时，应适当降低气井产能或采取必要的控水排水措施，以降低生产过程中的压力损失。

3.3 渗透率各向异性程度

图 8为水平段方向与K_y方向平行时，不同渗透率各向异性程度K_y/K_x和不同气井产能时的生产压差，可看出：①渗透率各向异性程度K_y/K_x < 1时，水平段垂直于渗透率主值K_x方向，此时各向异性程度越小，储层非均质性越强，相同产气量时气井所需生产压差越大，且气井产能越大，生产压差随各向异性程度的变化越明显；②渗透率各向异性程度K_y /K_x > 1时，水平段平行于渗透率主值K_y方向，此时随着渗透率各向异性程度的增加，气井生产压差略有减小。因此，从渗透率各向异性程度对压裂水平井储层压力损失的影响来考虑，布井时应尽可能使水平井水平段沿着渗透率主值方向，即河道方向分布。

3.4 井筒与渗透率主值方向夹角的影响

根据渗透率实验结果取渗透率各向异性程度K_y/K_x为1.5 ^[23-24]，分别设置水平段与渗透率主值Ky方向的夹角为0°，15°，30°，45°，60°，75°，90°，计算不同气井产能条件下的气井生产压差，得到图 9所示图版，可以看出：①水平段与渗透率主值方向的夹角θ < 30°时，相同产气量条件下的气井生产压差几乎不变；②当水平段与渗透率主值方向的夹角θ > 30°时，随着夹角θ的增加气井所需生产压差逐渐增大，且气井产能越大，压差增长越明显。因此，从降低压裂水平井储层压力损失的角度来考虑，布井时应尽可能将水平井水平段与渗透率主值方向的夹角控制在30°。

4 结论

（1）建立了多因素影响条件下的压裂水平井气液两相产能方程，并运用计算结果指导了气井生产初期的合理配产，实现了气井的稳产，可靠性高。

（2）泄气边界的敏感性分析表明，无因次泄气边界大于0.55时，受储层供气能力不足的影响，气井生产压差随配产呈下凹型快速增长，宜适当降低配产生产。

（3）生产水气比的敏感性分析表明，在保证气井携液能力的同时，宜适当降低气井产能或采取必要的控水排水措施，以降低生产过程中的压力损失。

（4）渗透率各向异性程度和水平段与渗透率主值方向夹角的敏感性分析表明，从降低压裂水平井储层压力损失的角度来考虑，布井时应尽可能将水平井水平段与渗透率主值方向的夹角控制在30°。