(19)国家知识产权局 (12)发明 专利申请 (10)申请公布号 (43)申请公布日 (21)申请 号 202210693931.9 (22)申请日 2022.06.19 (71)申请人 西南石油大 学 地址 610500 四川省成 都市新都区新都大 道8号 (72)发明人 张启辉　李海涛　王文星　王科　 陆宽　马欣　高素娟　聂松　 马寒松　代晶晶　 (51)Int.Cl. G06F 30/20(2020.01) G06F 17/16(2006.01) G01N 15/08(2006.01) G06F 119/14(2020.01) (54)发明名称 一种受限制游走模拟以计算多孔介质内流 动迂曲度的方法 (57)摘要 本发明公开一种受限制游走模拟以计算多 孔介质内流动迂曲度的方法， 包括基于分水岭算 法简化修正样品数字岩心； 包括结合数字岩心灰 度值表征的微观尺度连续相粒子在不同区域不 同运移方向的碰撞概率； 包括针对不同数字岩心 分辨率的连续相算子运移方向以及数字岩心镜 面拓展映射方法； 通过模拟在无限大镜像映射下 稳定扩散系数， 结合上文中真实情况下的算子运 移方式来量化不同状态下的迂曲度。 本发明基于 数字岩心， 通过对多孔介质内随机游走模拟的方 法获取的迂曲度数据， 相较于现有的koponen、 sousa以及barrande等理想数学模型的方法， 引 入了运移中复杂的孔隙结构以及受力影 响， 新方 法获取的迂曲度进行模拟计算， 相较于传统方法 的迂曲度， 结果符合 率提高了11％。 权利要求书3页 说明书5页 附图4页 CN 115081210 A 2022.09.20 CN 115081210 A 1.一种受限制游走模拟以计算多孔介质内流动迂曲度的方法， 其特点在于， 包括以下 步骤： Step1： 对岩 样进行360度扫描， 基于岩 样多角度图像进行三维修 正处理； Step2： 对处理后的岩样三维灰度数据通过图像阈值分割分离岩样孔隙与岩样骨架， 对 比矩阵统计孔隙与实际测试孔隙度调整阈值数据， 获取表征真实岩 样骨架与孔隙的数据； Step3： 结合现场实验， 使用分水岭算法分割 空隙相欧式距离矩阵， 将复杂矩阵转化为 较为简单的对应样品空间连接点 集合， 设定该矩阵区域 为随机游走中无碰撞区域； Step4： 通过不断去除岩石空隙相边界体素， 结合现场实验将对应获取对应岩样在对应 扫描分辨 率下单个体素的碰撞与反弹概 率； Step5： 选取样品数字岩心矩阵中孔隙内任意一点， 进行随机游走模拟， 游走过程中记 录对应步数对应游走距离， 对游走过程中矩阵边界进行镜像映射， 随机游走模拟坐标在接 触数字岩心骨架采取不同措施限制步数或游走距离； Step6： 通过随机游走步数与距离的关系表征样品储层迂曲度。 2.一种受限制游走模拟以计算多孔介质内流动迂曲度的方法， 其特点在于， 所述Step1 具体过程为： 使用计算机遍历 重构后的三 维数字岩心矩阵， 以灰度数据量为标准， 将灰度数 据划分为两个频段， 对频段划分后的数据， 在较难区分的部分基于对应灰度值与密度的重 构算法， 转 化灰度值 为密度数据， 并代替 其灰度值数据。 3.一种受限制游走模拟以计算多孔介质内流动迂曲度的方法， 其特点在于， 所述Step2 具体过程 为： Step201： 通过 液体饱和排液法或氦气法等岩 样孔隙度实验测定方法测定样品孔隙度； Step202： 基于修正后的三维数字岩心数据， 通过设定灰度或密度阈值划分骨架与孔 隙， 遍历划分后的二值化孔隙， 对矩阵内标记为孔隙的坐标进行bwlabel算法计算其独立 性， 标记并计算调整阈值后不连通孔隙所占体积， 对比Step201中实验所获的孔隙度φ1， 调 整阈值后统计孔隙度φ2以及不连通孔隙φu， 通过不断调 整阈值使得φ1＝φ2‑φu， 从而确 定能够反应真实岩心孔隙与骨架结构的数据。 4.一种受限制游走模拟以计算多孔介质内流动迂曲度的方法， 其特点在于， 所述Step3 具体过程为： 对数字岩心中根据Step2过程中划分为孔隙的区域中灰度值进行高斯平滑操 作， 抹除极小值后调整分水岭算法中孔隙灰度阈值， 根据数字岩心中灰度划分为多个相似 小区间并根据其灰度值大小进 行编号， 编号， 体素半径R以及其对应阈值Th储存在 对应数字 岩心对应[x， y， z]位置矩阵中。 5.一种受限制游走模拟以计算多孔介质内流动迂曲度的方法， 其特点在于， 所述Step4 具体过程 为： 基于波尔茨曼 方程中对于连续相速度对碰撞相的影响： 其中： f为无量纲外力， Xα为对应位置， ξα为无量纲粒子扩散 速度 考虑多孔介质中迂曲度计算忽略外力作用， 结合玻尔兹曼方程中对于速度场于碰撞相 的方程， 数字岩心中某一 孔隙点的扩散分布概 率可以等效为： 权　利　要　求　书 1/3 页 2 CN 115081210 A 2考虑数字岩心中孔隙网络可以等效为长圆管， 因此结合流动过程中泊肃叶方程以及 Step3过程中分水岭算法分隔的独立孔隙， 则岩样流动过程中各体素分布概率可以通过如 下方程表示： 化简， 则微观 尺度中岩心孔隙中单个 体素在进行自由扩散流动时各 方向概率为： 其中： ThN为数字岩心中灰度值； Ri为分水岭算法中计算所得等效半径； θ*为接触角；μ为 流体粘度。 6.一种受限制游走模拟以计算多孔介质内流动迂曲度的方法， 其特点在于， 所述Step5 具体过程 为： Step501:基于Step2步骤中获取的孔隙网络模型， 在孔隙网络中随机生成一算子， 读取 Step3过程中获取的数据 矩阵， 结合Step4过程计算该随机体素在下一迭代时间步长内的周 边N个体素的运移 概率， 若数字岩心分辨率较低， 即单个体素尺寸较大， 限于精度限制， 该模 型无法精细表达孔隙结构， 因此可能存在低于精度的孔隙结构， 考虑周边27个体素作为运 移方向； 若 数字岩心分辨率较高， 能够表征精细孔隙结构， 则取单个体素周边 14或6个方向， 对应体素运移概 率为： 当算子在随机游走过程中运移到数字岩心矩阵中骨架位置时， 则算子位置保持不动， 下一次迭代运移过程中忽略当前运移位置概 率； Step502： 经过Step501过程获取数字岩心矩阵中获取的不同方向概率， 调用编程软件 中自带随机数库， 使算子在数字岩心矩阵对应孔隙位置内进 行带有对应方向概率的随机游 走； 对数字岩心矩阵进 行横向以及纵向的镜像映射， 用以扩充边界， 对应边界镜面映射的过 程为: 其中： xi+1， yi+1代表映射前的算子所处矩阵内坐标； REF(x,y)代表映射后的数字岩心矩 阵坐标； b为矩阵边界大小； Step503： 使用编程软件重复迭代Step501 ‑Step502过程， 同时运行一次无限制空间下 的运移模拟， 即矩阵内全部为孔隙； 根据不同需要设计并记录不同时间步长下算子所处位权　利　要　求　书 2/3 页 3 CN 115081210 A 3