Arbitrary resolved-unresolved CFD-DEM coupling method and its application to seepage flow analysis in sandy soil
-
摘要: 基于欧拉-拉格朗日连续体与非连续体耦合理论进行岩土工程流固耦合问题分析是一种较新颖和盛行的方法。针对该理论下的全解流（Fully-resolved）耦合与非解流（Un-resolved）耦合方法各自的缺陷,在已有的半解流（Semi-resolved）流固耦合数值方法（将全解流与非解流联合）基础上,通过引入修正的高斯权函数,建立了新的任意解流流固耦合方法（Arbitrary Resolved-Unresolved CFD-DEM coupling method）。该任意解流流固耦合方法能够较好地解决全解流方法中由于对粗颗粒周围流场精细化所带来的计算量过大问题；同时,成功地解决了流体网格内细颗粒较大时无法获得局部平均化变量问题；因此,该方法能够对具有一定粒径级配砂土土体的流固耦合问题开展模拟分析。通过室内砂土向上渗流试验,对所建立的任意解流流固耦合方法的准确性和有效性进行了验证；进一步地,采用该任意解流流固耦合模型,从细观层面上分析和研究了砂土渗流过程中水力梯度、土体变形随渗流速度变化规律。建立的任意解流流固耦合方法能够为岩土工程土体渗流问题的研究提供新的方法和手段。Abstract: The Euler-Lagrange coupling scheme based on the continuous and discrete theories has been becoming increasingly popular in numerical analysis of fluid-particle interaction. In this study, by introducing the modified Gaussian weighting function, a new arbitrary resolved-unresolved CFD-DEM coupling method (ARU CFD-DEM) is proposed based on the authors’ previous developed semi-resolved coupling approach by combing the fully-resolved and un-resolved coupling methods. This ARU CFD-DEM method is powerful to relieve the overload in computation due to refining the flow field around the coarse particles in the fully-resolving method. At the same time, it is also able to solve the difficulty in computing the local averaging variables when fine particles with large diameter exist in fluid grids. By doing so, the ARU CFD-DEM is able to simulate the fluid-particle interaction in sand mass which contains a wide range of particle diameters. By comparing with the results of upward seepage flow tests in sand, the accuracy and effectiveness of ARU CFD-DEM model is verified. Furthermore, the hydraulic gradient-flow velocity relationship and soil deformation-flow velocity relationship in the upward seepage flow are analyzed on the particle-scale by the ARU CFD-DEM. The proposed ARU CFD-DEM model can provide a new tool for investigating the fluid-particle interaction in the seepage flow in geotechnical engineering.
-
Keywords:
- CFD-DEM coupling /
- seepage flow /
- sand particles /
- drag force /
- void fraction
-
0. 引言
近年来，随着中国城市化进程不断推进，一系列问题也日益突显，其中就包括交通拥堵问题。将地下交通隧道作为城市道路系统的延伸和补充，有利于城市可持续发展，引导城市交通系统向立体化扩展[1]，但地下交通隧道在施工期容易引起周边建筑、地层的扰动变形，必须实施严格的监测。
针对施工期隧道结构的内力及变形监测，已有较多的相关研究。在监测内容方面，张雪婵[2]对杭州地铁1号线湘湖站北2基坑坍塌事故进行分析后发现，实际监测方案比规范中少了周围地下关系位移、土体侧向变形和立柱沉降3项内容，这导致了相关单位对现场情况了解不足，对事故应对不及时；Brownjohn[3]对2004年新加坡某隧道坍塌事故的报告进行研究后发现，现场的离线监测仪器在两个月前已经记录到隧道有过大变形，但并未及时处理，这大大推动了后来无线自动监测技术的发展。
针对监测数据的处理，早期工程中主要利用参数监测与反分析方法结合[4]，依据跟踪实测数据，判别周围地层的稳定性；随着无线监测技术的发展，施工监测往往能获得大量的数据，对此，国内许多学者的监测试验结果均表明，监测数据多随时间呈离散变化，采用时间序列分析法能有效提高对监测数据变化趋势的预测精度[5]。
由此可见，施工期监测技术不论体系还是方法都已经较为成熟。然而，上述成果多数针对的是盾构法或明挖法施工的隧道。不同于盾构隧道，顶管隧道管片是移动的，这导致其在整个施工过程中的内力演化可能十分复杂。目前，对于这种演化规律的相关研究还较少，王承德[6]从理论计算上得出了顶管不直接承受土压力作用的结论，但监测数据资料较少，不足以反映此规律；魏纲等[7]通过现场试验，积累了一定的实测数据，但监测对象长度仅为22 m，规模较小。陆翔路-祁连山路贯通工程Ⅱ标段的顶管隧道一次顶进距离达445 m，本研究通过对该工程进行跟踪监测，使用高频率的自动采集工具获得了大量实测数据。基于这些监测成果，本文对浅覆土工况下施工期顶管隧道内力演化规律进行了分析研究，并与以往的研究结论进行对比和验证，得出了一些结论，能为顶管隧道的设计计算提供依据。
1. 工程概况
陆翔路—祁连山路贯通工程Ⅱ标段的地道工程埋深约4.5 m，包含两条相互平行的断面尺寸9.9 m×8.15 m、一次顶进距离445 m的矩形顶管隧道，因此兼具了大断面、浅覆土和超长距离的特点。顶管隧道依次下穿顾村公园、河浜段、S20外环高速，并且靠近上海地铁7号线和在建15号线，其平面位置如图 1所示。
场区所处的位置在正常沉积土层和古河道沉积的交替处，地层分布比较复杂，自上而下可分为9个大层，表层一般为杂填土，平均厚度在2 m左右。拟建线路沿线有3处明浜，浜底有0.7~1.6 m的含有机质淤泥，另有少数几处暗浜等不良地质条件，土质极差。
拟建场地地下水类型主要包括潜水和承压水，潜水主要赋存于浅部黏性土层和浅部粉性土层中，水位埋深在0.3~1.5 m，承压水包括分布在第5层的微承压水和中下部的承压水，对工程建设影响较大。
2. 监测方案
本项目的监测内容为管节的接触压力和内力，是通过土压力计和钢筋计来获取监测数据的。本次监测使用的振弦式柔性土压力计测量范围为0~1.0 MPa，分辨率≤0.06%FS。管节内力监测使用的是振弦式钢筋计，最大测量压应力为100 MPa，最大拉应力为200 MPa。
考虑监测数据量的要求以及现场工期协调问题，本项目选取了西线第16，50，130环管节用于埋设所需传感器。其中，土压力计均布置在4个侧面的中心位置；钢筋计的安装需要能反映多个位置的受力变化情况，同时又不致使传感器过多而影响浇筑，因此选择在第16环和第50环安装钢筋计各8支，并使其分布相对对称，兼顾管节的内侧和外侧，实际的安装位置示意如图 2所示。
![]()
图 2 第16环、50环钢筋计安装位置示意图Figure 2. Installation positions of reinforcement stress gauges on 16th and 50th pipe rings待上述特殊管节下放至工作井后，将各传感器的线缆接入采集箱内，施工前期设置的采集间隔为5 min，后期调整为10 min。
3. 监测结果分析
本次监测从第16环管节下放开始，持续整个施工期，在顶管工程贯通后，第50环和130环还继续监测了一段时间。尽管有部分传感器未能存活，本研究仍收集到了大量监测数据。实测的结果大致可以反映出一些规律，但也有部分数据存在问题，如第16环右侧的土压力计始终保持在一个固定值，或是某些钢筋计的测量值在某些时刻显著偏离正常范围甚至超出量程。本研究在后期分析时去掉了这部分数据，最终形成的监测结果按土压力和管节内力两个方面展开。
3.1 土压力监测结果分析
以往的理论分析认为[6-7]，顶管环向受到的接触压力主要由注浆压力和土压力两个因素决定，但主要承受的是注浆压力，而非土压力。现行的顶力计算规范中，使用的也是与土压力无关的公式。注浆能减小侧壁所受压力的水平，且对顶部和底部两侧的影响最为明显。管道底部的土体反力与顶部相同，大多数情况下这两个位置的受力都比较接近，且底部大一些。此次检测的结果大致符合这一规律，但也略有出入。
尽管从整体上看，监测数据是不断变化的，但在连续的一段时间内，其数据往往是较为稳定的。因此，本研究选择将测点各传感器监测数据按整体监测天数为区分，尽量选取各天相同监测时段内较稳定的监测值作为当天的代表值，绘制出趋势图，如图 3所示，图中各编号对应的实际测点位置见表 1。
表 1 编号与实际土压力测点位置对照表Table 1. Serial number and locations of measuring points for actual earth pressure编号 位置 #1 第16环底部 #2 第16环顶部 #3 第16环左侧 #4 第50环右侧 #5 第50环左侧 #6 第50环顶部 #7 第130环左侧 由图 3可见，在整个施工过程中，土压力值的水平都较低且波动明显，在后期施工结束但#4~#7测点监测仍在继续的阶段，其监测结果才较为稳定，压应力保持在小于0.05 MPa的状态。由于此次监测仅第16环底部土压力计存活，故顶部与底部的对照略有缺失，根据第16环的监测结果来看(即#1、#2)，顶部土压力稍大于底部土压力的受力水平，但二者的值十分接近，这与理论结果不符。
与顶部和底部相比，左右两侧的受力水平(即#3、#4、#5、#7)则要高一些，而这其中又以第16环最为明显，此处最大土压力为0.197 MPa，对比之下，此管节顶部最大土压力仅为0.022 MPa。以往的测试结果显示[7]，管节所受压力大小受覆土压力和注浆压力影响，其中，注浆压力又起到控制作用。通过注浆形成泥浆套，能够减小管壁与周围土体的接触压力，这种影响对顶部和底部尤为明显，而对左右两侧影响不大，本次监测结果也符合这一规律。
3.2 管节内力监测结果分析
顶管施工时，顶力的来源是位于始发井的液压千斤顶装置，顶力通过管节不断向前传递，一般认为，管节纵向钢筋受力包括顶力、开挖面阻力以及周围土体的摩阻力等，管节环向钢筋则主要受到自重、水土压力和注浆压力。由于顶管管道本身是弯曲的，在弯曲处管道接头容易产生应力集中现象，使钢筋受到拉应力的作用；在直线顶进段，或是更换受顶管节时，管节内部既存在受拉部分，也存在受压部分，且大体上呈对称分布[7]，本次监测结果也能反映这一规律。钢筋内力监测数据的处理方式与土压力数据相同，并将纵向和环向数据分开呈现，如图 4所示，各编号对应的实际测点位置见表 2。
表 2 编号与实际钢筋计测点位置对照表Table 2. Serial number and locations of measuring points for actual reinforcement stress编号 位置及安装方向 #1 第16环底部内侧、纵向 #2 第16环顶部外侧、纵向 #3 第16环右侧内部、纵向 #4 第50环顶部内侧、纵向 #5 第50环底部外侧、纵向 #6 第50环右侧内部、纵向 #7 第16环底部内侧、环向 #8 第16环底部外侧、环向 #9 第16环顶部外侧、环向 #10 第16环顶部内侧、环向 #11 第50环左侧内部、环向 #12 第50环左侧外部、环向 #13 第50环右侧内部、环向 #14 第50环右侧外部、环向 由图 4可见，受实际工况不断变化的影响，管节各部位的受力状态也在不断改变，总体受力水平保持在一定范围内。监测结果显示，纵向钢筋计所受最大拉应力为11.136 MPa，最大压应力为3.984 MPa；而环向钢筋计受到的最大拉应力为13.850 MPa，最大拉应力为5.390 MPa。相对而言，纵向各传感器的受力水平更为接近，环向受力则更为复杂多变。
进一步整理数据还发现，靠管节内侧的纵向钢筋计(#3、#4、#6)容易表现出受压的性质，靠外侧的钢筋计(#2、#5)则表现为受拉性质，这与以往的测试结果[7]相近，这表明，内侧钢筋受力主要受管节顶部内力影响，外侧钢筋受力则取决于左右两侧的内力。唯一不符合此规律的是#1测点，由于该测点位于管节底部，靠近机头堆载，与其他测点受力情况略有不同，故可认为该因素对钢筋的受力具有一定的影响。
对于环向钢筋的受力方面，各传感器在施工期间均表现出了受拉和受压两种状态，四个侧面受力性质分布比较均匀，但其数值变化情况则比较复杂，缺乏规律性。由此表明，环向钢筋受力水平受到地层性质、注浆参数、线路弯直等因素的影响，在实际施工过程中，这些因素都是随工程推进而不断变化的，多种因素叠加作用下导致了环向钢筋受力不断波动的现象。
4. 结论
通过对陆翔路顶管工程监测数据的处理分析，本研究总结的施工期顶管内力演化基本规律如下：
(1) 管节环向接触压力受注浆压力和土压力影响，但管节本身不直接承受土压力，而是注浆压力起控制作用。通过注浆能够有效减小管节受到的接触压力，这种作用对顶部和底部最为明显，而对左右两侧影响不大。顶部和底部的受力水平与实际工况有关，并无直接的大小关系。
(2) 实际工程中，纵筋受力性质受管节内力、内部堆载等因素的影响，变化规律比较复杂。一般情况下认为，管节内侧的纵向钢筋主要受顶部内力的影响，容易表现出受压的性质，外侧的纵筋主要受两侧内力的影响而更多地表现出受拉的性质。
(3) 环向钢筋在施工期可以表现出受拉、受压两种性质，且在管节上的分布区域较为均匀。环向钢筋的受力水平受到地层性质、注浆参数、线路弯直等多种因素的影响，这些因素的叠加作用会导致受力不断波动的现象，规律性较差。
-
![]()
图 1 高斯权函数支撑域示意图
Figure 1. Diagram of supporting domain of Gaussian-based weighting function
![]()
图 4 全解流与任意解流用时对比
Figure 4. Comparison in CPU time between full resolved and ARU method
![]()
图 7 水力梯度、滤层变形随表观渗流速度变化关系
Figure 7. Relationship between hydraulic gradient and coarse-particle layer deformation with superficial flow velocity
![]()
图 8 不同阶段粗颗粒滤层变形随表观渗流速度变化情况
Figure 8. Evolution of coarse-particle layer deformation with superficial flow velocity
![]()
图 9 土层孔隙率分布随表观渗流速度变化情况（粗细颗粒粒径比αs=6.7）
Figure 9. Evolution of porosity profile with superficial flow velocity (diameter ratio, αs=6.7)
![]()
图 10 土基颗粒力链分布随表观渗流速度变化情况（粗细颗粒.粒径比αs=6.7）
Figure 10. Evolution of force chain in fine-particle layer with superficial flow velocity (diameter ratio, αs=6.7)
表 1 基准渗流模拟所采用材料属性
Table 1 Material properties in the benchmark of seepage flow
颗粒属性 流体属性 颗粒密度ρs/(kg·m-3) 滑动摩擦系数μr 杨氏模量E/Pa 回弹系数e 泊松比 滚动摩擦系数μf 密度ρf/(kg·m-3) 动力黏度μf/(kg·m-1·s-1) 2650 0.84 2×107 0.9 0.2 0.26 1000 1×10-3 
下载: 导出CSV
表 2 基准渗流模拟相关设置
Table 2 Computational settings in the benchmark of seepage flow
计算设置 任意解流设置 全解流设置 DEM时间步长ΔtDEM/s CFD时间步长ΔtCFD/s 流体网格尺寸Lm/mm 高斯权函数带宽bw/mm 流体网格尺寸Lm/mm 2×10-5 2×10-4 1.32, 0.8, 0.5 8.0 0.5
下载: 导出CSV
表 3 滤层-土基数值模型中的颗粒属性
Table 3 Material properties in numerical simulations of seepage test
参数 数值 参数 数值 砂颗粒密度ρs/(kg·m-3) 2650 玻璃-玻璃滑动摩擦μs_gg 0.1545 玻璃珠密度ρg/(kg·m-3) 2450 砂-玻璃滑动摩擦μs_sg 0.3 砂颗粒杨氏模量Es/Pa 2×1010 砂-砂滚动摩擦μr_ss 0.26 玻璃珠杨氏模量Eg/Pa 5×1010 玻璃-玻璃滚动摩擦μr_gg 0.045 砂颗粒泊松比 0.2 砂-玻璃滚动摩擦μr_sg 0.1 玻璃珠泊松比 0.25 砂颗粒回弹系数es 0.9 砂-砂滑动摩擦μs_ss 0.84 玻璃珠回弹系数eg 0.9
下载: 导出CSV
表 4 滤层-土基数值模型中的流体属性
Table 4 Computational settings in numerical simulations of seepage test
密度ρf/(kg·m-3) 动力黏度μf/(kg·m-1·s-1) DEM时间步长ΔtDEM/s 1000 1×10-3 2×10-5
下载: 导出CSV
表 5 滤层-土基数值模型中的计算设置
Table 5 Computational settings in numerical simulations of seepage test
CFD时间步长ΔtCFD/s 流体网格尺寸Lm/mm 高斯权函数带宽bw/mm 2×10-4 1, 0.5 3.0
下载: 导出CSV
-
[1] 周健, 姚志雄, 张刚. 砂土渗流过程的细观数值模拟[J]. 岩土工程学报, 2007, 29(7): 977-981. doi: 10.3321/j.issn:1000-4548.2007.07.004 ZHOU Jian, YAO Zhi-xiong, ZHANG Gang. Mesomechanical simulation of seepage flow in sandy soil[J]. Chinese Journal of Geotechnical Engineering, 2007, 29(7): 977-981. (in Chinese) doi: 10.3321/j.issn:1000-4548.2007.07.004
[2] EL SHAMY U, ZEGHAL M, DOBRY R, et al. Micromechanical aspects of liquefaction-induced lateral spreading[J]. International Journal of Geomechanics, 2010, 10(5): 190-201. doi: 10.1061/(ASCE)GM.1943-5622.0000056
[3] ZHAO J D, SHAN T. Coupled CFD-DEM simulation of fluid-particle interaction in geomechanics[J]. Powder Technology, 2013, 239: 248-258. doi: 10.1016/j.powtec.2013.02.003
[4] 蒋明镜, 张望城. 一种考虑流体状态方程的土体CFD-DEM耦合数值方法[J]. 岩土工程学报, 2014, 36(5): 793-801. https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201405002.htm JIANG Ming-jing, ZHANG Wang-cheng. Coupled CFD-DEM method for soils incorporating equation of state for liquid[J]. Chinese Journal of Geotechnical Engineering, 2014, 36(5): 793-801. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTGC201405002.htm
[5] 王胤, 艾军, 杨庆. 考虑粒间滚动阻力的CFD-DEM流-固耦合数值模拟方法[J]. 岩土力学, 2017, 38(6): 1771-1780. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201706027.htm WANG Yin, AI Jun, YANG Qing. A CFD-DEM coupled method incorporating soil inter-particle rolling resistance[J]. Rock and Soil Mechanics, 2017, 38(6): 1771-1780. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201706027.htm
[6] HAGER A, KLOSS C, PIRKER S, et al. Parallel resolved open source CFD-DEM: method, validation and application[J]. The Journal of Computational Multiphase Flows, 2014, 6(1): 13-27. doi: 10.1260/1757-482X.6.1.13
[7] WANG Z L, FAN J R, LUO K. Combined multi-direct forcing and immersed boundary method for simulating flows with moving particles[J]. International Journal of Multiphase Flow, 2008, 34(3): 283-302. doi: 10.1016/j.ijmultiphaseflow.2007.10.004
[8] ANDERSON T B, JACKSON R. Fluid mechanical description of fluidized beds. equations of motion[J]. Industrial & Engineering Chemistry Fundamentals, 1967, 6(4): 527-539.
[9] TSUJI Y, KAWAGUCHI T, TANAKA T. Discrete particle simulation of two-dimensional fluidized bed[J]. Powder Technology, 1993, 77(1): 79-87. doi: 10.1016/0032-5910(93)85010-7
[10] CHENG K, WANG Y, YANG Q. A semi-resolved CFD-DEM model for seepage-induced fine particle migration in gap-graded soils[J]. Computers and Geotechnics, 2018, 100: 30-51. doi: 10.1016/j.compgeo.2018.04.004
[11] DI RENZO A, DI MAIO F P. Comparison of contact-force models for the simulation of collisions in DEM-based granular flow codes[J]. Chemical Engineering Science, 2004, 59(3): 525-541. doi: 10.1016/j.ces.2003.09.037
[12] DI FELICE R. The voidage function for fluid-particle interaction systems[J]. International Journal of Multiphase Flow, 1994, 20(1): 153-159. doi: 10.1016/0301-9322(94)90011-6
[13] WANG Y, ZHOU L X, WU Y, et al. New simple correlation formula for the drag coefficient of calcareous sand particles of highly irregular shape[J]. Powder Technology, 2018, 326: 379-392. doi: 10.1016/j.powtec.2017.12.004
[14] 王胤, 周令新, 杨庆. 基于不规则钙质砂颗粒发展的拖曳力系数模型及其在细观流固耦合数值模拟中应用[J]. 岩土力学, 2019, 40(5): 2009-2015. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201905046.htm WANG Yin, ZHOU Ling-xin, YANG Qing. New drag coefficient model for irregular calcareous sand particles and its application into fluid-particle coupling simulation[J]. Rock and Soil Mechanics, 2019, 40(5): 2009-2015. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201905046.htm
[15] YANG Q, CHENG K, WANG Y, et al. Improvement of semi-resolved CFD-DEM model for seepage-induced fine-particle migration: Eliminate limitation on mesh refinement[J]. Computers and Geotechnics, 2019, 110: 1-18. doi: 10.1016/j.compgeo.2019.02.002
[16] JASAK H. Error Analysis and Estimation for the Finite Volume Method with Applications to Fluid Flows[D]. London: Imperial College, University of London, 1996.
[17] ISSA R I. Solution of the implicitly discretised fluid flow equations by operator-splitting[J]. Journal of Computational Physics, 1986, 62(1): 40-65. doi: 10.1016/0021-9991(86)90099-9
[18] SHIRGAONKAR A A, MACIVER M A, PATANKAR N A. A new mathematical formulation and fast algorithm for fully resolved simulation of self-propulsion[J]. Journal of Computational Physics, 2009, 228(7): 2366-2390.
[19] KLOSS C, GONIVA C. LIGGGHTS - open source discrete element simulations of granular materials based on lammps[M]//Supplemental Proceedings. Hoboken: John Wiley & Sons, Inc., 2011: 781-788.
[20] WELLER H. OpenFOAM: The open source CFD toolbox, Programmer’s Guide[EB/OL]. UK: CFD Direct Ltd, 2009, Http://www.openfoam.com.
-
期刊类型引用(4)
1. 周兵权，李喜安，秦彪，赵文琪. 黄土持水曲线的CT分析及其在晋西洞穴发育特征中的应用. 农业工程学报. 2025(02): 155-162 . 
百度学术
2. 常洲，晏长根，安宁，兰恒星，石玉玲，包含，许江波. 干湿循环作用下原状黄土渗透性及其对土-水特征曲线的影响. 长江科学院院报. 2024(01): 143-150+166 . 百度学术
3. 李奇，陈文娟，石文豪，孙少波，张永根. 基于单一和集合土壤转换函数模型对土壤含水量的模拟性能分析. 土壤. 2023(03): 658-670 . 百度学术
4. 鄢士程，解广成，辛全明，王振霖. 多监测条件下地铁明挖车站周边地表沉降预测分析. 电子测量技术. 2022(02): 48-54 . 百度学术
其他类型引用(6)
计量
- 文章访问数: 344
- HTML全文浏览量: 39
- PDF下载量: 345
- 被引次数: 10
