简介:面向石油工程与试井分析专业人员的裂缝性气藏分支水平井试井模型完整复现资料。资源以1个docx文档呈现,约60KB,聚焦模型构建、Laplace变换求解、Stehfest数值反演及典型曲线绘制。内容从控制方程与无量纲定义出发,逐步给出Python完整可运行代码,涵盖裂缝储容比、窜流系数等参数敏感性分析,以及九阶段流动形态的自动识别与对比。读者可据此掌握鱼骨型分支水平井的渗流阶段划分逻辑,理解分支长度与间距对第二拟径向流段的影响规律。文档还结合工程场景提出分支参数优化、试井测试设计与数据解释建议,便于研究人员直接复现论文结果或用于现场方案论证。压缩包仅有docx说明文档,体积精简、模型主线清晰,适合已有一定渗流力学基础、需要代码级参考的开发与研究人员。已有73人学习/浏览此资源。
1. 裂缝性气藏分支水平井试井模型:为什么压力曲线和裂缝地质图总对不上
地质建模图上裂缝密度明明很大,分支水平井也穿过了好几条主缝,可压力恢复的双对数曲线看上去就是一个光滑的均质径向流形态,双孔介质的那个“凹槽”死活不出现。裂缝性气藏分支水平井试井模型就是用来解释这种矛盾的——把储层抽象成基质和裂缝两套介质,把分支井筒抽象成若干条线源,在拉普拉斯空间求解析解再用 Stehfest 反演回实空间,最终从实测压力曲线上反推出裂缝储存比、窜流系数和有效渗透率。适合正在做试井解释、气藏动态描述或地质模型校验的工程师。下面从模型假设、求解代码到典型曲线识别和现场踩坑,一条线讲完。
2. 模型怎么建:Warren-Root 双孔介质与分支井筒的三个关键假设
2.1 双孔介质假设:基质储气、裂缝导流的 ω 和 λ
裂缝性气藏最常见的简化是 Warren-Root 双孔介质模型。基质岩块孔隙度低但占据绝大部分储量,裂缝系统孔隙体积小但渗透率高出基质几个量级,压力波到达井筒时先沿裂缝快速传播,随后基质才开始向裂缝缓慢补给。这个“先快后慢”的补给过程,在压力导数曲线上表现为一个明显的下凹段。
两个无因次参数决定这个凹槽的形态。储存比 ω 定义为裂缝系统弹性储容占总弹性能的比例:
$$\omega = \frac{(\phi C_t)_f}{(\phi C_t)_f + (\phi C_t)_m}$$
ω 越小,说明裂缝储容占比越低,基质那部分弹性能需要更长时间释放,凹槽就越深。窜流系数 λ 表征基质向裂缝补给的强度,λ 越大补给越快,凹槽出现得越早。很多解释报告里只写“裂缝性气藏”而不标 ω 和 λ 的取值,实际上这两个数才是后续储量计算和产能预测的地质根。
双孔介质的无因次扩散方程在拉普拉斯变换后,会得到一个特征函数 f(s),所有后续求解都建立在它之上。气体不同于液体,严格说应该用拟压力 m(p) 代替压力 p,但方程结构和解的形式不变,拿到代码里只需要把输入压力换成拟压力即可。
2.2 分支水平井线源叠加:从分支坐标到无因次化
单支水平井的试井模型要沿井筒积分点源解,分支水平井则更复杂——每条分支都有自己的起点、长度和方位角,整口井等价于一个由多条线源组成的叠加系统。工程上多数人先把每条分支离散成若干段,把每段当成一个点源,再把所有点源的压力响应叠加起来,这就是线源叠加法。
实际操作时我会先建立全局坐标系,以主井筒的井底测试点为原点,X 轴取主支方向。每条分支用 [x0, y0, x1, y1] 四个数描述,其中 x1、y1 由分支长度 L 和方位角 θ 通过旋转矩阵算出。如果分支在钻探中有实测轨迹数据(LWD 或定向井轨迹),直接用轨迹坐标更准,不用硬套理想几何。所有这些坐标都要除以一个参考长度 L_ref 做无因次化,L_ref 我一般取主支长度或储层厚度,关键是时间无因次化必须用同一个 L_ref。
这里有一个新手最容易翻车的地方:点源解在井轴处是发散的,线源积分时如果观测点正好落在源点上,压力会趋于无穷。常见的处理办法是设一个无因次截止半径 rwD,积分时强制 r 不小于 rwD,相当于把线源当成半径有限的井筒。这个数值技巧不是可有可无,不设这个下限,Stehfest 反演出来的早期压力会完全失真。
2.3 内边界与外边界:井储、表皮和无限大地层的拉普拉斯写法
裂缝性气藏分支水平井试井模型的完整数学描述还包括两类边界条件。内边界是井筒处的井储效应和表皮效应。井储效应来自井筒和地面管线里的流体弹性能,压力计记录到的最早一段数据几乎全被井储主导,数学上用无因次井储系数 C_D 描述。表皮效应来自钻井液侵入和完井污染,用表皮系数 S 描述,S 为正表示污染、为负表示酸化压裂增产。这两个效应在拉普拉斯空间可以用 Van Everdingen-Hurst 公式一行代码叠加上去,不用回到实空间做卷积,这是解析试井模型最方便的地方。
外边界我默认取无限大地层,边界效应(封闭边界、定压边界)需要在点源解上叠加镜像源,代码量不大但会把文章带偏,现场试井早期的目标基本不涉及边界识别,先把无限大解跑通、把 ω 和 λ 拟合出来更重要。整个模型求解路径是:双重介质扩散方程做拉普拉斯变换 → 得到带 f(s) 的点源解 → 沿分支做线源积分 → 叠加井储表皮 → Stehfest 数值反演回实空间。下一步就直接进入代码。
3. 拉普拉斯空间求解与 Stehfest 反演:Python 代码从零跑通
3.1 双重介质点源解与多分支线源积分的代码实现
以下代码是整条求解链的核心,只依赖 numpy 和 scipy。先写双重介质的特征函数和多分支线源积分:
# 裂缝性气藏分支水平井试井模型 —— 拉普拉斯空间核心解 import numpy as np from math import factorial from scipy.special import kv def warren_root_f(s, omega, lamda): """双重介质Warren-Root特征函数 f(s),s 为拉普拉斯变量""" return s * (omega * (1 - omega) * s + lamda) / ((1 - omega) * s + lamda) def multibranch_pressure(s, omega, lamda, branch_geo, rwD=1e-4, n_seg=40): """ 多分支水平井的拉普拉斯空间无因次压力(均匀流量线源叠加) branch_geo: 列表,每个元素为 [x0, y0, x1, y1],坐标已用参考长度无因次化 rwD: 无因次截止半径,用于避开线源积分奇异性 """ f = warren_root_f(s, omega, lamda) L_total = 0.0 for x0, y0, x1, y1 in branch_geo: L_total += np.hypot(x1 - x0, y1 - y0) integral = 0.0 for x0, y0, x1, y1 in branch_geo: L_branch = np.hypot(x1 - x0, y1 - y0) xs = np.linspace(x0, x1, n_seg + 1) ys = np.linspace(y0, y1, n_seg + 1) for j in range(n_seg): xm = 0.5 * (xs[j] + xs[j + 1]) ym = 0.5 * (ys[j] + ys[j + 1]) r = max(np.hypot(xm, ym), rwD) integral += kv(0, r * np.sqrt(f)) * (L_branch / n_seg) return integral / (s * L_total)warren_root_f 返回的就是 2.1 节里那个特征函数,它在拉普拉斯空间替代了均质介质中的 s。multibranch_pressure 对每条分支均匀离散成 n_seg 段,取每段中点作为点源位置,用第二类零阶修正贝塞尔函数 kv(0, x) 计算点源响应,再按每段长度占总长度的比例加权。整个积分在 s 意义上除以 s,对应拉普拉斯变换中阶跃产量的积分性质。
参数说明:omega 和 lamda 的取值范围对裂缝性气藏来说,omega 通常在 0.01~0.3,lamda 在 1e-6~1e-2;rwD 取 1e-4 左右足够避开奇异性又不会引入明显误差;n_seg 取 40 时精度已经稳定,再多只是增加计算量。分支坐标必须预先除以 L_ref,如果直接传米制坐标,r 会变成几百,kv 函数的自变量过大直接溢出,这一步几乎人人踩过。
3.2 井储表皮叠加与 Stehfest 反演:参数选 8 还是 10
拉普拉斯空间算出的只是无井储、无表皮的裸地层响应,要叠加井储和表皮,用 Van Everdingen-Hurst 公式。同时把 Stehfest 反演的系数算出来:
def wellbore_storage_skin(p_d, s, S, C_D): """含井储C_D和表皮S的拉普拉斯空间井底压力公式""" A = s * p_d + S return A / (s + C_D * s * s * A) def stehfest_coeff(N): """Stehfest数值反演系数,N必须为偶数,常用8/10/12""" M = N // 2 V = np.zeros(N) for i in range(1, N + 1): acc = 0.0 for k in range((i + 1) // 2, min(i, M) + 1): num = k**M * factorial(2 * k) den = factorial(M - k) * factorial(k) * factorial(k - 1) den *= factorial(i - k) * factorial(2 * k - i) acc += num / den V[i - 1] = (-1)**(i + M) * acc return Vwellbore_storage_skin 里的逻辑可以这样理解:分子 sp_d+S 是地层压力加表皮压降,分母里的 C_DssA 则是井储效应在拉普拉斯空间的表达。C_D 越大,早期井储段持续越长。这个公式直接对拉普拉斯压力解作用,所以不需要在实空间做卷积,这也是整套方法计算上最高效的地方。
Stehfest 系数公式看着吓人,代码里就是两层循环加阶乘。N 的选择直接影响精度:N=8 对光滑曲线够用但偏粗,N=10 是多数试井程序的标准配置,N=12 精度更高但浮点舍入误差开始显现,N=14 以上在极端情况下会出现振荡。我的习惯是固定用 N=10,只有在压力导数早期有明显振荡时才降到 8。网上能搜到很多 Stehfest 的示例代码,但多数没有说明 N 的选取与数值稳定性关系,这里补一句:把 N 当成黑匣子直接拉到 16,大概率得到一条锯齿状曲线。
3.3 一整套可复现的算例:三支鱼骨井的压力响应曲线
最后是反演和主流程,用一个三支鱼骨井算例跑通全流程:
def stehfest_invert(tD, F_s, V): """Stehfest数值反演:由拉普拉斯空间函数值还原实空间解""" alpha = np.log(2.0) / tD total = 0.0 for i in range(1, len(V) + 1): total += V[i - 1] * F_s(i * alpha) return alpha * total def build_response(t_hours, q, k, h, phi, mu, Ct, rw, omega, lamda, S, C_D, branch_geo, L_ref): """计算压降阶段的Δp和对数压力导数,返回时间、压差、导数""" V = stehfest_coeff(10) rwD = rw / L_ref dp = np.zeros(len(t_hours)) for i, t in enumerate(t_hours): tD = 3.6 * k * t / (phi * mu * Ct * L_ref * L_ref) def F_s(s): pD = multibranch_pressure(s, omega, lamda, branch_geo, rwD) return wellbore_storage_skin(pD, s, S, C_D) dp[i] = stehfest_invert(tD, F_s, V) * (1.842e-3 * q * mu / (k * h)) deriv = t_hours * np.gradient(dp, t_hours) return t_hours, dp, deriv if __name__ == "__main__": L_ref = 500.0 branches = [ [0, 0, 400, 0], [150, 0, 150 + 200 * np.cos(np.pi / 4), 200 * np.sin(np.pi / 4)], [150, 0, 150 + 200 * np.cos(np.pi / 4), -200 * np.sin(np.pi / 4)], ] branches = [[v / L_ref for v in b] for b in branches] t = np.logspace(-3, 3, 60) t_h, dp, deriv = build_response( t_hours=t, q=1e5, k=2.5, h=30, phi=0.12, mu=0.025, Ct=0.005, rw=0.1, omega=0.08, lamda=1e-4, S=1.0, C_D=1e3, branch_geo=branches, L_ref=L_ref )build_response 里无因次时间用的是 3.6kt/(phimuCt*L_ref²),这是油藏工程常用单位下的换算,k 为 mD、t 为小时、Ct 为 1/MPa,和后面 1.842e-3 的压力换算系数配套。如果单位制不统一,结果会差几个量级,这是把油藏公式搬进 Python 最常见的血泪错误。分支几何里三条分支分别是 400m 主支和两条约 200m 的对称分支,在 x=150m 处分叉,这个几何形态模拟的是典型的鱼骨井,跑出来的曲线可以直接和实测压力对照。
跑完这段代码,会得到一组时间上对数均匀分布的压差和对数导数,这就是后续做典型曲线识别的原材料。代码里所有函数都是独立的,warren_root_f 和 stehfest_coeff 可以单独拿出来给别的模型复用,multibranch_pressure 换成单支时只需把 branch_geo 传一个线段即可。建议把 branch_geo 的构造单独提成一个函数,以后换井型时不用改主流程。
4. 典型曲线识别与参数敏感性:压力导数图上每一段在说什么
4.1 ω 和 λ 的敏感性:双孔介质“凹槽”有多深、多早出现
拿到压力和导数曲线后,第一件事不是拟合,而是分段认清每一段反映的是什么。典型曲线按时间顺序可以分四段:早期井储段、裂缝径向流段、双孔过渡段、系统径向流段。四段的识别特征和解释价值如下表:
| 阶段 | 压力导数双对数特征 | 可获取参数 |
|---|---|---|
| 早期井储段 | 斜率约等于 1 的直线 | 井储系数 C_D |
| 裂缝径向流段 | 第一条水平直线 | 裂缝系统 kh |
| 双孔过渡段 | 导数下凹再回升 | 储存比 ω、窜流系数 λ |
| 系统径向流段 | 第二条水平直线 | 裂缝+基质的系统 kh |
ω 的敏感性最直观地体现在凹槽深度上。ω 越小,凹槽越深,因为基质提供的弹性储容占比更大,压力波需要等待更长时间的补给才能进入系统径向流。λ 则控制凹槽出现的早晚和宽度,λ 越大凹槽越早出现、恢复也越快。现场解释中如果凹槽几乎没有,往往有两个可能:一是储层其实接近均质,裂缝不发育;二是井储段太长把过渡段完全掩盖了。这两种情况的处理方式完全相反,前者要换模型,后者要优化测试工艺,在第 5 章详说。
气体和液体在双孔介质过渡段上有个容易被忽略的差异:气体压缩性强,拟压力变换后在过渡段的表现比液体更缓和,同一条实测曲线用压力平方和用拟压力解释出来的 ω 可以差一倍。所以解释气井时不要拿液体的工作流直接套,至少要确认 PVT 计算里拟压力的积分处理是完整路径积分,而不是 Z 因子近似。
4.2 分支数、分支长度与方位角:哪些几何参数真正影响曲线
分支水平井比单支水平井多了一层不确定性:几何参数本身也是待求的。但试井曲线对分支几何的敏感度并不是均匀的。分支总长度(所有分支长度之和)对早期径向流出现的时间影响最大,总长度越大,等效井径越大,裂缝径向流段的水平直线出现得越晚、压降越平缓。分支数的影响其次,从单支变双支,曲线在前期会明显下移;从三支变四支,差异已经很小了,这就是边际效应递减。
方位角对压力响应的影响要看分支间的相互距离。两条分支夹角 90 度以上时,彼此干扰很小,线源叠加的结果接近独立分别计算后求和;夹角小于 30 度时,两条分支几乎平行,中间区域的地层压力波会重叠,表现在压力导数早期会有一个异常的“鼓包”,不是储层特征,而是几何引起的。现场往往拿不到精确的分支轨迹,只能靠随钻测量估算方位角,这个不确定性要直接放到拟合参数里去,而不是先固定死再调 ω 和 λ。
分支之间产量是否均匀也是一个隐性的强假设。代码里默认均匀流量假设,真实分支井各段产量受跟端和趾端压降影响,分布并不均匀。如果某条分支在实测产气剖面里明显高产出,而模型里还按平均分配处理,拟合出来的渗透率会偏高。遇到这种情况,我的做法是把分支评分段,按产气剖面给每段加权系数,而不是继续用均匀假设硬凑。
4.3 井储和表皮对早期段的干扰:如何判断数据是否可用
井储和表皮对典型曲线的干扰集中在早期,但会连锁影响整个拟合。C_D 越大,早期井储段越长,它会像一个窗帘一样逐渐拉上,把后面的裂缝径向流和双孔过渡段慢慢遮住。如果曲线在 24 小时之内还看不到第二条水平直线,基本可以放弃从这段数据里提取 ω 和 λ 的想法。表皮系数 S 不改变导数曲线形态,只使压差整体抬高(S>0)或降低(S<0),所以在拟合表皮时要用绝对压力段而不是导数段。
这里有一个判断数据可用性的经验法则:先看双对数导数曲线里井储段结束的时间点,再看双孔过渡段的起始时间,如果前者晚于后者,说明解释结果不可靠,不应该强行输出 ω 和 λ。很多试井软件会自动拟合,但不会警告你早期段被井储掩盖,输出一个置信区间很窄却不物理的结果。判断曲线是否“够用”,比拟合本身更重要。
5. 试井解释踩坑记录:四个最常见的翻车现场
5.1 压力导数出现假驼峰,别急着往里加封闭边界
现象:导数曲线中段出现一个明显的上升-下降驼峰,新手第一反应是“多分支干扰”或“封闭边界”,然后往模型里加各种参数去追。原因:绝大多数假驼峰来自井筒的相重新分布或液面降落,气体井尤其常见,关井后井筒内气液两相重新分层,井底压力计记录到的是井筒过程而不是地层过程。解决:先对比同一口井较早一次的测试记录,如果驼峰位置随测试时间漂移,基本可以断定是井筒效应。处理方法是把驼峰段标记为无效段,在自动拟合时排除,而不是加到模型里。
5.2 拟合太“完美”但地质不认,问题出在参数可辨识性
现象:自动拟合收敛很漂亮,压力导数和模型曲线几乎完全重合,但解释出的 ω=0.05,而地质建模给出的裂缝储容占比是 0.2,差了 4 倍。原因:分支井模型参数多(ω、λ、S、C_D、分支长度、分支数),彼此之间存在强相关,大量等效参数组合能产生几乎相同的压力曲线,这叫可辨识性不足,不是拟合质量差。解决:把分支几何参数固定下来,只用随钻轨迹和 LWD 解释结果作为输入;反演时固定 S 和 C_D 中至少一个,只让 ω、λ、k 自由;最终结果必须同时列地质约束下的参数置信区间,而不是给一个孤立的最优值。分支长度、井储系数这些本来就有独立测量手段的参数,交给试井去反演是最典型的资源浪费。
5.3 早期数据总是压不齐,先检查初值和井储段
现象:模型算出的早期段与实测点总是对不上,调整 C_D 和 S 也没有明显改善,自动拟合直接发散。原因:早期段主导因素是井储,如果初始 C_D 差了一个数量级,拟合算法的雅可比矩阵在早期段几乎是病态的;另外压力计分辨率和采样频率不足也会让早期点噪声过大。解决:先人工在双对数图上量斜率=1 直线段的长度,估算 C_D 初值再交给拟合;对早期点做降权处理,不让它们主导目标函数;如果数据噪声实在太大,直接删掉开井后前 0.5 小时的点再做拟合。先定 C_D、S,再定 k,最后定 ω 和 λ,不要一次全放开。
5.4 双孔过渡段被井储吃掉,换压力恢复数据或改测试工艺
现象:整条导数曲线光滑得像均质储层,没有任何凹槽,但岩心和测井都显示裂缝发育。原因:井储系数过大,关井后井筒持续向地层倒灌,双孔过渡段发生时间完全落在井储段内部,被掩盖了。解决:最有效的办法是改用压力恢复测试而不是压降测试,恢复测试用叠加原理把有效测试时间拉长;有条件的企业直接下井下关井工具,地面关井的井储效应通常比井下关井大一个量级以上。如果数据已经这样了,就不要硬套 Warren-Root 模型,试算一下双渗透模型或者均质等效模型,哪个能解释数据就先用哪个,模型选择服务于工程决策,不服务于“证明它是裂缝性气藏”。
6. 模型落地的三个验证习惯:从拟合曲线到可信报告
做完拟合只是第一步,真正让这套模型可信的是验证流程。第一个习惯是分段拟合而不是全曲线一键拟合。我一般先用斜率=1 的直线段定 C_D,用终期径向流段定系统渗透率,再用过渡段凹槽定 ω 和 λ,每段只调一两个参数,拟合结果回到整条曲线做一致性检查。这不是为了省事,而是因为全参数同时拟合在数学上能收敛,在地质上往往给出一组互相矛盾的值。
第二个习惯是至少做一次独立验证。解析模型里分支几何假设是理想化的,真实井筒有轨迹偏离、有跟端趾端压差,把这些因素一并交给解析解本身不现实。常见做法是建一个简单的数值试井模型(比如用有限元或有限差分网格),把解析模型解释出的 ω、λ、k 放进去,对比同样测试时间下的压力恢复曲线。两者偏差小于 5% 才认为解释可信,这个验证步骤在审查报告时几乎必看。
第三个习惯是报告里永远写参数的不确定性范围。试井曲线对参数组合的响应不是一一映射,同样的拟合质量可能对应 ω 在 0.03~0.12 之间变化。我会把等效拟合组合列出两到三组,标注各自对应的地质解释,再附上推荐值。盲目追求单一最优解,后续储量计算和产能预测就会把参数误差放大成开发决策风险。
我早年在一个裂缝性气藏项目里,为把曲线贴得完美,拼命加分支方位角参数,最后拟合出一个分支走向与随钻测井完全矛盾的结果,被地质师当场问住。自那以后我把“地质约束优先于数学拟合”写进自己的流程里,固定能独立的参数,只反演不可测的储层参数,再用数值模型做交叉验证。这套流程救了很多次场,希望帮到你。
本文还有配套的精品资源,点击获取