CN113191102A - 一种钻井液瞬态波动压力确定方法及系统 - Google Patents

一种钻井液瞬态波动压力确定方法及系统 Download PDF

Info

Publication number
CN113191102A
CN113191102A CN202110535627.7A CN202110535627A CN113191102A CN 113191102 A CN113191102 A CN 113191102A CN 202110535627 A CN202110535627 A CN 202110535627A CN 113191102 A CN113191102 A CN 113191102A
Authority
CN
China
Prior art keywords
model
pressure
node
characteristic
grid
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110535627.7A
Other languages
English (en)
Other versions
CN113191102B (zh
Inventor
刘书杰
黄熠
刘和兴
王成龙
马传华
郑金龙
田野
柳亚亚
吴艳辉
夏茜
蔡饶
陈鹏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
CNOOC China Ltd Zhanjiang Branch
CNOOC China Ltd Hainan Branch
Original Assignee
CNOOC China Ltd Zhanjiang Branch
CNOOC China Ltd Hainan Branch
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by CNOOC China Ltd Zhanjiang Branch, CNOOC China Ltd Hainan Branch filed Critical CNOOC China Ltd Zhanjiang Branch
Priority to CN202110535627.7A priority Critical patent/CN113191102B/zh
Publication of CN113191102A publication Critical patent/CN113191102A/zh
Application granted granted Critical
Publication of CN113191102B publication Critical patent/CN113191102B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/28Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/11Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13Differential equations
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00Details relating to the application field
    • G06F2113/08Fluids
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces

Abstract

本发明涉及一种钻井液瞬态波动压力确定方法及系统,包括以下步骤:对波动压力物理模型的流动区域进行按照设定的时间步长和距离步长进行网格划分,得到多个网格节点;获取网格节点原始控制模型的特征线模型,根据得到的特征线模型获得原始控制模型的常微分计算模型;获取特征线模型的在待求解时刻的前一时刻对应的特征节点的位置信息,根据与特征节点相邻的网格节点的速度和压力信息获取特征节点速度和压力,将获取的特征节点速度和压力带入常微分计算模型中,得到待计算时刻的网格节点的波动压力计算模型;利用波动压力计算模型得到不同网格节点在不同时刻的瞬态波动压力,采用本发明的方法降低了对网格划分的要求,精度高。

Description

一种钻井液瞬态波动压力确定方法及系统
技术领域
本发明涉及钻井技术领域,具体涉及一种钻井液瞬态波动压力确定方法及系统。
背景技术
这里的陈述仅提供与本发明相关的背景技术,而不必然地构成现有技术。
起下钻或下套管作业引起的波动压力是钻进安全钻井液密度窗口窄地层时井下压力控制的重要敏感因素,引起波动压力的主要因素为:钻井液胶凝强度、钻井液粘性和钻井液惯性。波动压力计算方法有稳态法和瞬态法。发明人发现,稳态法将钻井液看做不可压缩流体,不考虑流道弹性和惯性力,依据管柱运动引起的钻井液流量计算得到的摩阻压降作为波动压力,计算方法相对简单稳定,但计算结果比实测数据偏高30%-50%。瞬态法计算结果与实测数据吻合程度较高,但常规瞬态波动压力计算方法对时间步长和距离步长划分要求高,现场工况一般难以满足时空比要求,计算经常溢出、发散。
发明内容
本发明的目的是为克服现有技术的不足,提供一种钻井液瞬态波动压力确定方法,保证了计算过程的稳定收敛。
为实现上述目的,本发明采用如下技术方案:
第一方面,本发明的实施例提供了一种钻井液瞬态波动压力确定方法,包括以下步骤:
对波动压力物理模型的流动区域进行按照设定的时间步长和距离步长进行网格划分,得到多个网格节点;
获取网格节点原始控制模型的特征线模型,根据得到的特征线模型获得原始控制模型的常微分计算模型;
获取特征线模型的在待求解时刻的前一时刻对应的特征节点的位置信息,根据与特征节点相邻的网格节点的速度和压力获取特征节点速度和压力,将获取的特征节点速度和压力带入常微分计算模型中,得到待计算时刻的网格节点的波动压力计算模型;
利用波动压力计算模型得到不同网格节点在不同时刻的瞬态波动压力。
可选的,根据与特征节点相邻的网格节点的压力和速度信息采用插值算法得到特征节点的压力和速度信息。
可选的,计算边界条件,对流动区域施加边界条件,设定初始条件,对网格节点进行求解变量初始化并按照时间步长,利用波动压力计算模型推演得到不同时刻、不同网格节点处的瞬态波动压力。
可选的,所述流动区域划分为管柱内流动区域、环空流动区域及管柱下方流动区域,运动管柱流动区域为管柱模型内的区域,环空流动区域为管柱模型与井壁模型之间的区域,管柱下方流动区域为管柱模型下方的流动区域。
可选的,计算管柱内流动区域井口节点的压力、速度信息、计算环空流动区域井口节点的压力速度信息、计算管柱内流动区域、环空流动区域和管柱下方流动区域交汇处节点的压力和速度信息,将计算得到的压力和速度信息作为施加的边界条件。
可选的,所述原始控制模型包括质量守恒模型和动量守恒模型,相应的,分别根据质量守恒模型和动量守恒模型获取第一特征线模型和第二特征线模型,根据第一特征线模型得到第一常微分计算模型,根据第二特征线模型得到第二常微分计算模型,根据第一常微分计算模型和第二常微分计算模型得到波动压力计算模型。
可选的,时间步长和距离步长的确定方法为:
获取管柱模型的单管压力波速和流道弹性系数;
对单管距离步长进行初步计算;
对单管距离步长进行调节,根据调节后的单管距离步长得到时间步长;
根据调节后的单管距离部长和管柱的单管长度确定节点数;
根据确定的节点数得到距离步长;
可选的,根据管柱模型运动距离、管柱模型最大速度和加速度确定管柱模型运动时间,根据管柱模型运动时间和时间步长确定时间计算次数。
可选的,采用厚壁筒弹性理论得到管柱模型的流道弹性系数。
一种钻井液瞬态波动压力确定系统,包括:
网格划分模块:被配置为能够对波动压力物理模型的流动区域进行按照设定的时间步长和距离步长进行网格划分;
常微分计算模型获取模块:被配置为能够获取网格节点原始控制模型的特征线模型,根据得到的特征线模型获得原始控制模型的常微分计算模型;
波动压力计算模型获取模块:被配置为能够获取特征线模型的在待求解时刻的前一时刻对应的特征节点的位置信息,根据与特征节点相邻的网格节点的速度和压力信息获取特征节点速度和压力信息,将获取的特征节点速度和压力信息带入常微分计算模型中,得到待计算时刻的网格节点的波动压力计算模型。
波动压力计算模块:被配置为能够利用波动压力计算模型得到不同网格节点在不同时刻的瞬态波动压力。
本发明的有益效果:
1.本发明的方法,利用特征节点的速度和压力信息得到网格节点的波动压力计算模型,无需网格节点全部落在特征线上,允许特征节点的位置不必与网格节点的位置相同,保证了数值计算过程的稳定收敛,提高了计算方法的健壮性,降低了网格划分的要求,避免了数值计算崩溃溢出。
2.本发明的方法,计算边界条件时,三流域汇交部分是三个流动区域的耦合,距离步长和时间步长设定时考虑了流道弹性系数,相较于稳态法将钻井液看做不可压缩和将管柱看做完全刚性体,解决了稳态法计算结果比实测数据偏高30%-50%的问题,计算精度可达90%以上。
附图说明
构成本申请的一部分的说明书附图用来提供对本申请的进一步理解,本申请的示意性实施例及其说明用于解释本申请,并不构成对本申请的限定。
图1为本发明实施例1波动压力物理模型示意图;
图2为本发明实施例1管柱内流动区域网格节点压力和速度计算示意图;
图3为本发明实施例1环空流动区域网格节点压力和速度计算示意图;
图4为本发明实施例1管柱下方流动区域网格节点压力和速度计算示意图;
图5为本发明实施例1钻井系统流动示意图;
图6为本发明实施例1时间步长和距离步长确定流程图;
图7为本发明实施例1流道内节点布置示意图;
图8为本发明实施例1管柱井口节点压力和速度计算示意图;
图9为本发明实施例1管柱底部节点压力和速度计算示意图;
图10为本发明实施例1环空流动区域顶部节点压力和速度计算示意图;
图11为本发明实施例1环空流动区域底部节点压力和速度计算示意图;
图12为本发明实施例1管柱下方流动区域顶部节点压力和速度计算示意图;
图13为本发明实施例1管柱下方流动区域井底节点压力和速度计算示意图;
图14为本发明实施例1波动压力计算流程图;
图15为本发明实施例1变截面连接节点压力和速度计算示意图;
具体实施方式
实施例1
本实施例公开了一种钻井液瞬态波动压力确定方法,包括以下步骤:
步骤1:建立波动压力物理模型,如图1所示,所述波动压力物理模型包括钻井模型,位于钻井模型内的管柱模型及填充在管柱模型和钻井模型之间的钻井液模型,所述波动压力物理模型具有三个流动区域,分别为位于管柱模型内的管柱内流动区域,位于管柱模型下方的管柱下方流动区域及位于管柱和钻井模型井壁之间的环空流动区域。
以管柱底部为Z轴原点,建立两套坐标系。管柱内和环空流体Z轴向上,管柱下部流体Z轴向下
本实施例中,所述波动压力物理模型基于以下假设:
(1)井内水力系统各流道中,泥浆的流动均为一元流动;
(2)钻井液可压缩,钻井液物性受温度和压力影响而改变;
(3)流道和泥浆均为线弹性的,即应力与应变成正比;
(4)3个流动区域的横截面积均受固壁材料的弹性力学参数和流域内压力影响而变化;
(5)管柱底部运动速度受管柱轴向弹性和流体作用在运动管柱壁面上的作用力影响,使得运动管柱底部速度不同于井口速度;
(6)井壁弹性受地层弹性、套管弹性和水泥弹性复合影响;
(7)计算流道中稳定流摩阻压降的公式适用于瞬变流。
步骤2:对波动压力物理模型的流动区域案子设定的时间步长和距离步长划分网格,得到多个网格节点。
对获取波动压力物理模型中流动区域原始控制模型的特征线模型,根据得到的特征线模型获得原始控制模型的常微分计算模型,根据得到的常微分计算模型得到波动压力计算模型。
具体的,原始控制模型包括质量守恒模型和动量守恒模型:
针对管柱内流动区域:
质量守恒模型:
Figure BDA0003069507380000041
动量守恒模型:
Figure BDA0003069507380000042
针对环空流动区域:
质量守恒模型:
Figure BDA0003069507380000051
动量守恒模型:
Figure BDA0003069507380000052
针对管柱下方流动区域:
质量守恒模型:
Figure BDA0003069507380000053
动量守恒模型:
Figure BDA0003069507380000054
其中v为流速;p为总压力;A为流道横截面积;z为轴向坐标;ρ为钻井液密度;Δp为摩阻压降,是流体流速和管柱运动速度vp的函数;vp为管柱下入速度,是时间的函数;K为钻井液体积模量,下标1、2、3分别表示三个不同流动区域。
公式(1)-(6)构成了三对一阶拟线性双曲型二元偏微分方程组。
针对管柱内流动区域的原始控制模型:
对原始控制模型进行离散:
原始控制模型中,记
Figure BDA0003069507380000055
则控制方程展开为:
Figure BDA0003069507380000056
Figure BDA0003069507380000057
对于管柱内流动区域来说,假定p2已知,即
Figure BDA0003069507380000058
Figure BDA0003069507380000059
为常数。
应用双曲型偏微分方程特征方向分析方法,得该方程组的第一特征线模型和第二特征线模型为:
Figure BDA0003069507380000061
Figure BDA0003069507380000062
物理意义为压力波在钻井液中的传播速度。
现场起下管柱引起的管内钻井液运动速度一般为10m/s数量级,
Figure BDA0003069507380000063
一般为1000m/s数量级。故,计算特征方向时可忽略v1项。
将两个原始控制模型进行线性组合,沿特征线模型的特征线方向有;
dz-(v1+c1)dt=0 (10)
a1dp1+[a1ρ1(v1+c1)-a1ρ1v1]dv1+[K(v1+c1)-H]dt=0 (11)
dz-(v1-c)dt=0 (12)
a1dp1+[a1ρ1(v1-c1)-a1ρ1v1]dv1+[K(v1-c1)-H]dt=0 (13)
式中,K=a1[Δp(v1,vp)-ρ1gcosθ]
Figure BDA0003069507380000064
其中公式(11)和公式(13)分别为第一常微分计算模型和第二常微分计算模型。
将公式(10)-(13)进行简化有:
dz-(v1+c)dt=0 (15)
Figure BDA0003069507380000065
dz-(v1-c)dt=0 (17)
Figure BDA0003069507380000066
获取特征线模型的在待求解时刻的前一时刻对应的特征节点的位置信息,根据位置信息获取特征节点速度和压力信息,将获取的特征节点速度和压力信息带入常微分计算模型中,得到待计算时刻的网格节点的波动压力计算模型。
如图2所示,P点为网格节点,tN+1时刻为其待求解时刻,A、Q、B分别为tN时刻已知压力和速度的三个相邻的网格节点,设自P点的特征线C+和C-分别与上一时刻交于R和S点,即R和S点分别为特征线模型在待求解时刻的前一时刻对应的特征节点
Figure BDA0003069507380000071
Figure BDA0003069507380000072
由特征线模型得,得到R点和S点的z轴坐标即位置信息。
Figure BDA0003069507380000073
Figure BDA0003069507380000074
采用插值法获得R点和S点的压力p和速度v
Figure BDA0003069507380000075
Figure BDA0003069507380000076
Figure BDA0003069507380000077
Figure BDA0003069507380000078
式中,θ=ΔtΔz,为松弛时空比例因子。
则第一常微分计算模型和第二常微分计算模型一阶近似为:
Figure BDA0003069507380000079
Figure BDA00030695073800000710
参数a、b、c、ρ、v等参数下标i表示在节点i和节点i+1形成的流道的参数。
Figure BDA0003069507380000081
Figure BDA0003069507380000082
ui-1=a1,i-1ρ1,i-1c1,i-1
ui+1=a1,i+1ρ1i+1c1,i+1
则波动压力计算模型为:
Figure BDA0003069507380000083
Figure BDA0003069507380000084
Figure BDA0003069507380000085
针对环空流动区域:
质量守恒模型:
Figure BDA0003069507380000086
动量守恒模型:
Figure BDA0003069507380000087
对原始控制模型进行离散
原始控制模型中,记
Figure BDA0003069507380000088
则控制模型展开为:
Figure BDA0003069507380000089
Figure BDA00030695073800000810
对于环空流动区域来说,假定p1已知,即
Figure BDA00030695073800000811
Figure BDA00030695073800000812
为常数。
应用双曲型偏微分方程特征方向分析方法,得到该方程组的第一特征线模型和第二特征线模型为:
Figure BDA0003069507380000091
Figure BDA0003069507380000092
物理意义为压力波在钻井液中的传播速度。
现场下管柱引起的管内钻井液运动速度一般为10m/s数量级,
Figure BDA0003069507380000093
一般为1000m/s数量级,因此,计算特征线模型时可忽略v2项。
将原始控制模型进行线性组合,沿特征线模型的特征线方向有:
dz-(v2+c2)dt=0 (37)
a2dp2+[a2ρ2(v2+c2)-a2ρ2v2]dv2+[K(v2+c2)-H]dt=0 (38)
dz-(v2-c2)dt=0 (39)
a2dp2+[a2ρ2(v2-c2)-a2ρ2ν2]dv2+[K(v2-c2)-H]dt=0 (40)
式中,K=a2[Δp(v2,vp)-ρ2gcosθ]
Figure BDA0003069507380000094
公式(38)和公式(40)为第一常微分计算模型和第二常微分计算模型。
公式(37)-(40)简化后:
dz-(v2+c2)dt=0 (41)
Figure BDA0003069507380000095
dz-(v2-c2)dt=0 (43)
Figure BDA0003069507380000096
如图3所示,P点为网格节点,tN+1时刻为其待求解时刻,A、Q、B分别为tN时刻已知压力和速度的三个相邻的网格节点,设自P点的特征线C+和C-分别与上一时刻交于R和S点,即R和S点分别为特征线模型在待求解时刻的前一时刻对应的特征节点
Figure BDA0003069507380000101
Figure BDA0003069507380000102
由特征线模型得,得到R点和S点的z轴坐标即位置信息。
Figure BDA0003069507380000103
Figure BDA0003069507380000104
采用插值法获得R点和S点的压力p和速度v
Figure BDA0003069507380000105
Figure BDA0003069507380000106
Figure BDA0003069507380000107
Figure BDA0003069507380000108
式中,θ=ΔtΔz,为环空松弛时空比例因子。
则第一常微分计算模型和第二常微分计算模型一阶近似为:
Figure BDA0003069507380000109
Figure BDA00030695073800001010
参数a、b、c、ρ、v等参数下标i表示在节点i和节点i+1形成的流道的参数。
Figure BDA0003069507380000111
Figure BDA0003069507380000112
ui-1=a2,i-1ρ2,i-1c2,i-1
ui+1=a2,i+1ρ2,i+1c2,i+1
则波动压力计算模型为:
Figure BDA0003069507380000113
Figure BDA0003069507380000114
Figure BDA0003069507380000115
针对管柱下方流动区域
原始控制模型的离散:
该区域内节点原始控制模型的离散形式与管柱内流动区域的离散形式相同,唯一不同之处在于该区域无b1相关项,因此在此不进行重复叙述。
如图4所示,P、A、R、Q、S、B节点的设置方式与管柱内流动区域的设置方式相同。
由特征线模型得,得到R点和S点的z轴坐标即位置信息。
Figure BDA0003069507380000116
Figure BDA0003069507380000117
则采用插值法获取R点和S点的压力速度信息。
Figure BDA0003069507380000121
Figure BDA0003069507380000122
Figure BDA0003069507380000123
Figure BDA0003069507380000124
式中θ=ΔtΔz为松弛时空比例因子。
方程组简化为:
Figure BDA0003069507380000125
p3,0-p3,S3,1c3,1(v3,0-v3,S)+[-c3,1Δp(v3,1)+c3,1ρ3,1gcosθ]Δt=0 (65)
dp33c3dv3+[c3Δp(v3,vp)-c3ρ3gcosθ]dt=0 (66)
dp33c3dv3+[-c3Δp(v3,vp)+c3ρ3gcosθ]dt=0 (67)
Figure BDA0003069507380000126
Figure BDA0003069507380000127
ui-1=p3,i-1c3,i-1
ui+1=ρ3,i+1c3,i+1
则波动压力计算模型为:
Figure BDA0003069507380000128
pP=pS+ui+1vP-ui+1vS-n (69)
划分网格的具体方法为:如图5所示,根据井身结构和运动管柱组合,分析三个流动区域内的流道组成情况。根据套管下入深度在井中的位置和井身结构类型,确定每个流道由几个单管组成,记录单管的外径、内径和长度等信息。如图所示,点划线表示运动管柱下入深度,裸眼井较简单,只有一种情况;套管+裸眼井有三种情况;套管衬管+裸眼井身结构最复杂,共5种情况。以Courant准则(ΔZi/Δt≤1)为约束条件,基于用户设定的时间步长,利用系统流道组合分析确定的数据,根据下述算法确定单管i的节点数、距离步长和调整后的时间步长。
如图6所示,所述设定的时间步长和距离步长的确定方法包括以下步骤:
获取管柱模型的单管压力波速ci和流道弹性系数;
对单管距离步长进行初步计算ΔZi=(ci+50)Δt;
若Li<2ΔZi,则对单管距离步长进行调节ΔZi=Li/2,根据调节后的单管距离步长得到时间步长:Δt=ΔZi/(ci+50);否则不进行调整。
根据调节后的单管距离部长和管柱的单管长度确定节点数Ni=(int)Li/ΔZi
根据确定的节点数得到距离步长:ΔZi=Li/Ni
其中,Li为管柱模型的单管i长度,Ni为单管节点数,ΔZi为单管i距离步长。
基于管柱运动距离、管柱运动最大速度和加速度确定管柱运动时间,基于时间步长确定时间计算次数即需要计算多少个时间部长。
Figure BDA0003069507380000131
Nt=(int)(tt+3)/Δt (143)
式中,L为管柱运动距离,vmax为管柱运动最大速度,a为管柱运动加速度。
本实施例中,离散方程
Figure BDA0003069507380000132
为泥浆的压缩系数,
Figure BDA0003069507380000133
为管道的弹性系数,泥浆的压缩系数与泥浆性质,温度,压力有关,一般应通过实测得特定条件下的压缩系数。研究表明,压力在0~49MPa范围内变化时,泥浆压缩系数对波动压力的影响很小,因此,为方便应用,实际计算环空波动压力时,以水在50℃,50MPa下的压缩系数代替泥浆的实际压缩系数:
k=0,39×10-9(Pa)-1 (137)
流道内节点布置如图7所示,对于流道弹性系数,采用厚壁筒弹性理论进行计算,
1.管柱内流道的弹性系数,仅考虑波动压力引起的流道变形,而管内外静压平衡:
Figure BDA0003069507380000134
式中R=D2/D1
2.裸眼井流道弹性系数:
Figure BDA0003069507380000141
3.套管与套管间环空流道的弹性系数,假设套管材料相同:
Figure BDA0003069507380000142
4.裸眼井壁与套管间环空流道的弹性系数:
Figure BDA0003069507380000143
式中ES,μS,Ef,μf分别为套管和地层的弹性模量和泊松比,计算中取Es=0.2068×1012Pa,μs=0.3,Ef=0.17237×1011,μf=0.28。
本实施例中,考虑了钻井液的可压缩性、管柱和井壁的弹性,相对于稳态法提高了计算精度。
步骤3:计算边界条件,对流动区域施加边界条件。
具体的,计算管柱内流动区域井口节点的压力、速度信息、计算环空流动区域井口节点的压力速度信息、计算管柱内流动区域、环空流动区域和管柱下方流动区域交汇处节点的压力和速度信息,将计算得到的压力和速度信息作为施加的边界条件。
管柱井口节点:
如图8所示,井口节点Z方向索引号为is,该节点为求解区域的下游边界,C+为过P点的特征线,可根据该节点原始控制模型得到,根据特征线得到R点的z轴坐标:
Figure BDA0003069507380000144
进而利用A点和Q点已知的压力和速度信息采用插值得到R点的压力和速度信息。
Figure BDA0003069507380000145
Figure BDA0003069507380000146
Figure BDA0003069507380000151
井口压力为大气压,表压为0,以表压表示压力,则井口节点的压力和速度即管柱井口的边界条件为:
Figure BDA0003069507380000152
Figure BDA0003069507380000153
管柱底部节点
管柱底部节点编号为0,为求解区域的上游边界点,求解节点如图9所示:
C-为过P点的特征线,根据P点的原始控制模型求得,特征节点S的在轴坐标为:
Figure BDA0003069507380000154
则利用Q点和B点的压力和速度信息采用插值法得到S点的压力和速度信息。
Figure BDA0003069507380000155
Figure BDA0003069507380000156
Figure BDA0003069507380000157
该节点为运动管柱内流动区域、环空流动区域和运动管柱下部流动区域的汇交点,需要应用汇交点边界条件求解该节点出的压力和流速。
环空流动区域井口节点
如图10所示,C+为过P点的特征线,由P点的原始控制模型得到,特征节点R点的Z轴坐标为:
Figure BDA0003069507380000158
采用插值法获得R点的压力和速度:
Figure BDA0003069507380000161
Figure BDA0003069507380000162
Figure BDA0003069507380000163
井口压力为大气压,表压为0,以表压表示压力,则井口节点的边界条件为:
Figure BDA0003069507380000164
Figure BDA0003069507380000165
环空流动区域底部节点
环空底部节点Z方向索引号为0,节点如图11所示。
过P点的特征线为C-,由P点的原始控制模型获得,特征节点S的z轴坐标为:
Figure BDA0003069507380000166
采用插值法获得S点的速度和压力:
Figure BDA0003069507380000167
Figure BDA0003069507380000168
Figure BDA0003069507380000169
该节点为运动管柱内流动区域、环空流动区域和运动管柱下部流动区域的汇交点,需要应用汇交点边界条件求解该节点出的压力和流速。
管柱下方流动区域顶部节点
运动管柱下方区域的顶部也是运动管柱的底部,设该位置出的节点Z方向索引号为0,
如图12所示,过P点的特征线为C-,由P点的原始控制模型求得,特征节点S的z轴坐标为:
Figure BDA0003069507380000171
采用插值法得到S点的压力和速度:
Figure BDA0003069507380000172
Figure BDA0003069507380000173
Figure BDA0003069507380000174
管柱下方流动区域井底节点:
井底节点Z方向的索引号为ib,如图13所示,
过P点的特征线为C+,由P点的原始控制模型得到,则特征节点R的z轴坐标为:
Figure BDA0003069507380000175
采用插值法得到R点的压力和速度:
Figure BDA0003069507380000176
Figure BDA0003069507380000177
Figure BDA0003069507380000178
井底处流速vp=0(vp即为
Figure BDA0003069507380000179
),得
Figure BDA00030695073800001710
管柱底部三流动区域汇交区域
管柱流动区域底部节点边界条件满足:
Figure BDA00030695073800001711
环空流动区域底部节点边界条件满足:
Figure BDA0003069507380000181
管柱下方流动区域顶部节点边界条件满足:
p3,0-p3,s3,1c3,1(v3,0-v3,S)+[-c3,1Δp(v3,1)+c3,1ρ3,1gcos|Δt=0
井底汇交点处三个流动区域的流量之和为运动管柱单位时间排开的钻井液体积,即
v1,0A1,0+v2,0A2,0+v3,0A3,0=vp(Ap-Ae) (101)
环空流动区域底部节点压力等于管柱下方区域顶部节点压力减去钻头和井眼形成的窄环空压降,即
Figure BDA0003069507380000182
运动管柱内区域底部节点压力等于运动管柱下方区域顶部节点压力减去钻头喷嘴压降,即
Figure BDA0003069507380000183
式中,Aa是钻头与井眼之间的间隙面积,Aε为钻头喷嘴截面面积。
公式(99)-公式(103)的六个方程含有六个未知量,p1,0、p2,0、p3,0、v1,0、v2,0、v3,0,且摩阻压降也是未知量的函数,为降低计算难度,用上一时间节点的速度计算摩阻压降梯度,令:
Figure BDA0003069507380000191
Figure BDA0003069507380000192
K3=[-c3,1Δp(v3,1)+c3,1ρ3,1gcosθ]Δt (106)
Figure BDA0003069507380000193
Figure BDA0003069507380000194
则上述六个方程组简化为:
a1,1(p1,0-p1,s)-a1,1ρ1,1c1,1(v1,0-v1,s)+K1+L1v1,0=0 (109)
a2,1(p2,0-p2,s)-a2,1ρ2,1c2,1(v2,0-v2,s)+K2+L2v2,0=0 (110)
p3,0-p3,s-u3,1(v3,0-v3,s)+K3=0 (111)
v1,0A1,0+v2,0A2,0+v3,0A3,0=vp(Ap-Ae) (112)
Figure BDA0003069507380000195
Figure BDA0003069507380000196
令B1=-u1,1+L1,B2=-u2,1+L2,B3=-u3,1,C1=u1,1v1,s+K1-a1,1p1,s,C2=u2,i1v2,s+K2-a2,1p2,s
Figure BDA0003069507380000197
Qt=vp(Ap-Ae)则方程(109)-方程(114)简化为:
Figure BDA0003069507380000201
Figure BDA0003069507380000202
Figure BDA0003069507380000203
v1,0A1,0+v2,0A2,0+v3,0A3,0=vp(Ap-Ae) (118)
Figure BDA0003069507380000204
Figure BDA0003069507380000205
取v2,0最小值vmin为0,最大值vmax为vp(Ap-Ae)/A2,0,用二分法求解上述6个代数方程组,步骤如下:
①v2,0=(vmin+vmax)/2;
②将v2,0代入式(116)求得p2,0
③将p2,0代入式(119),求得p3,0
④将p3,0代入式(117),求得v3,0
⑤将v2,0和v3,0代入式(118)求得v1,0
⑥将v1,0代入式(115)求得
Figure BDA0003069507380000206
⑦将p3,0和v1,0代入式(120)
Figure BDA0003069507380000207
⑧比较
Figure BDA0003069507380000208
Figure BDA0003069507380000209
Figure BDA00030695073800002010
则vmin=ν2,0,否则,vmax=v2,0,重复步骤①~⑧,直至
Figure BDA00030695073800002011
步骤4:如图14所示,输入井身数据(包括各井段套管顶深、套管下深、套管外径、套管内径、井眼尺寸和井眼深度,如果是裸眼段则套管外径输入“-1”。)、输入管柱数据(运动管柱组合外径、内径、长度、运动最大速度Vmax和加速度a。)、输入钻井液数据(钻井液密度,钻井液流变性:K、n,钻井液等温压缩系数(缺省取0.39×10-9Pa-1)。)、输入材料参数(地层弹性模量(缺省取Ef=0.17237×1011Pa),地层泊松比(缺省取0.28),钢弹性模量(缺省取Es=0.2068×1012Pa),钢泊松比(缺省取0.3))。
初始化初始条件,节点压力为依据垂深计算的静液柱压力,节点速度为0,根据划分的网格和施加的边界条件,利用步骤2得到的三个流动区域的波动压力计算模型循环迭代求解各个网格节点在不同时刻的速度和压力,压力即为波动压力。
本实施例中,如图15所示,管柱内流动区域、环空流动区域和管柱下方流动区域存在变截面。
管柱内流动区域和环空流动区域变截面处连接点的离散方程相同,如图11所示,节点i位于流道变截面处,该节点处存在两个不同特征线确定的流速vpu和vpd以及1个压力pp,沿特征线得差分方程:
Figure BDA0003069507380000211
Figure BDA0003069507380000212
AuvPv=AdvPd-(Au-Ad)vp (123)
解得波动压力计算模型为:
Figure BDA0003069507380000213
Figure BDA0003069507380000214
Figure BDA0003069507380000215
管柱下方流动区域连接点离散方程:
Figure BDA0003069507380000221
Figure BDA0003069507380000222
AuvPu=AdvPd(129)
Figure BDA0003069507380000223
Figure BDA0003069507380000224
pP=pS+ui+1vPd-ui+1vS-n (132)
Figure BDA0003069507380000225
Figure BDA0003069507380000226
ui-1=ρ3,i-1c3,i-1 (135)
ui+1=ρ3,i+1c3,i+1 (136)
实施例2:
本实施例公开了一种钻井液瞬态波动压力确定系统,包括:
网格划分模块:被配置为能够对波动压力物理模型的流动区域进行按照设定的时间步长和距离步长进行网格划分;
常微分计算模型获取模块:被配置为能够获取网格节点原始控制模型的特征线模型,根据得到的特征线模型获得原始控制模型的常微分计算模型;
波动压力计算模型获取模块:被配置为能够获取特征线模型的在待求解时刻的前一时刻对应的特征节点的位置信息,根据与特征节点相邻的网格节点的速度和压力信息获取特征节点速度和压力信息,将获取的特征节点速度和压力信息带入常微分计算模型中,得到待计算时刻的网格节点的波动压力计算模型;
波动压力计算模块:被配置为能够利用波动压力计算模型得到不同网格节点在不同时刻的瞬态波动压力。
上述虽然结合附图对本发明的具体实施方式进行了描述,但并非对本发明保护范围的限制,所属领域技术人员应该明白,在本发明的技术方案的基础上,本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。

Claims (10)

1.一种钻井液瞬态波动压力确定方法,其特征在于,包括以下步骤:
对波动压力物理模型的流动区域进行按照设定的时间步长和距离步长进行网格划分,得到多个网格节点;
获取网格节点原始控制模型的特征线模型,根据得到的特征线模型获得原始控制模型的常微分计算模型;
获取特征线模型的在待求解时刻的前一时刻对应的特征节点的位置信息,根据与特征节点相邻的网格节点的速度和压力信息获取特征节点速度和压力,将获取的特征节点速度和压力带入常微分计算模型中,得到待计算时刻的网格节点的波动压力计算模型;
利用波动压力计算模型得到不同网格节点在不同时刻的瞬态波动压力。
2.如权利要求1所述的一种钻井液瞬态波动压力确定方法,其特征在于,根据与特征节点相邻的网格节点的压力和速度信息采用插值算法得到特征节点的压力和速度信息。
3.如权利要求1所述的一种钻井液瞬态波动压力确定方法,其特征在于,计算边界条件,对流动区域施加边界条件,设定初始条件,对网格节点进行求解变量初始化并按照时间步长,利用波动压力计算模型推演得到不同时刻、不同网格节点处的瞬态波动压力。
4.如权利要求3所述的一种钻井液瞬态波动压力确定方法,其特征在于,所述流动区域划分为管柱内流动区域、环空流动区域及管柱下方流动区域,运动管柱流动区域为管柱模型内的区域,环空流动区域为管柱模型与井壁模型之间的区域,管柱下方流动区域为管柱模型下方的流动区域。
5.如权利要求4所述的一种钻井液瞬态波动压力确定方法,其特征在于,计算管柱内流动区域井口节点的压力、速度信息、计算环空流动区域井口节点的压力速度信息、计算管柱内流动区域、环空流动区域和管柱下方流动区域交汇处节点的压力和速度信息,将计算得到的压力和速度信息作为施加的边界条件。
6.如权利要求1所述的一种钻井液瞬态波动压力确定方法,其特征在于,所述原始控制模型包括质量守恒模型和动量守恒模型,相应的,分别根据质量守恒模型和动量守恒模型获取第一特征线模型和第二特征线模型,根据第一特征线模型得到第一常微分计算模型,根据第二特征线模型得到第二常微分计算模型,根据第一常微分计算模型和第二常微分计算模型得到波动压力计算模型。
7.如权利要求1所述的一种钻井液瞬态波动压力确定方法,其特征在于,时间步长和距离步长的确定方法为:
获取管柱模型的单管压力波速和流道弹性系数;
对单管距离步长进行初步计算;
对单管距离步长进行调节,根据调节后的单管距离步长得到时间步长;
根据调节后的单管距离部长和管柱的单管长度确定节点数;
根据确定的节点数得到距离步长;
8.如权利要求7所述的一种钻井液瞬态波动压力确定方法,其特征在于,根据管柱模型运动距离、管柱模型最大速度和加速度确定管柱模型运动时间,根据管柱模型运动时间和时间步长确定时间计算次数。
9.如权利要求7所述的一种钻井液瞬态波动压力确定方法,其特征在于,采用厚壁筒弹性理论得到管柱模型的流道弹性系数。
10.一种钻井液瞬态波动压力确定系统,其特征在于,包括:
网格划分模块:被配置为能够对波动压力物理模型的流动区域进行按照设定的时间步长和距离步长进行网格划分;
常微分计算模型获取模块:被配置为能够获取网格节点原始控制模型的特征线模型,根据得到的特征线模型获得原始控制模型的常微分计算模型;
波动压力计算模型获取模块:被配置为能够获取特征线模型的在待求解时刻的前一时刻对应的特征节点的位置信息,根据与特征节点相邻的网格节点的速度和压力信息获取特征节点速度和压力信息,将获取的特征节点速度和压力信息带入常微分计算模型中,得到待计算时刻的网格节点的波动压力计算模型;
波动压力计算模块:被配置为能够利用波动压力计算模型得到不同网格节点在不同时刻的瞬态波动压力。
CN202110535627.7A 2021-05-17 2021-05-17 一种钻井液瞬态波动压力确定方法及系统 Active CN113191102B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110535627.7A CN113191102B (zh) 2021-05-17 2021-05-17 一种钻井液瞬态波动压力确定方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110535627.7A CN113191102B (zh) 2021-05-17 2021-05-17 一种钻井液瞬态波动压力确定方法及系统

Publications (2)

Publication Number Publication Date
CN113191102A true CN113191102A (zh) 2021-07-30
CN113191102B CN113191102B (zh) 2024-04-09

Family

ID=76982170

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110535627.7A Active CN113191102B (zh) 2021-05-17 2021-05-17 一种钻井液瞬态波动压力确定方法及系统

Country Status (1)

Country Link
CN (1) CN113191102B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114893174A (zh) * 2022-04-07 2022-08-12 中海石油(中国)有限公司海南分公司 一种基于多因素耦合的砂岩储层可压裂性评价方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5408638A (en) * 1990-12-21 1995-04-18 Hitachi, Ltd. Method of generating partial differential equations for simulation, simulation method, and method of generating simulation programs
CN101710468A (zh) * 2009-12-16 2010-05-19 西南石油大学 钻井模拟器压力控制模拟方法
US20150066455A1 (en) * 2013-08-27 2015-03-05 Halliburton Energy Services, Inc. Proppant Transport Model for Well System Fluid Flow Simulations
CN107035327A (zh) * 2017-05-09 2017-08-11 中国石油大学(北京) 确定启停泵过程中瞬时波动压力的方法和装置
US20180371902A1 (en) * 2017-06-21 2018-12-27 Schlumberger Technology Corporation Downhole characterization of formation pressure

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5408638A (en) * 1990-12-21 1995-04-18 Hitachi, Ltd. Method of generating partial differential equations for simulation, simulation method, and method of generating simulation programs
CN101710468A (zh) * 2009-12-16 2010-05-19 西南石油大学 钻井模拟器压力控制模拟方法
US20150066455A1 (en) * 2013-08-27 2015-03-05 Halliburton Energy Services, Inc. Proppant Transport Model for Well System Fluid Flow Simulations
CN107035327A (zh) * 2017-05-09 2017-08-11 中国石油大学(北京) 确定启停泵过程中瞬时波动压力的方法和装置
US20180371902A1 (en) * 2017-06-21 2018-12-27 Schlumberger Technology Corporation Downhole characterization of formation pressure

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
周开吉,钟兵,袁其骥,柳萍: "拟瞬态井内波动压力预测模型", 西南石油学院学报, no. 04, 30 November 1995 (1995-11-30), pages 58 - 64 *
樊洪海,褚元林,刘希圣: "起下钻时井眼内动态波动压力的预测", 石油大学学报(自然科学版), no. 05, 5 October 1995 (1995-10-05), pages 35 - 41 *
王竹林;陈恭洋;: "钻井泵瞬态波动压力研究", 油气地球物理, no. 02, 26 April 2009 (2009-04-26), pages 39 - 42 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114893174A (zh) * 2022-04-07 2022-08-12 中海石油(中国)有限公司海南分公司 一种基于多因素耦合的砂岩储层可压裂性评价方法
CN114893174B (zh) * 2022-04-07 2022-12-16 中海石油(中国)有限公司海南分公司 一种基于多因素耦合的砂岩储层可压裂性评价方法

Also Published As

Publication number Publication date
CN113191102B (zh) 2024-04-09

Similar Documents

Publication Publication Date Title
CN111581786B (zh) 用于分析缝洞串联模式双孔复合储层的试井解释模型的试井解释方法
Ouyang et al. A general single-phase wellbore/reservoir coupling model for multilateral wells
CN113191102B (zh) 一种钻井液瞬态波动压力确定方法及系统
CN104727811A (zh) 鱼骨状水平井分段耦合的产能预测方法
CN110847894B (zh) 一种井下节流气井流压的确定方法
Bhalla Coiled tubing extended reach technology
Kong et al. A new model for predicting dynamic surge pressure in gas and drilling mud two-phase flow during tripping operations
Yang et al. Dynamic response mechanism of borehole breathing in fractured formations
Liu et al. Study on the coupling model of wellbore temperature and pressure during the production of high temperature and high pressure gas well
CN109812236B (zh) 一种确定异形井眼中的井眼清洁效果的方法
CN108104800A (zh) 一种井下节流气井的节流器入口压力计算方法和装置
Yuan et al. An experimental and analytical study of single-phase liquid flow in a horizontal well
CN114254520A (zh) 一种超深井油套环空液面高度的确定方法
Haige et al. Study on steady surge pressure for yield-pseudoplastic fluid in a concentric annulus
Wang et al. Flow simulation of a horizontal well with two types of completions in the frame of a wellbore–annulus–reservoir model
White et al. An analytical concept of the static and dynamic parameters of intermittent gas lift
Hu et al. Mechanism and prevention method of drill string uplift during shut-in after overflow in an ultra-deep well
CN112177598B (zh) 一种考虑压裂液压缩性的地层起裂压力预测方法
CN111520132B (zh) 一种确定地层中洞距离的方法及系统
CN109812237B (zh) 一种满足异形井眼清洁的钻井液排量确定方法
Haciislamoglu et al. Effect of pipe eccentricity on surge pressures
Tang et al. Transient dynamic characteristics of the gas-lift unloading process
Ding et al. Prediction and analysis of sustained annular pressure of gas wells considering the influence of production system and liquid thermodynamic properties
Owusu et al. BUCKLEY-LEVERETT DISPLACEMENT THEORY FOR WATERFLOODING PERFORMANCE IN STRATIFIED RESERVOIR.
Hatletvedt et al. Modelling of a hydro-pneumatic system for heave compensation

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant