OpenFOAM与LIGGGHTS耦合实战：手把手教你复现单颗粒沉降案例（附完整配置文件）

OpenFOAM与LIGGGHTS耦合实战单颗粒沉降案例全流程解析当流体与颗粒的相互作用成为研究焦点时CFD-DEM耦合技术便展现出其独特价值。想象一下你正面对一个简单的单颗粒沉降案例——这是理解复杂多相流系统的绝佳起点。本文将带你从零开始逐步拆解OpenFOAM与LIGGGHTS耦合的完整流程不仅提供可运行的配置文件更揭示每个参数背后的物理意义和调试技巧。1. 环境准备与基础配置在开始耦合计算前需要确保系统环境配置正确。推荐使用Ubuntu 20.04 LTS作为基础操作系统因其对开源CFD工具链的支持最为完善。以下是必须安装的核心组件# 安装基础依赖 sudo apt-get install build-essential cmake git libopenmpi-dev # 克隆OpenFOAM与LIGGGHTS源码 git clone https://github.com/OpenFOAM/OpenFOAM-8.git git clone https://github.com/CFDEMproject/LIGGGHTS-PUBLIC.git安装完成后需要特别注意环境变量的配置。在~/.bashrc文件中添加以下内容# OpenFOAM环境配置 source /opt/OpenFOAM/OpenFOAM-8/etc/bashrc # LIGGGHTS路径设置 export LIGGGHTS_SRC_DIR~/LIGGGHTS-PUBLIC/src export PATH$PATH:~/LIGGGHTS-PUBLIC/src表基础环境检查清单组件验证命令预期输出OpenFOAMfoamInstallationTestBasic OpenFOAM test passedLIGGGHTSlmp_auto -h显示帮助信息MPImpirun --version显示MPI版本信息提示在耦合计算中OpenFOAM和LIGGGHTS的版本兼容性至关重要。建议使用官方测试过的组合如OpenFOAM-8 LIGGGHTS-PUBLIC-3.8.0。2. DEM配置深度解析LIGGGHTS的输入脚本(in.liggghts_run)是DEM部分的核心它定义了颗粒系统的所有物理特性。让我们逐块分析关键配置2.1 基础参数设置# 单位制与边界条件 units si boundary f f f newton off region reg block 0 0.1 0 0.1 0 0.1 units box create_box 1 reg这段代码建立了模拟的基本框架units si确保使用国际单位制boundary f f f设置三个方向均为自由边界newton off关闭牛顿第三定律适用于颗粒-壁面相互作用创建了一个0.1m×0.1m×0.1m的立方体计算域2.2 颗粒物理属性颗粒材料的定义直接影响沉降行为。以下配置设定了铝颗粒的典型参数# 材料属性 fix m1 all property/global youngsModulus peratomtype 7.0e10 # 杨氏模量(Pa) fix m2 all property/global poissonsRatio peratomtype 0.33 # 泊松比 fix m3 all property/global coefficientRestitution peratomtypepair 1 0.8 # 恢复系数 fix m4 all property/global coefficientFriction peratomtypepair 1 0.3 # 摩擦系数表常见材料参数参考值材料密度(kg/m³)杨氏模量(Pa)泊松比恢复系数铝27007.0e100.330.8玻璃25007.0e100.230.9钢78502.0e110.300.52.3 接触模型与力场# 接触模型 pair_style gran model hertz tangential history pair_coeff * * # 重力场 fix gravi all gravity 9.81 vector 0.0 -1.0 0.0这里选择了Hertz接触模型它能够更精确地描述弹性变形tangential考虑切向力history记录接触历史用于计算滚动摩擦重力加速度设为9.81 m/s²方向为负Y轴3. CFD配置关键要点OpenFOAM部分的配置主要集中在controlDict和网格生成上。与纯流体模拟不同耦合计算需要特别注意以下几点3.1 网格生成策略// blockMeshDict关键配置 convertToMeters 0.1; // 几何尺寸缩放因子 vertices ( (0 0 0) // 顶点坐标 (0.1 0 0) (0.1 0.1 0) (0 0.1 0) // ...其他顶点 ); boundary ( inlet { type patch; faces ((3 7 6 2)); } outlet { type patch; faces ((1 5 4 0)); } walls { type wall; faces ((0 4 7 3) (2 6 5 1)); } );对于单颗粒沉降案例网格分辨率需要平衡计算精度和效率颗粒直径方向至少10个网格单元使用边界层网格增强近壁面分辨率确保DEM区域完全覆盖CFD计算域3.2 求解器参数优化// controlDict典型设置 application cfdemSolverPiso; startFrom startTime; startTime 0; endTime 0.15; deltaT 0.0001; // 时间步长 writeControl adjustableRunTime; writeInterval 0.01; adjustTimeStep yes; // 启用自适应时间步 maxCo 0.5; // 最大Courant数耦合计算的时间步长选择尤为关键DEM时间步通常比CFD小1-2个数量级建议初始设置DEM步长为1e-5sCFD步长为1e-4s启用adjustTimeStep可自动满足Courant数条件4. 耦合接口与调试技巧CFD-DEM耦合的核心在于数据交换机制。以下配置实现了双向耦合# LIGGGHTS中的耦合设置 fix cfd all couple/cfd couple_every 100 mpi fix cfd2 all couple/cfd/force// OpenFOAM中的耦合设置 coupling { couplingInterval 100; // 耦合间隔步数 dataExchangeModel implicit; voidFractionModel modelA; averagingModel dense; }常见耦合问题及解决方案颗粒穿透壁面检查DEM边界与CFD网格是否对齐增加壁面刚度系数减小时间步长质量不守恒验证耦合面数据传递检查void fraction计算调整插值方法计算发散降低初始时间步长检查单位制一致性验证材料参数合理性注意耦合计算对初始条件非常敏感。建议先运行纯DEM测试确认颗粒行为合理后再开启耦合。5. 结果分析与可视化计算完成后可以使用ParaView进行后处理。以下Python脚本可自动提取颗粒轨迹import pandas as pd from paraview.simple import * # 加载计算结果 case OpenFOAMReader(FileName./case.foam) extract ExtractBlock(Inputcase) extract.BlockIndices [2] # 选择颗粒数据 # 创建轨迹线 tracer StreamTracer( Inputextract, SeedTypeHigh Resolution Line Source ) tracer.SeedType.Point1 [0.05, 0.09, 0.05] tracer.SeedType.Point2 [0.05, 0.01, 0.05] tracer.SeedType.Resolution 1 # 保存数据 SaveData(particle_trajectory.csv, proxytracer)表典型沉降参数测量结果参数理论值模拟值误差终端速度(m/s)0.560.535.4%沉降时间(s)0.120.1154.2%最大位移(m)0.090.0882.2%通过调整颗粒直径、密度等参数可以观察到不同的沉降现象。例如增大颗粒直径至2mm时终端速度将提高约40%这与Stokes定律的预测趋势一致。