[发明专利]一种多孔材料的联合测试方法

权利要求书

1.一种多孔材料的联合测试方法，其特征在于：包括如下步骤：

步骤(1)、搭建试验系统；

本试验系统包括气压加载装置、围压加载装置和试验腔体；圆柱形试样放置于圆形试验腔体中部，上下端分别与出气孔和进气孔相连，出气孔直接通向大气，进气孔连接气压加载装置，并在连接管路上安装阀门；围压加载装置连接到试验腔体，给试样施加相当于最大进气压2～5倍的围压；

步骤(2)、实现非稳态渗流过程；

步骤(3)、实现稳态渗流过程；

步骤(4)、构建非稳态气体渗流控制方程；

步骤(5)、利用全隐式有限体积格式求解非稳态气体渗流控制方程；

步骤(6)、基于稳态渗流结果计算绝对渗透率k_∞与滑移因子b的对应关系；

步骤(7)、基于最小二乘原理计算滑移因子、绝对渗透率和有效孔隙率；

所述步骤(2)的具体步骤为，关闭气压加载装置和进气孔之间的阀门，在气压加载装置中形成气压P₀，即P(x＝0,t＝0)＝P₀，形成的储气空间体积为V(t)≡V₀；

之后打开阀门，在压力差作用下气体由试样进气端渗流向出气端；同时，储气空间中的气体摩尔数由于渗流而逐渐减小，因此，储气空间中的气压逐渐下降；记录气压下降过程，形成P(0,t)v.s.t曲线；

当储气空间中的气压降到P₁时，结束非稳定渗流过程，将气压加载装置调整为稳压加载模式；

所述P₁为小于初始压力P₀的非稳态渗流结束控制压力，取值为0.8～0.9P₀；

所述步骤(3)的具体步骤为：进气端气压P₁保持不变，记录稳压加载过程中的储气空间体积变化过程，形成V(t)v.s.t曲线；实时计算曲线斜率，当两次计算相对差异不大于ζ时，取当前斜率的绝对值作为P₁压力条件下的气体流量Q；

所述ζ取值为1％；

所述步骤(4)的具体步骤为：

多孔材料渗流质量守恒方程：

多孔材料渗流达西定律：

多孔材料渗流气体滑移方程：

其中，φ为多孔材料孔隙率；t为时间；P为多孔材料孔隙中的气体压力，为空间、时间以及温度的函数；v为等效流速，此流速并不是气体在孔隙中的真实流动速度，而是以试样整个截面面积作为流动截面的等效流速；Q为气体流量；A为试样垂直于渗流方向上的截面面积；k为试样有效气体渗透率，为空间和时间的函数；μ为气体粘性系数，看作恒定常数；k_∞为绝对渗透率，仅与孔隙结构相关；b为滑移因子，与气体类型和孔隙结构相关；

根据公式(1)、(2)、(3)，并定义U(x,t)＝P(x,t)+b，获得气体渗流控制方程：

所述步骤(5)的具体步骤为：

非稳态渗流过程的初始状态表述为：

非稳态渗流过程的出气端边界条件表述为：

U(H,t)＝P_atm+b (6)

非稳态渗流过程的进气端边界条件，通过物质守恒原理获得：

其中，P_atm为大气压，取为1×10⁵Pa；

步骤(4)中获得控制方程以及上述初始状态和边界条件方程，组成一个二维非线性抛物线型偏微分方程定解问题；对于此类非线性问题，采用有限体积法中的全隐式Picard-Newton迭代格式进行数值求解；Picard-Newton迭代格式具有良好的二阶收敛特性；为构建Picard-Newton迭代格式，对控制方程左边时间项进行向后欧拉时间离散；对控制方程右边的k_∞U/μ项进行Taylor级数展开，继而完成Newton线性化，获得如下线性迭代格式：

式中，j为当前时间步标号，取为1,2,…,n-1；n为总时间步数；j＝0对应于初始状态；(s)和(s+1)分别代表第s次和第s+1次迭代计算结果，s取为0,1,2,…；Δt为时间步长；

采用全隐Picard空间离散格式，则上式最终变换为：

式中：i为当前空间网格点编号，取为1,2,…,m-1；i＝0为进气端；i＝m为出气端；为第i个空间网格点在j-1时间步上的最终解，表达式中未带迭代计算次数(s)的U值均为最终解；

为启动公式(9)的计算，需要定义所有网格点的初始压力值；

代表初始条件的方程离散为：

公式(9)形成含有m+1个未知量的m-1个方程，需要代入上下边界的两个边界条件方程，才能完成计算；代表边界条件的方程离散为：

其中

同时，需要假其中定第j时间步的初始值为第j-1时间步的最终解，即：

其中i＝1,2,...,m-1，j＝0,1,2,...,n，根据上述离散的控制方程公式(9)、初始条件方程公式(10)和11)和边界条件方程公式(12)和(13)，形成m+1个线性方程，继而求解此线性方程组m+1个未知量以此类推求解获得当两次迭代计算结果差值的最大值小于精度要求时，即认为获得第j时间步的最终解；其中，ξ为气压最小分辨率，取为10^-8P₀；

所述步骤(6)的具体步骤为：

对于稳定渗流，步骤(4)控制方程左侧的时间域偏微分项为0，控制方程具有解析解：

将上式以及步骤(2)获得的P₁压力条件下的流量Q代入步骤(4)所示的渗流达西定律，获得绝对渗透率k_∞与滑移因子b的换算关系：

其中，H为试样高度(m)；

所述步骤(7)的具体步骤为：

假定滑移因子为b、有效孔隙率为φ，并根据公式(16)计算绝对渗透率k_∞；利用步骤(5)中公式(9)、(10)、(11)、(12)和(13)形成线性方程组，计算获得的第j时间步的储气空间压力为采用最小二乘原理设定目标函数如下：

式中，为步骤(2)获得的P(0,t)v.s.t曲线在j时间步处的实测压力值；采用不需要直接计算梯度的修正Powell方向加速法进行计算上述目标函数的最小值，继而获得绝对渗透率k_∞、滑移因子b和有效孔隙率φ最终结果；

为求解上述目标函数最小值及其对应的绝对渗透率k_∞、滑移因子b和有效孔隙率φ，采用不需要直接计算梯度的修正Powell方向加速法进行计算；主要计算步骤为：

步骤7.1：

定义算法起始点v₀⁽¹⁾＝[lg(b_r),lg(φ_r)]、初始搜索方向s₁⁽⁰⁾＝[1,0]和s₂⁽⁰⁾＝[0,1]、收敛精度ε＝10^-8；其中，b_r和φ_r分别为设定的滑移因子和有效孔隙率的初始值，根据测试材料的经验值进行确定，对于花岗岩分别取为1.0×10^-19m²和0.005；

步骤7.2：

对于第j轮计算j＝1,2,…，沿方向s₁^(j)和s₂^(j)逐次进行一维搜索，即

步骤7.3：

分别进行以下计算

判断下述两个条件是否同时满足

f₃＜f₁ (27)

步骤7.4：

当上述两个条件同时满足时，进行以下线性搜索

令

如果停止计算，最终参数值为

否则，重复步骤7.2和7.3进行第j+1轮计算；

步骤7.5：

当步骤7.3中两个条件不能同时满足时，则取

如果停止计算，最终参数值为

否则，重复步骤7.2和7.3进行第j+1轮计算；

步骤7.6：

根据步骤7.4或7.5获得的滑移因子b，利用公式(16)计算绝对渗透率k_∞。