一种基于格子‑玻尔兹曼模型的油层流体模拟方法技术

本发明专利技术公开了一种基于格子‑玻尔兹曼模型的流体模拟方法。所述模拟方法包括，将多孔介质的图像网格化，用fi(x,y,t)表示网格点I(x,y)处，运动速度为ci的粒子所对应的粒子分布；判断粒子运动方向ci是否朝向固壁边界，是则令粒子分布fi(x,y,t)执行反向函数，否则对fi(x,y,t)执行格子‑玻尔兹曼模型的碰撞函数，然后根据演化后的粒子分布fi(x,y,t)获得流体密度ρ’(x,y)和流体速度u’(x,y)；直至满足演化结束条件。本发明专利技术通过将多孔介质的图像网格化，并将网格中的流体离散化为运动的粒子，从而根据粒子分布fi(x,y,t)获得流体密度和速度，提高了模拟的运行效率。

【技术实现步骤摘要】

本专利技术属于流体力学和计算机高性能计算领域，更具体地，涉及一种基于格子-玻尔兹曼模型的流体模拟方法。
技术介绍
岩石和土壤皆为天然的多孔介质，其孔隙直径大约在1μm～500μm之间。多孔介质中流体流动的模拟问题广泛存在于石油开采、环境保护、生物芯片和化学工程等领域，例如，在石油开采过程中，压力影响着油层的渗透率和采收率，研究压力驱动对流体流动的影响可对开采方法起到指导作用。因此，利用数值方法和计算机技术进行模拟对解决上述领域的实际问题有着十分重要的意义。在实际应用中，由于多孔介质的结构复杂、边界处理困难，同时使用者对大规模、高效和快速模拟的要求不断提高，这使得模拟多孔介质中流体流动面临重要挑战。专利文献CN 103776739A公开了一种罗伯逊－斯蒂夫流体在多孔介质中的启动压力梯度的预测方法。然而该方法在计算中需要涉及多孔介质的最小孔隙半径、最大孔隙半径以及毛细管的直线长度等多个参数，使得计算步骤较为繁琐；其次，该方法仅能模拟启动压力，无法预测流体的速度以及压力对流体流动后续的影响，预测能力有限。
技术实现思路
针对现有技术的以上缺陷或改进需求，本专利技术提供了一种基于格子-玻尔兹曼模型的流体模拟方法，其目的在于通过将多孔介质的图像网格化，并将网格中的流体离散化为运动的粒子，从而根据粒子分布fi(x,y,t)获得流体密度和速度，提高模拟的运行效率。按照本专利技术的一方面，提供了一种流体的模拟方法，包括以下步骤：S1.对多孔介质的图像进行网格化处理后获得模拟网格I(x,y)，并获得所述模拟网格I(x,y)对应的无量纲化的流体粘度ν、流体密度ρ(x,y,t)、流体速度u(x,y,t)以及固壁边界；其中，x、y分别表示模拟网格I(x,y)上网格点的横坐标和纵坐标，x为1～Nx，y为1～Ny；S2.令时刻t＝0，利用格子-玻尔兹曼模型建立流体的离散密度函数、离散速度函数以及平衡态分布函数；所述离散密度函数为所述离散速度函数为所述平衡态分布函数为fi(x,y,t)表示网格点(x,y)处，运动速度为ci的粒子所对应的粒子分布，i表示粒子速度的方向编号；表示粒子分布fi(x,y,t)的平衡态，ωi为权系数；S3.判断粒子速度ci的运动方向是否朝向固壁边界，是则通过格子-玻尔兹曼模型的反向函数，令速度ci对应的粒子分布fi(x,y,t)在时刻t+δt的取值为与i反向的粒子分布，否则执行格子-玻尔兹曼模型的碰撞函数其中，时间步长δx为网格点的边长，cs为声速，θi为方向i沿x轴的正方向逆时针旋转的夹角，松弛时间S4.通过演化后的离散密度函数以及离散速度函数，获得时刻t+δt时的流体密度ρ’(x,y)和流体速度u’(x,y)；所述演化后的离散密度函数为演化后的离散速度函数为S5.判断是否满足演化结束条件，是则模拟结束，否则令时刻t＝t+δt，返回步骤S3。优选地，在所述步骤S1中，无量纲化的流体粘度为流体的运动学粘性系数，L0为流体的特征长度，u0为流体的特征速度，且u0＝L0/1秒；在所述步骤S5之后，还包括，获得流体的实际速度作为进一步优选地，所述流体的特征长度或为模拟网格I(x,y)的总宽度，为模拟网格I(x,y)的总高度。优选地，在所述步骤S1中，对多孔介质的图像进行网格化处理的方法具体为：将多孔介质的图像转化为网格化的灰度图像，二值化所述灰度图像，并去除噪点，获得所述模拟网格I(x,y)。作为进一步优选地，所述二值化的方法为双峰法、Otsu法、Bernsen法或循环阈值法。作为进一步优选地，所述去除噪点的方法为连通域法。优选地，在所述步骤S2中，i＝0时，ci＝(0,0)；i＝1～8时，θi分别为0、π/2、π、3π/2、π/4、3π/4、5π/4以及7π/4，对应粒子分布fi(x,y,t)中的粒子在时间步长δt内运动的距离ciδt为从网格点I(x,y)到其八邻域内的另一网格点的距离；在所述步骤S3中，所述反向函数为：fi(x,y,t)＝fi+2(x,y,t)，i＝1、2、5或6；fi(x,y,t)＝fi－2(x,y,t)，i＝3、4、7或8。作为进一步优选地，在所述步骤S2中，权系数ωi满足： ω i = 4 / 9 i = 0 , ω i = 1 / 9 i = 1 ~ 4 , ω i = 1 / 36 i = 5 ~ 8. ]]>优选地，所述步骤S5中的所述演化结束条件为：时刻t是否≥时间增量T，或者流体速度u’(x,y)是否满足err为速度收敛条件。作为进一步优选地，所述步骤S5中的所述时间增量T为10000δt～1000000δt。作为进一步优选地，所述步骤S5中的所述速度收敛条件err为0～0.1。优选地，在所述步骤S5之后，还包括将流体速度u’(x,y)可视化。本专利技术具有以下有益效果：1、本专利技术通过将多孔介质的图像网格化，以多孔介质所在的区域作为固壁边界，实现了对于在多孔介质中的复杂计算工况的设定；2、本专利技术将网格中的流体离散化为运动的粒子，使得粒子间的碰撞以及粒子与固壁边界的作用能够并行处理，从而提高了模拟的计算效率和计算速度；3、由于本专利技术中不同网格点的粒子的计算可同时运行，本专利技术可利用OpenACC等方法进行并行优化，适于在多种加速器平台并行化执行。附图说明图1为本专利技术组成流体的粒子的运动方向示意图；图2为本专利技术组成流体的粒子运动朝向固壁边界并反向的示意图；图3为本专利技术组成流体的粒子向相邻网格点运动的示意图；图4为本专利技术实施例1流程图；图5为本发本文档来自技高网...

【技术保护点】
一种流体模拟方法，其特征在于，包括以下步骤：S1.对多孔介质的图像进行网格化处理后获得模拟网格I(x,y)，并获得所述模拟网格I(x,y)对应的无量纲化的流体粘度ν、流体密度ρ(x,y,t)、流体速度u(x,y,t)以及固壁边界；其中，x、y分别表示模拟网格I(x,y)上网格点的横坐标和纵坐标，x为1～Nx，y为1～Ny；S2.令时刻t＝0，利用格子‑玻尔兹曼模型建立流体的离散密度函数、离散速度函数以及平衡态分布函数；所述离散密度函数为所述离散速度函数为所述平衡态分布函数为fi(x,y,t)表示网格点(x,y)处，运动速度为ci的粒子所对应的粒子分布，i表示粒子速度的方向编号；表示粒子分布fi(x,y,t)的平衡态，ωi为权系数；S3.判断粒子速度ci的运动方向是否朝向固壁边界，是则通过格子‑玻尔兹曼模型的反向函数，令速度ci对应的粒子分布fi(x,y,t)在时刻t+δt的取值为与i反向的粒子分布，否则执行格子‑玻尔兹曼模型的碰撞函数其中，时间步长δx为网格点的边长，cs为声速，θi为方向i沿x轴的正方向逆时针旋转的夹角，松弛时间S4.通过演化后的离散密度函数以及离散速度函数，获得时刻t+δt时的流体密度ρ’(x,y)和流体速度u’(x,y)；所述演化后的离散密度函数为演化后的离散速度函数为S5.判断是否满足演化结束条件，是则模拟结束，否则令时刻t＝t+δt，返回步骤S3。...

【技术特征摘要】
1.一种流体模拟方法，其特征在于，包括以下步骤：S1.对多孔介质的图像进行网格化处理后获得模拟网格I(x,y)，并获得所述模拟网格I(x,y)对应的无量纲化的流体粘度ν、流体密度ρ(x,y,t)、流体速度u(x,y,t)以及固壁边界；其中，x、y分别表示模拟网格I(x,y)上网格点的横坐标和纵坐标，x为1～Nx，y为1～Ny；S2.令时刻t＝0，利用格子-玻尔兹曼模型建立流体的离散密度函数、离散速度函数以及平衡态分布函数；所述离散密度函数为所述离散速度函数为所述平衡态分布函数为fi(x,y,t)表示网格点(x,y)处，运动速度为ci的粒子所对应的粒子分布，i表示粒子速度的方向编号；表示粒子分布fi(x,y,t)的平衡态，ωi为权系数；S3.判断粒子速度ci的运动方向是否朝向固壁边界，是则通过格子-玻尔兹曼模型的反向函数，令速度ci对应的粒子分布fi(x,y,t)在时刻t+δt的取值为与i反向的粒子分布，否则执行格子-玻尔兹曼模型的碰撞函数其中，时间步长δx为网格点的边长，cs为声速，θi为方向i沿x轴的正方向逆时针旋转的夹角，松弛时间S4.通过演化后的离散密度函数以及离散速度函数，获得时刻t+δt时的流体密度ρ’(x,y)和流体速度u’(x,y)；所述演化后的离散密度函数为演化后的离散速度函数为S5.判断是否满足演化结束条件，是则模拟结束，否则令时刻t＝t+δt，返回步骤S3。2.如权利要求1所述的模流体模拟方法，其特征在于，在所述步骤S1中，无量纲化的流体粘度为流体的运动学粘性系数，L0为流体的特征长度，u0为流体的特征速度，且u0＝L0/1秒；在所述步骤S5之后，还包...

【专利技术属性】
技术研发人员：武琛，施保昌，柴振华，赵勇，黄昌盛，汪垒，唐冲，
申请(专利权)人：华中科技大学，
类型：发明
国别省市：湖北;42