简介:面向石油工程、试井分析领域研究人员与现场工程师的裂缝性气藏分支水平井试井模型解析文档,完整呈现从模型建立、Laplace变换求解到Stehfest数值反演的全流程。文档内嵌可运行的Python代码及逐段解释,覆盖典型曲线绘制、九个渗流阶段划分(含纯井储段、裂缝系统线性流段、第二拟径向流段、窜流段等)、分支参数敏感性分析及流动阶段自动识别,并给出分支参数优化、试井测试设计和数据解释流程等工程建议。此外还特别分析了分支长度与分支间距对压力动态曲线的影响,为理解复杂井型渗流机制提供直观依据。包体为单个docx文件,大小约60KB,结构紧凑,便于直接阅读与代码复用。已有73人学习下载,适合从事油气田开发、试井解释的技术人员作为理论推导与编程实现的参考。通过这份文档,读者可掌握鱼骨型分支水平井复杂渗流机制的压力动态特征,并可直接修改代码参数用于自家井况的初步模拟。
1. 裂缝性气藏分支水平井试井模型:为什么一口井的压降曲线,能在软件里算成两种答案
裂缝性气藏分支水平井试井模型,解决的是这样一个现实问题:储层是“基质储气、裂缝导气”的双重介质,井筒是带多个分支的复杂结构,两件事叠加之后,常规单孔介质、单分支水平井的解析解已经不适用,同一个压降数据用不同模型解释,渗透率和表皮能差出几倍。围绕裂缝性气藏分支水平井试井模型的构建、求解及应用,本文会从双重介质渗流方程讲起,给出完整的 Python 求解代码和典型曲线拟合流程。适合要做试井解释的气藏工程师、研究分支水平井压力动态的研究生,以及想给自编解释程序补上双孔介质模块的开发者。
2. 双孔介质与分支水平井的物理图景:压力波在三套系统里的传播路径
2.1 为什么裂缝性气藏必须用双孔介质模型,而不是把裂缝“平均”掉
裂缝性气藏最典型的矛盾是:基质岩块有着绝大部分储集空间,渗透率极低;裂缝系统只占一小部分孔隙度,却承担了几乎全部导流能力。如果按常规单孔介质处理,把基质和裂缝的物性加权平均,得到的是一个“看起来合理”的渗透率,但它无法解释压力导数曲线上那个标志性的下凹段——裂缝先供液、基质后补给的滞后响应。Warren-Root 在 1963 年提出的双孔介质模型,用两个重叠的连续介质来描述这种结构:裂缝系统作为主流动通道,基质系统作为分散的储集体,两者之间通过拟稳态窜流发生质量交换。
这套模型在试井里的核心价值是引入了两个无因次参数:储容比(裂缝系统储容占总储容的比例)和窜流系数(基质与裂缝间流动能力对比)。裂缝性气藏分支水平井试井模型之所以难做,正是因为它要在双孔介质的基础上再叠加水平井的分支结构——压力波从井筒出发,先经过近井裂缝网络,再在基质与裂缝之间反复交换,最后才进入拟径向流。任何一个环节处理不对,拟合出来的参数都没有物理意义。
2.2 分支水平井压力响应的三个流态阶段:早期线性流、窜流过渡段、晚期拟径向流
分支水平井的压力动态可以粗分成三个阶段。早期,压力波在井筒周围裂缝系统内传播,主分支和各分支各自形成线性流,压力导数曲线上表现为与单位斜率平行的线段;产量沿分支分布不匀时,这一段的形态会比较杂乱。中期,基质开始向裂缝供气,压力导数出现典型的下凹槽——这个凹槽的深度和宽度直接对应储容比和窜流系数,是裂缝性气藏解释中信息量最大的一段。晚期,各分支的压力波相互叠加成一个整体,储层进入拟径向流,导数趋近于 0.5 水平线。
分支水平井与单分支水平井在流态上的本质差别在于:不同分支长度不一、角度各异,压力波到达边界的时间错开,导致过渡段被拉长、凹槽被抹平。如果分支间距小于压力波在某个时刻的探测半径,分支之间会发生干扰,表现为压力导数在过渡段出现多个小幅波动。这是解释中必须靠几何参数去约束的部分,绝不能只靠调储容比和窜流系数来硬凑。
2.3 建模的边界条件选择:均匀流量、无限导流与分支间干扰的处理
分支水平井试井模型的数学解,依赖于井筒边界条件的假设。常用的是两种:均匀流量假设认为沿井筒单位长度进气量相同,数学处理简单,对中低渗气藏、打开程度不高的井近似度好;无限导流假设认为井筒压力处处相等,适用于高渗气藏或井筒压降不可忽略的情况。对于裂缝性气藏,我一般先按均匀流量建模,因为这个假设对双重介质的适应性更好,而且它给出的早期线性流斜率与实测更接近。
多分支之间的空间干扰通过叠加原理处理:每一分支作为独立线源,对任一点产生的压力响应为该分支单位流量贡献之和。分支越密、间距越小,叠加效应越强,表现在压力导数上就是过渡段凹槽变浅、变宽。这里有一个容易忽略的点:分支角度影响的是各分支到观测点的等效距离,而不是简单地改变分支长度。代码里如果只用总长度代替分支几何,就会把分支井拟合成本质上的直井+裂缝,解释完全失真。
3. 从渗流方程到可运行的 Python 模型:拉普拉斯空间求解与 Stehfest 反演
3.1 无因次化:把气藏参数压缩成四个控制数
构建模型的第一步是无因次化。对裂缝性气藏,无因次压力定义为:
[ p_D = \frac{2\pi k_f h}{q \mu B} \Delta p ]
无因次时间定义为:
[ t_D = \frac{k_f t}{\phi_f c_{tf} \mu r_w^2} ]
加上前面提到的储容比 (\omega) 和窜流系数 (\lambda),整个裂缝性气藏分支水平井试井模型的控制参数就是这四组。无因次化的好处是让计算结果与具体井的产量、厚度、粘度无关,典型曲线可以直接用于拟合。对气井,压差要用拟压力 (\psi(p) = 2 \int_{p_0}^{p} \frac{p}{\mu z} dp) 代替,这部分我在第 4 章展开。
3.2 拉普拉斯空间的点源解与线源积分:Warren-Root 特征函数怎么进入解
双孔介质的核心技巧在拉普拉斯空间:基质与裂缝的窜流关系被吸收进一个特征函数 (f(s))。对于拟稳态窜流(Warren-Root 模型),特征函数为:
[ f(s) = \frac{\omega(1-\omega)s + \lambda}{(1-\omega)s + \lambda} ]
当无量纲时间极短((s \rightarrow \infty))时,(f(s) \rightarrow \omega),对应裂缝系统单独响应;当无量纲时间很长((s \rightarrow 0))时,(f(s) \rightarrow 1),整个系统趋于均质。这一性质在 Python 里实现非常直接:无限大地层三维点源在拉普拉斯空间的压力解为 ( \frac{1}{s} \frac{e^{-r\sqrt{s f(s)}}}{r} ),分支水平井就是沿着每个分支段对点源解做线积分,再按分支数做流量加权叠加。
以下代码实现了 f(s)、单分支线源积分和多分支叠加的拉普拉斯空间解:
import math import numpy as np from scipy.integrate import quad def f_dual(s, omega, lam): """ Warren-Root 拟稳态窜流特征函数。 s: 拉普拉斯变量 omega: 裂缝储容比,0~1 lam: 窜流系数(无因次) 极限检查:s->0 时 f->1(整体均质),s->inf 时 f->omega(裂缝独自响应) """ return (omega * (1.0 - omega) * s + lam) / ((1.0 - omega) * s + lam) def segment_pressure(s, x0, y0, x1, y1, x_obs, y_obs, omega, lam): """ 单个直线分支在拉普拉斯空间的均匀流量线源解。 分支从 (x0,y0) 到 (x1,y1),观测点 (x_obs,y_obs),坐标为无因次长度。 返回单位流量下的无因次压力 p_D(s)。 """ seg_len = math.hypot(x1 - x0, y1 - y0) if seg_len < 1e-12: return 0.0 Ld = seg_len / 2.0 alpha = math.sqrt(s * f_dual(s, omega, lam)) # 对分支参数 t 从 0 到 1 积分,dl = seg_len * dt def integrand(t): xp = x0 + t * (x1 - x0) yp = y0 + t * (y1 - y0) r = math.hypot(x_obs - xp, y_obs - yp) if r < 1e-6: r = 1e-6 return math.exp(-alpha * r) / r val, _ = quad(integrand, 0.0, 1.0, limit=200) # 线源归一化:总流量除以分支全长的均匀流量假设 return (val * seg_len) / (2.0 * Ld * s) def branch_well_pressure(s, branches, x_obs, y_obs, omega, lam): """ 多分支水平井叠加:每个分支按等流量分配,共 N 条分支。 branches: list of [(x0,y0,x1,y1), ...] 返回无因次压力 p_D(s) """ n_branch = len(branches) total = 0.0 for (x0, y0, x1, y1) in branches: total += segment_pressure(s, x0, y0, x1, y1, x_obs, y_obs, omega, lam) / n_branch return total这里的f_dual是整套模型的灵魂。omega 和 lam 不直接出现在流动方程的微分项里,而是通过这个特征函数改变点源解的衰减系数。分段积分时,每个分支在观测点处的距离不同,衰减快慢也不同,所以分支长、角度开、离观测点远的分支,贡献被指数项自然压低。quad积分的limit=200是因为分支长、alpha 大时积分函数振荡剧烈,默认上限可能不够;如果压力随时间出现锯齿,优先把这个值加大到 400。
要注意无因次坐标的缩放:这段代码里分支长度直接用无因次数,比如主分支从 -10 到 10 表示半长 10。实际换算时,无因次长度 (L_D = \frac{L}{r_w}\sqrt{\frac{k_z}{k_h}}),各向异性通过这个根号项进入。代码没有显式写各向异性,因为它在无因次坐标里已经被吸收掉了,不要重复乘。
3.3 Stehfest 数值反演与典型曲线绘制:把拉普拉斯空间的解拉回时间域
拉普拉斯空间的解无法直接用于拟合,必须做数值反演。Stehfest 方法是试井解释里最常用的反演算法,核心是用一组加权系数对拉普拉斯空间的函数值做加权求和。N 取 12 是工程经验的折中点:N 太小,曲线过于光滑,早期细节丢失;N 太大,舍入误差放大,压力出现震荡。
def stehfest_coeff(n): """计算 Stehfest 反演系数 V_i,N=12 为工程推荐值""" v = np.zeros(n + 1) m = n // 2 for i in range(1, n + 1): s = 0.0 k_lo = (i + 1) // 2 k_hi = min(i, m) for k in range(k_lo, k_hi + 1): a = k ** m * math.factorial(2 * k) b = math.factorial(m - k) * math.factorial(k - 1) c = math.factorial(i - k) * math.factorial(2 * k - i) s += a / (b * c) v[i] = (-1) ** (i + m) * s return v def stehfest_invert(F, t, n=12): """对拉普拉斯空间函数 F(s) 反演到实空间 t""" v = stehfest_coeff(n) ln2 = math.log(2.0) acc = 0.0 for i in range(1, n + 1): s = i * ln2 / t acc += v[i] * F(s) return ln2 / t * acc # 分支几何:主分支沿 x 轴 + 两条斜分支 branches = [ (-10.0, 0.0, 10.0, 0.0), # 主分支 (0.0, 0.0, 12.0, 9.0), # 右斜分支 (0.0, 0.0, 12.0, -9.0) # 左斜分支 ] omega = 0.1 lam = 1e-5 tD_arr = np.logspace(-4, 4, 120) def F(s): return branch_well_pressure(s, branches, 0.02, 0.0, omega, lam) pD_arr = np.array([stehfest_invert(F, t, 12) for t in tD_arr]) deriv = np.gradient(pD_arr, np.log(tD_arr))反演函数接受一个 Python 函数F作为参数,这个设计是为了拟合时能复用同一套反演流程,只替换参数。观测点取 (0.02, 0.0) 而不是 (0, 0),是为了避免积分时出现距离为 0 的奇异点,这个偏移量对应无因次井半径量级。np.gradient(pD_arr, np.log(tD_arr))得到的是对 ln(t_D) 的导数,也就是压力导数曲线,这一段反映了裂缝供气的快速压降和基质补给的滞后。
绘制典型曲线时,压力导数比压力曲线更值得关注。导数凹槽的最低点位置对应窜流发生的开始时刻,凹槽最低点的纵坐标和储容比成反比关系,横坐标则主要受窜流系数控制。如果画出来导数曲线在过渡段没有凹槽,优先检查是 lambda 设置过大还是分支数太多了。
4. 从典型曲线到储层参数反演:储容比与窜流系数是怎么被“读”出来的
4.1 储容比和窜流系数对曲线形态的控制作用
储容比 (\omega) 和窜流系数 (\lambda) 是裂缝性气藏解释的两个输出核心。储容比决定了凹槽的深度:(\omega) 越小,说明裂缝储容占比越低,过渡段压力下跌越明显,凹槽越深。窜流系数决定了凹槽出现的早晚:(\lambda) 越大,基质向裂缝补气越及时,凹槽越靠近早期并与早期线性流重叠,这时看起来反而像一个均质储层。两个参数一个管“深度”,一个管“横向位置”,形态上耦合并不高,但分支几何会同时影响这两个维度的判断。
典型曲线拟合的工程步骤是:先靠晚期拟径向流段确定裂缝总渗透率 (k_f),再用压力导数凹槽的深度和位置确定 (\omega) 和 (\lambda),最后用早期线性流段确定分支的有效长度。这个顺序不能乱。如果一开始就去拟合 (\omega) 和 (\lambda),晚期段的轻微偏移就会把它们带偏,给出虚假的双孔响应。
4.2 用 scipy 做自动拟合:目标函数怎么写才不容易翻车
自动拟合的本质是最小化模型曲线与实测曲线的对数域距离。直接对压力拟合会遇到早期数据权重过大的问题;更稳的做法是对压力导数拟合,或者在目标函数里同时给压力和导数各分配 0.5 的权重。以下代码演示了如何拟合 (\omega) 和 (\lambda):
from scipy.optimize import minimize # 模拟一组“实测”压力导数,作为拟合目标 # 真实模型参数: omega=0.08, lam=2e-5 omega_true = 0.08 lam_true = 2e-5 def F_true(s): return branch_well_pressure(s, branches, 0.02, 0.0, omega_true, lam_true) pD_true = np.array([stehfest_invert(F_true, t, 12) for t in tD_arr]) deriv_true = np.gradient(pD_true, np.log(tD_arr)) # 拟合目标:导数曲线的对数误差 def objective(params): w, l = params if not (0.01 <= w <= 0.5) or not (1e-7 <= l <= 1e-2): return 1e6 def F_model(s): return branch_well_pressure(s, branches, 0.02, 0.0, w, l) pD_model = np.array([stehfest_invert(F_model, t, 12) for t in tD_arr]) deriv_model = np.gradient(pD_model, np.log(tD_arr)) return np.mean((np.log10(deriv_model + 1e-8) - np.log10(deriv_true + 1e-8)) ** 2) res = minimize(objective, x0=[0.15, 1e-4], method='Nelder-Mead', options={'xatol': 1e-6, 'fatol': 1e-8, 'maxiter': 200}) print("拟合结果 omega=%.4f, lambda=%.2e" % (res.x[0], res.x[1]))这里把目标函数定义为对数域的导数误差,是为了让小数值的 (\lambda) 和大数值的 (\omega) 在同一尺度下比较。+1e-8是防止导数过零时取对数报错。初值选择上,(\omega) 取 0.15、(\lambda) 取 1e-4 是常见默认值,如果拟合结果落到边界上,说明分支几何与曲线形态不匹配,而不是参数问题。
让我跑一下思路验证:真实参数是 (0.08, 2e-5),初值 (0.15, 1e-4),Nelder-Mead 在这个光滑的二维目标函数上一般几秒内收敛。注意这段代码里我没有把分支几何放进拟合变量,实际项目中分支长度和角度如果有测井或地质约束,应该固定;如果完全未知,建议先用手工拟合确定几何,再放开自动拟合。
拟合完成后要检查残差曲线:如果导数拟合得很好但压力拟合差,大概率是井储效应或表皮没处理;如果压力拟合好但导数在凹槽处偏差大,说明 (\omega) 方向错了,而不是 (\lambda)。
4.3 气藏与油藏试井公式的差异:拟压力和多产量效应
前面所有代码用的都是油藏形式的无因次压力,换到气藏时必须在输入数据阶段把实测压力转换成拟压力。拟压力 (\psi(p) = 2\int_{p_0}^{p} \frac{p}{\mu(p) z(p)} dp) 把气体物性随压力的变化吸收掉,之后模型公式与油藏完全一致。这个变换必须在进入模型前完成,而不是在模型里乘一个系数。常见的错误是直接用压力平方差 (p^2) 代替拟压力,这只在低压范围近似成立,对高压气藏会造成系统性偏差。
气井还有一点特殊:高速流动造成的非达西效应会以附加表皮的形式进入模型。解释时如果发现表皮系数随产量明显增大,就要怀疑非达西项。处理方式是做多产量测试,回归出速率相关表皮。裂缝性气藏里,非达西效应主要在裂缝内发生,影响的是早期段,所以不要因为在早期拟合不好就盲目加大分支数量。
5. 裂缝性气藏分支水平井试井解释的5个避坑点:从数据采集到模型选型
5.1 早期井储效应淹没裂缝响应,把双孔曲线误判成均质
现象:压力双对数曲线最早期斜率为 1,导数也贴着单位斜率上升,看不到任何分支线性流的影子;强行拟合时发现无论怎么调 (\omega) 都拟合不出凹槽。
原因:井筒储集效应在关井初期占主导,井筒里气体的压缩和液面升降掩盖了地层响应。裂缝性气藏往往伴随较大井筒体积,这个问题更严重。如果井储系数未知,直接拿原始数据拟合,模型会把井储吸收掉,导致解释结果完全偏离。
解决:先做井储诊断。在双对数图上把压力导数早期的单位斜率段单独拟合出井储系数 (C),然后对压力数据做去井储处理,或者使用带井储边界条件的完整模型让 (C) 作为独立拟合参数。我习惯是把 (C) 和表皮 (S) 放进拟合变量,而不是手动扣掉。分支几何在早期数据里反映很弱,拟合时先固定分支参数,让 (C) 去吸收早期误差。
5.2 分支角度和间距设置不当,模型解不唯一
现象:两组完全不同的分支几何和参数组合,能拟合出几乎重合的压降曲线;解释结果里渗透率和窜流系数相差悬殊,但拟合误差一样小。
原因:分支水平井对角度和间距的敏感性远低于对总长度的敏感性。远离主分支的细小分支对压力响应的贡献被指数衰减抹平,本质上等同于一个较短的主分支。这是裂缝性气藏分支水平井试井模型的固有不适定问题。
解决:不要指望单靠试井曲线确定分支几何。必须用现场的分段流量测试、随钻测量数据或成像测井来固定分支数、长度和夹角。如果这些数据都没有,就做“简约化处理”:把多个细分支合并成一个等效分支,用一个等效半长去表征,再在解释报告里明确说明这个近似。强行拟合分支角度属于自欺欺人。
5.3 窜流系数预估偏大,双孔凹槽消失被解释成均质储层
现象:压力导数曲线在过渡段没有明显的下凹,整条曲线和单孔介质模型一致;用单孔模型拟合,得到的渗透率介于裂缝渗透率和基质渗透率之间,看似合理但无法预测生产。
原因:(\lambda) 偏大时,基质向裂缝的补给非常及时,窜流发生在早期线性流还没有完全展开的时候,凹槽被早期段“吃掉”了。这在裂缝发育好、基质渗透率不算太低的储层里很容易出现。注意气藏与油藏不同,气体粘度随压力变化会让过渡段拉长,如果气藏压力高且 (\lambda) 较大,凹槽可能非常浅。
解决:把拟合窗口延长到中晚期,重点检查拟径向流段的形状。如果晚期导数高于 0.5 但压力呈持续上升趋势,考虑地层边界;如果晚期导数正常而只是早期没有凹槽,不要急着排除双孔模型,先让 (\lambda) 从 1e-7 左右的小值扫一遍,看拟合残差是否有明显改善。扫参比人工判断可靠得多。
5.4 高速非达西流动让表皮系数虚高,早期数据报废
现象:实测压力导数早期段明显高于模型曲线,拟合出来的表皮系数达到几十甚至上百,与压裂改造预期不符;产量越高,虚高越明显。
原因:气井在裂缝内流速过高,达西定律失效,附加压降与流量平方成正比。裂缝性气藏中裂缝导流能力有限,这个问题在高产井里很常见。把非达西效应硬塞进表皮系数,解释结果无法用于预测不同产量下的压力动态。
解决:至少做两个产量下的测试,分别拟合得到表皮 (S_1) 和 (S_2),用 (S = S_0 + D q) 回归出速率相关系数 (D)。在模型里把非达西项显式加入边界条件,而不是并入表皮。这个步骤看起来费事,但分支水平井本就参数多,不拆开非达西项,其他参数的置信区间全是假的。
5.5 Stehfest 反演的 N 值选择:早期数据震荡是数值问题还是地质响应
现象:反演后的压力曲线在早期出现锯齿状波动,特别是 (t_D < 0.01) 时;改变 Stehfest 的 (N) 值,曲线形态变化明显。
原因:Stehfest 系数随 (N) 增大而急剧增大,舍入误差随之放大。(N=12) 是多数文献推荐值,但它对函数光滑性有要求;当分支数多、积分分段多时,拉普拉斯空间函数本身带数值噪声,反演会把噪声放大。这不是油藏响应,是纯数值问题。
解决:先用 (N=12) 计算,如果早期震荡,把分支积分精度提高(增大 quad 的limit和epsabs),再尝试 (N=10)。对比两组结果的差异:只有早期段震荡且中晚期重合,可以确认是数值问题;如果整个曲线形态都变,说明是模型本身的问题。另外,观测点离分支太近(如小于 1e-6)也会导致积分奇异,检查一下无因次井径设置。
6. 进阶:用叠加原理验证多分支模型的物理一致性,一个值得养成的习惯
在把模型交给拟合之前,先做一次叠加原理验证,这个习惯帮我抓出过不止一次代码错误。原理很简单:多分支井的总压降,应当等于每个分支单独生产时在观测点产生的压降之和。具体做法是设计两组计算——第一组直接用我们的branch_well_pressure计算三分支完整模型的压力;第二组分别计算三个分支单独存在时的压力,再按 1/3 权重叠加。如果两者在全部时间范围内重合且误差小于千分之一,说明分支叠加实现正确。
这个验证极度敏感,能同时暴露三类问题:一是分支积分方向写反导致符号错误;二是流量权重分配不对;三是观测点与分支距离太近时奇异点污染积分结果。代码层面对比很简单,不需要引入额外库——在branch_well_pressure函数里临时传入只含一个分支的列表,然后把三分支单独结果的平均值与三分支整体模型的输出做对数域误差比对。误差来源纯粹是数值积分的精度差异,工程上小于 1e-3 就像。
跑通叠加验证只是第一步,我更习惯的做法是顺手做一个“守恒检查”:把观测点放在主分支井筒表面与放在远离所有分支 100 倍半径处,分别计算压力响应。前者应反映井筒压力,后者应趋向于点源解。如果远端响应与点源解偏差超过百分之几,说明分支长度和观测点坐标的无因次换算有错。
最后一个收尾习惯:所有拟合结果必须输出到一张包含压力残差和导数残差的对照图上,并按 (\omega)、(\lambda)、(k_f h)、(S) 的顺序记录置信区间。这样报告里的曲线才经得起审阅。我在实际项目里见过太多只给“拟合好”的漂亮曲线、却拿不出参数不确定性范围的解释报告,那对气藏开发决策几乎没参考价值。希望这篇拆解能帮你把裂缝性气藏分支水平井试井模型的代码跑通,也把解释中的那些暗坑提前亮出来——希望帮到你。
本文还有配套的精品资源,点击获取