一种考虑动态毛管力的数值模拟方法与流程

本发明涉及油藏数值模拟技术应用领域，尤其涉及一种考虑动态毛管力的数值模拟方法。

背景技术：

毛管力是毛细管压力的简称，是指毛细管中弯液面两侧两种流体(非湿相流体与湿相流体)的压力差。油气藏开发过程中，储层中的油、水、气是始终流动的，由于流体流动必然会造成动润湿滞后(流体流动速度大于三相周界移动速度引起的润湿角改变)，从而影响流体流动过程中毛管力的数值。渗流力学和油藏数值模拟过程中涉及到的毛管力，应当是与流体流动速度相关的动态毛管力。但是，目前毛管力的测试实验所应用的岩心夹持器都是一端进入高压流体，另一端封闭，测试过程中流体的流动方式与油藏条件下的流体流动方式截然不同，整个实验过程的压力变化体现的是孔喉半径的差异，其实质是“静态毛管力”，并非油藏条件下流体流动的毛管力。

在水动力学领域，动态毛管力的研究成果较多，但都是测试的气水两相在填砂模型中的动态毛管力，从camps-roach测试得到的气水两相的实验结果可知静态毛管力与动态毛管力的存在明显差异。kalaydjian利用经过改装的测试装置，测试得到了油水两相的动态毛管力的实验结果，从测试结果分析，静态毛管力与动态毛管力之间依然存在较大差异。很多学者利用物质平衡方法计算得到了动态毛管力数值，并对动态毛管力的影响因素开展了分析，认为动态毛管力不仅受控于孔隙结构、渗透率、流体性质等静态因素，而且受控于流体流动速度等动态因素。

hassanizadehandgray提出动态毛管力与静态毛管力之间的差异可用如下的数学方程进行描述：

式(1)中，pdyn是动态毛管力；pstat是静态毛管力(处于平衡状态的毛管力，即τ是动态系数，用于描述饱和度变化对毛管力的影响程度。已有大量学者研究了动态毛管力在不同条件下的变化规律和取值方法。hassanizadeh利用室内实验结果，提出τ的取值范围是5×10⁴～6×10⁴pa·s。stauffer对大量室内实验的结果进行总结，提出了τ的具体计算公式：

式中，α是无因次系数，通常取0.1；φ是孔隙度，单位为f；k是渗透率，单位为md；μ是湿相流体粘度，单位为mpa·s；g是重力加速度，单位为n/kg；ρ是湿相流体密度，单位为g/m³；p^d是孔隙进入压力，单位为mpa；λ是粒径分布系数，单位为f。

然而，目前应用较为广泛的油藏数值模拟方法在处理毛管力时，都认为毛管力仅仅是含水饱和度的函数，数值模拟所用的毛管力在取值时，参考的是压汞法、离心法或者半渗透隔板法所测试得到的毛管力，考虑到润湿滞后效应，将毛管力分为吸吮过程的毛管力和驱替过程的毛管力。基本渗流方程建立过程和基本渗流方程差分求解过程中，没有考虑流体流动速度和饱和度变化速度对毛管力的影响。sandoval-torres等利用有限元解法，将动态毛管力引入到计算木头真空干燥的数值模拟方程中，将数值模拟的计算结果与实验结果进行了对比，认为建立的数值模拟模型准确可靠。shu等同样利用有限元解法，将动态毛管力引入到计算油水两相渗流的数值模拟方程中，计算得到了一维条件下考虑动态毛管力的数值模拟结果，分析了一维条件下动态毛管力对油水运移规律的影响。但是，油藏数值模拟器更多的选用有限差分方法，且油藏数值模拟针对的是宏观的三维空间的油水运移规律。上述两种模型并不能很好的用于油藏条件下的数值模拟计算。tian等将动态毛管力引入到油藏数值模拟基本渗流方程中，并用有限差分方法进行了求解，但是求解方法选用的是impes(显式饱和度，隐式压力)方法，计算的稳定性和准确性不如全隐式求解，并且，计算结果仅与室内实验结果进行了对比，并未与实际生产数据进行对比。

现有技术中在油藏数值模拟中不考虑动态毛管力的方式影响模拟结果的精确性。

技术实现要素：

为了解决上述技术问题，本发明一种考虑动态毛管力的数值模拟方法。

本发明提供的考虑动态毛管力的数值模拟方法，包括：

步骤1，对岩心设置多个测点，通过动态毛管力测试方法获得各测点下不同驱替速度下的动态毛管力曲线；

步骤2，利用渗流力学方法推导得到动态毛管力与含水饱和度、流速之间的函数关系，使用此函数关系对实验测试得到的动态毛管力曲线进行拟合，得到函数关系中的待定系数；

步骤3，根据所述函数关系对基本油相方程和基本水相方程进行差分，得到油相方程的差分离散方程，并且得到水相方程的与动态毛管力相关的差分离散方程，对所述油相方程的差分离散方程和所述水相方程的与动态毛管力相关的差分离散方程进行全隐式迭代求解，得到线性化的系数矩阵。

上述考虑动态毛管力的数值模拟方法还具有以下特点：

所述步骤1中动态毛管力曲线的横坐标为含水饱和度sw，纵坐标为动态毛管力数值vw表示水相驱替流速；

所述步骤2中，所述函数关系为：

所述函数关系中的待定系数为k1和k2；

其中，

指静态毛管力，αs是实验常数，μw是润湿相相粘度，是孔隙度，ct是压缩系数，k是绝对渗透率，λ和pe是brook–corey模型中的毛管压力因子，ρw是水密度，g是重力加速度。

上述考虑动态毛管力的数值模拟方法还具有以下特点：

所述步骤3中水相方程的与动态毛管力相关的差分离散方程为:

是二阶差商算子，是第n+1个时间步的传导率，是第n+1个时间步的势，是指第n+1时刻的水相注入速度，vijk是指坐标为(i，j，k)网格的孔隙体积，δtⁿ⁺¹第n+1个时间步的时间间隔，φ是孔隙度，ρw是指润湿相密度。

上述考虑动态毛管力的数值模拟方法还具有以下特点：

步骤3中对所述油相方程的差分离散方程和所述水相方程的与动态毛管力相关的差分离散方程进行全隐式迭代求解包括：

设置参数势为

是指第n时刻的非润湿相压力，是指第n个时间步的毛管力，x是指坐标x的方向，p是指压力，t是指时间，γw是指水的重度，d是指垂向高度；

根据所述参数势对所述油相方程的差分离散方程和所述水相方程的与动态毛管力相关的差分离散方程进行全隐式迭代求解，得到线性化的系数矩阵。

上述考虑动态毛管力的数值模拟方法还具有以下特点：

所述动态毛管力测试方法包括：

步骤一，在岩心不同测点上设置油相压力传感器和水相压力传感器，在岩心的油输入口连接油容器，在岩心的水输入口连接水容器，在实验用油容器和实验用水容器的输入端连接驱替泵；

步骤二，确定油水比例，确定驱替泵的驱替速度；

步骤三，分别在实验用油容器和实验用水容器放置符合所述油水比例的实验用油和实验用水；

步骤四，将所述油容器的实验用油和所述水容器的实验用水驱替入所述岩心，在驱替稳定后，通过所述油相压力传感器和所述水相压力传感器实时检测并记录油相压力和水相压力；

步骤五，改变油水比例和驱替速度，重复执行步骤二和步骤三，获得不同驱替速度条件下的油相压力、水相压力和饱和度变化，从而获得各测点下不同驱替速度下的动态毛管力曲线离散参数。

本发明解决了现有技术中无法获得线性化的系数矩阵的问题。本发明将动态毛管力引入到油藏数值模拟基本渗流方程中，并选用全隐式的方法对饱和度和压力进行了差分求解，提高了模拟结果的精确性。申请人将计算结果分别与室内实验结果和实际油田生产数据进行了对比，验证了本方法的可靠性和准确性。

附图说明

图1是本发明中考虑动态毛管力的数值模拟方法的流程图。

具体实施例

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。需要说明的是，在不冲突的情况下，本申请中的实施例及实施例中的特征可以相互任意组合。

图1是本发明中考虑动态毛管力的数值模拟方法的流程图。此考虑动态毛管力的数值模拟方法包括：

步骤1，对岩心设置多个测点，通过动态毛管力测试方法获得各测点下不同驱替速度下的动态毛管力曲线；

步骤2，利用渗流力学方法推导得到动态毛管力与含水饱和度、流速之间的函数关系，使用此函数关系对实验测试得到的动态毛管力曲线进行拟合，得到函数关系中的待定系数；

步骤3，根据所述函数关系对基本油相方程和基本水相方程进行差分，得到油相方程的差分离散方程，并且得到水相方程的与动态毛管力相关的差分离散方程，对所述油相方程的差分离散方程和所述水相方程的与动态毛管力相关的差分离散方程进行全隐式迭代求解，得到线性化的系数矩阵。

其中，

步骤1中动态毛管力曲线的横坐标为含水饱和度sw，纵坐标为动态毛管力数值vw表示水相驱替流速。

步骤1中，所述动态毛管力测试方法包括：

步骤一，在岩心不同测点上设置油相压力传感器和水相压力传感器，在岩心的油输入口连接油容器，在岩心的水输入口连接水容器，在实验用油容器和实验用水容器的输入端连接驱替泵；

步骤二，确定油水比例，确定驱替泵的驱替速度；

步骤三，分别在实验用油容器和实验用水容器放置符合所述油水比例的实验用油和实验用水；

步骤四，将所述油容器的实验用油和所述水容器的实验用水驱替入所述岩心，在驱替稳定后，通过所述油相压力传感器和所述水相压力传感器实时检测并记录油相压力和水相压力；

步骤五，改变油水比例和驱替速度，重复执行步骤二和步骤三，获得不同驱替速度条件下的油相压力、水相压力和饱和度变化，从而获得各测点下不同驱替速度下的动态毛管力曲线离散参数。

从不同驱替速度的动态毛管力曲线可以看出，驱替速度对动态毛管力曲线的影响较为敏感，驱替速度越大，动态毛管力曲线数值越大，在束缚水饱和度附近和残余油饱和度附近，不同驱替速度下的动态毛管力差别较小，在两相共渗区的中间段，不同驱替速度下的动态毛管力差别较大。

步骤2中此函数关系为：

所述函数关系中的待定系数为k1和k2；

其中，

指静态毛管力，αs是实验常数，μw是润湿相相粘度，是孔隙度，ct是压缩系数，k是绝对渗透率，λ和pe是brook–corey模型中的毛管压力因子，ρw是水密度，g是重力加速度。

步骤3中水相方程的与动态毛管力相关的差分离散方程为

是二阶差商算子，是第n+1个时间步的传导率，是第n+1个时间步的势，是指第n+1时刻的水相注入速度，vijk是指坐标为(i，j，k)网格的孔隙体积，δtⁿ⁺¹第n+1个时间步的时间间隔，φ是孔隙度，ρw是指润湿相密度。

其中，

步骤3中对所述油相方程的差分离散方程和所述水相方程的与动态毛管力相关的差分离散方程进行全隐式迭代求解包括：

设置参数势为

是指第n时刻的非润湿相压力，是指第n个时间步的毛管力，x是指坐标x的方向，p是指压力，t是指时间，γw是指水的重度，d是指垂向高度。

根据所述参数势对所述油相方程的差分离散方程和所述水相方程的与动态毛管力相关的差分离散方程进行全隐式迭代求解，得到线性化的系数矩阵。

本发明解决了现有技术中无法获得线性化的系数矩阵的问题。本发明将动态毛管力引入到油藏数值模拟基本渗流方程中，并选用全隐式的方法对饱和度和压力进行了差分求解，提高了模拟结果的精确性。申请人将计算结果分别与室内实验结果和实际油田生产数据进行了对比，验证了本方法的可靠性和准确性。

上面描述的内容可以单独地或者以各种方式组合起来实施，而这些变型方式都在本发明的保护范围之内。

需要说明的是，在本文中，术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含，从而使得包括一系列要素的物品或者设备不仅包括那些要素，而且还包括没有明确列出的其他要素，或者是还包括为这种物品或者设备所固有的要素。在没有更多限制的情况下，由语句“包括……”限定的要素，并不排除在包括所述要素的物品或者设备中还存在另外的相同要素。

以上实施例仅用以说明本发明的技术方案而非限制，仅仅参照较佳实施例对本发明进行了详细说明。本领域的普通技术人员应当理解，可以对本发明的技术方案进行修改或者等同替换，而不脱离本发明技术方案的精神和范围，均应涵盖在本发明的权利要求范围当中。