一种采用微观周期性孔隙结构模型获取胶体沉积率的方法

本发明属于地下水污染迁移和修复领域，尤其涉及一种采用微观周期性孔隙结构模型获取胶体沉积率的方法。

背景技术：

1、胶体颗粒的迁移和沉积现象在自然界和工程中均广泛存在，比如将特殊的胶体颗粒注入含水层以修复地下水污染，受胶体的沉积、堵塞特性影响而导致修复效果十分有限。对于多孔介质中胶体沉积行为的理解，有助于人们制定合理的污染防控策略，从微观机理层面研究胶体的迁移和沉积行为具有重要意义。

2、过去的研究将多孔介质概化为单一孔隙结构，如经典happel模型、收缩管模型等，然后采用理论分析或数值模拟研究胶体在孔隙内的沉积率。实际上，受介质颗粒的接触和堆叠效应的影响，即便介质颗粒的大小统一，也会产生孔隙结构空间不均匀。过去的研究大多数忽视了复杂孔隙结构对胶体沉积的影响，导致了模拟结果与实验之间的部分偏差。

技术实现思路

1、发明目的：本发明的目的在于提供一种采用微观周期性孔隙结构模型获取胶体沉积率的方法。根据规则球形介质颗粒的堆积情况，提出三种周期性孔隙结构，并根据多孔介质的孔隙率调整模型结构，采用多种周期结构组合代表多孔介质的方法。三种孔隙结构模拟结果以及理论预测值的对比情况，证明了本发明的有效性和可靠性。

2、技术方案：本发明的一种采用微观周期性孔隙结构模型获取胶体沉积率的方法，包括如下步骤：

3、步骤1、构建多孔介质周期孔隙模型，包括立方体内球体结构模型i、体心结构模型ii和面心结构模型iii；

4、步骤2、基于多孔介质周期孔隙模型，通过stokes方程求解流场，通过eikonal方程求解多孔介质孔隙内与介质表面的距离场；

5、步骤3、基于拉格朗日方法模拟胶体的迁移运动轨迹；

6、步骤4、计算胶体的沉积率并进行验证。

7、进一步的，步骤1具体为：假定多孔介质颗粒尺寸相等，且按照周期性重复的规律堆积，根据多孔介质的堆积结构，将孔隙形态划分为三种，分别为立方体内球体结构模型i、体心结构模型ii和面心结构模型iii；根据多孔介质的孔隙率确定介质颗粒的位置和间距，以确保模型的孔隙率一致；

8、其中，模型i真实孔隙率为：

9、

10、模型ii和模型iii的真实孔隙率数学表达式为：

11、

12、式中：ε为模型或者多孔介质的孔隙率；rc为介质颗粒半径；δd为连接两个相邻介质颗粒中心，颗粒表面之间的间距；对于模型i，根据公式(1)计算得到的孔隙率恒为0.476，因此将立方体收缩为圆柱从而保证多孔介质的孔隙率可任意变化；对于模型ii和iii，根据公式(2)反过来计算δd；

13、首先建立介质颗粒堆积实体模型，然后将介质颗粒实体挖掉，最后剩下孔隙空间；因模型入口截面大小存在差异，为保证流场的一致性，根据入口总流量计算横截面流速；其中模型i采用达西流速u∞，而模型ii和iii采用孔隙间平均流速，且采用充分发展流动入口边界；根据流量计算模型ii和iii的孔隙间流速，即：

14、

15、式中：u∞为达西流速；u为孔隙间平均流速。

16、进一步的，步骤2中，所述通过stokes方程求解流场具体为：

17、stokes方程求解流场公式如下：

18、

19、式中：p为流体压强，μ为流体粘滞系数；u为流速；模型上边界设定恒定流速，下边界为零压强。

20、进一步的，步骤2中，所述通过eikonal方程求解多孔介质孔隙内与介质表面的距离场具体为：

21、令胶体表面与介质表面间距为h，定义一个最小间距σ＝0.5nm，当h<σ时冻结胶体，为了在数值计算中准确计算h，需要获取孔隙内与介质表面的距离场；考虑到孔隙结构内会存在多个介质颗粒表面，采用eikonal方程获取距离场，令介质表面上距离δw＝0，其他边界上为修正eikonal方程形式：

22、

23、式中：为nabla算子符号；δw为壁面距离；为壁面距离倒数；为最近壁面指向胶体颗粒的法向矢量；γs为光滑因子；

24、则最终胶体与介质表面的分离距离为：

25、h＝δw,p-rp (9)

26、式中：δw,p为胶体中心坐标与壁面的距离。

27、进一步的，步骤3具体为：基于拉格朗日方法模拟胶体的运动：

28、

29、颗粒受力的具体表达式为：

30、fd＝6πμrp(u-v) (11)

31、

32、式中：m为胶体质量；m*为流体中胶体的虚拟质量；v为颗粒速度矢量；f为施加在胶体上的作用力，fd为曳力，fg为有效重力，fl,saffman为升力，fb为随机布朗力，fp为压强梯度力；rp为胶体半径；ρp和ρf分别为胶体和流体密度；r(0,1)为平均值等于0、方差等于1的高斯随机数矢量；kb＝1.38×10-23j/k，为boltzmann常数；t为绝对温度；δt为数值计算的时间步长；公式(11)采用stokes曳力定律，适用于相对雷诺数rer＝dpρf||u-v||/μ＜＜1的条件，即流速较慢且胶体尺寸较小的条件。

33、进一步的，步骤4具体为：

34、采用以下公式计算胶体沉积率：

35、η＝ndep/ntot (16)

36、式中：η为胶体沉积率；ndep为最终沉积的粒子数量；ntot为从上游入口释放的粒子数量；

37、针对三种孔隙结构的模拟结果，采用解析解和非线性回归预测方程两种预测结果进行对比验证。

38、进一步的，所述解析解的表达式为：

39、

40、所述非线性回归预测方程表达式为：

41、

42、式中：γ＝(1-ε)1/3和as为孔隙率相关变量；na为吸引数；ng为重力数；npe为佩克莱数；nr为拦截数，均为无量纲变量，具体表达式如下：

43、

44、式中：a123为成对可加性的hamaker常数；kb＝1.38×10-23j/k为boltzmann常数；t为绝对温度；u∞为入口平均流速；d＝kbt/(6πμrp)为胶体扩散系数；ρp为胶体密度；ρf为流体密度；μ为流体动力粘滞系数。

45、有益效果：与现有技术相比，本发明具有如下显著优点：

46、考虑了多孔介质中孔隙结构的非均一性对胶体沉积的影响。提出采用三种代表性周期孔隙模型模拟胶体沉积，对模拟结果进行加权平均以表征多孔介质的复杂孔隙条件，解决了过去采用单一孔隙模型代表多孔介质的问题。

47、能够同时考虑扩散、重力沉积和拦截的耦合机制。令胶体颗粒在给定的分离距离处停滞，以模拟胶体颗粒的沉积状态。这种方法能够同时满足布朗(粒径小于2μm)和非布朗(粒径大于2μm)胶体沉积模拟，在方法上更为统一，克服了过去布朗胶体只考虑扩散沉积而非布朗胶体只考虑重力沉降和拦截的假设。

48、适用性和拓展性很强。本发明采用解析方法和修正的eikonal方程计算胶体与介质表面的距离，拓展了本方法的适用性，可进一步应用到介质表面粗糙、存在静电力或者表面物理化学异质性的条件。

49、综上，本方法适用于更接近真实地下水渗流和胶体沉积的条件，实现准确的模拟结果，突破简化解析解和室内试验的局限，从而辅助和指导工程实践。