CN113868978A - 一种页岩基质内气体非稳态微观传输特征的确定方法 - Google Patents

一种页岩基质内气体非稳态微观传输特征的确定方法 Download PDF

Info

Publication number
CN113868978A
CN113868978A CN202111149220.7A CN202111149220A CN113868978A CN 113868978 A CN113868978 A CN 113868978A CN 202111149220 A CN202111149220 A CN 202111149220A CN 113868978 A CN113868978 A CN 113868978A
Authority
CN
China
Prior art keywords
gas
transmission
adsorption
density
calculation model
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
CN202111149220.7A
Other languages
English (en)
Other versions
CN113868978B (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.)
China University of Petroleum Beijing
Original Assignee
China University of Petroleum Beijing
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 China University of Petroleum Beijing filed Critical China University of Petroleum Beijing
Priority to CN202111149220.7A priority Critical patent/CN113868978B/zh
Publication of CN113868978A publication Critical patent/CN113868978A/zh
Application granted granted Critical
Publication of CN113868978B publication Critical patent/CN113868978B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00—Details relating to the application field
    • G06F2113/08—Fluids
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Fluid Mechanics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Computing Systems (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Solid-Sorbent Or Filter-Aiding Compositions (AREA)

Abstract

本发明提供了一种页岩基质内气体非稳态微观传输特征的确定方法。该方法包括：获取目标页岩基质无吸附作用和/或吸附作用气体传输参数；获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型，和/或获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型；获取孔隙内无吸附作用和/或吸附作用气体密度关于压力的计算模型；在视渗透率计算模型和气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附和/或吸附气体传输模型；基于气体传输模型，结合气体传输参数，确定目标页岩基质内气体微观传输特征。

Description

一种页岩基质内气体非稳态微观传输特征的确定方法
技术领域
本发明属于页岩气开发技术领域，特别涉及一种页岩基质内气体非稳态微观传输特征的确定方法。
背景技术
面对中国乃至世界强劲的能源需求，页岩气作为重要的非常规油气资源，愈加受到重视，其开发规模逐渐扩大。页岩气储层与非常规储层存在明显差异，微纳尺度孔隙发育，呈现低孔、低渗特征。由于微尺度效应的存在，储层条件下页岩基质内气体传输机理复杂。现阶段，对页岩气微观传输特征的认识有待加深。
为研究页岩气微观传输特征，现有研究者开展了相应研究。页岩基质渗透率是描述页岩传输特征的重要参数。Javadpour等人基于滑脱传输、努森扩散等传输机理，建了页岩纳米孔气体传输模型，首次提出了视渗透率的概念(Javadpour F.Nanopores andapparent permeability of gas flow in mudrocks(shales and siltstone).Journalof Canadian Petroleum Technology2009；48(08):16-21)。Wu等人拓展了视渗透率理论，通过权重系数整合了滑脱传输和努森扩散传输(Wu K,Chen Z,Li X,Guo C,Wei M.A modelfor multiple transport mechanisms through nanopores of shale gas reservoirswith real gas effect–adsorption-mechanic coupling.International Journal ofHeat and Mass Transfer 2016；93:408-26)，并获得了考虑滑脱传输、努森扩散传输和表面扩散传输的视渗透率模型(吴克柳,李相方,陈掌星.页岩气纳米孔气体传输模型.石油学报2015；36(7):837-48)，实现了多种传输模式的高度统一。
现有CN108710723A公开了一种页岩储层多孔介质表观渗透率的计算方法，通过叠加不同尺寸毛细管内滑脱传输、努森扩散传输和表面扩散传输渗透率，获得了岩心尺度表观渗透率；在该发明专利公开文本中，计算参数需预先给定，且未曾涉及真实页岩基质内非稳态传输的表征。
现有CN110472348A公开了一种页岩气藏非稳态渗流模型的建立办法，给出了页岩基质和裂缝的非稳态渗流模型；该模型依托达西渗流规律建立，未涉及努森扩散传输等微观传输规律，且未考虑页岩微纳米孔隙内吸附相和体相的协同传输特征。在页岩孔隙中，体相和吸附相气体同时传输，且相互影响。孔隙内气体传输是孔隙和气体本身共同完成的一个过程。所以，为准确表征气体传输，除了孔隙特征，需考虑真实气体效应，采用对应的压缩系数、粘度和偏差因子等真实气体参数。
现有基质传输模型未能在与实际页岩样品气体非稳态微观传输相对应的基础上，综合考虑滑脱流、努森扩散和表面扩散等传输规律及压缩系数、粘度和偏差因子等真实气体参数，准确表征真实页岩基质内微纳米孔内吸附相和体相的协同传输特征。因此，不能真实反映页岩气储层基质内无吸附作用和有吸附作用的真实气体的传输机理，影响对气体传输预测的准确性。
发明内容
本发明的目的在于提供一种考虑滑脱流、努森扩散和表面扩散等传输规律及气体传输参数，准确确定真实页岩基质内微纳米孔内吸附相和体相的协同传输特征的方法。
为了实现上述目的，本发明提供了一种页岩基质内气体非稳态微观传输特征的确定方法，其中，该方法包括：
获取目标页岩基质无吸附作用气体传输参数；和/或；获取目标页岩基质吸附作用气体传输参数；
获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型；和/或；获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型；
获取孔隙内无吸附作用气体密度关于压力的计算模型；和/或；获取孔隙内吸附作用气体密度关于压力的计算模型；
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型；和/或；在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型；
基于无吸附气体传输模型，结合目标页岩基质无吸附作用气体传输参数，确定目标页岩基质内无吸附气体微观传输特征；和/或；基于吸附气体传输模型，结合目标页岩基质吸附作用气体传输参数，确定目标页岩基质内吸附气体微观传输特征。
本发明提供的技术方案，建立页岩基质内考虑综合考虑滑脱流、努森扩散和表面扩散等传输规律及气体参数(在优选技术方案中考虑了压缩系数、粘度和偏差因子等)，准确表征真实页岩基质内微纳米孔内气体传输特征的非稳态气体传输模型，进而确定真实页岩基质内微纳米孔内吸附相和体相的协同传输特征，为页岩基质内气体的传输机理、页岩储层的产能预测提供理论支撑。与现有技术相比，本发明提供的技术方案具备如下优势：
(1)在气体传输模型中综合考虑了滑脱流、努森扩散和表面扩散等传输规律，提高了确定得到的页岩基质内微观气体传输特征的准确性。
(2)本发明提供的技术方案以气体密度代替压力作为传输模型变量，实现了滑脱流、努森扩散和表面扩散等微观传输模型在宏观气体传输模型中的精准应用，可以实现页岩基质内体相和吸附相的协同传输模拟。
(3)本发明提供的技术方案结合气体传输参数，实现了页岩基质内无吸附作用和有吸附作用的气体的传输特征的确定，有助于构建微观传输理论与实际页岩气储层气体渗流间的联系，为实际页岩储层的产能预测提供理论支撑。
附图说明
图1为本发明一实施例提供的页岩基质内气体非稳态微观传输特征的确定方法的流程示意图。
图2为本发明一实施例孔隙内吸附相和体相气体分布示意图。
图3A为实施例1中3MPa、9MPa压力下样品1氦气传输压降曲线图。
图3B为实施例1中4MPa、8MPa压力下样品1氦气传输验证压降曲线图。
图3C为实施例1中3MPa、9MPa压力下样品2氦气传输压降曲线图。
图3D为实施例1中4MPa、8MPa压力下样品2氦气传输验证压降曲线图。
图3E为实施例1中3MPa、9MPa压力下样品3氦气传输压降曲线图。
图3F为实施例1中4MPa、8MPa压力下样品3氦气传输验证压降曲线图。
图4A为实施例1中3MPa压力下样品1、2、3甲烷传输压降曲线图。
图4B为实施例1中5MPa压力下样品1、2、3甲烷传输压降曲线图。
图4C为实施例1中7MPa压力下样品1、2、3甲烷传输压降曲线图。
图4D为实施例1中8MPa压力下样品1、2、3甲烷传输压降曲线图。
图5A为实施例1中不同时刻基质颗粒内氦气压力和努森数分布图。
图5B为实施例1中不同时刻基质颗粒内氦气渗透率分布图。
图6A为实施例1中不同时刻基质颗粒内甲烷体相压力、总密度分布特征图。
图6B为实施例1中不同时刻总密度、吸附相密度和体相密度分布特征图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚完整的描述。显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其它实施例，都属于本发明的保护范围。
针对砂砾岩储层评价中，为了区分储层类型来建立可靠的储层参数和含油饱和度评价模型，以提高储层参数评价精度，在全岩矿物分析获取的浊沸石含量基础上，分析浊沸石胶结砂砾岩储层和不含浊沸石胶结砂砾岩储层的测井响应特征，在此基础上，提取对区分储层是否含浊沸石较敏感的测井参数，以建立利用测井资料准确识别浊沸石胶结砂砾岩储层的方法，进而提出了本发明用以实现浊沸石胶结砂砾岩储层测井识别的方法。
参见图1，本发明一具体实施方式提供了一种页岩基质内气体非稳态微观传输特征的确定方法，其中，该方法包括：
步骤S1：获取目标页岩基质无吸附作用气体传输参数；和/或；获取目标页岩基质吸附作用气体传输参数；
步骤S2：获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型；和/或；获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型；
步骤S3：获取孔隙内无吸附作用气体密度关于压力的计算模型；和/或；获取孔隙内吸附作用气体密度关于压力的计算模型；
步骤S4：在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型；和/或；在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型；
步骤S5：基于无吸附气体传输模型，结合目标页岩基质无吸附作用气体传输参数，确定目标页岩基质内无吸附气体微观传输特征；和/或；基于吸附气体传输模型，结合目标页岩基质吸附作用气体传输参数，确定目标页岩基质内吸附气体微观传输特征。
在一具体实施方式中，页岩基质内气体非稳态微观传输特征的确定方法包括：
获取目标页岩基质无吸附作用气体传输参数；
获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型；
获取孔隙内无吸附作用气体密度关于压力的计算模型；
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型；
基于无吸附气体传输模型，结合目标页岩基质无吸附作用气体传输参数，确定目标页岩基质内无吸附气体微观传输特征。
在一具体实施方式中，页岩基质内气体非稳态微观传输特征的确定方法包括：
获取目标页岩基质吸附作用气体传输参数；
获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型；
获取孔隙内吸附作用气体密度关于压力的计算模型；
在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型；
基于吸附气体传输模型，结合目标页岩基质吸附作用气体传输参数，确定目标页岩基质内吸附气体微观传输特征
在一具体实施方式中，页岩基质内气体非稳态微观传输特征的确定方法包括：
获取目标页岩基质无吸附作用气体传输参数和吸附作用气体传输参数；
获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型；获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型；
获取孔隙内无吸附作用气体密度关于压力的计算模型和吸附作用气体密度关于压力的计算模型；
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型；在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型；
基于无吸附气体传输模型，结合目标页岩基质无吸附作用气体传输参数，确定目标页岩基质内无吸附气体微观传输特征；基于吸附气体传输模型，结合目标页岩基质吸附作用气体传输参数，确定目标页岩基质内吸附气体微观传输特征。
进一步地，无吸附作用气体传输参数包括孔隙度、孔径和迂曲度。
进一步地，所述吸附作用气体传输参数包括孔隙度、孔径、迂曲度、兰氏压力、兰氏吸附量和表面扩散系数。
在一具体实施方式中，结合室内模拟实验结果，通过拟合的方式获取无吸附作用气体传输参数和/或吸附作用气体传输参数；例如参照CN112362552A中公开的页岩基质渗透率测定装置以及的页岩基质渗透率测定方法进行室内模拟实验获取压降曲线，进而通过拟合的方式获取无吸附作用气体传输参数和/或吸附作用气体传输参数。
进一步地，获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型包括：
获取滑脱流视渗透率计算模型；
获取努森扩散视渗透率计算模型；
获取单组分无吸附作用气体总视渗透率计算模型。
进一步地，获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型包括：
获取滑脱流视渗透率计算模型；
获取努森扩散视渗透率计算模型；
获取表面扩散视渗透率计算模型；
获取单组分吸附作用气体总视渗透率计算模型。
进一步地，滑脱流视渗透率计算模型为滑脱流视渗透率关于孔隙度、孔径、迂曲度、努森数以及气体滑脱常数的计算模型；
更进一步地，滑脱流视渗透率计算模型为：
Figure BDA0003286587200000061
式中，kv为滑脱流视渗透率，m2；kn为努森数，无因次；τ为迂曲度，无因次；bc为气体滑脱常数，无因次；φ为样品孔隙度，无因次；rt为孔半径，m；
其中，气体滑脱常数bc可以利用本领域常规方式进行确定，例如基于二阶滑脱边界条件进行确定(参见吴克柳,李相方,陈掌星,李靖,梁羽丰,吴修华.页岩气复杂孔裂隙真实气体传输机理和数学模型.中国科学:技术科学2016；46(8):851-863)；通常bc可以取值为-1。
进一步地，努森扩散视渗透率计算模型为努森扩散视渗透率关于气体参数、孔隙度、孔径和迂曲度的计算模型；
更进一步地，气体参数包括单气体黏度、气体压缩系数、气体分子直径、偏差因子和摩尔质量；
再进一步地，努森扩散视渗透率计算模型为：
Figure BDA0003286587200000071
式中，kk为努森扩散视渗透率，m2；τ为迂曲度，无因次；bc为气体滑脱常数，无因次；φ为样品孔隙度，无因次；rt为孔半径，m；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；δ为气体分子直径与孔径比值，无因次；Df为表面粗糙度，无因次；Z为偏差因子，无因次；R为气体常数，通常取8.314472m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol；
其中，表面粗糙度Df的取值按照本领域常规方式进行即可，在一具体实施方式中取值2.5。
进一步地，单组分无吸附作用气体总视渗透率计算模型为总视渗透率关于努森扩散视渗透率和滑脱流视渗透率的计算模型；
更进一步地，单组分无吸附作用气体总视渗透率计算模型为：
k＝kv·ωv+kk·ωk，其中，
Figure BDA0003286587200000072
式中，k为单组分无吸附作用气体总视渗透率，m2；kv为滑脱流视渗透率，m2；kk为努森扩散视渗透率，m2；ωv为滑脱流权重系数，无因次；ωk努森扩散权重系数，无因次；kn为努森数，无因次。
进一步地，表面扩散视渗透率计算模型为表面扩散视渗透率关于气体参数、孔隙度、孔径和迂曲度的计算模型；
更进一步地，气体参数包括单气体黏度、吸附相气体浓度、游离态吸附气传输对应的孔隙半径、偏差因子和摩尔质量；
再进一步地，表面扩散视渗透率计算模型为：
Figure BDA0003286587200000073
式中，ks为表面扩散视渗透率，m2；Ds为吸附气表面扩散系数，m2/s；τ为迂曲度，无因次；p为压力，Pa；M为摩尔质量，kg/mol；φ为样品孔隙度，无因次；μ为气体黏度，Pa·s；Ccon为吸附相气体浓度，kg/m3；Z为偏差因子，无因次；R为气体常数，通常取8.314472m3·Pa/(K·mol)；T为温度，K；rt为孔半径，m；reff为游离态吸附气传输对应的孔隙半径，m；
其中，吸附相气体浓度Ccon可以利用本领域常规方式进行确定，在一具体实施方式中，吸附相气体浓度Ccon可表示为：
Figure BDA0003286587200000081
式中，Ccon为吸附相气体浓度，kg/m3；θ为吸附气壁面覆盖度，无因次；dm为分子直径，m；NA为阿伏伽德罗常数；
其中，游离态吸附气传输对应的孔隙半径reff可以利用本领域常规方式进行确定，在一具体实施方式中，游离态吸附气传输对应的孔隙半径reff可表示为：
reff＝(rt-dmθ)
式中，rt为孔半径，m；reff为游离态吸附气传输对应的孔隙半径，m；dm为分子直径，m；θ为吸附气壁面覆盖度，无因次。
进一步地，单组分吸附作用气体总视渗透率计算模型为总视渗透率关于努森扩散视渗透率、滑脱流视渗透率和表面扩散视渗透率的计算模型；
更进一步地，单组分吸附作用气体总视渗透率计算模型为：
k＝kv·ωv+kk·ωk+ks，其中，
Figure BDA0003286587200000082
式中，k为单组分吸附作用气体总视渗透率，m2；kv为滑脱流视渗透率，m2；kk为努森扩散视渗透率，m2；ωv为滑脱流权重系数，无因次；ωk努森扩散权重系数，无因次；ks为表面扩散视渗透率，m2；kn为努森数，无因次。
在一具体实施方式中，各类式渗透率的计算模型可以参考下述文献获得：Wu K,Chen Z,Li X,Guo C,Wei M.A model for multiple transport mechanisms throughnanopores of shale gas reservoirs with real gas effect–adsorption-mechaniccoupling.International Journal of Heat and Mass Transfer 2016；93:408-26。
其中，孔隙内无吸附作用气体密度关于压力的计算模型可以基于赋存状态分析结果进行确定；
进一步地，孔隙内无吸附作用气体密度关于压力的计算模型基于气体状态方程确定；
更进一步地，孔隙内无吸附作用气体密度关于压力的计算模型为：
Figure BDA0003286587200000083
式中，ρ为气体密度，kg/m3；p为压力，Pa；M为摩尔质量，kg/mol；Z为偏差因子，无因次；R为气体常数，通常取8.314472m3·Pa/(K·mol)；T为温度，K。
其中，孔隙内吸附作用气体密度关于压力的计算模型可以基于赋存状态分析结果进行确定；
进一步地，孔隙内吸附作用气体密度关于压力的计算模型结合孔隙内吸附相和体相气体赋存特征综合确定(孔隙内吸附相和体相气体分布状态参见图2)；
更进一步地，孔隙内吸附作用气体密度关于压力的计算模型为：
Figure BDA0003286587200000091
其中，
Figure BDA0003286587200000092
式中，ρ为基质中游离态和吸附态气体的综合密度，kg/m3；p为压力，Pa；pL为兰氏压力，Pa；M为摩尔质量，kg/mol；Z为偏差因子，无因次；R为气体常数，通常取8.314472m3·Pa/(K·mol)；T为温度，K；rt为孔半径，m；nL为兰氏吸附量，mol/kg；θ为吸附气壁面覆盖度，无因次；ρs为目标页岩样品密度，g/mL；φ为样品孔隙度，无因次；dm为分子直径，m。
在一具体实施方式中，努森数可表示为：
Figure BDA0003286587200000093
其中，
Figure BDA0003286587200000094
式中，λ为分子自由程，m；μc为气体在1.01325×105Pa和423K条件下的黏度，Pa·s；Z为偏差因子，无因次；R为气体常数，m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol；p为压力，Pa。
在一具体实施方式中，吸附气壁面覆盖度θ可表示为：
Figure BDA0003286587200000095
其中，pL为兰氏压力，Pa；Z为偏差因子，无因次；p为压力，Pa；θ为吸附气壁面覆盖度，无因次。
进一步地，在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型包括：
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输解析模型；
在无吸附气体传输解析模型基础上，确定无吸附气体传输数值模型作为无吸附气体传输模型；
更进一步地，在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输解析模型包括：
假设基质颗粒为刚性球体，基于质量守恒定理，构建球形空间控制方程；其中，所述球形空间控制方程为气体密度关于时间和运移距离的偏微分方程且球形空间控制方程包含渗透率参数；
获取内边界条件、外边条件和初始条件；
在球形空间控制方程的基础上，联合内边界条件、外边条件、初始条件以及单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型，形成以气体密度为传递变量的无吸附气体传输解析模型；
再进一步地，在无吸附气体传输解析模型基础上，确定无吸附气体传输数值模型作为无吸附气体传输模型包括：
将球形空间控制方程进行离散处理，得到离散后的控制方程；
将内边界条件、外边条件和初始条件分别进行离散处理，得到离散后内边界条件、离散后外边条件和离散后初始条件；
将孔隙内无吸附作用气体密度关于压力的计算模型进行离散处理，得到离散后孔隙内无吸附作用气体密度关于压力的计算模型；
在离散后的控制方程的基础上，联合离散后内边界条件、离散后外边条件、离散后初始条件以及单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和离散后孔隙内无吸附作用气体密度关于压力的计算模型，形成无吸附气体传输数值模型。
进一步地，在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型包括：
在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输解析模型；
在吸附气体传输解析模型基础上，确定吸附气体传输数值模型作为吸附气体传输模型；
更进一步地，在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输解析模型包括：
假设基质颗粒为刚性球体，基于质量守恒定理，构建球形空间控制方程；其中，所述球形空间控制方程为气体密度关于时间和运移距离的偏微分方程且球形空间控制方程包含渗透率参数；
获取内边界条件、外边条件和初始条件；
在球形空间控制方程的基础上，联合内边界条件、外边条件、初始条件以及单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型，形成以气体密度为传递变量的吸附气体传输解析模型；
再进一步地，在吸附气体传输解析模型基础上，确定吸附气体传输数值模型作为吸附气体传输模型包括：
将球形空间控制方程进行离散处理，得到离散后的控制方程；
将内边界条件、外边条件和初始条件分别进行离散处理，得到离散后内边界条件、离散后外边条件和离散后初始条件；
将孔隙内吸附作用气体密度关于压力的计算模型进行离散处理，得到离散后孔隙内吸附作用气体密度关于压力的计算模型；
在离散后的控制方程的基础上，联合离散后内边界条件、离散后外边条件、离散后初始条件以及单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和离散后孔隙内吸附作用气体密度关于压力的计算模型，形成吸附气体传输数值模型。
在一具体实施方式中，球形空间控制方程为：
Figure BDA0003286587200000111
式中，ρ为气体密度，kg/m3；φ为样品孔隙度，无因次；t为时间，s；r为流体运移距离，m；k为气体总视渗透率，m2；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa。
在一具体实施方式中，内边界条件为：
Figure BDA0003286587200000112
式中，p为压力，Pa；r为流体运移距离，m。
在一具体实施方式中，外边条件为：
Figure BDA0003286587200000113
式中，Vt为基质颗粒所处外界环境的体积(在一具体实施方式中指排除基质颗粒颗粒本身体积之外的样品室以及与样品室连接的中间室的总体积)，m3；St为基质颗粒总外表面积，m2；ρout为基质颗粒外部气体密度，kg/m3；R为基质颗粒半径，m；t为时间，s；r为流体运移距离，m；k为气体总视渗透率，m2；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa。
在一具体实施方式中，初始条件为
Figure BDA0003286587200000121
式中，p0为初始时刻的压力，Pa；R为基质颗粒半径，m；t为时间，s；r为流体运移距离，m；pout-ini为初始阶段页岩颗粒外部的压力，Pa。
在一具体实施方式中，将球形空间控制方程进行离散处理通过将时间差分、空间差分的方式进行；
更进一步地，离散后的控制方程为：
Figure BDA0003286587200000122
其中，
Figure BDA0003286587200000123
Figure BDA0003286587200000124
Figure BDA0003286587200000125
式中，
Figure BDA0003286587200000126
为第n+1步第i节点的气体密度，kg/m3；
Figure BDA0003286587200000127
为第n步第i+1节点的气体密度，kg/m3；
Figure BDA0003286587200000128
为第n步第i节点的气体密度，kg/m3；
Figure BDA0003286587200000129
为第n步第i-1节点的气体密度，kg/m3；φ为样品孔隙度，无因次；i为第i节点；Δt为各步间隔时间，s；Δr为各节点间隔流体运移距离，m；
Figure BDA00032865872000001210
为第n步第i+1节点的气体总视渗透率，m2；
Figure BDA00032865872000001211
为第n步第i节点的气体总视渗透率，m2；
Figure BDA00032865872000001212
为第n步第i-1节点的气体总视渗透率，m2；
Figure BDA00032865872000001213
为第n步第i+1节点的气体黏度，Pa·s；
Figure BDA00032865872000001214
为第n步第i节点的气体黏度，Pa·s；
Figure BDA00032865872000001215
为第n步第i-1节点的气体黏度，Pa·s；
Figure BDA00032865872000001216
为第n步第i+1节点的气体压缩系数，1/Pa；
Figure BDA00032865872000001217
为第n步第i节点的气体压缩系数，1/Pa；
Figure BDA00032865872000001218
为第n步第i-1节点的气体压缩系数，1/Pa。
在一具体实施方式中，离散后内边界条件为：
Figure BDA0003286587200000131
式中，
Figure BDA0003286587200000132
为第n步第1节点的压力，Pa；
Figure BDA0003286587200000133
为第n步第2节点的压力，Pa。
在一具体实施方式中，离散后外边条件为：
Figure BDA0003286587200000134
式中，Vt为基质颗粒所处外界环境的体积(在一具体实施方式中指排除基质颗粒颗粒本身体积之外的样品室以及与样品室连接的中间室的总体积)，m3；St为基质颗粒总外表面积，m2；Δt为各步间隔时间，s；Δr为各节点间隔流体运移距离，m；
Figure BDA0003286587200000135
为第n+1步第m节点的基质颗粒外部气体密度，kg/m3；
Figure BDA0003286587200000136
为第n步第m节点的基质颗粒外部气体密度，kg/m3；
Figure BDA0003286587200000137
为第n步第m节点的气体总视渗透率，m2；
Figure BDA0003286587200000138
为第n步第m节点的气体黏度，Pa·s；
Figure BDA0003286587200000139
为第n步第m节点的气体压缩系数，1/Pa；
Figure BDA00032865872000001310
为第n步第m节点的气体密度，kg/m3；
Figure BDA00032865872000001311
为第n步第m-1节点的气体密度，kg/m3。
在一具体实施方式中，离散后初始条件为：
Figure BDA00032865872000001312
式中，
Figure BDA00032865872000001313
为第1步第i+1节点的压力，Pa；
Figure BDA00032865872000001314
为第1步第m节点的压力，Pa；pout-ini为初始阶段页岩颗粒外部的压力，Pa。
在一具体实施方式中，离散后孔隙内无吸附作用气体密度关于压力的计算模型为：
Figure BDA00032865872000001315
式中，
Figure BDA00032865872000001316
为第n+1步第i节点的压力，Pa；
Figure BDA00032865872000001317
为第n+1步第i节点的气体密度，kg/m3；
Figure BDA00032865872000001318
为第n步第i节点的偏差因子，无因次；R为气体常数，通常取8.314472m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol。
在一具体实施方式中，离散后孔隙内吸附作用气体密度关于压力的计算模型为：
Figure BDA00032865872000001319
式中，
Figure BDA0003286587200000141
为第n+1步第i节点的气体密度，kg/m3；p为压力，Pa；
Figure BDA0003286587200000142
为第n+1步第i节点的压力，Pa；pL为兰氏压力，Pa；M为摩尔质量，kg/mol；
Figure BDA0003286587200000143
为第n步第i节点的偏差因子，无因次；R为气体常数，通常取8.314472m3·Pa/(K·mol)；T为温度，K；rt为孔半径，m；nL为兰氏吸附量，mol/kg；
Figure BDA0003286587200000144
为第n+1步第i节点的吸附气壁面覆盖度，无因次；ρs为目标页岩样品密度，g/mL；φ为样品孔隙度，无因次。
进一步地，目标页岩基质内无吸附气体微观传输特征包括不同时刻压力和/或渗透率分布。
进一步地，目标页岩基质内吸附气体微观传输特征包括不同时刻密度、压力和/或渗透率分布；更优选地，所述渗透率包括总渗透率、吸附相密度和/或体相密度。
其中，气体参数例如粘度、压缩系数、偏差因子等气体物性参数可通过美国国家标准和技术会(National Institute of Standards and Technology，NIST)所推荐的数学模型直接获得。氦气真实气体物性参数计算采用的模型参见文献1-文献3；甲烷的真实气体物性参数计算采用的模型参见文献4-文献6；
其中，文献1为Ortiz-Vega D,Hall K,Holste J,Arp V,Harvey A,LemmonE.Helmholtz equation of state for helium.J Phys Chem Ref Data 2020；
文献2为Hands B,Arp V.A correlation of thermal conductivity data forhelium.Cryogenics1981；21(12):697-703；
文献3为Arp VD,McCarty RD.Thermophysical Properties of Helium-4from0.8to 1500Kwith Pressures to 2000MPa.1989；
文献4为Setzmann U,Wagner W.A new equation of state and tables ofthermodynamic properties for methane covering the range from the melting lineto 625K at pressures up to 100MPa.Journal of Physical and Chemical referencedata 1991；20(6):1061-1155；
文献5为Friend DG,Ely JF,Ingham H.Thermophysical properties ofmethane.Journal of Physical and Chemical Reference Data 1989；18(2):583-638；
文献6为Vogel E.Reference viscosities of gaseous methane and hydrogensulfide at low density in the temperature range from(292to 682)K.Journal ofChemical&Engineering Data2011；56(7):3265-3272。
实施例1
本实施例提供了一种页岩基质内气体非稳态微观传输特征的确定方法，其中，该方法包括：
1、获取目标页岩基质无吸附作用气体传输参数和吸附作用气体传输参数；其中，无吸附作用气体传输参数包括孔隙度、孔径和迂曲度；所述吸附作用气体传输参数包括孔隙度、孔径、迂曲度、兰氏压力、兰氏吸附量和表面扩散系数。
获取目标页岩基质实验用样品，进行室内模拟实验获取压降曲线(压降曲线获取步骤可以参见CN112362552A中压降曲线的获取方式)，通过调整无吸附作用气体孔隙度、孔径和迂曲度等传输参数，及吸附作用气体孔隙度、孔径、迂曲度、兰氏压力、兰氏吸附量及表面扩散系数等传输参数，视理论曲线和实验数据的近似重合。此时，调整后的无吸附作用气体孔隙度、孔径和迂曲度等传输参数，及吸附作用气体孔隙度、孔径、迂曲度、兰氏压力、兰氏吸附量及表面扩散系数等传输参数即可视为拟合所得实际样品的传输参数。
实验用样品参数如表1所示，通过无吸附作用气体氦气和吸附作用气体甲烷实验数据拟合，获得相应传输参数(分别如表2和表3所示)；压降曲线如图3A-图4D所示。
表1实验用样品参数
Figure BDA0003286587200000151
表2氦气传输拟合参数
样品 φ r<sub>t</sub>/nm τ
1 0.03 1.2 40
2 0.021 1 42
3 0.038 1.5 37
表3甲烷传输实验拟合参数
Figure BDA0003286587200000152
2、获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型；获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型；
单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型为：
k＝kv·ωv+kk·ωk，其中
Figure BDA0003286587200000161
Figure BDA0003286587200000162
其中，
Figure BDA0003286587200000163
式中，k为单组分无吸附作用气体总视渗透率，m2；kv为滑脱流视渗透率，m2；kk为努森扩散视渗透率，m2；ωv为滑脱流权重系数，无因次；ωk努森扩散权重系数，无因次；kn为努森数，无因次；τ为迂曲度，无因次；bc为气体滑脱常数，无因次，取值-1；φ为样品孔隙度，无因次；rt为孔半径，m；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；δ为气体分子直径与孔径比值，无因次；Df为表面粗糙度，无因次，取值2.5；Z为偏差因子，无因次；R为气体常数，取8.314472m3·Pa/(K·mol)；T为温度，K(本实施例中为303.15K)；M为摩尔质量，kg/mol；氦气作为无吸附作用气体分子直径为0.26nm。
单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型为：
k＝kv·ωv+kk·ωk+ks
Figure BDA0003286587200000164
Figure BDA0003286587200000165
Figure BDA0003286587200000166
其中，
Figure BDA0003286587200000167
reff＝(rt-dmθ)
式中，k为单组分吸附作用气体总视渗透率，m2；kv为滑脱流视渗透率，m2；kk为努森扩散视渗透率，m2；ωv为滑脱流权重系数，无因次；ωk努森扩散权重系数，无因次；ks为表面扩散视渗透率，m2；kn为努森数，无因次；τ为迂曲度，无因次；bc为气体滑脱常数，无因次，取值-1；φ为样品孔隙度，无因次；rt为孔半径，m；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；δ为气体分子直径与孔径比值，无因次；Df为表面粗糙度，无因次，取值2.5；Z为偏差因子，无因次；R为气体常数，取8.314472m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol；Ds为吸附气表面扩散系数，m2/s；p为压力，Pa；Ccon为吸附相气体浓度，kg/m3；reff为游离态吸附气传输对应的孔隙半径，m；θ为吸附气壁面覆盖度，无因次；dm为分子直径，m；NA为阿伏伽德罗常数；甲烷作为吸附作用气体分子直径为0.3758nm。
其中，努森数表示为：
Figure BDA0003286587200000171
Figure BDA0003286587200000172
式中，λ为分子自由程，m；μc为气体在1.01325×105Pa和423K条件下的黏度，Pa·s；Z为偏差因子，无因次；R为气体常数，m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol；p为压力，Pa。
其中，吸附气壁面覆盖度θ为
Figure BDA0003286587200000173
其中，pL为兰氏压力，Pa；Z为偏差因子，无因次；p为压力，Pa；θ为吸附气壁面覆盖度，无因次。
3、获取孔隙内无吸附作用气体密度关于压力的计算模型和吸附作用气体密度关于压力的计算模型；
通过赋存状态分析，确定孔隙内无吸附作用气体和吸附作用气体密度；其中，
对于无吸附作用气体，通过气体状态方程简单表征气体密度，得到的孔隙内无吸附作用气体密度关于压力的计算模型为：
Figure BDA0003286587200000174
式中，ρ为气体密度，kg/m3；p为压力，Pa；M为摩尔质量，kg/mol；Z为偏差因子，无因次；R为气体常数，取8.314472m3·Pa/(K·mol)；T为温度，K；
对于吸附作用气体，结合吸附相和体相气体赋存特征获得气体综合密度，得到的孔隙内吸附作用气体密度关于压力的计算模型为：
Figure BDA0003286587200000181
其中，
Figure BDA0003286587200000182
式中，ρ为基质中游离态和吸附态气体的综合密度，kg/m3；p为压力，Pa；pL为兰氏压力，Pa；M为摩尔质量，kg/mol；Z为偏差因子，无因次；R为气体常数，取8.314472m3·Pa/(K·mol)；T为温度，K；rt为孔半径，m；nL为兰氏吸附量，mol/kg；θ为吸附气壁面覆盖度，无因次；ρs为目标页岩样品密度，g/mL，见表4；φ为样品孔隙度，无因次；dm为分子直径，m。
4、在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型；在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型；
基于球形页岩微粒，构建宏观气体传输模型；对于基质孔隙内吸附作用气体，不仅存在体相传输，还存在吸附相的表面扩散传输；所以，不同于以压力为传递变量的常规模型，需以包含体相和吸附相的综合密度为新的传递变量，构建新的宏观气体传输模型：
4.1、假设基质颗粒为刚性球体，基于质量守恒定理，构建球形空间控制方程；其中，所述球形空间控制方程为气体密度关于时间和运移距离的偏微分方程且球形空间控制方程包含渗透率参数；
Figure BDA0003286587200000183
式中，ρ为气体密度，kg/m3；φ为样品孔隙度，无因次；t为时间，s；r为流体运移距离，m；k为气体总视渗透率，m2；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；
4.2、获取内边界条件、外边条件和初始条件；
内边界条件为：
Figure BDA0003286587200000191
式中，p为压力，Pa；r为流体运移距离，m；
外边条件为：
Figure BDA0003286587200000192
式中，Vt为基质颗粒所处外界环境的体积，m3；St为基质颗粒总外表面积，m2；ρout为基质颗粒外部气体密度，kg/m3，基于颗粒外压力由气体状态方程获得；R为基质颗粒半径，m，取值0.3625mm；t为时间，s；r为流体运移距离，m；k为气体总视渗透率，m2；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；初始条件为
Figure BDA0003286587200000193
式中，p0为初始时刻的压力，Pa；R为基质颗粒半径，m；t为时间，s；r为流体运移距离，m；pout-ini为初始阶段页岩颗粒外部的压力，Pa；
4.3、在球形空间控制方程的基础上，联合内边界条件、外边条件、初始条件以及单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型，形成以气体密度为传递变量的无吸附气体传输解析模型；
4.4、在球形空间控制方程的基础上，联合内边界条件、外边条件、初始条件以及单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型，形成以气体密度为传递变量的吸附气体传输解析模型；
4.5、在无吸附气体传输解析模型基础上，确定无吸附气体传输数值模型作为无吸附气体传输模型；在吸附气体传输解析模型基础上，确定吸附气体传输数值模型作为吸附气体传输模型；具体包括：
4.5.1、将球形空间控制方程进行离散处理，得到离散后的控制方程；
离散后的控制方程为：
Figure BDA0003286587200000194
其中，
Figure BDA0003286587200000201
Figure BDA0003286587200000202
Figure BDA0003286587200000203
式中，
Figure BDA0003286587200000204
为第n+1步第i节点的气体密度，kg/m3；
Figure BDA0003286587200000205
为第n步第i+1节点的气体密度，kg/m3；
Figure BDA0003286587200000206
为第n步第i节点的气体密度，kg/m3；
Figure BDA0003286587200000207
为第n步第i-1节点的气体密度，kg/m3；φ为样品孔隙度，无因次；i为第i节点；Δt为各步间隔时间，s；Δr为各节点间隔流体运移距离，m；
Figure BDA0003286587200000208
为第n步第i+1节点的气体总视渗透率，m2；
Figure BDA0003286587200000209
为第n步第i节点的气体总视渗透率，m2；
Figure BDA00032865872000002010
为第n步第i-1节点的气体总视渗透率，m2；
Figure BDA00032865872000002011
为第n步第i+1节点的气体黏度，Pa·s；
Figure BDA00032865872000002012
为第n步第i节点的气体黏度，Pa·s；
Figure BDA00032865872000002013
为第n步第i-1节点的气体黏度，Pa·s；
Figure BDA00032865872000002014
为第n步第i+1节点的气体压缩系数，1/Pa；
Figure BDA00032865872000002015
为第n步第i节点的气体压缩系数，1/Pa；
Figure BDA00032865872000002016
为第n步第i-1节点的气体压缩系数，1/Pa；
4.5.2、将内边界条件、外边条件和初始条件分别进行离散处理，得到离散后内边界条件、离散后外边条件和离散后初始条件；
离散后内边界条件为：
Figure BDA00032865872000002017
式中，
Figure BDA00032865872000002018
为第n步第1节点的压力，Pa；
Figure BDA00032865872000002019
为第n步第2节点的压力，Pa；
离散后外边条件为：
Figure BDA00032865872000002020
式中，Vt为基质颗粒所处外界环境的体积，m3；St为基质颗粒总外表面积，m2；Δt为各步间隔时间，s；Δr为各节点间隔流体运移距离，m；
Figure BDA00032865872000002021
为第n+1步第m节点的基质颗粒外部气体密度，kg/m3；
Figure BDA00032865872000002022
为第n步第m节点的基质颗粒外部气体密度，kg/m3；
Figure BDA00032865872000002023
为第n步第m节点的气体总视渗透率，m2；
Figure BDA00032865872000002024
为第n步第m节点的气体黏度，Pa·s；
Figure BDA00032865872000002025
为第n步第m节点的气体压缩系数，1/Pa；
Figure BDA00032865872000002026
为第n步第m节点的气体密度，kg/m3；
Figure BDA00032865872000002027
为第n步第m-1节点的气体密度，kg/m3；
离散后初始条件为：
Figure BDA0003286587200000211
式中，
Figure BDA0003286587200000212
为第1步第i+1节点的压力，Pa；
Figure BDA0003286587200000213
为第1步第m节点的压力，Pa；pout-ini为初始阶段页岩颗粒外部的压力，Pa；
4.5.3、将孔隙内无吸附作用气体密度关于压力的计算模型进行离散处理，得到离散后孔隙内无吸附作用气体密度关于压力的计算模型；
离散后孔隙内无吸附作用气体密度关于压力的计算模型为：
Figure BDA0003286587200000214
式中，
Figure BDA0003286587200000215
为第n+1步第i节点的压力，Pa；
Figure BDA0003286587200000216
为第n+1步第i节点的气体密度，kg/m3；
Figure BDA0003286587200000217
为第n步第i节点的偏差因子，无因次；R为气体常数，取8.314472m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol。
4.5.4、将孔隙内吸附作用气体密度关于压力的计算模型进行离散处理，得到离散后孔隙内吸附作用气体密度关于压力的计算模型；
离散后孔隙内吸附作用气体密度关于压力的计算模型为：
Figure BDA0003286587200000218
式中，
Figure BDA0003286587200000219
为第n+1步第i节点的气体密度，kg/m3；p为压力，Pa；
Figure BDA00032865872000002110
为第n+1步第i节点的压力，Pa；pL为兰氏压力，Pa；M为摩尔质量，kg/mol；
Figure BDA00032865872000002111
为第n步第i节点的偏差因子，无因次；R为气体常数，取8.314472m3·Pa/(K·mol)；T为温度，K；rt为孔半径，m；nL为兰氏吸附量，mol/kg；
Figure BDA00032865872000002112
为第n+1步第i节点的吸附气壁面覆盖度，无因次；ρs为目标页岩样品密度，g/mL；φ为样品孔隙度，无因次；
4.5.5、在离散后的控制方程的基础上，联合离散后内边界条件、离散后外边条件、离散后初始条件以及单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和离散后孔隙内无吸附作用气体密度关于压力的计算模型，形成无吸附气体传输数值模型；
4.5.6、在离散后的控制方程的基础上，联合离散后内边界条件、离散后外边条件、离散后初始条件以及单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和离散后孔隙内吸附作用气体密度关于压力的计算模型，形成吸附气体传输数值模型。
5、基于无吸附气体传输模型，结合目标页岩基质无吸附作用气体传输参数，确定目标页岩基质内无吸附气体微观传输特征；基于吸附气体传输模型，结合目标页岩基质吸附作用气体传输参数，确定目标页岩基质内吸附气体微观传输特征；
基于气体传输数值模型，带入无吸附作用气体和吸附作用气体传输参数，可准确获得真实页岩基质内气体微观传输特征，为实际页岩气储层开发提供理论指导。目标页岩基质内无吸附气体微观传输特征如图5A、图5B所示，为真实页岩基质颗粒不同时刻氦气努森数、压力和渗透率分布，用以分析真实页岩基质内无吸附作用气体非稳态传输特征。目标页岩基质内吸附气体微观传输特征如图6A、图6B所示，为真实页岩基质颗粒不同时刻甲烷压力、总密度、体相密度和吸附相密度分布，用以分析真实页岩基质内吸附作用气体非稳态传输特征。
在上述实施例中，采用由美国国家标准和技术会(National Institute ofStandards and Technology，NIST)所推荐的关于气体粘度μ、气体压缩系数Cg、偏差因子Z的数学模型，获得氦气的真实气体物性参数和甲烷的真实气体物性参数；
氦气的真实气体物性参数和甲烷的真实气体物性参数利用NIST所推荐的数学模型进行获取即可。
以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种页岩基质内气体非稳态微观传输特征的确定方法，其中，该方法包括：
获取目标页岩基质无吸附作用气体传输参数；和/或；获取目标页岩基质吸附作用气体传输参数；
获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型；和/或；获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型；
获取孔隙内无吸附作用气体密度关于压力的计算模型；和/或；获取孔隙内吸附作用气体密度关于压力的计算模型；
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型；和/或；在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型；
基于无吸附气体传输模型，结合目标页岩基质无吸附作用气体传输参数，确定目标页岩基质内无吸附气体微观传输特征；和/或；基于吸附气体传输模型，结合目标页岩基质吸附作用气体传输参数，确定目标页岩基质内吸附气体微观传输特征。
2.根据权利要求1所述的方法，其中，
获取单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型包括：
获取滑脱流视渗透率计算模型；
获取努森扩散视渗透率计算模型；
获取单组分无吸附作用气体总视渗透率计算模型；
获取单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型包括：
获取滑脱流视渗透率计算模型；
获取努森扩散视渗透率计算模型；
获取表面扩散视渗透率计算模型；
获取单组分吸附作用气体总视渗透率计算模型。
3.根据权利要求2所述的方法，其中，
滑脱流视渗透率计算模型为滑脱流视渗透率关于孔隙度、孔径、迂曲度、努森数以及气体滑脱常数的计算模型；优选地，滑脱流视渗透率计算模型为：
Figure FDA0003286587190000021
式中，kv为滑脱流视渗透率，m2；kn为努森数，无因次；τ为迂曲度，无因次；bc为气体滑脱常数，无因次；φ为样品孔隙度，无因次；rt为孔半径，m；
努森扩散视渗透率计算模型为努森扩散视渗透率关于气体参数、孔隙度、孔径和迂曲度的计算模型；优选地，气体参数包括单气体黏度、气体压缩系数、气体分子直径、偏差因子和摩尔质量；更优选地，努森扩散视渗透率计算模型为：
Figure FDA0003286587190000022
式中，kk为努森扩散视渗透率，m2；τ为迂曲度，无因次；bc为气体滑脱常数，无因次；φ为样品孔隙度，无因次；rt为孔半径，m；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；δ为气体分子直径与孔径比值，无因次；Df为表面粗糙度，无因次；Z为偏差因子，无因次；R为气体常数，m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol；
单组分无吸附作用气体总视渗透率计算模型为总视渗透率关于努森扩散视渗透率和滑脱流视渗透率的计算模型；优选地，单组分无吸附作用气体总视渗透率计算模型为：
k＝kv·ωv+kk·ωk，其中，
Figure FDA0003286587190000023
式中，k为单组分无吸附作用气体总视渗透率，m2；kv为滑脱流视渗透率，m2；kk为努森扩散视渗透率，m2；ωv为滑脱流权重系数，无因次；ωk努森扩散权重系数，无因次；kn为努森数，无因次；
表面扩散视渗透率计算模型为表面扩散视渗透率关于气体参数、孔隙度、孔径和迂曲度的计算模型；优选地，气体参数包括单气体黏度、吸附相气体浓度、游离态吸附气传输对应的孔隙半径、偏差因子和摩尔质量；更优选地，表面扩散视渗透率计算模型为：
Figure FDA0003286587190000024
式中，ks为表面扩散视渗透率，m2；Ds为吸附气表面扩散系数，m2/s；τ为迂曲度，无因次；p为压力，Pa；M为摩尔质量，kg/mol；φ为样品孔隙度，无因次；μ为气体黏度，Pa·s；Ccon为吸附相气体浓度，kg/m3；Z为偏差因子，无因次；R为气体常数，m3·Pa/(K·mol)；T为温度，K；rt为孔半径，m；reff为游离态吸附气传输对应的孔隙半径，m；
单组分吸附作用气体总视渗透率计算模型为总视渗透率关于努森扩散视渗透率、滑脱流视渗透率和表面扩散视渗透率的计算模型；优选地，单组分吸附作用气体总视渗透率计算模型为：
k＝kv·ωv+kk·ωk+ks其中，
Figure FDA0003286587190000031
式中，k为单组分吸附作用气体总视渗透率，m2；kv为滑脱流视渗透率，m2；kk为努森扩散视渗透率，m2；ωv为滑脱流权重系数，无因次；ωk努森扩散权重系数，无因次；ks为表面扩散视渗透率，m2；kn为努森数，无因次。
4.根据权利要求1所述的方法，其中，
孔隙内无吸附作用气体密度关于压力的计算模型基于气体状态方程确定；优选地，孔隙内无吸附作用气体密度关于压力的计算模型为：
Figure FDA0003286587190000032
式中，ρ为气体密度，kg/m3；p为压力，Pa；M为摩尔质量，kg/mol；Z为偏差因子，无因次；R为气体常数，m3·Pa/(K·mol)；T为温度，K；
孔隙内吸附作用气体密度关于压力的计算模型结合孔隙内吸附相和体相气体赋存特征综合确定；优选地，孔隙内吸附作用气体密度关于压力的计算模型为：
Figure FDA0003286587190000033
其中，
Figure FDA0003286587190000034
式中，ρ为基质中游离态和吸附态气体的综合密度，kg/m3；p为压力，Pa；pL为兰氏压力，Pa；M为摩尔质量，kg/mol；Z为偏差因子，无因次；R为气体常数，Pa/(K·mol)；T为温度，K；rt为孔半径，m；nL为兰氏吸附量，mol/kg；θ为吸附气壁面覆盖度，无因次；ρs为目标页岩样品密度，g/mL；φ为样品孔隙度，无因次；dm为分子直径，m。
5.根据权利要求1所述的方法，其中，
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输模型包括：
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输解析模型；
在无吸附气体传输解析模型基础上，确定无吸附气体传输数值模型作为无吸附气体传输模型；
在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输模型包括：
在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输解析模型；
在吸附气体传输解析模型基础上，确定吸附气体传输数值模型作为吸附气体传输模型。
6.根据权利要求5所述的方法，其中，
在单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建无吸附气体传输解析模型包括：
假设基质颗粒为刚性球体，基于质量守恒定理，构建球形空间控制方程；其中，所述球形空间控制方程为气体密度关于时间和运移距离的偏微分方程且球形空间控制方程包含渗透率参数；
获取内边界条件、外边条件和初始条件；
在球形空间控制方程的基础上，联合内边界条件、外边条件、初始条件以及单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和孔隙内无吸附作用气体密度关于压力的计算模型，形成以气体密度为传递变量的无吸附气体传输解析模型；
在单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型的基础上，以气体密度为传递变量构建吸附气体传输解析模型包括：
假设基质颗粒为刚性球体，基于质量守恒定理，构建球形空间控制方程；其中，所述球形空间控制方程为气体密度关于时间和运移距离的偏微分方程且球形空间控制方程包含渗透率参数；
获取内边界条件、外边条件和初始条件；
在球形空间控制方程的基础上，联合内边界条件、外边条件、初始条件以及单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和孔隙内吸附作用气体密度关于压力的计算模型，形成以气体密度为传递变量的吸附气体传输解析模型。
7.根据权利要求6所述的方法，其中，
在无吸附气体传输解析模型基础上，确定无吸附气体传输数值模型作为无吸附气体传输模型包括：
将球形空间控制方程进行离散处理，得到离散后的控制方程；
将内边界条件、外边条件和初始条件分别进行离散处理，得到离散后内边界条件、离散后外边条件和离散后初始条件；
将孔隙内无吸附作用气体密度关于压力的计算模型进行离散处理，得到离散后孔隙内无吸附作用气体密度关于压力的计算模型；
在离散后的控制方程的基础上，联合离散后内边界条件、离散后外边条件、离散后初始条件以及单组分无吸附作用气体滑脱传输、努森扩散传输视渗透率计算模型和离散后孔隙内无吸附作用气体密度关于压力的计算模型，形成无吸附气体传输数值模型；
在吸附气体传输解析模型基础上，确定吸附气体传输数值模型作为吸附气体传输模型包括：
将球形空间控制方程进行离散处理，得到离散后的控制方程；
将内边界条件、外边条件和初始条件分别进行离散处理，得到离散后内边界条件、离散后外边条件和离散后初始条件；
将孔隙内吸附作用气体密度关于压力的计算模型进行离散处理，得到离散后孔隙内吸附作用气体密度关于压力的计算模型；
在离散后的控制方程的基础上，联合离散后内边界条件、离散后外边条件、离散后初始条件以及单组分吸附作用气体滑脱传输、努森扩散传输、表面扩散传输视渗透率计算模型和离散后孔隙内吸附作用气体密度关于压力的计算模型，形成吸附气体传输数值模型。
8.根据权利要求6或7所述的方法，其中，
球形空间控制方程为：
Figure FDA0003286587190000061
式中，ρ为气体密度，kg/m3；φ为样品孔隙度，无因次；t为时间，s；r为流体运移距离，m；k为气体总视渗透率，m2；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；
内边界条件为：
Figure FDA0003286587190000062
式中，p为压力，Pa；r为流体运移距离，m；
外边条件为：
Figure FDA0003286587190000063
式中，Vt为基质颗粒所处外界环境的体积，m3；St为基质颗粒总外表面积，m2；ρout为基质颗粒外部气体密度，kg/m3；R为基质颗粒半径，m；t为时间，s；r为流体运移距离，m；k为气体总视渗透率，m2；μ为气体黏度，Pa·s；Cg为气体压缩系数，1/Pa；
初始条件为
Figure FDA0003286587190000064
式中，p0为初始时刻的压力，Pa；R为基质颗粒半径，m；t为时间，s；r为流体运移距离，m；pout-ini为初始阶段页岩颗粒外部的压力，Pa。
9.根据权利要求7所述的方法，其中，
离散后的控制方程为：
Figure FDA0003286587190000065
其中，
Figure FDA0003286587190000066
Figure FDA0003286587190000071
Figure FDA0003286587190000072
式中，
Figure FDA0003286587190000073
为第n+1步第i节点的气体密度，kg/m3；
Figure FDA0003286587190000074
为第n步第i+1节点的气体密度，kg/m3；
Figure FDA0003286587190000075
为第n步第i节点的气体密度，kg/m3；
Figure FDA0003286587190000076
为第n步第i-1节点的气体密度，kg/m3；φ为样品孔隙度，无因次；i为第i节点；Δt为各步间隔时间，s；Δr为各节点间隔流体运移距离，m；
Figure FDA0003286587190000077
为第n步第i+1节点的气体总视渗透率，m2；
Figure FDA0003286587190000078
为第n步第i节点的气体总视渗透率，m2；
Figure FDA0003286587190000079
为第n步第i-1节点的气体总视渗透率，m2；
Figure FDA00032865871900000710
为第n步第i+1节点的气体黏度，Pa·s；
Figure FDA00032865871900000711
为第n步第i节点的气体黏度，Pa·s；
Figure FDA00032865871900000712
为第n步第i-1节点的气体黏度，Pa·s；
Figure FDA00032865871900000713
为第n步第i+1节点的气体压缩系数，1/Pa；
Figure FDA00032865871900000714
为第n步第i节点的气体压缩系数，1/Pa；
Figure FDA00032865871900000715
为第n步第i-1节点的气体压缩系数，1/Pa；
离散后内边界条件为：
Figure FDA00032865871900000716
式中，
Figure FDA00032865871900000717
为第n步第1节点的压力，Pa；
Figure FDA00032865871900000718
为第n步第2节点的压力，Pa；
离散后外边条件为：
Figure FDA00032865871900000719
式中，Vt为基质颗粒所处外界环境的体积(在一具体实施方式中指排除基质颗粒颗粒本身体积之外的样品室以及与样品室连接的中间室的总体积)，m3；St为基质颗粒总外表面积，m2；Δt为各步间隔时间，s；Δr为各节点间隔流体运移距离，m；
Figure FDA00032865871900000720
为第n+1步第m节点的基质颗粒外部气体密度，kg/m3；
Figure FDA00032865871900000721
为第n步第m节点的基质颗粒外部气体密度，kg/m3；
Figure FDA00032865871900000722
为第n步第m节点的气体总视渗透率，m2；
Figure FDA00032865871900000723
为第n步第m节点的气体黏度，Pa·s；
Figure FDA00032865871900000724
为第n步第m节点的气体压缩系数，1/Pa；
Figure FDA00032865871900000725
为第n步第m节点的气体密度，kg/m3；
Figure FDA00032865871900000726
为第n步第m-1节点的气体密度，kg/m3；
离散后初始条件为：
Figure FDA00032865871900000727
式中，
Figure FDA00032865871900000728
为第1步第i+1节点的压力，Pa；
Figure FDA00032865871900000729
为第1步第m节点的压力，Pa；pout-ini为初始阶段页岩颗粒外部的压力，Pa；
离散后孔隙内无吸附作用气体密度关于压力的计算模型为：
Figure FDA0003286587190000081
式中，
Figure FDA0003286587190000082
为第n+1步第i节点的压力，Pa；
Figure FDA0003286587190000083
为第n+1步第i节点的气体密度，kg/m3；
Figure FDA0003286587190000084
为第n步第i节点的偏差因子，无因次；R为气体常数，m3·Pa/(K·mol)；T为温度，K；M为摩尔质量，kg/mol；
离散后孔隙内吸附作用气体密度关于压力的计算模型为：
Figure FDA0003286587190000085
式中，
Figure FDA0003286587190000086
为第n+1步第i节点的气体密度，kg/m3；p为压力，Pa；
Figure FDA0003286587190000087
为第n+1步第i节点的压力，Pa；pL为兰氏压力，Pa；M为摩尔质量，kg/mol；
Figure FDA0003286587190000088
为第n步第i节点的偏差因子，无因次；R为气体常数，m3·Pa/(K·mol)；T为温度，K；rt为孔半径，m；nL为兰氏吸附量，mol/kg；
Figure FDA0003286587190000089
为第n+1步第i节点的吸附气壁面覆盖度，无因次；ρs为目标页岩样品密度，g/mL；φ为样品孔隙度，无因次。
10.根据权利要求1所述的方法，其中，
无吸附作用气体传输参数包括孔隙度、孔径和迂曲度；
吸附作用气体传输参数包括孔隙度、孔径、迂曲度、兰氏压力、兰氏吸附量和表面扩散系数。
CN202111149220.7A 2021-09-29 2021-09-29 一种页岩基质内气体非稳态微观传输特征的确定方法 Active CN113868978B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202111149220.7A CN113868978B (zh) 2021-09-29 2021-09-29 一种页岩基质内气体非稳态微观传输特征的确定方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202111149220.7A CN113868978B (zh) 2021-09-29 2021-09-29 一种页岩基质内气体非稳态微观传输特征的确定方法

Publications (2)

Publication Number Publication Date
CN113868978A true CN113868978A (zh) 2021-12-31
CN113868978B CN113868978B (zh) 2022-09-02

Family

ID=78992478

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111149220.7A Active CN113868978B (zh) 2021-09-29 2021-09-29 一种页岩基质内气体非稳态微观传输特征的确定方法

Country Status (1)

Country Link
CN (1) CN113868978B (zh)

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107461192A (zh) * 2017-06-01 2017-12-12 西南石油大学 储层条件下页岩动态表观渗透率的计算方法
CN107480316A (zh) * 2017-06-01 2017-12-15 西南石油大学 页岩纳米级孔隙中气体传输流量的计算方法
CN108343433A (zh) * 2018-02-28 2018-07-31 西南石油大学 用于计算页岩微裂缝缝宽变化下气体质量传输参数的方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107461192A (zh) * 2017-06-01 2017-12-12 西南石油大学 储层条件下页岩动态表观渗透率的计算方法
CN107480316A (zh) * 2017-06-01 2017-12-15 西南石油大学 页岩纳米级孔隙中气体传输流量的计算方法
CN108343433A (zh) * 2018-02-28 2018-07-31 西南石油大学 用于计算页岩微裂缝缝宽变化下气体质量传输参数的方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
吴克柳等: "页岩纳米孔吸附气表面扩散机理和数学模型", 《中国科学:技术科学》 *
成景烨等: "页岩气储层基质纳米孔隙渗流模型研究进展", 《常州大学学报(自然科学版)》 *

Also Published As

Publication number Publication date
CN113868978B (zh) 2022-09-02

Similar Documents

Publication Publication Date Title
Cai et al. Shale gas transport model in 3D fractal porous media with variable pore sizes
Yu et al. Multiscale gas transport behavior in heterogeneous shale matrix consisting of organic and inorganic nanopores
Yuan et al. An analytical model of apparent gas permeability for tight porous media
Sun et al. Transport capacity of gas confined in nanoporous ultra-tight gas reservoirs with real gas effect and water storage mechanisms coupling
Wang et al. Apparent permeability prediction of organic shale with generalized lattice Boltzmann model considering surface diffusion effect
CN107461192B (zh) 储层条件下页岩动态表观渗透率的计算方法
Kazemi et al. An analytical model for shale gas permeability
CN105910971B (zh) 富有机质致密岩心气体渗透率和扩散系数的联测方法
Firouzi et al. Klinkenberg effect on predicting and measuring helium permeability in gas shales
CN108959789B (zh) 一种考虑应力敏感效应的页岩表观渗透率计算方法
CN108229048B (zh) 一种考虑页岩基质差异的多场耦合渗透率计算方法
CN108710723B (zh) 一种页岩储层多孔介质表观渗透率的计算方法
Alnoaimi et al. Characterization and measurement of multiscale gas transport in shale-core samples
CN106932323B (zh) 一种页岩气藏气体有效孔隙度反演方法
CN108801879B (zh) 一种页岩基质颗粒孔隙率和渗透率一体化测量系统及方法
CN106525684A (zh) 一种基于孔喉结构的致密砂岩克氏渗透率的校正方法
Li et al. Improved mathematical model of apparent permeability: A focused study on free and multilayer adsorptive phase flow
Pang et al. An innovative method to characterize sorption-induced kerogen swelling in organic-rich shales
ZUO et al. Limitations of Lattice Boltzmann Modeling of Micro‐Flows in Complex Nanopores
CN113825996A (zh) 用于在脉冲衰减实验中确定岩心渗透率的方法和系统
CN113128040A (zh) 一种基于孔隙结构特征的页岩表观渗透率预测方法
CN113868978B (zh) 一种页岩基质内气体非稳态微观传输特征的确定方法
Civan et al. Theoretical fundamentals, critical issues, and adequate formulation of effective shale gas and condensate reservoir simulation
JP2008157826A (ja) 透過細孔を識別する方法
CN113624847B (zh) 一种页岩水化损伤系数的预测模型建立方法及预测方法

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