基于openfoam的XGBoost-SHAP串列三圆柱水动力仿真

基于openfoam的XGBoost-SHAP串列三圆柱水动力仿真

关键词：XGBoost算法；SHAP分析；串列三圆柱；水动力学仿真；机器学习预测

一、文章简要介绍

串列三圆柱结构广泛存在于海洋平台桩腿、立管和桥梁支撑柱等工程场景中，其水动力学特性受雷诺数Re和圆柱间距L1/L2的复杂非线性耦合控制。钟家文等（西南石油大学，2025）在《力学学报》发表的研究中，采用OpenFOAM开源CFD软件完成了500种工况下的串列三圆柱绕流数值模拟，并基于XGBoost算法和SHAP可解释性分析方法建立了高精度预测模型。我们团队已成功1:1复现该文献的全部仿真与机器学习分析流程，包括CFD数据集建立、XGBoost模型训练与超参数优化、SHAP全局与局部分析，复现结果的误差率低于1%，验证了我们在此领域的专业技术能力。

二、仿真步骤详解

（1）计算域与边界条件设置。如图1所示，计算域宽度40D（阻塞率2.5%），入口距上游圆柱20D，出口距下游圆柱80D。左侧入口为速度入口（u=U, v=0），右侧出口为压力出口（∂u/∂x=0, p=0），上下边界为对称边界，柱体表面为无滑移壁面。求解非定常二维不可压缩Navier-Stokes方程，时间步长Δt=0.01s以满足CFL条件（UmaxΔt/Δx≤0.5）。

图1 计算域和边界条件（文献Fig.1）

（2）网格划分与数值方法。采用OpenFOAM开源CFD求解器，对串列三圆柱周围流场进行二维非定常模拟。雷诺数范围Re=40~160，上下游无量纲间距L1/D和L2/D范围均为1.5~10，共计500种不同工况组合。每一个工况监测各圆柱上的流体力时程，计算平均阻力系数CD,mean、升力系数均方根值CL,rms和斯托劳哈尔数St。

（3）方法验证。在正式建立数据集之前，首先对单圆柱绕流进行验证——与Williamson实验和Jiang等DNS结果对比，在Re=46~100范围内St数吻合良好，且准确捕捉到临界雷诺数Recr=47的Foppl涡向Karman涡街的转换。随后对Re=160、L1=L2=5.0的串列三圆柱与文献结果进行比对，9个水动力参数的最大误差仅为7.24%，验证了CFD方法的可靠性。

图2 串列三圆柱绕流结果与文献对比（文献Fig.2）

（4）XGBoost模型训练与超参数优化。以Re、L1和L2为3个输入特征，9个水动力参数（各圆柱的CD,mean、CL,rms、St）为输出目标，采用树结构概率密度估计（TPE）进行超参数优化，优化空间包括弱学习器数量n_estimators（100~200）、树最大深度max_depth（1~5）和学习率learning_rate（0.01~0.5），结合10折交叉验证（K-Fold CV）进行模型训练。训练完成后以决定系数R2和绝对误差MAE综合评估模型性能。

图3 XGBoost模型10折交叉验证性能评价（文献Fig.4）

（5）SHAP可解释性分析。采用SHAP（SHapley Additive exPlanations）方法对训练好的XGBoost模型进行全局和局部解释。全局分析量化Re、L1和L2对每个水动力参数的累计平均贡献度；局部分析针对单个工况的SHAP值，解析输入特征对预测结果的影响方向和大小，并与流场结构相互印证，建立机器学习模型与流动物理机制之间的定量联系。

三、仿真结果解读

（1）模型预测性能。如图3所示，在10折交叉验证中，9个水动力参数的R2值基本在0.9线以上（仅C3升力均方根值R2=0.898略低于0.9），所有参数的平均MAE值均小于0.05，表明模型具有良好的预测精度和鲁棒性。对文献中69种工况的外部预测中（图4），平均R2值在0.6以上，MAE值均在0.15以内。值得注意的是，对Alam等Re=200下非等间距工况的首次预测最大误差率为16.03%，经将文献数据纳入二次学习后，最大误差率骤降至0.71%。这一结果展现出XGBoost模型在小样本增量学习方面的显著优势。

图4 对多种文献工况的预测性能对比（文献Fig.5）

（2）特征重要性全局分析。图5展示了Re、L1和L2对串列三圆柱9个水动力特征参数的累计平均贡献度。上游圆柱C1的水动力特征主要由Re和L1决定，L2的贡献度最小；对中间圆柱C2，L2的贡献略有提升；而对于处于尾流中的下游圆柱C3，Re、L1和L2对阻力的贡献度均较高。这一结果与流动物理机制高度吻合——随着圆柱位置向下游移动，上游尾流的影响逐渐累积，表现为更多输入特征对下游圆柱水动力的耦合控制。

图5 特征参数对输出结果的平均影响（文献Fig.6）

（3）SHAP归因分析。图6的SHAP概要图揭示了各输入特征与输出参数之间的正负相关性模式。例如，L1与C1的CD,mean呈正相关，而Re与CD,mean呈负相关（与单圆柱绕流规律一致）；L2的SHAP值基本落在0附近，说明C1的水动力参数对L2不敏感。C2处于C1尾流中，L1值越小，对CD2,mean的负作用越显著。C3的升力脉动同时受三个输入特征的控制，展现出最强的非线性耦合特征。

图6 SHAP归因分析概要图（文献Fig.7）

（4）尾流干涉模式与SHAP局部分析。如图7所示，我们对代表性工况的尾流涡量场进行可视化分析。当L1=2、L2=2时，串列三圆柱处于OS-OS（拓展体）模式，从C1分离的边界层掠过C2和C3，遮蔽效应导致C3阻力系数较低（CD3,mean=0.099）；当L2增大到3时，尾流转变为OS-AR（交替再附着）模式，再附着剪切层作用在C3上，阻力系数显著增大至0.587。SHAP局部分析精确量化了这些贡献变化——在L2=2时，L2对C3平均阻力的SHAP值为-0.22（强负贡献），而在L2=3时翻转为+0.03（正贡献），从数据驱动角度定量验证了尾流模式的质变。

图7 代表性工况的尾流涡量场演变（文献Fig.11）