流速敏感性是指因流体流动速度变化引起储层岩石中微粒运移从而堵塞喉道，导致储层岩石渗透率变化的现象.通过速敏评价实验可以判断开始发生速敏现象的临界流速，其大小对确定油井合理产能及注入速度具有重要意义.目前，岩心速敏实验主要参考SYT 5358—2010《储层敏感性流动实验评价方法》(以下简称标准)，因为标准主要适用于气测渗透率大于1×10^-3 μm²的碎屑岩储层岩样，对气测渗透率小于1×10^-3 μm²的岩样没有说明，所以小于1×10^-3 μm²的低渗透岩心进行速敏实验时，有的参考标准执行(标准中注入流量：0.1，0.25，0.5，0.75，1，2，3，4，5，6 mL/min)，有的按照自己的理解设定注入流量，如吴亚红等^[1]流量设置为0.05，0.07，0.09，0.12，0.15，0.18，0.21，0.25，0.3 mL/min，吴小斌等^[2]流量设置0.099 2，0.245，0.492，0.506 mL/min，谢娜等^[3]流量设置为0.019，0.036，0.068，0.098，0.115，0.131 mL/min.不同的人注入流量不一样，得到的结论也不同，影响对储层客观的认识和评价，因此超低渗透岩心的速度敏感性实验的注入流量设置研究十分必要.

国内外许多学者通过实验^[4-6]或者数值模拟方法^[7-9]对流体微粒运移进行了广泛研究.早期数值模拟方法主要建立在渗透率与孔隙度的关系式基础上.最新研究成果表明，低渗透储层渗流能力与孔隙半径关联性较差，主要受平均喉道半径及其分布形态的控制^[10-11].因此本文根据岩心平均喉道半径及其分布建立新的不等径毛细管数学模型，分析特低、超低渗透岩心速敏伤害规律及注入流量对临界流量判断的影响.

1 速敏伤害机理

速敏实验一般用地层水作为注入水，因此忽略黏土膨胀对渗透率的影响，渗透率的伤害主要由于岩心内部微粒的释放、表面沉积和堵塞造成.微粒从孔隙壁面启动前，必须克服微粒和周围介质之间阻力和阻力矩^[12-16].一旦微粒发生运移，大部分微粒会重新滞留在岩心内部，少量微粒因流体携带作用流出岩心.微粒滞留主要有以下几种机理^{[7, 17-19]}.一是表面沉积，释放的微粒在重力、静电力等作用下重新在孔喉壁面沉积，如果释放的微粒大部分随流体流出岩心，则毛细管半径变大，渗透率变大，如果微粒大部分沉积在岩心内部，会在沉积位置使局部毛细管半径减小，导致渗透率变差.二是单孔堵塞，如果释放的微粒半径相当于孔喉半径或大于孔喉半径将会直接堵塞孔喉，使网络连通性变差.三是形成桥堵，当流体中微粒达到一定浓度时，多个微粒通过架桥方式在孔喉处发生堵塞.根据BarKman和Davidson^[20]对悬浮物在多孔介质渗流过程的研究结果，对于直径为d的悬浮颗粒和直径为D的孔喉，当10d<D时颗粒可以自由通过孔喉，3d<D<10d时颗粒在孔喉处可形成桥堵，3d>D时颗粒可形成卡堵.

2 不等径毛细管模型

由于低渗透岩心渗流能力主要受平均喉道半径及其分布形态控制^[10-11]，岩心内部孔喉结构用半径不相等的一系列毛细管束表示^{[12, 21-22]}，单根毛细管内流体的流动遵循泊肃叶定律，毛细管管径尺寸与实际孔径尺寸大小相当并服从一定的分布规律^[21]，毛细管束模型和真实岩心孔隙体积相等，外观尺寸相等.

2.1 单根毛细管的有效半径

假设半径为r_i的毛细管内任意时刻微粒浓度为C_i，则微粒遵守连续性方程

$ - v\frac{{\partial {C_i}}}{{\partial x}} + D\frac{{{\partial ^2}{C_i}}}{{\partial {x^2}}} + R = \frac{{\partial {C_i}}}{{\partial t}}, $

(1)

式中，v为单根毛细管内平均流速, cm/s，R为单位时间内微粒净释放量，等于微粒释放的质量浓度C_{i, r}和微粒沉积的质量浓度C_{i, t}之差

$ R = \frac{{\partial {C_{i,r}}}}{{\partial t}} - \frac{{\partial {C_{i,t}}}}{{\partial t}}. $

(2)

岩心实际内部空间主要由孔隙和孔喉两部分组成，且微粒流动通道弯曲度较高，当流体中微粒浓度越高时，流速越大，微粒与孔隙表面或孔喉表面碰撞概率越大，越容易在表面沉积，因此表面沉淀速率与流动悬浮液中颗粒质量通量vC成正比，则表面沉淀速率方程

$ \frac{{\partial {C_i}}}{{\partial t}} = {K_t}{C_i}\left( {\alpha + v} \right), $

(3)

$ {C_i}{|_{t = 0}} = {C_m}, $

(4)

式中，α为静止条件下沉淀的静止沉淀系数；K_t为滞留系数，1/cm；C_m为初始状态下介质壁面微粒浓度，g/cm³.

毛细管壁面微粒的释放速率与孔壁剪切力(τ-τ_c)和孔隙表面处沉积微粒浓度成正比^[23]，则

$ \frac{{\partial {C_{i,r}}}}{{\partial t}} = {K_r}\left( {{C_m} - {C_i}} \right)\left( {\tau - {\tau _c}} \right), $

(5)

式中，K_r为微粒释放系数，1/(N·s)，当τ>τ_c时，K_r不等于0，否则K_r等于0；τ为毛细管壁面剪切力；τ_c为微粒平均临界剪切力，即所有微粒的临界剪切力的平均值.

速敏实验时，一旦岩心内部微粒被水动力启动，对流成为影响微粒启动、沉积和堵塞的主要因素，此时扩散项可以忽略，将式(2)—式(5)代入式(1)得

$ \begin{array}{l} - v\frac{{\partial {C_i}}}{{\partial x}} + {K_r}\left( {{C_m} - {C_i}} \right)\left( {\tau - {\tau _c}} \right) - \\ \;\;{K_t}{C_i}\left( {\alpha + v} \right) = \frac{{\partial {C_i}}}{{\partial t}}. \end{array} $

(6)

半径r_t的毛细管连续性方程的初始条件和边界条件满足

$ \left\{ {\begin{array}{*{20}{c}} {{C_i}{|_{t = 0}} = 0}\\ {{C_i}{|_{x = 0}} = 0}\\ {{C_i}{|_{x = L}} = {w_i}} \end{array},} \right. $

式中，w_i为岩心出口端流出的微粒浓度，g/cm³;L为岩心外观长度，m.

由牛顿内摩擦定律和泊肃叶方程，半径r_i的毛细管壁面剪切应力为

$ \tau = \mu \frac{{{\rm{d}}v}}{{{\rm{d}}r}} = \frac{{{P_i} - {P_o}}}{{2L}}{r_{i,\max }}, $

(7)

式中，P_i为入口注入压力；P_o为出口压力.

孔隙中的微粒主要受水流动力、举升力、重力和范德华力作用^{[12, 15-18]}，当微粒的力矩不平衡时，微粒发生迁移，微粒启动的条件为^{[12, 15-16]}

$ {F_{\rm{d}}}{L_{\rm{d}}} + {F_{\rm{l}}}{L_{\rm{n}}} \ge \left( {{F_{\rm{e}}} + {F_{\rm{g}}}} \right){L_{\rm{n}}}, $

(8)

式中，F_d为水流动力，当式(8)取等号时，水流动力即为微粒启动的临界剪切力；L_d为牵力臂；L_n为法向力臂；F_l举升力^[16]；F_e范德华力；F_g微粒受到的重力.

岩心内部微粒质量和平均粒径分布服从下面统计模型^[24]

$ {\rm{f}}\left( {{r_{\rm{d}}} < \overline {{r_{{\rm{d}}i}}} } \right) = \frac{{m\left( {{r_{\rm{d}}} < \overline {{r_{{\rm{d}}i}}} } \right)}}{{{m_T}}} = {\left( {\frac{{\overline {{r_{{\rm{d}}i}}} }}{{{r_{{\rm{d}}\max }}}}} \right)^{\left( {3 - \overline D } \right)}}, $

(9)

式中，r_d为微粒半径，μm；r_dmax为最大可运移微粒半径，μm；m_T为可运移微粒总质量，g；D为分布分维，介于2.270 8～2.948 3之间.

假设岩心在给定注入流量条件下进行流动实验，经过时间T达到平衡状态，在单位体积半径为r_i的毛细管内部微粒净沉积质量为$\sum\nolimits_{n = 0}^{n = T} {\left( {C_{i,t}^n - C_{i,r}^n} \right)} $，净沉积量引起孔喉半径变化δt，根据质量守恒

$ \sum\nolimits_{n = 0}^{n = T} {\left( {C_{i,t}^n - C_{i,r}^n} \right){\rm{\cdot}}\pi {\rm{\cdot}}r_i^2} = {\rho _{\rm{P}}}{\rm{\cdot}}\pi \left( {{r_i} - \delta t} \right){\rm{\cdot}}\delta t. $

(10)

因为δt非常小，忽略平方项，得

$ \delta {t_i} = \frac{{\sum\nolimits_{n = 0}^{n = T} {\left( {C_{i,t}^n - C_{i,r}^n} \right){\rm{\cdot}}{r_i}} }}{{{\rho _{\rm{P}}}}}, $

(11)

$ {r_{i,{\rm{eff}}}}{\rm{ = }}{r_i} - \delta t, $

(12)

式中，C_{i, t}ⁿ为n时刻岩心单位体积内沉积质量，C_{i, r}ⁿ为n时刻岩心单位体积内微粒释放质量; ρ_p为黏土密度，g/cm³；r_{i, eff}为毛细管有效半径, um.

2.2 未堵塞毛细管数量 2.2.1 微粒桥堵引起毛细管数量减少

当微粒体积浓度达到一定值时，微粒在孔喉处可能会发生架桥堵塞，为引入微粒桥堵对渗透率的影响，假设发生桥堵概率P与流动悬浮液中微粒浓度成正比，当微粒浓度为零时，单根毛细管不发生堵塞，当微粒浓度达到最大时，单根毛细管一定发生堵塞，半径r_i毛细管内部发生桥堵概率P表示为

$ P = {{\rm{e}}^{\left( {{C_i} - {C_m}} \right)}} - \frac{{\left( {{C_i} - {C_m}} \right)\left( {{C_i} - {{\rm{e}}^{ - {C_m}}}} \right)}}{{{C_m}}}, $

(13)

半径r_i毛细管堵塞的毛细管数量为

$ {N_{i,{\rm{bp}}}} = {N_i}{\rm{\cdot}}P, $

(14)

式中，N_{i, bp}为半径r_i毛细管由于微粒桥堵引起的毛细管数量的减少，N_i为半径r_i为毛细管初始状态下毛细管数量.

2.2.2 单孔堵塞引起毛细管数量的减少

当启动的颗粒尺寸与孔隙喉道大小相当或大于孔隙喉道时，将发生单孔堵塞.堵塞毛细管数量N_{i, sp}用堵塞孔喉颗粒的总体积与临界尺寸的单一颗粒体积之比来估算^[25-26].

$ {N_{i,{\rm{sp}}}}{\rm{ = }}\left( {{q_i}{f_d}\int\limits_0^t {\frac{{{C_i}{\rm{dt}}}}{{{\rho _{\rm{P}}}}}} } \right)/\left( {\frac{{\pi {d^3}}}{6}} \right){\rm{\cdot}}{10^4}, $

(15)

式中，q_i为毛细管内流量, cm³/s；d为发生卡堵的临界微粒尺寸，大小相当于毛细管直径；f_d为流动相中微粒的体积分数，无量纲；C_i为半径r_i毛细管内流体中的微粒浓度，g/cm³.

2.2.3 未堵塞毛细管数量

由于把岩心看成由一系列不同大小管径组成的平行管束模型^{[12, 20-21]}，管径尺寸与实际孔喉大小相当，实际储层喉道尺寸和数量分布符合一定的分布规律，管径数量的概率密度函数采用对数正态分布

$ f\left( {{r_i}} \right) = \frac{1}{{{r_i}\sigma \sqrt {2\pi } }}{{\rm{e}}^{\frac{{{{\left( {\ln {r_i} - \overline \mu } \right)}^2}}}{{2{\sigma ^2}}}}} $

(16)

式中，μ为统计量均值，σ为统计量方差.

根据不等径毛细管孔隙体积和岩心孔隙体积相等，得到约束条件

$ AL\phi = \sum {\pi r_i^2{\rm{\cdot}}{L_{\rm{e}}}{\rm{\cdot}}N\left( {{r_i}} \right){\rm{\cdot}}{{10}^{ - 8}}} , $

(17)

式中，L为岩心长度，cm; L_e为毛细管实际长度，cm; N(r_i)为半径r_i对应的毛细管数量；ø为岩心孔隙体积，小数；A为岩心横截面积，cm²; r_i为毛细管半径，μm.

根据式(16)和式(17), 得不等径毛细管半径总的毛管数为

$ N = \sum {N\left( {{r_i}} \right)} = {10^8}{\rm{\cdot}}\frac{{A\phi }}{{\sum {\pi r_i^2f\left( {{r_i}} \right)\varsigma } }}, $

(18)

式中，ζ为迂曲度，等于孔隙实际长度和岩心长度的比值.

半径为r_i的毛细管初始毛细管数量为N_i=N·f(r_i)，N_{i, np}表示半径为r_i毛细管未堵塞的毛细管数量，则发生桥堵和卡堵后，未堵塞的毛细管数量为

$ {N_{i,{\rm{np}}}}{\rm{ = }}{N_i} - {N_{i,{\rm{bp}}}} - {N_{i,{\rm{sp}}}}. $

(19)

2.3 渗透率与有效半径和未堵塞毛细管数量关系

根据泊肃叶定律，单根毛细管内流量为

$ {q_i} = {10^{ - 8}} \times \frac{{\pi r_{i,{\rm{eff}}}^4\Delta P}}{{8\mu {L_{\rm{e}}}}}, $

(20)

式中，ΔP为毛细管两端压差，10^-1 MPa; r_i，eff毛细管有效孔喉半径，μm；μ流体黏度，mPa·s.

毛细管内总流量为所有半径r_i毛细管内流量和对应的毛细管数量的乘积

$ Q{\rm{ = }}\sum {{q_i}{\rm{\cdot}}{N_{i,{\rm{np}}}}} . $

(21)

根据达西定律，岩心渗透率为

$ K = {10^3}{\rm{\cdot}}\frac{{\mu LQ}}{{A\Delta P}}, $

(22)

式中，K为渗透率，10^-3μm².

联立式(19)—式(22)，渗透率与有效半径和连通毛细管数量关系式

$ K = \frac{\pi }{{8A\varsigma }}\sum\nolimits_i^N {r_{i,{\rm{eff}}}^4{\rm{\cdot}}{N_{i,{\rm{np}}}}} . $

(23)

3 模型求解

模型求解过程首先根据式(6)计算半径r_i的毛细管内可运移微粒在一定速度下达到稳定时间T的微粒浓度，然后分别根据式(12)、式(19)计算微粒释放和沉积后的有效毛细管的半径r_i，eff和毛细管未堵塞的毛细管数量N_{i, np}，最后通过式(23)求得岩心渗透率.

对半径r_i的毛细管内微粒连续方程，式(6)是一阶拟线性偏微分方程，采用Lax-Friedrichs格式进行差分

$ \begin{array}{l} C_j^{n + 1} = \left( {0.5 + \lambda v/2} \right)C_{j + 1}^n + \left( {0.5 - \lambda v/2} \right)C_{j - 1}^n - \\ \;\;\;\;\;\;\;\;\;\Delta t\left[ {{K_r}\left( {\tau - {\tau _{\rm{c}}}} \right) + {K_t}\left( {\alpha + v} \right)} \right]C_j^n + \\ \;\;\;\;\;\;\;\;\;\Delta t{K_r}\left( {\tau - {\tau _{\rm{c}}}} \right){C_m}. \end{array} $

(24)

当vλ≤1时，格式稳定.微粒沉积方程(3)和微粒释放方程(5)采用中心差分格式：

$ \begin{array}{l} C_{r,j}^{n + 1} = C_{r,j}^n + {K_r}\Delta t\left( {\tau - {\tau _{\rm{c}}}} \right)\left( {{C_m} - C_j^n} \right),\\ C_{t,j}^{n + 1} = C_{t,j}^n + {K_t}C_j^n\left( {\alpha + v} \right)\Delta t. \end{array} $

4 实例模拟分析

由流速敏感性引起的渗透率伤害率按下式计算

$ {D_{vn}} = \frac{{\left| {{K_n} - {K_i}} \right|}}{{{K_i}}} \times 100\% , $

(25)

式中，D_vn不同流速下所对应的岩样渗透伤害率；K_n实验中不同流速下所对应的岩心渗透；K_i实验中最小流速所对应的岩样渗透率.首先采用标准中注入流量，研究微粒释放、表面沉积和堵塞等3种因素对3种不同低渗透率岩心(K₀=0.52×10^-3μm²;5×10^-3μm²;25×10^-3μm²)的渗透率的影响，模拟参数见表 1.

4.1 微粒释放对渗透率的影响

令表面释放系数K_r=2.5，表面沉积系数K_t=0，堵塞管子数N_p=0，模拟微粒释放对岩心渗透率的影响，结果如图 1所示.从图 1可以看出，渗透率伤害率与注入流量呈较好的线性关系，注入流量越大，渗透率伤害率越大；岩心渗透率越低，渗透率伤害率越大.这是因为内部微粒脱离孔喉表面，孔喉有效半径变大，岩心渗透能力增强，渗透率变大.在注入流量为6 mL/min时，0.52×10^-3μm²岩心的渗透率增加40%, 5×10^-3μm²岩心的渗透率增加3.7%，25×10^-3μm²岩心的渗透率仅增加0.7%.

4.2 表面沉积对渗透率的影响

令表面释放系数K_r=2.5，堵塞管子数N_p=0，沉积系数K_t=0.01，为了方便观察表面沉积对渗透率的影响，将模拟结果和图 1数据进行差值计算，结果如图 2所示.从图 2可以看出，注入流量越大，渗透率伤害率越大，岩心渗透率降低.当注入流量为6 mL/min时，0.52×10^-3μm²岩心的渗透率降低41.07%, 5×10^-3μm²岩心的渗透率降3.72%，25×10^-3μm²岩心的渗透率降低0.73%.

4.3 微粒堵塞对渗透率的影响

令表面释放系数K_r=2.5，沉积系数K_t=0.01，为了方便观察堵塞对渗透率的影响，将模拟结果和图 1的差值作图，结果如图 3所示.从图 3可以看出，注入流量越大，堵塞作用对岩心渗透率伤害越大，渗透率减小.当注入流量为6 mL/min时，0.52×10^-3μm²岩心的渗透率降低17.04%，5×10^-3μm² 岩心的渗透率降低1.44%，25×10^-3μm²岩心的渗透率降低0.27%.

总结微粒释放、表面沉积和堵塞对渗透率伤害率的影响，结果如表 2所示.结果表明，表面沉积、微粒释放和堵塞对渗透率的影响程度为：表面沉积>微粒释放>堵塞，其中表面沉积和微粒释放是渗透率伤害的主要因素.渗透率越低，微粒运移对渗透率伤害越明显.这是因为渗透率越小，毛细管半径越小，微粒运移对毛细管半径影响越大.值得注意的是，模拟中的堵塞是指微粒完全堵塞喉道，导致模型中对应毛细管半径的毛细管数量的减少.

4.4 注入流量设置对临界流量的影响

为研究低渗透岩心合理注入流量，模拟7种不同低渗透率岩心，分别用标准注入流量(记为方案1)和新的注入流量(0.025, 0.05, 0.075, 0.1, 0.25, 0.5, 0.75, 1.0, 1.5, 2.0, 3.0, 4.0 mL/min，记为方案2)，进行速敏实验.取释放系数K_r=2.5;滞留系数K_t=0.01.当注入流量为方案1时，结果如图 4所示；当注入流量为方案2时，结果如图 5所示.

根据渗透率变化超过20%时，前一个注入流量为临界流量，7种低渗透率岩心临界流量，判断结果如表 3所示.模拟结果表明：方案2和方案1相比，当渗透率小于等于 1×10^-3μm²时，岩心编号2, 3, 4, 5的临界流量判断结果偏小，且平均减小1个间隔；而渗透率为5×10^-3μm²和2.5×10^-3μm²的岩心，临界流量判断结果一致；但是渗透率为0.05×10^-3μm²时，临界流量判断结果偏差较大.因此，当渗透率大于1×10^-3μm²时，模拟结果与标准推荐注入流量一致，当渗透率为(0.1~1)×10^-3μm²时，不建议继续参考标准注入流量，建议在标准中最小流量0.1 mL/min之前添加1或2个点进行实验.

5 结论

本文基于毛细管半径分布、毛细管数量分布及微粒分布，综合考虑微粒释放、表面沉积及堵塞等3种因素对渗透率的影响，建立新的不等经毛细管模型，给出渗透率伤害率的数值计算方法并进行实例仿真，分析3种因素对不同渗透率岩心渗透率伤害的影响程度以及注入流量设置对低渗透率岩心的临界流量判断的影响，得到以下两点认识：

1)模拟结果表明，表面沉积、微粒释放和堵塞对低渗透岩心渗透率伤害率影响程度为：表面沉积>微粒释放>堵塞，其中微粒释放和表面沉积是渗透率伤害的主要因素；渗透率越低，微粒运移对渗透率伤害越明显.本模型中堵塞是指孔隙完全封堵，其结果导致毛细管数量减少.

2)1×10^-3μm²以上的低渗透率岩心，参考标准注入流量时，临界流量判断没有影响； 小于等于1×10^-3μm²渗透率岩心，继续参考标准中注入流量进行速敏实验时，临界流量判断结果偏大，超低渗透率岩心(渗透率为(0.1~1)×10^-3μm²)建议在0.1 mL/min以下添加1或2个注入点.