推土应用中的土壤变形与力的耦合预测模型

《Journal of Terramechanics》：A coupled soil deformation and force prediction model for bulldozing applications

　　

## 摘要

本研究提出了一种耦合建模方法，将连续土壤切割过程中的土壤变形预测与推土力估算相结合。本研究开发了两种土壤质量转移模型：一个半经验模型和一个基于离散元法（DEM）仿真数据训练的神经网络模型。半经验模型提供了一种具有物理可解释性的公式，其中土壤质量转移仅取决于切割深度和推土距离。相比之下，神经网络通过引入五个输入特征来捕捉复杂的变形行为。基于预测的质量转移，可以估算不同位置的土壤变形。半经验模型适用于内聚力土壤和非内聚力土壤，而当前的神经网络模型作为初步测试仅限于非内聚力土壤。本研究深入探讨了影响土壤质量转移的多个因素，包括内聚力、推土速度、内摩擦等。本文提出了一种耦合方法，将土壤质量转移预测模型与半经验推土力模型相结合。与DEM仿真结果的验证表明，两种模型均能准确预测土壤质量转移和表面变形。该耦合框架在广泛的切割条件和不同重力水平下，均能提供准确且一致的推土力预测。

### 1. 背景

越野车辆的性能评估至关重要，因为其工作环境始终复杂多变。已有各种类型的模型用于预测车辆-土壤相互作用，从半经验模型到基于离散元法（DEM）或有限元法（FEM）的表达方式不等。虽然DEM和FEM模型能提供准确的结果，但计算成本较高；特别是DEM计算密集，而FEM相对更高效。另一方面，半经验模型简单快速，但缺乏土壤变形效应的纳入。因此，对于需要快速且精确模拟的场景，如果半经验模型能有效地纳入土壤变形模型，则将是理想选择。为沿此方向提高预测能力，德国宇航中心（DLR）开发了土壤接触模型（SCM），引入了一种可与改进的半经验公式配合使用的土壤变形模型（Krenn和Hirzinger, 2009）。该模型使用高度图表示土壤几何形状，并动态更新该图以反映土壤变形。然而，原始模型和后续模型中的土壤变形算法均缺乏坚实的理论基础（Krenn和Hirzinger, 2009, Buse等, 2016）。

为提高SCM类模型中土壤变形建模的可靠性，本研究探讨了由车辆的切向运动（如转向和整平）引起的土壤质量转移。这些运动可视为连续的土壤切割（推土）动作，导致切割工具前方形成土堆，例如车轮或挖掘机铲斗。在本研究中，使用SCM类二维高度图表示土壤表面。土壤变形通过估算地图上跨越不同区域的质量转移来表征，在土壤密度均匀的假设下，通过更新高度场来计算形状变化。此外，本研究观察到并随后证明，准确的土壤切割力预测需要可靠的土壤质量转移模型。相应地，本研究提出了一种新模型来表征切向车辆运动在土壤网格中引起的质量转移和土壤变形。为此，本研究引入了两种不同的方法。

第一种方法基于太沙基楔体理论（Terzaghi's wedge theory）。在该方法中，将对数螺旋滑移线近似为直线以简化力的求解。因此，剪切破坏区被转化为三角形。与基于神经网络的方法相比，该模型采用显著简化的形式，并具有更强的可解释性。由于其半经验性质（无需训练数据），它在不同类型的土壤中具有更广泛的适用性。然而，作为权衡，它无法捕捉切割过程中某些土壤变形特征。

作为模型开发的一部分，选择了若干参数来考察其对土壤质量转移的影响，包括内聚力、内摩擦、推土速度和密度。

第二种方法基于机器学习神经网络。使用'EDEM'软件生成的DEM仿真结果用于神经网络的训练和验证。用于训练的数据集包含100,000个样本，本研究提出了一种浅层多层感知器（MLP）神经网络来执行土壤质量预测任务。神经网络的详细架构和超参数选择将在本文中呈现。

通过输入表征土壤切割过程的当前状态，神经网络可估算转移到各个网格位置的土壤质量。在土壤密度空间均匀的假设下，可以简单地计算任何目标节点的高度变化。值得注意的是，该神经网络目前仅考虑非内聚力土壤作为初步起点，因为数据收集工作量呈乘积式增长，需要大量时间。将模型扩展到内聚力土壤将是未来工作有趣的方向。

现有的推土力模型，如太沙基楔体理论和Hegedus推土力模型（Terzaghi, 1941, Hegedus, 1960），未考虑土壤切割过程中的土壤变形或质量转移，这可能导致不准确的力预测。本研究提出了一种耦合方法，将所提出的土壤质量转移预测模型与半经验推土力模型相结合。该耦合方法实现了在一系列切割角度和深度下准确估算推土力。此外，耦合方法的固有特性使其在理论上适用于不同重力环境下的推土力预测。因此，本研究还考察了不同重力水平下的推土力。耦合方法在不同重力条件下的性能可作为交叉验证形式。如果该方法在不同重力条件下表现良好，则表明耦合方法、土壤变形模型和半经验推土力模型各自均适用于不同重力水平。此外，此类性能本身也是对该耦合策略的一种概念验证。

本研究中使用的所有离散元法（DEM）实验均使用商业DEM仿真软件'EDEM'进行。地形力学实验由于其实验环境的复杂性，往往难以实施并提取准确结果。准备特定的实验条件，如火星土壤或减重力条件，可能特别具有挑战性。因此，DEM虚拟实验已被广泛采用作为物理实验的替代方案（Karpman等, 2020, Smith和Peng, 2013, Johnson等, 2015）。遵循此先例，本研究也采用DEM仿真用于数据集生成和模型验证。

总之，本研究的主要目的是提出土壤质量转移预测模型，并提出一种将太沙基楔体理论与所提出的土壤质量转移模型耦合以预测推土力的方法。此外，还介绍了与两种质量转移模型兼容的土壤变形算法。

### 2. 离散元法实验设置

DEM虚拟实验用于模型验证和生成神经网络训练数据。

虚拟实验装置包括一个土壤箱和一个宽度为0.08米的矩形刀片，如图1所示。土壤流由放置在工具两侧的两个挡板在横向方向（垂直于推土方向）上进行约束。切割刀片长度从0.2米到0.3米不等（取决于切割深度），厚度为0.001米。在虚拟实验中，工具与侧挡板之间的摩擦设为零。此外，颗粒与工具/切割刀片材料之间的摩擦系数 $f_t$ 也设为零。本研究中考察了这一摩擦效应，结果在"概念证明.2：参数化研究"一节中呈现。结果表明，土壤质量转移几乎不受此效应影响。

由于测试旨在模拟土壤变形区（楔体）的剖面，箱体（两个侧挡板）与土壤之间的摩擦 $f_b$，代表作用于楔体两个侧面的摩擦，设为零。此摩擦效应仅在全面规模的土壤变形中才有意义。本研究仍对此进行了考察，"概念证明.2"中的结果表明土壤质量转移几乎不受影响。因此，为简便起见将其设为零。

图中不同颜色表示不同的颗粒速度，红色、绿色和蓝色分别表示较高、中等和较低速度。以红色三角形突出显示了一个由绿色颗粒组成的明显太沙基楔体（剪切区）。

共使用60,000个颗粒，颗粒密度为3500 kg/m3，选择该值以再现2000 kg/m3的干容重，颗粒半径遵循正态分布，范围为0.003米到0.004米。恢复系数设为0.3，泊松比为0.3，颗粒杨氏模量为2.5×10? Pa。颗粒接触采用Hertz–Mindlin接触模型建模，内聚力通过Johnson–Kendall–Roberts（JKR）模型引入，通过调整颗粒表面能实现（Johnson等, 1971, Mindlin, 1949）。更高的表面能对应更大的土壤内聚力。在本研究中，表面能主要在0 J/m2到5 J/m2范围内变化，以生成不同水平的内聚力进行实验，而在概念证明案例中使用了一些更极端的值。为简便起见，本文其余部分将单位J/m2缩写为J。

尽管杨氏模量和颗粒形状等其他参数也与土壤内摩擦有关，但摩擦系数与其直接相关。在本研究中，内摩擦角主要通过调整颗粒摩擦系数来控制；因此，在大多数情况下，内摩擦效应使用不同的颗粒摩擦系数来描述。本研究中使用的摩擦系数范围为0.1到1.5，最高值1.5对应42.4°的内摩擦角。神经网络训练数据覆盖的系数范围为0.2到1.0，而验证案例在两端范围之外进行了扩展。

本研究中关注的土壤特性是内摩擦角和内聚力；它们分别通过最直接的接触参数——颗粒间滑动摩擦系数和JKR表面能——来控制。其余接触参数（杨氏模量、泊松比、恢复系数、滚动摩擦）也会影响土壤响应，但无需调节它们。固定这些参数保留了每个校准参数与其目标特性之间的清晰映射，避免了同时调节可能引入的混淆相互作用。因此，这些参数保持在EDEM提供的参考值。颗粒杨氏模量2.5×10? Pa对应约10? Pa的剪切模量，这是一个常用的降低值，可以在不影响整体响应的情况下允许更大的稳定时间步长（Hlosta等, 2020）。泊松比设为0.3，这是粒状材料的标准值，其对整体响应的影响可忽略不计（Hlosta等, 2020）。3–4毫米的颗粒半径在计算成本和物理保真度之间取得了平衡，与其他DEM地形力学仿真中使用的数值相当（Karpman等, 2020, Smith和Peng, 2013）。

仿真过程中，刀片以恒定速度向前推进并保持恒定的切割深度。'EDEM'软件中的质量传感器功能可用于监测特定感兴趣区域内的土壤质量。如图1所示，由黑色方框表示的质量传感器测量地表面上方的总土壤质量。该测量用于支持概念证明和半经验土壤质量转移预测模型的验证。

下载：下载高分辨率图片（466KB）
下载：下载全尺寸图片
图1.

图1展示了DEM仿真设置示意图。（关于图中颜色标注的解释，读者请参见本文的网络版本。）用于神经网络的测量方法与半经验方法有显著差异，将在后续章节单独说明。3. 方法3.1. 半经验方法如背景部分所述，在Terzaghi楔体理论中，滑移线被假定为直线，剪切破坏区被表示为楔形，如图2所示。角度α对应土壤切割角，角度β代表楔角，随土壤类型不同而变化。在该理论中，土壤变形被假定为沿工具宽度方向（图2中的y轴）完全均匀，这一实用假设得到了实验验证的支持。该方法也常用于其他推土力模型中，力通常按单位宽度计算。遵循这一惯例，本研究也采用单位宽度公式，土壤质量传输模型也相应地基于单位宽度建立。本研究提出了一种受Terzaghi楔体理论启发的方法，用于预测连续推土过程中的土壤质量传输。在该方法中，土壤剪切被视为主要的变形模式。下载：下载高分辨率图片（140KB）下载：下载全尺寸图片图2. Terzaghi楔体理论示意图（McKyes, 1985）。下载：下载高分辨率图片（199KB）下载：下载全尺寸图片图3. 半经验方法示意图。（关于图中颜色标注的解释，读者请参见本文的网络版本。）本研究提出的推土变形过程可概括如下：•当土壤楔体（剪切破坏区）被工具推动时，它沿滑移线发生剪切，同时沿刀片表面滑动，如图3(a)所示。虚线勾勒的三角形代表楔体的原始位置，而橙色楔体表示其在工具推动一个无穷小时间步后的位置。•在每个后续无穷小时间步开始时，一个新的土壤楔体被工具推动，如图3(b)所示。红色线条勾勒的三角形代表新形成的楔体。•在每个无穷小时间步期间，重复前述两步所述的剪切过程，被剪切楔体的上层土壤被推到地面表面。随着剪切持续发生，连续楔体上层的反复积累导致土壤堆的形成。推论如下：由于土壤层在过程中保持相同的长度，当切割角为α=90°时，地面上积累的土壤预计呈近似矩形。如果刀片不是垂直切割的，表面几何形状类似于平行四边形，假设土壤具有非常大的安息角。在一般情况中考虑侵蚀时，土壤自然向下流动并积累成丘状形态，如图1所示。本方法受楔体理论启发，描述了推土过程背后的一种可能的土壤变形机制。概念验证研究见概念验证部分，用于提供初步验证。3.2. 概念验证部分一从累积土壤形状的视角进行了验证。它既验证了形状相关的推论，也验证了基于该推论所提出的方法的概念。通过调节土壤颗粒的表面能，可以调节土壤的凝聚力，从而改变安息角。图4所示案例中所用的土壤样本具有较高的表面能，为30 J，因此具有较大的安息角。结果是安息角违例极少，从而有效地消除了自然土壤向下流动对推土后土壤剖面形态的影响。最终形成的形状仅代表由推土过程本身塑造的未扰动剖面。如图所示，在整个虚拟实验过程中，地面上累积的土壤保持矩形形状，与之前的推论一致。图中观察到的矩形剖面进一步支持了每一个无穷小时间步都剪切新楔体并将其顶层沉积在地面表面的假设。图4中的土壤剖面是用高凝聚力土壤获得的。随着凝聚力降低，累积土壤开始形成丘状堆，因为土壤颗粒因安息角违例而向下流动。图5和图6所示的情况分别用5 J和0 J的表面能值进行模拟。随着凝聚力降低，所得土壤堆的坡度变平缓。这一现象很重要，因为它导致土壤的自然流动（非直接由推土作用引起），并影响整体土壤几何形状。这些效应在后续章节中进行了详细讨论。下载：下载高分辨率图片（569KB）下载：下载全尺寸图片图4. 累积土壤的矩形剖面。在推土过程中，土壤剪切被视为主要机制，压实作用被认为作用较小。因此，每个无穷小时间步的增量土壤体积等于工具刮取的土壤体积。这意味着每个无穷小时间步地面上土壤质量的增量等于被推土工具"切割"的土壤量。单位宽度下每个时间步的增量体积可估算为切割深度与该时间步推土距离的乘积。换句话说，如果假设切割深度恒定且土壤密度均匀，则单位宽度的总质量传输与总推土距离成正比。下载：下载高分辨率图片（439KB）下载：下载全尺寸图片图5. 表面能为5 J的累积土壤剖面。下载：下载高分辨率图片（372KB）下载：下载全尺寸图片图6. 无表面能的累积土壤剖面。这一等价关系构成了土壤质量传输预测模型的基础。该等价关系的推导及预测模型将在下节提供。3.2.1. 半经验土壤质量传输预测模型及其推导如图7所示，在每个时间步中，地面上方的红色矩形用于表示在仿真中一个时间步内传输到表面的土壤质量。在图中，H表示切割深度，D=vΔt是工具在当前时间步的推土距离，a是楔体上边的长度，h是楔体在当前时间步的高度变化。由于剪切是主导机制，被刮取或移除的土壤体积应等于累积在地面表面的体积。换言之，以下关系成立：H·D=a·h。以下用图7所示的具体土壤变形实例简要推导该关系。对于图中所示的实例，工具处于垂直切割姿态（α=90°）。当前时间步单位切割宽度下的地表土壤体积增量可计算为：ΔVp=h·a。切割工具在当前时间步移除的土壤体积为：ΔVr=D·H。同时有：h=D·tan(β)，a=H/tan(β)。将这两个等式代入ΔVp：ΔVp=h·a=D·H。因此得出结论：ΔVp=ΔVr。如果不考虑膨胀，则：Vp=Vr，其中Vp和Vr分别代表位移到表面的土壤总体积和工具刮取（移除）的土壤总体积。体积变化可能源于地表上方土壤约束力的丧失。例如，当累积土壤堆超过安息角时，土壤自然向下流动，这可能导致空隙比和总体积增加。因此，总体积的等价性并非严格成立。通过为正向切割情况（α<90°）和负向切割情况（α>90°）进行类似的推导，观察到ΔVp与ΔVr之间相同的等量关系。下载：下载高分辨率图片（172KB）下载：下载全尺寸图片图7. 半经验方法示意图。就土壤质量而言，可以得到以下关系：Δmp=ρ·ΔVp=ρ·ΔVr=Δmr。这表明每个时间步土壤堆的质量增量等于工具刮取的质量。自然地，对于地面总质量和工具移除的质量，同样的等价关系也适用：mp=mr。同时，从方程可以看出，地表累积土壤质量取决于当前时间步的切割深度H和推土距离。半经验土壤质量传输预测模型可表示如下：•地表的土壤质量或体积累积速率与切割工具移除速率一致。总累积量和总移除量的等价关系同样成立。•地表土壤的增量体积可计算为：ΔVp=ΔVr=v·Δt·H，其中v是推土速度，Δt是时间步长。地面上土壤累积高度的增量可计算为：h=v·Δt·H/a。值得注意的是，在该半经验方法中，楔体上边的长度a被假设为常数。•假设土壤密度均匀，总土壤质量传输（累积在地表）可计算为：mp=mr=ρ·Lb·H，其中Lb代表总推土距离。如果推土过程中切割深度保持恒定，则总土壤质量传输与推土距离成线性关系。•在该模型中，沉积在地表的增量和总土壤质量仅取决于切割深度和推土距离。根据实验观察，忽略了土壤参数、切割角、推土速度、工具-土壤摩擦力和楔体几何形状的影响。这一假设将在下节中结合实验结果进行讨论和论证。为了验证被移除土壤质量与累积土壤质量之间的理论等价性，并评估累积土壤质量与土壤类型、推土速度等因素的无关性，下一节（参数研究）提供了概念验证。3.3. 概念验证部分二：参数研究DEM实验设置与图1所示相同。3.3.1. 土壤参数通过改变摩擦系数和表面能来区分土壤类型，它们与土壤的内摩擦角和凝聚力密切相关。图8所示的图线呈现了使用不同土壤样本进行土壤切割过程中收集到的累积土壤质量数据，其内摩擦角范围为33.2°至40.5°（无凝聚力）。这四个案例中使用的摩擦系数范围为0.3至1.2。在后续章节中，粒子摩擦系数用于代表内摩擦效应，因为在本研究中它们直接相关。内摩擦角是通过常规的虚拟剪切盒试验确定的，使用两个16 cm × 16 cm × 8 cm的盒子。x轴代表推土距离，y轴表示地面上累积的土壤质量。总土壤质量传输的预测使用提出的半经验模型：mp=ρ·Vp=ρ·Vr=ρ·Lb·H。如图所示，土壤质量传输随内摩擦角增大呈轻微上升趋势，但增长幅度非常小。总土壤质量传输仍与推土距离呈线性关系，且不同内摩擦角下相应的斜率几乎不变。半经验模型的预测很好地捕捉了DEM数据。图9中的图线展示了使用表面能值从1 J到5 J不等的土壤样本获得的累积土壤质量。所有案例中摩擦系数固定为0.5。图中所示的凝聚力值同样通过虚拟剪切盒试验确定，范围为766 Pa至3010 Pa。与无凝聚力案例类似，所有曲线均与推土距离呈线性关系。

尽管总土壤质量迁移量随凝聚力增大而增加，但增幅很小，且未观察到（推土距离与总土壤质量迁移量之间的线性关系）斜率有任何明显变化。半经验模型略微低估了数据，但总体吻合度仍然良好。

下载：下载高分辨率图像（170KB）
下载：下载全尺寸图像
图8. 不同摩擦条件下迁移至地表的总土壤质量

下载：下载高分辨率图像（207KB）
下载：下载全尺寸图像
图9. 不同凝聚力条件下迁移至地表的总土壤质量

在密度方面，其对变形行为的影响也很小。如图10所示，所采用的案例使用了比例为1、2和4的密度。相应的质量数据按此比例进行了归一化处理（例如，以2000 kg/m3为参考密度，则4000 kg/m3获得的质量除以2）。如图所示，半经验模型的预测结果与离散元（DEM）模拟结果高度吻合，且不同密度下的归一化质量数据几乎完全一致。

下载：下载高分辨率图像（178KB）
下载：下载全尺寸图像
图10. 不同土壤密度下迁移至地表的归一化土壤质量

下载：下载高分辨率图像（162KB）
下载：下载全尺寸图像
图11. 不同切削角度下迁移至地表的累积土壤质量

3.3.2. 切削角度

图11中的曲线是使用60°至120°范围内的切削角度、针对不同土壤类型获得的。如图所示，土壤质量迁移几乎不受切削角度的影响。

在黏性土和非黏性土两种情况下的预测值相同，因为其仅取决于推土距离和切削深度。这表明，使用所有十种土壤类型所获得的土壤质量几乎完全一致，无论摩擦系数或凝聚力如何，从而验证了地表上的增量和累积土壤质量与土壤参数无关。此外，如图11所示，切削角度对累积土壤质量没有影响，这与累积土壤质量独立于切削角度的结论一致。

3.3.3. 推土速度

推土速度也是可能影响土壤质量迁移的参数之一。本研究系统地研究了切削速度对土壤质量迁移和推土力的影响，涵盖了广泛的工况范围。离散元（DEM）模拟在0.05 m/s至1.0 m/s的切削速度范围内，对不同非黏性土和黏性土试样进行了仿真。

如图12和图13所示，对于具有两种不同内摩擦角（颗粒摩擦系数分别为0.4和0.6）的非黏性土，推土速度的变化不会影响总土壤质量迁移的线性趋势或量级；两种情况下的结果几乎完全一致。

下载：下载高分辨率图像（169KB）
下载：下载全尺寸图像
图12. 不同推土速度下的总土壤质量迁移——0.4摩擦系数

下载：下载高分辨率图像（136KB）
下载：下载全尺寸图像
图13. 不同推土速度下的总土壤质量迁移——0.6摩擦系数

对具有不同凝聚力的土壤也进行了类似测试，因为凝聚力可能会放大速度效应。结果如图14、图15和图16所示。对于三种不同凝聚力水平（表面能分别为1 J、3 J和5 J，摩擦系数为0.5）的土壤，推土速度的变化不会改变总土壤质量迁移的趋势或量级；所有情况下的测试结果同样几乎完全一致。

下载：下载高分辨率图像（115KB）
下载：下载全尺寸图像
图14. 不同推土速度下的总土壤质量迁移——1 J凝聚力

下载：下载高分辨率图像（143KB）
下载：下载全尺寸图像
图15. 不同推土速度下的总土壤质量迁移——3 J凝聚力

下载：下载高分辨率图像（149KB）
下载：下载全尺寸图像
图16. 不同推土速度下的总土壤质量迁移——5 J凝聚力

综上所述，在不同推土速度下，所有考虑的土壤类型的土壤质量迁移均保持一致，表明推土速度对推土过程中土壤质量迁移的影响可以忽略不计。

3.3.4. 补充测试：推土速度对推土力的影响

本补充测试旨在确定推土速度是否对推土力产生不可忽略的影响，因为本研究中后文提出的耦合方法旨在用于推土力预测。

与之前的测试集类似，在四种不同土壤条件下、以不同推土速度测得的推土力结果如下所示：两组零凝聚力但颗粒摩擦系数不同的测试集，以及两组颗粒摩擦系数为0.5但凝聚力不同的测试集。

如图17至图20所示，推土速度的增加导致测得的推土力波动更加明显。然而，尽管较高速度的力振荡水平更高，推土力随推土距离的整体线性趋势仍保持不变。所有测试速度下的力斜率和平均力量级均一致，表明推土速度主要影响力的稳定性，而非平均力量级或其增长速率。

下载：下载高分辨率图像（313KB）
下载：下载全尺寸图像
图17. 不同推土速度下的推土力——0.4颗粒摩擦系数

下载：下载高分辨率图像（269KB）
下载：下载全尺寸图像
图18. 不同推土速度下的推土力——0.6颗粒摩擦系数

下载：下载高分辨率图像（259KB）
下载：下载全尺寸图像
图19. 不同推土速度下的推土力——1 J凝聚力

下载：下载高分辨率图像（266KB）
下载：下载全尺寸图像
图20. 不同推土速度下的推土力——3 J凝聚力

下载：下载高分辨率图像（164KB）
下载：下载全尺寸图像
图21. 不同侧摩擦效应下的总土壤质量迁移

下载：下载高分辨率图像（143KB）
下载：下载全尺寸图像
图22. 不同工具-土壤摩擦效应下的总土壤质量迁移

3.3.5. 试楔侧摩擦和工具-土壤摩擦的影响

图21和图22展示了在不同侧摩擦系数和工具-土壤摩擦系数条件下的总土壤质量迁移。如图所示，楔体的侧摩擦和工具-土壤摩擦对总土壤质量迁移的影响均可忽略不计。控制土壤质量迁移的主要因素是切削深度和土壤失效机制的几何形状。根据本研究中提出的土壤运动机制，工具-土壤摩擦不应影响总土壤质量迁移。侧摩擦同样不影响土壤质量迁移的根本原因可能需要进一步研究。

3.3.6. 测试总结

综上所述，在所有案例中，预测的（刮除/移除）土壤质量与地表的总累积土壤质量相匹配，从而验证了移除土壤质量 $m_r$ 与总（累积）土壤质量 $m_p$ 之间的等价关系。

对于推土行为由楔体理论控制的材料，内摩擦角、密度、土壤切削角度和推土速度对所提方法的影响极小。凝聚力的影响更为显著，但仍不足以被视为土壤质量迁移的主要贡献因素。总质量迁移主要由推土距离和切削深度控制。

楔体侧摩擦和工具-土壤摩擦对总土壤质量迁移的影响同样发现可以忽略不计。

关于推土力，在切削深度恒定的情况下，推土力随推土距离线性增加。较高的推土速度主要引入更大的力测量波动，而整体力趋势和平均量级基本不变。

总体而言，验证结果符合预期。

3.4. 土壤变形算法——半经验法

如前所述，使用带有正方形网格的土壤高度图来表示土壤，如图23所示。对应于单个土柱的网格节点在图中标记为 $N_1$ 至 $N_n$。本节中的节点表示相对于刀片切削点的位置，以简化方法展示。然而，在实际模拟中，随着刀片的推进，楔体范围内的土壤网格节点在每个时间步内会重新对齐。土壤变形的几何演化通过土壤网格节点高度的变化来表示。

在图示案例中，当土壤与工具相互作用时，楔体区域内每个土柱的高度在每个时间步根据土壤变形算法进行更新。

每个时间步的土堆高度增量 $h$ 可按前述方法计算。地表上的土壤高度增量是由楔体的向上运动引起的（即形成新楔体顶层的土壤）。因此，在每个时间步，所有土柱 $N_1$ 至 $N_n$ 的高度增量可估计为 $h = v \cdot \Delta t \cdot H_a$。

下载：下载高分辨率图像（123KB）
下载：下载全尺寸图像
图23. 土壤网格高度图示意图

如前所述，楔体沿滑动线运动，同时沿刀片向上移动。由于土壤变形区域（或楔体范围）假定为与刀片固定，这种复合运动导致在切削角度不为90°时，楔体顶层相对于土壤变形区域发生偏移。

图24中的三个图简化展示了无侵蚀情况（即无因重力引起的自然土壤流动）下三种不同切削角度场景的土壤变形模式。楔体区域表示土壤变形区域，红色矩形表示每个时间步中每个剪切楔体的顶层。

下载：下载高分辨率图像（175KB）
下载：下载全尺寸图像
图24. 土壤变形的简化示意（无侵蚀）

当切削角度恰好为90°时，楔体顶层与土壤变形区域保持对齐，不发生相对偏移，如图24(a)所示。

当切削角度大于90°时，楔体沿倾斜刀片运动，导致每个剪切楔体的顶层（由红色矩形表示）相对于变形区域（楔体区域）向前偏移，如图24(b)所示。其基本机制与垂直切削情况相同：在每个时间步中，一层土壤被新剪切楔体中的一层向上推升。唯一的区别在于切削角度。这种运动导致土堆形成平行四边形形状，其侧面与刀片表面平行。每个时间步中每层楔体的偏移距离可估计为 $d = h \cdot |\tan(\alpha)|$，其中 $h$ 为高度增量。

另一方面，当切削角度小于90°时，楔体沿前倾刀片运动，导致楔体层在每个时间步相对于变形区域向后偏移距离 $d$。这种运动形成平行四边形的土堆，如图24(c)所示。

3.4.1. 模拟流程

在本案例的实际模拟中，红色矩形层未被使用，以保持简洁性。仍然只使用二维土壤高度图。在垂直方向上没有分层。

对于垂直切削情况，如图23所示，楔体范围内（$N_1$ 至 $N_n$）所有节点的高度统一增加一个高度增量 $h$，该增量取决于给定时间步内的切削深度和推土距离。随后进行侵蚀处理，以考虑因重力引起的自然土壤流动，并消除任何局部休止角违规行为。侵蚀机制在前文已介绍，用于模拟土壤变形的详细实现将在后续章节中呈现。

对于切削角度小于90°的情况，其过程与垂直切削情况类似，不同之处在于每个时间步需要将要处理的土壤向刀片方向偏移距离 $d$。高度增量的施加方式与垂直切削情况相同。之后，高度值被重新分配给位于距刀片距离 $d$ 更近的节点（例如，图25中节点 $N_2$ 至 $N_{n+1}$ 的高度值被重新分配给节点 $N_1$ 至 $N_n$）。切削刀片上积累的土壤由节点 $N_0$ 表示，该节点不对应实际的土壤网格节点。

最后，对所有相关的土壤网格节点应用侵蚀程序，整个过程在下一个时间步重复。对于切削角大于90°的情况，如图25所示，节点高度在每个时间步理想上应向前移动距离$d$，然后按照前面讨论的方式增量$h$。然而，为了在实际应用中简化操作，本研究提出了一种更为实用的方法。在该方法中，高度增量直接应用于楔形范围内的节点，如同之前切削情况中的做法。然后，如前所述，由于楔形层沿刀片滑动，任何其高度超过刀片在其水平位置高度的节点都会被降低至匹配刀片高度。总的减少体积由图25中的绿色三角形表示。另一方面，楔形范围外的红色三角形所示的土壤区域是由于前面描述的土壤层的水平移动产生的。该区域的体积与减少量相等（共享相同的高度和底边）。在每个时间步中，这部分被重新分配到紧邻楔形范围外部并与其相邻的节点上，因为平行四边形形状的土壤堆积在实际中很难实现（除非没有侵蚀发生）。以图25中的情况为例，重新分配节点对应节点$N+1$。此步骤完成后，应用侵蚀程序，整个过程在下一个时间步重复。如图所示，侵蚀程序对所有情况都是必要的。它通过检查每个土壤网格节点处安息角违规情况来运作。安息角定义为相对于水平面，颗粒材料堆积在不发生滑动的情况下能够保持稳定的最陡角度（Al-Hashemi和Al-Amoudi，2018）。侵蚀过程将在下一节中详细描述。

下载：下载高分辨率图片（170KB）
下载：下载全尺寸图片
图25. 切削角$α>90°$情况的模拟方法。

3.4.2 侵蚀

侵蚀程序的演示如图26所示。如图所示，在两种情况下，$N$、$N_c$和$N_i$是三个相邻的土壤网格节点。在此背景下，与安息角相比的对应角度为$θ_c$和$θ_i$。这两个角度可以计算为：
$$θ_c = \arctan\left(\frac{h_c}{l_{grid}}\right)$$
$$θ_i = \arctan\left(\frac{h_i}{l_{grid}}\right)$$
其中$l_{grid}$表示网格尺寸，$h_c$和$h_i$分别是节点$N_c$和$N_i$相对于节点$N$的高度差。

在两种情况下，$h_c$是导致角度$θ_c$等于土壤安息角时的临界高度差。因此，如果$h_i$小于$h_c$，如图左侧情况所示，则不发生侵蚀。相反，如果$h_i$超过$h_c$，则触发侵蚀过程。在这种情况下，土壤从节点$N$转移到节点$N_i$，直到高度差减少到$h_c$。

由于本研究使用的是2D高度图，该程序看起来像是将高度从一个土壤网格节点转移到另一个节点。然而，其物理本质是土壤质量的转移。侵蚀程序完成后，半经验土壤变形过程对于当前时间步即被认为完成。

下载：下载高分辨率图片（96KB）
下载：下载全尺寸图片
图26. 土壤网格高度图说明。

通过半经验模型和侵蚀模型的结合，土壤变形的几何演化可以得到合理的近似。

3.4.3 半经验模型的局限性

然而，某些非线性特征可能无法被捕捉到。一个突出的例子是在推土过程中观察到的Terzaghi楔形的几何变换。如图27所示，这里展示的三种DEM模拟使用了不同的土壤性质和切削角。前两种情况在颗粒摩擦系数上有所不同，这反映在表面堆积的形状中。第三种情况使用60°的切削角而不是90°。

如前所述，图中颗粒的颜色代表其速度，并且在所有情况下都可以观察到明显的楔形形状，由红色三角形标示。在所有情况下，楔形区域随着推土过程的推进而逐渐增大，如图27所示。这种行为主要归因于表面累积土壤造成的附加荷载积累，进而增加了颗粒间摩擦效应。

下载：下载高分辨率图片（1MB）
下载：下载全尺寸图片
图27. Terzaghi楔形的形状变化。（对于本图图例中颜色的解读，读者请参考本文章的网页版本。）

如前所述，在提出的半经验方法中，累积土壤质量仅取决于切削深度和推土距离。因此，尽管楔形形状在过程中发生变化（例如通过扩大或楔形角变化），增量和累积土壤质量预测仍然有效且不变。如图8、图9、图10、图11、图12、图13、图14、图15、图16所示，在推土过程中观察到的楔形形状变化并不影响提出的模型与DEM数据之间的一致性。预测在整个过程中与基于DEM的累积土壤质量保持一致，尽管楔形几何形状在不断演变。这为所提出模型的有效性提供了额外支持。

然而，尽管楔形几何形状不影响累积土壤的总体积或质量，但其形状变化确实影响土壤变形模式。由于地面上的土壤堆积是由连续楔形层的顶层在每个时间步剪切后积累产生的，楔形大小的任何变化都会直接影响土壤变形的范围和速率。

在半经验方法中，假设楔形形状是固定的，因此无法捕捉其几何变化对土壤变形和力预测的影响。为了应对这一及其他半经验方法的局限性，并探索数据驱动方法是否能够捕捉半经验方法未能充分表示的土壤变形模式，本文提出了基于神经网络的土壤质量转移模型。神经网络的发展将在下一节中介绍。

3.5 神经网络方法

由于土壤变形的复杂性，采用多层感知器（MLP）模型来预测质量分布，因为它具有近似高维输入特征与连续输出变量之间复杂非线性关系的能力。其分层架构能够进行分层特征提取和非线性变换，使其非常适合本工作。

在本实现中，神经网络（NN）被设计为预测沿土壤切削方向各个表面位置（即土壤网格节点）累积的土壤质量。神经网络的输入，如图28所示，包括土壤颗粒摩擦系数$μ$、切削深度$H$、切削角$α$、工具-节点距离$S$以及工具行进距离$D$。值得注意的是，工具-节点距离$S$指的是沿工具运动方向测量的从楔形顶面起始点（或图中突出的工具切削点）到节点的线性距离，如图28所示。$S_1$和$S_n$分别代表节点1和节点$n$的工具（切削点）-节点距离。网络产生一个单一输出：在输入距离$S$指定的目标节点处累积的土壤质量$Δm$。完整的沉积剖面$m_1, m_2, ..., m_n$通过对每个节点$S_1, S_2, ..., S_n$评估一次网络来获得。

值得一提的是，对于切削角$α<90°$的情况，在楔形顶面左边缘之前，一些土壤位于工具表面上。这部分总土壤质量转移，如图29中橙色突出所示，被分配至位置0 cm处（相对于切削点）。对于切削角$α>90°$的情况，位置0 cm处的质量设为零，因为土壤仅堆积在楔形区域上。

如文献中所记载并在本研究中验证的，土壤的内聚力和内摩擦系数都影响楔形的大小，如图30所示。三张图片展示了在不同摩擦系数$μ$和表面能$γ$（内聚力）下测试开始时获得的楔形。左侧的参考楔形是在切削深度为10 cm、$μ=0.2$且$γ=0$的情况下获得的；另外两个楔形在相同切削深度下，分别使用$μ=0.8$、$γ=0$和$μ=0.4$、$γ=5\text{J}$获得。如图所示，楔形大小随着摩擦系数或内聚力的增加而增大，从11.2 cm增大到16.4 cm和14 cm。

下载：下载高分辨率图片（122KB）
下载：下载全尺寸图片
图28. DEM设置中模型变量的定义。

下载：下载高分辨率图片（474KB）
下载：下载全尺寸图片
图29. 切削角$α<90°$的推土情况。（对于本图图例中颜色的解读，读者请参考本文章的网页版本。）

由于神经网络方法直接预测在楔形几何形状定义的变形范围内各个位置的土壤质量转移，因此土壤性质起着关键作用。

下载：下载高分辨率图片（762KB）
下载：下载全尺寸图片
图30. 摩擦和内聚力对楔形几何形状影响的说明。

准确预测土壤质量转移是神经网络的主要目标之一；同时，该网络还被设计为调查其捕捉超出半经验模型范围的土壤变形行为的能力。将额外的土壤参数作为输入特征会增加输入空间的维度，使网络架构复杂化，降低可解释性，并且需要大得多的训练数据集。因此，为了保持模型的简洁和可解释性，仅使用颗粒摩擦系数作为输入，而内聚力在此初始阶段被排除在外，将在未来工作中加以处理。内摩擦系数未被用作输入，因为它是一个通过实验标定确定的土壤性质，因此存在额外不确定性。

另一方面，密度对楔形几何形状没有显著影响，如图A.1附录中"密度对楔形几何形状的影响"所示。因此，密度作为输入特征在模型中被排除。

用于数据挖掘和网络验证的实验与之前介绍的实验设置相同。然而，不是测量总土壤质量，而是测量在距离工具水平距离$S$不同位置处累积的土壤质量，如图31所示。

为了获得训练神经网络所需的质量数据，如图31中黑色矩形框所示的质量传感器被放置在扰动区域的各个位置，以测量切削过程中每个位置的土壤质量变化。每个框上方显示的数字代表框内包含的总颗粒质量。根据情况的不同，最多使用十个长度各为2 cm的质量传感器来覆盖土壤变形范围。

下载：下载高分辨率图片（439KB）
下载：下载全尺寸图片
图31. DEM模拟设置说明。

在本研究中，假设每个传感器内的质量均匀分布。因此，如果模拟中使用较小的网格尺寸（例如1 cm），可以通过简单地将MLP输出除以质量传感器框长度（2 cm）与用户定义网格尺寸之间的比率来获得预测值。

神经网络架构如图32所示。隐藏层数和每层神经元数将在后续部分确定。

神经网络的训练数据通过DEM模拟获取。数据集包括五种不同的切削角$α$，范围从120°到60°，每15°递增；五种切削深度$H$，范围从2 cm到10 cm，每2 cm递增；80个推土距离$D$，从0.2 cm到16 cm，每0.2 cm递增；十个节点到工具的距离$S$，范围从0 cm到18 cm，每2 cm递增；以及五种不同的土壤颗粒摩擦系数，范围从0.2到1，每0.2递增。

总共使用约100,000个数据点进行神经网络训练。神经网络的输入参数总结于表1。下载：下载高分辨率图像（275KB）下载：下载全尺寸图像图32。神经网络架构及其输入和输出。

输出表示不同位置处原始地面以上的土壤质量增加量。该质量使用长度为2 cm的质量传感器进行测量。因此，如果在模拟中使用较小的网格尺寸（例如1 cm），可通过将MLP输出除以传感器长度与网格尺寸之比来获得土壤质量增加量。

表1。神经网络的输入。

| 量 | 符号 | 范围 | 增量 |
|---|---|---|---|
| 切削深度 | ?? | 2–10 cm | 2 cm |
| 切削角度 | ?? | 60°–120° | 15° |
| 颗粒摩擦系数 | ?? | 0.2–1.0 | 0.2 |
| 工具–节点距离 | ?? | 0–18 cm | 2 cm |
| 推土距离 | ?? | 0.2–16 cm | 0.2 cm |

相同的缩放方法也适用于土壤密度。如果模拟中使用的土壤密度与训练数据中的不同，可通过将MLP输出除以训练密度与模拟密度之比来获得土壤质量增加量。

本研究中的超参数调优分两个阶段进行。在第一阶段，通过比较不同隐藏层数量和神经元数量组合的性能来选择网络架构，同时保持与训练相关的超参数固定。一旦确定了合适的架构，则对其他超参数进行调整以进一步改善性能。这种分阶段方法有助于将网络容量与训练动力学的影响分离开来，并与若干现有研究中使用的方法一致（Gorishniy等，2021；Shavitt和Segal，2018）。

神经网络模型使用MATLAB中的深度学习工具箱实现。所有网络架构和训练超参数，包括隐藏层数量、每层神经元数量、丢弃率和相关训练参数，均进行了显式调整。具体过程在以下章节中描述。

**3.5.1 超参数：神经网络的复杂度**

神经网络的复杂度是最重要的超参数之一。选择隐藏层数量和神经元数量是MLP模型设计的关键方面，因为它直接影响模型表达能力与泛化能力之间的权衡。虽然具有更多神经元的深层网络可以学习更复杂的行为且通常泛化能力更好，但它们也会增加过拟合的风险（Bengio，2012；Goodfellow等，2016）。因此，为确定隐藏层数量与每层神经元数量的最优组合，进行了网格搜索。

值得一提的是，每次搜索的训练时间相对较短；使用Intel i9-13900K处理器，每次训练运行大约需要20–30秒。

**层数**

先前的研究表明，1到4个隐藏层的MLP通常足以处理结构化回归任务。Gorishniy等人在11个具有6到136个特征、样本量从6000到超过700,000的公开数据集上评估了深度学习模型，得出结论：具有最多三个隐藏层的浅层MLP与更深的基于Transformer的架构表现相当（Gorishniy等，2021）。Shavitt和Segal（2018）研究了具有最多200个输入特征、样本量从10,000到超过1,000,000的数据集，表明具有1–3个隐藏层的正则化MLP始终取得了较强的性能，而更深的架构往往导致过拟合。尽管如此，对于没有序列依赖关系的结构化数据，浅层网络通常提供足够的容量（Goodfellow等，2016）。

因此，在本工作中，探索了具有1到3个隐藏层的MLP架构搜索范围，以平衡模型表达能力和过拟合风险，同时考虑了输入特征的结构化特性和可用数据集的大小（100,000个样本）。

**每层神经元数量**

每层隐藏层的神经元数量决定了模型捕捉非线性关系的能力。本工作中使用8到20个神经元/层的搜索范围。下限8的选择基于一种直观且由经验驱动的启发式方法，该方法建议将隐藏层大小设为输入特征数量的2到4倍作为合理起点（Heaton，2008）。由于有5个输入特征，下限设为8个神经元，以确保搜索范围足够宽。同时，上限20得到了表格数据研究先前发现的支撑，在这些研究中，超过该大小的更宽层显示出收益递减，并且在低维问题中过拟合风险增加（Gorishniy等，2021）。该范围定义了调整模型容量的搜索空间，同时在当前结构化数据集上保持表达能力与泛化能力之间的平衡。

因此，在本工作中，隐藏层数量在1到3之间变化，每层神经元数量在8到20的范围内探索，依据先前研究的发现和本研究中使用的数据集大小。

为评估每种层数–神经元组合的性能，使用均方根误差（RMSE）作为主要指标。在不同网络配置下比较训练和验证RMSE，并使用热图进行可视化，其中RMSE值以隐藏层数量和每层神经元数量为坐标进行绘制。

如前所述，在架构搜索阶段，所有与训练相关的超参数均保持固定，以便独立评估隐藏层数量和每层神经元数量的影响。数据集分为75%用于训练和25%用于验证，初始超参数设置列于表2。这些超参数及其选择值的理由在下一节中讨论。

训练和验证RMSE的热图如图33、图34所示。

如图34所示，在所有测试的架构组合中，具有1个隐藏层和12个神经元的配置取得了最低的验证RMSE值0.0096，表明泛化性能最佳。相对于总体质量范围（最高可达2.4 kg），该RMSE非常小。虽然该配置未实现最低训练RMSE，但其0.0078的值仍然较低且与其他配置相当，如图33所示，表明模型仍能有效拟合训练数据。相比之下，更深或更宽的网络有时实现了略低的训练RMSE值，但代价是验证误差增加，反映了过拟合现象。值得注意的是，随着网络深度增加，训练RMSE持续下降，而验证RMSE趋于上升，表明泛化能力降低。相反，随着神经元数量增加未观察到明显趋势，表明在所测试范围内，模型宽度对性能的影响不那么显著或更为复杂。所选配置中训练RMSE与验证RMSE之间的微小差异表明该模型达到了良好平衡，能够可靠泛化且不过拟合或欠拟合。具有2个隐藏层和14个神经元的配置达到了可比较的低验证RMSE；但单层结构因其以更低计算成本实现同等精度而更受青睐。基于此分析，选择1个隐藏层和12个神经元的架构作为最终模型。

下载：下载高分辨率图像（151KB）下载：下载全尺寸图像图33。训练RMSE。

下载：下载高分辨率图像（155KB）下载：下载全尺寸图像图34。验证RMSE。

表2。架构搜索期间使用的初始训练超参数。

| 超参数 | 值 |
|---|---|
| 优化器 | Adam |
| 学习率 | 0.001 |
| 批大小 | 64 |
| 丢弃率 | 0.1 |
| 权重衰减（L2） | 0.0001 |
| 早停 | 耐心值 = 10 |

**3.5.2 其他超参数**

本节将介绍表2中列出的超参数，并提供其选择值的理由。鉴于这些初始值已取得强大的基线性能且训练与验证RMSE之间差距相对较小，未进行广泛的超参数调优，因为这可能导致过度正则化并降低模型性能。相反，仅对一个关键正则化超参数——权重衰减——进行了有限探索，以考察其对泛化的影响。仅在讨论中解释为何仅调整该参数而保持其他参数固定。

优化器在训练过程中更新模型权重以最小化损失函数。本研究中使用了Adam优化器，因其在训练MLP时的高效性和稳定性。Adam结合了动量和自适应学习率，能够根据过去的梯度对每个参数进行独立调整（Kingma和Ba，2015）。这使其非常适合结构化数据，并减少了对大量手动调优的需求。

学习率决定了每个训练步骤中权重更新的大小，是影响优化器行为的重要组成部分，直接影响收敛速度和训练过程的稳定性。本研究中使用了0.001的值，因为它被广泛采用为Adam优化器的默认值，并且在各种深度学习应用中展现了可靠的性能（Goodfellow等，2016；Kingma和Ba，2015）。对于本研究中使用的数据集，该值被认为合适且可靠，因此保持不做进一步调优。

批大小同时影响训练稳定性和计算效率，并可通过梯度更新中的噪声水平对泛化产生适度影响。然而，对于样本量少于100,000的结构化回归任务，当使用128以下的值时，批大小对最终模型精度和泛化的影响通常较小，而更大的批大小甚至可能降低性能（Masters和Luschi，2018）。因此，本研究中将批大小固定为64，因为它提供了平衡的权衡并避免了调优过程中不必要的复杂性。

丢弃是一种通过训练期间随机停用部分神经元来改善泛化的正则化技术，从而减少协同适应并帮助模型学习更鲁棒的特征（Srivastava等，2014）。它在冗余表示更可能形成的更宽网络中尤为有效。然而，在狭窄网络中其益处减弱，因为每个神经元在表示输入–输出映射中扮演更关键的角色。本研究中每层仅包含12个神经元；因此，即使停用少量单元也可能导致欠拟合。鉴于验证RMSE已经很低，应用了较小的丢弃率'0.1'并在调优过程中保持固定。

权重衰减，也称为L2正则化，是一种常用于防止过拟合的技术，通过对较大的模型权重进行惩罚来实现。这是通过增强原始损失函数，加入与模型权重欧几里得范数平方成正比的惩罚项来实现的。所得目标函数定义为：

(1) ??_total = ??_original + ??∑??_i2

其中??_original是原始损失（例如均方误差），??_i是模型权重，??是权重衰减系数。此正则化鼓励模型学习更平滑、更不复杂的函数，从而对未见数据具有更好的泛化能力（Goodfellow等，2016）。

本研究中使用的初始权重衰减值为0.0001，这是使用Adam优化器训练模型的常用基线（Loshchilov和Hutter，2019）。为探索其对模型性能的影响，在[1e?5, 1e?4, 3e?4, 5e?4, 1e?3]范围内进行了网格搜索，涵盖弱正则化和相对强正则化。调优范围的选择基于先前研究中小型到中型多层感知器中常用的权重衰减值（Goodfellow等，2016）。

图35展示了权重衰减变化对训练和验证RMSE的影响。随着权重衰减增加，训练RMSE从0.00232略微上升至0.0028，反映了模型拟合训练数据能力的约束增加。验证RMSE最初从0.00326下降到3×10??处的0.00308，表明较小的权重衰减值允许更好的训练性能，但可能以过拟合为代价。超过该点后，验证RMSE在1×10?3处增加到0.0041，表明过度正则化导致欠拟合。由于将权重衰减从1×10??增加到3×10??所获得的改善并不显著，本研究中保留了1×10??的值。

下载：下载高分辨率图像（181KB）下载：下载全尺寸图像图35。权重衰减RMSE。

为防止过拟合并减少不必要的计算，训练期间采用了早停。如果验证RMSE在连续10个epoch内未改善，则终止训练过程。这种方法确保了稳定收敛，且无需额外调优。

总之，神经网络架构通过两阶段调优过程确定。

在第一阶段，通过网格搜索选择隐藏层的数量和每层的神经元数量，最终配置为1个隐藏层，包含12个神经元。在第二阶段，根据现有文献确定关键训练超参数。经过有针对性的网格搜索，使用了1×10?4的权重衰减，而其他参数如批量大小和丢弃率则固定不变，以避免不必要的复杂性。对神经网络进行了验证，相应结果将在后文讨论。

3.5.3. 基于神经网络的土壤变形
与半经验方法相比，土壤变形过程被简化了。通过使用神经网络预测每个节点的质量，在假设密度均匀的情况下，可以直接估计每个位置的土壤体积变化。然后，可以通过将体积变化除以网格面积来计算每个土壤网格位置的高度变化。

3.6. 耦合方法：推土力计算
3.6.1. 力计算
力的估计基于特萨吉（Terzaghi）的试楔模型。如图2所示，假设作用在楔体上的力处于平衡状态，可以通过求解以下方程来计算土壤切削力（McKyes, 1985）。
(2)∑????=???sin?(??+??)+????????cos???????sin?(??+??)??????1?cos???=0
(3)∑????=????cos?(??+??)+????????sin???????cos?(??+??)+?????1?sin???+??+??=0
由于该模型旨在捕捉单位宽度的土壤切削力，因此假设该力沿y轴均匀分布。在方程中，L和L1表示楔体的两侧。ca和c分别表示粘聚力和内聚力。Q表示超载力，而W表示楔体的重量。P表示工具与楔体之间的法向反作用力和摩擦力的总和。此外，R表示土壤与（土壤）楔体之间的法向反作用力和内摩擦力的总和。由于工具与土壤之间的摩擦和粘附力与土壤的内摩擦力和剪力相比相对较小，本研究忽略了这些力，方程可简化为：
(4)∑????=???sin?(??)????sin?(??+??)??????1?cos???=0
(5)∑????=????cos?(??)????cos?(??+??)+?????1?sin???+??+??=0
通过求解这两个平衡方程，可以确定P和R的值，并计算出作用在工具上的水平合力，即推土力。

3.6.2. 土壤变形与力计算之间的耦合
在本研究中，为了将土壤变形模型与半经验推土力模型耦合，将转移到原始地表面以上的土壤质量视为超载Q的一部分。
转移到地表面以上的土壤不会对原始楔体提供的剪切和内摩擦效应做出贡献。这是因为，与楔体不同，地面上的土壤不受任何周围材料的约束，因此在工具向前移动时，它不会像楔体的L1侧面那样产生阻力（剪切力和内摩擦力）。因此，它不被视为楔体的一部分，并且在这种耦合模型中，切削深度保持不变。
剪切到表面的土壤总重量（质量）可以使用提出的两种土壤质量转移模型进行预测。通过将这个总重量视为超载Q的一部分，因此可以估算推土力。
所提出耦合方法的验证在结果部分进行了展示。用这种方法预测的推土力与推土力神经网络预测的力进行了比较，并且都与离散元（DEM）虚拟实验的结果进行了评估。
值得注意的是，耦合方法依赖于转移的土壤总质量；因此，半经验和神经网络土壤质量转移模型都可以与这种方法结合使用。

3.6.3. 重力对耦合推土力预测的影响
按照这种耦合方法，由于土壤变形对推土力的影响是通过将堆积土壤（表面以上）的重量视为超载来捕捉的，因此预计耦合方法也能在不同的重力环境下预测推土力。如果假设楔体不受重力影响，那么理论上所提出的耦合方法可以不作修改地应用。
从另一个角度来看，变化的重力场对土壤重量的影响类似于改变密度。由于密度不影响土壤变形机制或地表面累积土壤质量的预测，预计重力的影响也是极小的。
因此，进行了离散元（DEM）推土实验，以评估耦合模型在不同重力条件下的力预测性能。重力对土壤质量转移的影响也在结果部分进行了研究和展示。在这些条件下的成功预测将证明耦合方法（包括土壤变形模型和半经验推土力模型）在一系列重力水平下的适用性。这也作为耦合策略本身的概念性交叉验证。

4. 结果
对于质量预测验证，数据是使用前面描述的相同的DEM虚拟实验装置获得的。由EDEM质量传感器直接测量的土壤质量作为参考（真实值），用于与半经验和神经网络模型的预测进行比较，其一致性通过均方根误差（RMSE）进行量化。对于神经网络，没有任何验证案例包含在训练集中，它们涵盖了插值和外推。进行了两种类型的验证：转移到表面的土壤总质量（在下面的小节中展示），以及在目标节点处累积的土壤质量（在后面的土壤变形比较中展示）。

4.1. 累积土壤质量预测
使用九种具有不同切削角和切削深度组合的案例进行验证。如表3所示，案例1至5的参数值在神经网络训练数据范围内，但未包含在训练集中，而其余四个案例的参数在训练范围之外。
值得注意的是，该神经网络被设计用于预测转移到不同位置的土壤质量。因此，土壤总质量转移是通过求和每个位置的预测质量获得的。

表3. 九个验证案例使用的参数组合。
案例 | 切削角 (°) | 切削深度 (cm) | 摩擦系数
案例 1 | 6 | 5 | 10 | 0.5
案例 2 | 9 | 5 | 10 | 0.5
案例 3 | 11 | 5 | 10 | 0.5
案例 4 | 9 | 0 | 3 | 0.9
案例 5 | 9 | 0 | 5 | 0.9
案例 6 | 4 | 5 | 10 | 0.3
案例 7 | 13 | 0 | 5 | 0.3
案例 8 | 9 | 0 | 1.5 | 1.5
案例 9 | 9 | 0 | 15 | 1.5

对于参数在神经网络训练范围内的案例，验证结果如图36所示。地表面上的总累积土壤质量与土壤特性、切削角和楔体几何形状的变化无关。具有相同切削深度的案例表现出非常相似的总质量累积，表明累积土壤质量主要取决于切削深度和推土距离。
在图中，x轴表示工具沿推土方向行驶的距离，y轴表示土壤堆的总质量。
由于具有相同切削深度的案例的预测几乎完全相同，因此每种方法仅展示一个代表性结果。这两种方法都与DEM模拟结果吻合良好。基于半经验方法的预测由于方法本身的性质呈线性，并且在初期略微高估了质量累积。神经网络预测更接近DEM模拟结果，特别是在早期阶段。
然而，半经验模型预测与虚拟实验结果之间的差异可能归因于“EDEM”中使用的质量估计方法。在本研究使用的DEM模拟中，土壤质量是根据重心落在采样范围内的颗粒数量来估计的。因此，在测试初期部分进入该区域的颗粒被排除在外，这可能导致早期阶段的土壤质量被低估。因此，在实践中很难确定哪种方法提供更准确的估计。尽管如此，两种方法都表现出非常令人满意的表现。

对于切削角和切削深度在训练范围之外的案例，验证结果如图37、图38所示。案例6和7的切削角在训练范围之外，而案例8和9的切削深度在范围之外。
如图所示，两种方法的预测仍然令人满意。然而，对于切削深度为15厘米的案例，与DEM模拟结果相比，神经网络略微高估了质量。这种偏差在切削角在训练范围之外的两个案例中更为明显。特别是，对于α=45°的案例，在过程开始时，神经网络预测明显高于DEM数据。这种差异是可以预期的，因为神经网络在应用于训练域之外的输入时，通常表现出有限的泛化能力。尽管如此，预测仍然相当准确。另一方面，半经验方法在这些超出范围的案例中给出了非常好的预测，特别是由于其解析结构和较少依赖于训练数据。

下载：下载高分辨率图片 (282KB) 下载：下载完整大小图片 图36. 预测累积土壤质量的验证（案例1至案例5）。
下载：下载高分辨率图片 (233KB) 下载：下载完整大小图片 图37. 验证案例6和案例7（切削角）。
下载：下载高分辨率图片 (201KB) 下载：下载完整大小图片 图38. 验证案例8和案例9（切削深度）。

土壤质量转移预测的RMSE值总结于表4中。总体而言，在训练范围之外，半经验模型的表现优于神经网络模型，在大多数情况下表现出更低的RMSE。网络仅在1.5厘米的案例中达到了相当的可比精度，但在其他情况下表现不佳，特别是在较大深度或斜角切削条件下。尽管如此，两种模型与真实值（归一化RMSE或百分比误差范围为2.38%至8.23%）总体吻合良好。

表4. 训练范围外案例的RMSE（土壤质量）。
条件 | 半模型RMSE (kg) | 网络模型RMSE (kg)
90° - 1.5 cm | 0.0229 | 0.0161
90° - 15 cm | 0.1579 | 0.2550
45° - 10 cm | 0.0595 | 0.1478
130° - 5 cm | 0.0286 | 0.0760

4.2. 土壤切削力预测
本节中推土力的预测基于一种耦合方法，该方法将土壤总质量转移作为超载与特萨吉（Terzaghi）楔体理论相结合。本节使用的验证案例与前一节中介绍的案例相同。
工具（切削刀片）与土壤之间的摩擦并没有明确包含在推土力验证中，因为它已经在原始的特萨吉楔体理论中被考虑在内，因此可以被轻易捕捉。内聚力也是如此。
如图2所示，定义近似滑移线的角度β是从虚拟实验初始阶段形成的楔体上测量得到的。由于该角度在实验初始阶段会在一定范围内波动，因此在该范围内选择一个代表性值以最好地拟合DEM模拟数据。
由于本节的验证侧重于耦合方法，因此使用半经验模型预测的土壤总质量转移进行验证。在训练范围内，两种方法的预测与DEM模拟结果非常吻合，而在训练范围之外，半经验模型保持了良好的预测准确性并优于神经网络；因此，选择它进行此项验证。
表3中列出的案例1-5首先用于验证，因为它们是更具代表性的案例。

表5. 训练范围内案例的RMSE（推土力）。
条件 | RMSE (N)
65° - 10 cm | 4.5
95° - 10 cm | 8.2
115° - 10 cm | 15.6
90° - 5 cm | 2.4
90° - 3 cm | 4.2

下载：下载高分辨率图片 (350KB) 下载：下载完整大小图片 图39. 预测推土力的验证（不同的切削角）。
下载：下载高分辨率图片 (343KB) 下载：下载完整大小图片 图40. 预测推土力的验证（不同的切削深度）。

耦合方法的预测性能如图39、图40所示。第一张图展示了在摩擦系数为0.5时三个不同切削角α案例的结果（案例1、2、3），而第二张图展示了在摩擦系数为0.9时不同切削深度的结果（案例4、5）。
两张图中的验证结果表明，耦合方法在连续土壤切削过程中对推土力提供了相当准确的预测。除了5厘米切削深度的案例外，所有案例都与参考结果吻合良好；在该案例中，预测梯度略有偏差。
这些普通案例的RMSE值总结于表5中。

相应的归一化均方根误差范围约为10.2%（工况1）至14.7%（工况4），对于推土力的预测而言，这一误差水平是可以接受的。为了进一步验证耦合方法，本文考虑了另外四个测试工况。前两个工况对应不同的切削深度。第一个工况采用切削深度1.5 cm、切削角度90°和摩擦系数1.5。第二个工况采用切削深度20 cm、切削角度90°和摩擦系数0.5。后两个工况则关注不同的切削角度。其中一个工况采用切削深度10 cm、切削角度45°和摩擦系数0.3。最后一个工况采用切削深度10 cm、切削角度140°和摩擦系数0.4。

前两个工况的结果如图41和图42所示。如图所示，两种方法的总体预测均与离散元模拟数据表现出合理的吻合。与之前的推土力验证工况类似，这两个工况的均方根误差百分比约为10%。

具有不同切削角度的工况的力预测验证结果如图43和图44所示。

（图41：切削深度1.5 cm时的力预测验证。图42：切削深度20 cm时的力预测验证。）

如图所示，耦合半经验方法的总体预测结果仍然令人满意，尽管估计值始终略高于离散元模拟结果。对于切削角度为45°的工况，预测值超过了离散元数据，但仍保持合理的准确性。对于切削角度为140°的工况，预测值与离散元模拟结果之间存在明显的偏差。尽管这一偏差不大，但它代表了使用耦合方法进行的预测中观察到的最大差异。这两个工况的归一化均方根误差值分别为14.2%和18.3%。这些是所有验证工况中表现最差的，但预测精度仍在可接受范围内。

总体而言，在所有工况中，耦合方法与离散元模拟结果表现出令人满意的吻合。此外，由于离散元结果是在不同摩擦系数下获得的，该模型还对多种土壤类型表现出良好的适用性。

（图43：切削角度45°时的力预测验证。图44：切削角度140°时的力预测验证。）

**4.3 不同重力效应下的土壤切削力预测**

如前所述，耦合方法可以直接应用于不同重力条件下的推土力预测，无需对其公式做任何修改。在本节中，使用在0.166g（月球）至1.5g之间六个重力等级下获得的推土力进行验证。所有模拟均采用切削角度90°和固定切削深度10 cm。

如图45所示，重力对推土力的影响清晰可见。在每个重力等级下，耦合方法预测的推土力均紧密跟随离散元数据中观察到的上升趋势。在较低重力等级（如0.1666g和0.4g）下，预测精度尤其高，预测曲线与离散元结果紧密吻合。最差的工况出现在1.5g重力等级下，均方根误差为22.01 N，相对于离散元模拟数据的归一化误差为9.37%。

（图45：不同重力效应下的推土力。图46：不同重力效应下的归一化推土力。）

另一方面，如图46所示，重力对总累积土壤质量几乎没有影响，这与预期一致。

总之，耦合方法在不同重力等级下预测推土力表现出可靠的性能。

**4.4 基于不同位置质量转移的土壤变形比较**

在本节中，将半经验方法和基于网络的方法预测的土壤质量分布与离散元模拟结果进行比较。由于土壤的几何形状由不同位置的高度决定，而高度与质量分布直接相关，因此在本分析中使用高度作为比较指标。

如图47所示，比较中总共考虑了四个不同位置。这些位置沿切削方向距工具的距离分别为：1 cm、3 cm、12 cm和16 cm（覆盖土壤变形范围）。

选取了三种具有不同摩擦系数的土壤类型进行验证，因为摩擦系数影响楔形形状，进而影响最终土堆的几何形状。所选的三个摩擦系数为：0.9（处于神经网络的训练范围内），以及0.1和1.5（处于训练范围之外）。这三个工况均采用90°切削角度。

此外，还使用两个不同切削角度的数据集进行验证：45°和125°（处于训练范围之外）。这两个工况的摩擦系数设置为0.5，以确保比较的一致性。

如图48所示，对于摩擦系数为0.9的工况，神经网络在预测所有位置的土壤高度增加方面表现出优异的性能。半经验方法的准确性较低，但仍产生令人满意的结果。此外，在切削角度为90°的工况下，节点1和节点2处的土壤累积量非常相似。

（图47：不同位置的土壤质量估计。图48：摩擦系数=0.9时的高度变化。）

对于摩擦系数在训练范围之外的工况，如图49和图50所示，两种方法的预测性能与前一工况相当。神经网络通常提供更准确的预测，而半经验方法的准确性略低，但仍产生令人满意的结果。

（图49：摩擦系数=0.1时的高度变化。图50：摩擦系数=1.5时的高度变化。）

值得注意的是，在摩擦系数为1.5和0.9的工况中，半经验方法的预测未能捕捉到节点4处的土壤累积现象。这表明仅依靠安息角违反条件不足以准确捕捉初始土壤变形范围（初始楔形）之外发生的土壤累积，因为土壤楔形的几何形状在整个过程中会不断扩大。相比之下，神经网络在捕捉这些复杂的土壤变形行为方面具有明显优势。

使用两个切削角度在训练范围之外的工况来验证土壤高度变化，如下述结果所示。对于切削角度为45°的工况，节点0处的土壤高度被纳入验证范围（即堆积在工具上的土壤）。如前所述，节点0在几何位置上位于切削点之前，用于表示切削角度小于90°时工具上的土壤累积。该工况的拟合结果如图51所示（节点0至节点2），图52所示（节点3和节点4）。

如图所示，基于神经网络预测的质量转移的高度预测开始与离散元模拟出现一些偏差。特别是，节点2处的预测不如半经验模型准确。这一偏差是可以接受的，因为45°远低于训练范围的下限（下限为60°）。然而，在所有验证工况中，这是神经网络预测出现明显偏差（在节点2处）的唯一实例。尽管如此，预测仍能合理地捕捉到这些节点处离散元曲线的趋势，整体拟合仍可接受。

（图51：切削角度=45°时的高度变化：节点0至节点2。图52：切削角度=45°时的高度变化：节点3和节点4。）

另一方面，半经验方法的预测表现出良好的一致性和令人满意的性能。尽管在某些工况中其预测不如神经网络准确，但它们从未与离散元模拟数据出现大幅偏差。特别是，节点0处的预测在半经验方法下明显更准确。值得注意的是，节点0至节点2处的土壤堆高度几乎相同，半经验模型对这三个节点的预测也完全相同。这是因为，根据半经验模型，在没有侵蚀的情况下，土壤变形范围内每个节点的高度增量应该是相同的。

切削角度为125°的工况拟合结果如图53和图54所示。

如图53所示，转移到节点1和节点2的土壤质量最初增加，随后由于切削工具形状带来的几何约束而趋于稳定。半经验模型获得的预测与离散元模拟数据吻合良好。同时，神经网络在预测这两个节点处的高度变化方面也表现出良好的一致性。

（图53：切削角度=125°时的高度变化：节点1和节点2。图54：切削角度=125°时的高度变化：节点3和节点4。）

如图54所示，神经网络的预测与离散元数据几乎完美吻合。半经验模型显示出一些明显偏差，明显不如神经网络准确，但其性能仍然令人满意且一致。

**5. 结论**

本研究提出了一种预测推土过程中土壤质量转移的半经验模型。该模型假设在每一个微小时间步内，新的土壤楔形沿Terzaghi楔形理论定义的滑移线发生剪切。土壤质量增量等于工具或刀片在每个时间步移除的土壤量，总土壤质量转移仅取决于切削深度和推土距离。如果切削深度恒定，则总质量与推土距离之间的关系是线性的。

在模型开发过程中，考察了推土速度、凝聚力、密度以及各种摩擦系数等因素的影响。研究发现，在大多数情况下，土壤参数对总质量转移的影响可以忽略不计。与其余参数相比，凝聚力具有更明显的影响，但其效应仍然有限。未来可以进一步研究这一效应。

借助本研究引入的土壤变形算法和相容的侵蚀过程，推土过程中的土壤变形几何形状也可以被确定。

然而，在半经验方法中，楔形形状被假设为固定的，因此几何变化对土壤变形和力预测的影响未被考虑。为了克服这一局限性，并研究数据驱动方法是否能够捕捉半经验模型范围之外的土壤变形模式，本文基于离散元模拟数据开发了一个神经网络，用于预测楔形区域内不同位置的土壤质量转移。

土壤变形验证分为两部分：总土壤质量转移和不同位置（节点）的土壤高度变化。对于总土壤质量，半经验模型在所有九个具有不同切削角度和切削深度的工况中均提供准确且一致的预测。相比之下，神经网络在训练范围内的工况中与离散元数据吻合更好，但在训练范围之外的工况中其准确性有所下降。

对于土壤高度变化，半经验模型的预测与离散元结果存在一些偏差。在某些情况下，特别是在距切削点较远的节点处，模型未能准确捕捉土壤累积的起始时刻。尽管如此，整体趋势仍保持一致，土壤变形预测被认为是令人满意的。

另一方面，神经网络在大多数情况下预测更为准确。

这表明神经网络能够捕捉半经验模型未能表征的土壤变形特征，例如土楔的扩大。此外，离散元（DEM）推土模拟所包含的信息比本文所利用的更为丰富；特别是识别变形区域内的弹性区、塑性区和过渡区，并将这些信息纳入神经网络，可以进一步增强预测能力，这将在未来工作中进行研究。对于推土力预测，提出了一种耦合方法，将土壤质量转移预测模型与半经验推土力模型相结合。该方法通过将转移到地表的土壤质量视为附加荷载来实现。在一般参数范围内的常规情况下，耦合方法与离散元数据吻合良好。对于超出此范围的情况，预测结果仍然合理且一致，尽管精度有所降低。与此同时，耦合方法在广泛的引力条件下表现出稳健性，保持了较高的预测精度，尤其是在较低的重力水平下。仅在1.5g时观察到明显的偏差，而在所有其他条件下的预测结果仍保持一致且可靠。除了推土场景之外，所提出的土壤质量转移模型、土壤变形公式和耦合方法还可以扩展到轮土相互作用问题。在这种情况下，可以通过横向离散化车轮，分别计算每个横向段的土壤变形和相应的切削力，来估算车轮产生的侧向推土力，同时考虑由于车轮曲面几何形状导致的切削深度变化。这些结果表明，半经验模型和神经网络土壤质量转移预测模型各有独特优势：前者提供了具有物理可解释性且可泛化的框架，适用于外推场景；而后者则在其训练域内提供高精度的预测。两者结合，为捕捉复杂的土壤变形提供了一种稳健且灵活的建模策略。耦合方法为将土壤质量转移预测模型与半经验公式相结合提供了一个框架，由于其本质上属于半经验方法，因此能在广泛条件下实现准确预测。总体而言，所有提出的模型都通过准确且一致地预测土壤质量转移、土壤变形和推土力，提高了推土模拟的保真度。CRediT作者贡献声明：Yang Jiao：撰写—审阅与编辑，撰写—原始草稿，验证，方法论，调查研究，形式分析，数据整理。Jozsef K?vecses：监督指导。Marek Teichmann：监督指导。Yang Jiao|Jozsef K?vecses|Marek Teichmann