使用耦合的离散元方法（DEM）- 有限差分方法（FDM）模拟对砂土中进行圆锥贯入试验的数值研究

《Computers and Geotechnics》：Numerical investigation of cone penetration test in sand using a coupled DEM-FDM simulation

　　**摘要**
本研究开发了一种结合离散元方法（DEM）和有限差分方法（FDM）的耦合数值模拟技术，用于在虚拟校准室（VCC）中模拟锥形贯入试验（CPT），旨在研究影响锥端阻力qc、靴体摩擦力fs以及摩擦比Rf的因素。通过将模拟的宏观应力-应变行为与排水三轴试验的实验结果进行匹配，对Fontainebleau砂的微观参数进行了校准。基于校准后的模型，系统地评估了VCC尺寸、边界条件及应力状态对qc和fs的影响。此外，还分析了水平围压σ3对Rf的影响。训练了一种人工神经网络（ANN）模型，以预测qc、fs和Rf作为相对密度Dr、围压p、颗粒杨氏模量E和峰值内摩擦角?peak的函数。结果表明，DEM-FDM模拟的计算时间显著低于单独使用DEM的计算时间。当DEM和FDM的径向尺寸分别设置为锥体直径的8倍和12倍时，耦合 DEM-FDM模拟中的尺寸效应可以忽略不计。模拟中的位移和速度等值线显示出了良好的连续性。通过模拟得到的径向应力σr与基于腔体膨胀方法的预测结果之间的高度吻合，验证了该耦合方法的可靠性。最终，训练好的ANN模型表现出令人满意的预测能力，其预测值与实测值在qc、fs和Rf方面具有合理的吻合度。

**引言**
快速的城市化和大规模基础设施的发展增加了对高效岩土工程勘察方法的需求（Patel, 2019）。因此，锥形贯入试验（CPT）作为一种原位测试方法得到了广泛应用，因其简便性、高效率以及能够提供随深度变化的连续测量数据（Robertson, 2009; Mo et al., 2020; Ouyang and Mayne, 2020）。在测试过程中，锥形探头以恒定速率穿透土壤，并同时记录锥端阻力qc和靴体摩擦力fs等数据。CPT并不直接测量土壤参数，而是通常利用已建立的经验关系式从CPT测量数据中推断土壤参数（Jamiolkowski et al., 2003; Mayne, 2007; Chen and Abousleiman, 2012; Kouretzis et al., 2014）。尽管有这些优势，但从CPT测量数据估计土壤参数仍存在不确定性，因为经验关系式可能无法完全考虑土壤状态、应力历史、排水条件和材料属性的变化。为了改进CPT的解释，许多研究通过将CPT测量结果与原位测试、实验室测试、经验关系式及数值模拟得到的土壤参数进行关联，开发了基于现场的数据相关性（Huang and Hsu, 2005; Damavandi-Monfared and Sadrekarimi, 2015; Ahmadi and Dariani, 2017; Yang et al., 2022; Hu and Wang, 2024）。这些经验关系式是基于原位CPT测量结果建立的，然后将其与其他测试方法（如实验室实验或现场测试）得到的土壤参数进行关联。Tonni和Simonini（2013）利用原位测试数据研究了qc与剪切波速度之间的关系，改进了适用于中等密实度砂的经验公式。Saye等人（2017）提出了qc与fs之间的线性关系，以加强对不同压实程度砂的分类，并建立了其与颗粒大小和Atterberg界限等指标属性的联系。同样，Tang等人（2017）在未饱和残余土中进行了CPT测试，并结合基质吸力和有效应力开发了一种估计有效摩擦角的方法，提高了摩擦角估算的准确性。尽管这些研究推进了基于CPT的土壤特性分析，但由于初始应力状态、密度、土壤结构和排水条件无法独立控制，因此基于现场的数据相关性往往具有现场特异性。因此，校准室测试为建立CPT测量结果与土壤属性之间的关联提供了更为可控的环境，因为可以在试验中更明确地规定初始密度、应力状态和边界条件。

**校准室测试**
校准室测试已被广泛用于在受控应力及边界条件下研究CPT测量结果与土壤属性之间的关系。早期研究（Been et al., 1986; Mayne and Kulhawy, 1991; Salgado et al., 1998）探讨了室尺寸和边界条件对qc的影响，并建立了将QC与土壤状态参数相关联的经验关系式。最近，Ghali等人（2020）开发了一种专门用于模拟现场条件的校准室，证明了其在非塑性淤泥和颗粒材料中捕捉应力-应变行为及体积变化的能力。Sadrekarimi和Hashemi（2025）研究了实验室校准室测试中边界条件对qc的影响，研究指出刚性边界会导致密实样品中的qc值升高，并提出了相应的修正因子。然而，实验室校准室测试仍存在一些局限性，尤其是在制备具有预定相对密度的砂样及保证高重复性方面（Parkin and Lunne, 1982）。此外，锥形贯入过程中的应力状态和位移场难以直接测量。因此，数值模拟提供了一种有用的补充方法，因为它可以明确控制应力状态和边界条件，同时能够详细研究土壤在贯入过程中的响应。

**数值模拟**
数值模拟已成为研究CPT响应的有效方法，其中有限元方法（FEM）和离散元方法（DEM）是主要应用的技术（Arroyo et al., 2011; Falagush et al., 2015; Chen et al., 2020; Sharif et al., 2021; Chouhan and Chavda, 2023）。由于其计算效率高，FEM被广泛用于模拟锥形贯入问题。Ahmadi和Dariani（2017）使用Drucker-Prager本构模型对砂中的CPT进行了FEM模拟，分析了探头周围应力状态的演变，并开发了一种解析解，用以从qc和fs估算原位水平应力及相对密度。然而，FEM无法直接解析锥体与土壤颗粒之间的相互作用。相比之下，DEM能够提供颗粒级别的微观力学行为见解。Khosravi等人（2020）在虚拟校准室（VCC）中进行了CPT的DEM模拟，量化了不同应力状态和室直径下的室边界效应。结果显示，qc和fs随室直径的增加而系统性增加，在各向同性应力状态下的值明显高于k0应力状态下的值，这反映了它们对施加的侧向应力状态和室壁边界效应的敏感性。尽管在有限的情况下研究了孔隙比，但尚未系统地研究砂的相对密度对qc的影响。在此基础上，Song等人（2025）在VCC中利用DEM模拟研究了不同相对密度砂的qc受模型尺寸和应力状态的影响。然而，DEM的计算效率相对较低。此外，以往的DEM研究主要集中在qc上，而fs的演变及k0应力状态的影响尚未得到系统研究。为了在保持锥体附近颗粒级描述的同时提高计算效率，本研究开发了一种耦合DEM-FDM框架用于CPT模拟。该框架通过将k0应力条件下的模拟径向应力分布与基于非塑性本构模型的腔体膨胀解的预测结果进行比较来验证（Salgado et al., 1998; Cudmani and Osinov, 2001）。耦合方法的改进效率使得能够进行广泛的参数化模拟，从而生成了涵盖不同密度、应力状态和模型尺寸的数据库。随后利用该数据库训练了人工神经网络（ANN）模型，以数据驱动的方式预测CPT响应。

**本研究**
本研究开发了一种耦合DEM-FDM框架，用于在虚拟校准室中模拟砂中的CPT。DEM域用于捕捉锥体周围的颗粒级响应，而周围的FDM域代表远场土壤，以降低计算成本。Fontainebleau砂的模型参数是根据排水三轴试验结果进行校准的。通过与完整DEM模拟和腔体膨胀方法的比较，验证了该耦合框架的有效性。基于验证后的模型，系统研究了室尺寸、边界条件、相对密度Dr和应力状态对锥端阻力qc、靴体摩擦力fs以及摩擦比Rf的影响。最后，开发了一种人工神经网络（ANN）模型，用于根据关键土壤参数（包括Dr、围压p、杨氏模量E和峰值内摩擦角?peak）预测CPT响应。

**部分内容摘录**
**校准室模型**
耦合校准室模型和CPT探头如图1所示。该校准室为圆柱形样品，直径为Dc，高度为356毫米（见图1(a)）。它包含一个直径为Dd的中央颗粒域和一个尺寸为Df的外围区域域。中央域使用PFC3D（Itasca Consulting Group, 2015）中的DEM进行模拟，包含大约50,000个球形颗粒。外围域使用FLAC3D（Itasca）中的FDM进行建模。

**校准室水平尺寸的影响**
图3展示了Rdf对测试序列G-1中模拟结果的影响。对于每次测试，从2Dp贯入深度到贯入结束所测得的锥端阻力qc和靴体摩擦力fs用于计算平均锥端阻力qcs和平均靴体摩擦力fss。在松散砂中（见图3(a)），当Rdf为0.43时，qcs达到3.2 MPa，高于其他Rdf值。随着Rdf的增加，qcs稳定在约2.1 MPa。如图3(b)所示，在中等密实度砂中……

**腔体膨胀法预测qcs**
如公式（4）-（8）所示，腔体膨胀法可用于建立径向膨胀应力σr与qcs之间的关系（Salgado et al., 1997; Cudmani and Osinov, 2001）。在本研究中，采用了一个初始应力状态为k0的圆柱形腔体膨胀解，以考虑侧向土压力的影响。使用的非塑性本构模型是Wolffersdorff（1996）提出的，该模型基于Gudehus（1996）和Bauer的公式……

**锥形贯入试验的预测模型**
由于Young氏模量E、峰值内摩擦角?peak、相对密度Dr和围压p等岩土工程参数的相互作用，预测锥形贯入试验中的土壤响应（如qc、fss和Rf）面临挑战。图20显示了qcs、fss和Rf与物理参数（E、?peak、Dr和p）的关系图，并给出了皮尔逊相关系数ρ。从视觉上看，qcs、fss和Rf的值在统计上与这些参数相关……

**结论**
本研究采用耦合离散元方法（DEM）-有限差分方法（FDM）技术，在虚拟校准室中模拟了锥形贯入试验（CPT），研究了室尺寸和围压p对锥端阻力qc和靴体摩擦力fs的影响。通过将模拟结果与腔体膨胀方法的预测结果进行比较，评估了该耦合方法的有效性。此外，还开发了一种人工神经网络（ANN）模型……

**作者贡献声明**
Yue Song：撰写——原始草稿、可视化、方法论、调查、形式分析、概念化。Xiaoqiang Gu：撰写——审稿与编辑、监督、方法论、资金获取、概念化。Jing Hu：撰写——审稿与编辑、监督、方法论。

**利益冲突声明**
作者声明他们没有已知的可能会影响本文所述工作的财务利益或个人关系。

**致谢**
本项工作得到了国家自然科学基金（项目编号42361164615）和中央高校基本科研业务费的支持。对此表示衷心感谢。

**作者信息**
Yue Song | Xiaoqiang Gu | Jing Hu
中国上海同济大学岩土工程系，200092