非线性数值模拟进展评价

1非线性数值模拟进展评价通过对程林松、姜瑞忠、黄延章、何顺利、姚军等人的相关文章调研发现，目前对于低速非达西模型的描述主要为启动压力梯度和毛管力动态效应两种表述方式。1.单相启动压力梯度大量的实验表明，关于低渗透性油层中渗流特征和渗流规律的表述可归纳为：具有启动压力梯度的非线性律，如图1所示。其基本渗流特征为：(1)当压力梯度小于某个值时，流体不流动。这个压力梯度称为最小启动压力梯度(真实启动压力梯度)，以𝜆𝑎表示.(2)当压力梯度在比较低的范围时，渗流速度的增加呈下凹型曲线。(3)当压力梯度比较高时，渗流速度呈直线性增加。(4)该直线段的延伸与压力梯度轴交于某点而不通过坐标原点.此点称为平均启动压力梯度(拟启动压力梯度)，以𝜆𝑏表示,𝜆𝑐为最大启动压力梯度。图1渗流速度与压力梯度的关系目前常用的单相启动压力梯度描述模型主要分为拟启动压力梯度模型、分段模型和连续模型。1.1拟启动压力梯度模型一些研究人员基于拟启动压力梯度模型，针对拟启动压力梯度的求解方法进行了实验和理论研究，实验流体从单一流体到多种流体，考虑的影响因素也从单2一的渗透率变为渗透率、粘度等多种因素的综合作用。李忠兴，程林松[1-3]等人（2004）使用抽真空煤油进行室内驱替实验，依据考虑拟启动压力梯度的达西定律：Q=𝐾𝜇A(∆𝑃∆𝐿−𝜆)对两边取常用对数得到：lg(∆𝑃∆𝐿−𝜆)=lgQ−lg(𝐾𝜇A)该式中流量为零时驱替压力梯度即为启动压力梯度，通过回归实验数据，得出幂律关系式：𝜆=𝛼𝐾−𝑛最后研究人员绘制了确定低渗透油田启动压力梯度的理论图版，如图2所示。该求解方法仅考虑了渗透率和启动压力梯度的关系。郝斐，程林松[4]等人（2006）使用模拟油、地层水、注入水和蒸馏水进行物理模拟，实验结果表明：双对数坐标中启动压力梯度与渗透率呈线性关系且斜率为-1(图3所示)，这与理论研究相符合：lgG=−lgk+lg(𝑏×𝜇)验证了压力梯度与渗透率的幂函数关系，并且研究了不同粘度对启动压力梯度的影响，结果表明粘度越大，启动压力梯度越大。许建红，程林松[5-6]等人（2007）使用抽真空煤油进行驱替实验，拟合得到驱动压力梯度与渗流速度的关系式：v=a(∇𝑝)2+𝑏(∇𝑝)+𝑐(|∇𝑝|≥𝐺)再与考虑拟启动压力梯度的达西定律进行联立求解，可以得到下式：𝐺𝑖=−[𝑎(𝛻𝑃)2+(𝑏−𝐾𝑖𝜇𝑖)(𝛻𝑃)+𝑐]𝐾𝑖𝜇𝑖⁄该求解方法不仅考虑了渗透率与启动压力梯度的关系，还考虑了流体粘度、驱动压力梯度对启动压力梯度的影响。3图2确定低渗透地层启动压力梯度大小理论图版图3不同流体启动压力梯度与渗透率的关系也有一些研究人员讨论了考虑拟启动压力梯度的数学模型，并进行了数值模拟或者寻找简化方法求得解析解。韩洪宝，程林松[1-3,6-9]等人（2004）引入了边界层理论[5]，将油相粘度用表观粘度来代替：𝜇𝐵=𝐴𝜇1+(1−𝐴)𝜇2同时考虑了拟启动压力梯度（如图4点C），将结果应用到黑油模型中，得到了修正后的三维三相黑油模型：{油方程：∇∙⌈𝑘𝐾𝑟𝑜𝐵𝑜𝜇𝐵(∇𝑝𝑜−𝜌𝑜𝑔∇D−𝐺𝑜)⌉+𝑞𝑜𝜌𝑜𝑠𝑐=𝜕𝜕𝑡(𝜑𝑆𝑜𝐵𝑜)水方程：∇∙⌈𝑘𝐾𝑟𝑤𝐵𝑤𝜇𝑤(∇𝑝𝑤−𝜌𝑤𝑔∇D)⌉+𝑞𝑤𝜌𝑤𝑠𝑐=𝜕𝜕𝑡(𝜑𝑆𝑤𝐵𝑤)气方程：∇∙⌈𝑘𝐾𝑟𝑜𝑅𝑠𝑜𝐵𝑜𝜇𝐵(∇𝑝𝑜−𝜌𝑜𝑔∇D−𝐺𝑜)+𝑘𝐾𝑟𝑤𝑅𝑠𝑤𝐵𝑤𝜇𝑤(∇𝑝𝑤−𝜌𝑤𝑔∇D)+𝑘𝐾𝑟𝑔𝐵𝑔𝜇𝑔(∇𝑝𝑔−𝜌𝑔𝑔∇D)⌉+𝑞𝑔𝜌𝑔𝑠𝑐=𝜕𝜕𝑡[𝜑(𝑆𝑜𝑅𝑠𝑜𝐵𝑜+𝑆𝑤𝑅𝑠𝑤𝐵𝑤+𝑆𝑔𝐵𝑔)]4图4典型的非达西渗流关系曲线李莉，董平川等人（2006）基于拟启动压力梯度模型，建立了三维油水两相非达西渗流数学模型，利用差分法得到数值模型，并设计了正方形注采单元进行数值模拟，通过改变注采井距，的饿到了见水时地层压力、压力梯度等参数的空间分布，最后得到了有效驱动注采井距。一些研究人员也以拟启动压力梯度为基础建立了低渗透油藏非达西渗流模型，使用不同的方法探索模型的精确解。郝斐，程林松[10-12]等人（2006）建立了考虑拟启动压力梯度低渗透油藏不稳定渗流模型，并在内边界条件中考虑了动边界：{(𝜕𝑝𝜕𝑟−𝐺)|𝑟=𝑅(𝑡)=0𝑝=𝑝𝑒𝑟≥𝑅(𝑡)利用油藏压力近似表达式：p=𝑎0𝑙𝑛𝑟𝑅(𝑡)+𝑎1+𝑎2𝑟𝑅(𝑡),𝑟𝑤≤𝑟≤𝑅(𝑡)求得近似解和半解析解，分析了启动压力梯度对地层能量传播和流体运移的影响。两种求解结果有较好的一致性，解决了数值解法计算量大的弊端。李松泉，程林松[13]等人（2008）在前人研究的基础上，建立了考虑拟启动压力梯度和介质变形的稳定、不稳定渗流数学模型，利用近似求解方法，给出了定产、变产和定压条件下的解。同时也低渗透油藏产能递减规律、井距确定和含水率变化进行了分析。含有动边界条件的低渗透油藏非达西不稳定渗流数学模型具有较强的非线性，难以求出精确解析解。刘文超，姚军[14-17]等人（2013）建立了能反映动边界移动规律的低渗透油藏不稳定非达西渗流数学模型，模型中不含动边界条件，经过无因次化，模型如下：5∂𝑝𝐷∂𝑡𝐷=1𝑟𝐷[𝑣𝐷(∂𝑝𝐷∂𝑟𝐷)−𝑟𝐷𝑣𝐷′(∂𝑝𝐷∂𝑟𝐷)×𝜕2𝑝𝐷𝜕𝑟𝐷2]𝑝𝐷|𝑡𝐷=0=0𝑟𝐷∂𝑝𝐷∂𝑟𝐷|𝑟𝐷=1=−(1+𝜆𝐷)𝑝𝐷|𝑟𝐷=𝑟𝑒𝐷=0𝑣𝐷(∂𝑝𝐷∂𝑟𝐷)={00≤−∂𝑝𝐷∂𝑟𝐷≤𝜆𝐷∂𝑝𝐷∂𝑟𝐷+𝜆𝐷−∂𝑝𝐷∂𝑟𝐷𝜆𝐷𝑣𝐷′(∂𝑝𝐷∂𝑟𝐷)={00≤−∂𝑝𝐷∂𝑟𝐷≤𝜆𝐷−1−∂𝑝𝐷∂𝑟𝐷𝜆𝐷通过采用Douglas-Jones预估校正有限差分法，求得了该模型瞬时压力的精确数值解。计算结果表明井底瞬时压力的双对数曲线存在拐点，拐点前的曲线反映动边界到达定压边界前的状态，拐点后的曲线反映动边界到达定压边界后的状态，如图5所示。图5瞬时井底压力的双对数曲线上述关于拟启动压力梯度问题的探讨存在一定的局限性：井-网格方程未考虑拟启动压力梯度或因处理太简单而不合理，流体在网格、网格之间的流动与井、网格之间的流动会产生矛盾，前者可使用模型直接离散化实现，后者需要对井附近的流动进行解析研究，建立井-网格方程，从而修正了内边界条件。赵国忠[18,19]（2006）建立了考虑拟启动压力梯度的三维三相渗流数值模型，通过等效半径（如图6所示）6𝑟𝑒=(∆𝑥∆𝑦/𝜋)12将井与其所在网格关联起来。图6单位时间步内的井-网格稳定流动示意图王建忠，姚军等人（2013）建立了三维两相数值模型，利用Peaceman等效半径将网格压力和井底压力联系起来：𝑝𝑖,𝑗,𝑘−𝑝𝑤𝑓=𝑄𝜇2𝜋(𝐾𝑥𝐾𝑦)12⁄ℎ𝑙𝑛𝑟𝑜𝑟𝑤+(𝑟𝑜−𝑟𝑤)𝐺进而推导了考虑启动压力梯度影响的定井底压力、定产油量和定产液量条件下的内边界条件。现有的低渗透油藏非线性流动数值模拟，大都采用传统的有限差分法或控制体积有限差分，此类方法基于相邻单元中心点的物理量来构造数值计算个事，求取压力梯度，进而进行求解。而非均质低渗油藏中不同区域甚至不同网格单元的渗流模式不尽相同，必须获取各个网格单元的压力梯度，并选择相应渗流模式，因此传统的数值模拟方法无法满足要求。姚军，黄涛等人（2014）基于拟启动压力梯度模型，结合动态渗透率概念，建立了单相数值模型。在控制体积有限差分法的基础上引入界面中心点压力和速度，提出一种新的数值算法，该算法基于单个网格来构造数值计算格式，课精确获得每个网格的压力梯度，单元间通过流量连续条件联系起来。总结：拟启动压力梯度模型较为成熟，但不能体现出非线性渗流的特点，且7当压力梯度小于拟启动压力梯度时油藏流体无法动用，大大缩小了特低渗油藏的流动范围。1.2分段模型分段模型起步较早，黄延章[22]（1998）总结出目前主要的三种分段模型为：直线逼近法、幂律函数逼近法和二次函数逼近法，如下式所示。直线逼近法{𝑣=0∆𝑃𝐿≤𝑎𝑣=(𝐾𝜇)1(𝛥𝑃𝐿−𝑎)𝑎∆𝑃𝐿𝑏𝑣=(𝐾𝜇)2(𝛥𝑃𝐿−𝑐)𝑏≤∆𝑃𝐿幂律函数逼近法{𝑣=0∆𝑃𝐿≤𝑎𝑣=𝑎1𝐾𝜇(𝛥𝑃𝐿−𝐺𝑖)𝑎2𝑎∆𝑃𝐿𝑏𝑣=(𝐾𝜇)2(𝛥𝑃𝐿−𝑐)𝑏≤∆𝑃𝐿二次函数逼近法{𝑣=0∆𝑃𝐿≤𝑎𝑣=𝐴(∆𝑃𝐿)2+𝐵(∆𝑃𝐿)+𝑐𝑎∆𝑃𝐿𝑏𝑣=(𝐾𝜇)2(𝛥𝑃𝐿−𝑐)𝑏≤∆𝑃𝐿尹芝灵，孙文静等人（2011）以分段模型为基础，结合了姚军教授提出的动态渗透率概念：K=𝑓1(∆𝑝)𝐾𝑒(𝐺𝐴∆𝑝𝑚𝐺𝑐)𝑓1(∆𝑝)=2𝑎∆𝑝+𝑏2𝑎𝐺𝑐+𝑏建立了涵盖非线性段和拟线性段的油藏数值模拟模型，其运动方程为：𝑣1={0,∆𝑝𝑙𝐺𝐴𝑓1(∆𝑝𝑙)𝐾𝑚𝑎𝑥𝐾𝑟𝑙𝜇𝑙(∆𝑝𝑙−𝐺𝐴),(𝐺𝐴∆𝑝𝑙𝐺𝑐)𝐾𝑚𝑎𝑥𝐾𝑟𝑙𝜇𝑙(∆𝑝𝑙−𝐺𝐵),(∆𝑝𝑙𝐺𝑐)通过全隐式差分格式建立了数值模型，采用预处理共轭梯度法进行求解。结果表明：启动压力梯度对压力场有显著地影响，且用此非线性渗流模型计算的压力变化比拟线性情况平缓，更符合油藏开发实际（如图7所示）8图7不同渗流模型条件下注采井间的压力分布鉴于现有分段模型的连续性和确定性问题，难以直接应用于建模，刘文超，姚军[24]等人（2012）从低渗透多孔介质渗透率渐变理论出发，将幂律分段模型进行了公式化，得到了非线性运动方程实数域内的表示：{𝑣=0,|𝑑𝑃𝑑𝑟|≤𝜆𝐶𝑣=−(𝑑𝑃𝑑𝑟|𝑑𝑃𝑑𝑟|⁄)∙𝑘𝜇∙𝜆𝐴−𝜆𝐵(𝜆𝐴−𝜆𝐶)𝑛(|𝑑𝑃𝑑𝑟|−𝜆𝐶)𝑛,𝜆𝐶|𝑑𝑃𝑑𝑟|≤𝜆𝐴𝑣=−(𝑑𝑃𝑑𝑟|𝑑𝑃𝑑𝑟|⁄)∙𝑘𝜇∙(|𝑑𝑃𝑑𝑟|−𝜆𝐵),|𝑑𝑃𝑑𝑟|≥𝜆𝐴该表达式中渗流速度关于压力梯度的一阶导数是连续的，得到的是渗透率和拟启动压力梯度也是连续变化的，符合低渗介质渗透率连续变化的客观规律：𝑘𝑏=−𝑣′∙𝜇𝐺𝑏=𝑑𝑃𝑑𝑟−𝑣𝑣′之后，通过建立单相非线性径向流动的数值模型，采用Douglas-Jones预估校正发求得了精确数值解，结果表明非线性渗流模型世界语拟线性和达西模型之间的一种理想模型。总结：分段模型比较精确，描述了不同阶段的渗流规律，对于精确的科学研究可以使用。但是应用较为麻烦，并且对于各渗流阶段的临界点判断较难。1.3连续模型邓英尔，刘慈群[25]（2001）首次提出了连续函数模型，经典的三参数模型，为非线性渗流定量分析奠定了基础：𝑣(𝑎1+𝑎21+𝑏∙𝑣)=−∇𝑃该模型的缺点是无法反应渗流存在最小启动压力梯度的现象。9杨清立，杨正明[26]等人（2007）提出了一个新的模型，既反映了最小启动压力梯度现象，又很好地描述了非线性段渗流特征，该模型如下：𝑣=𝐾𝜇𝛻𝑃(1−1𝑎+𝑏|𝛻𝑃|)时宇，杨正明[27]（2008）根据大庆油田不同区块低渗透岩心的实验结果，特出了反应真实启动压力梯度和非线性渗流段的两参数模型：𝑣=𝐾𝜇𝑔𝑟𝑎𝑑𝑝(1−𝑏1+𝑎𝑔𝑟𝑎𝑑𝑝)以上模型为早期的连续函数模型，这些模型主要利用唯像学方法，根据非线性渗流