参考这篇文献,我需要做一个关于熔盐液滴滴落在低温平面凝固的ansys fluent仿真,请给出从模型...

作成日: 2026年7月29日

使用モデル GPT-5.6 Thinking by Chat01

質問

关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf

参考这篇文献,我需要做一个关于熔盐液滴滴落在低温平面凝固的ansys fluent仿真,请给出从模型建立到制作仿真动画的全过程和详细操作步骤

思考

Configuring wall coupling

I’m setting the fluid–solid interface as thermally coupled and adding contact resistance where available. If direct resistance is unavailable, I’ll represent it with an equivalent thin layer, then verify the option in the installed version.

Worked for 2m 2s

回答

下面给出一套适合 “高温熔盐液滴垂直撞击低温水平平面并凝固” 的 ANSYS Fluent 全流程。建模思路沿用文献中的 VOF 气—液界面捕捉 + Solidification/Melting 焓—多孔介质凝固模型,但把文献中的椭球液态锡、球形曲面改成熔盐液滴和水平平板。文献使用二维轴对称计算域、压力基瞬态求解器、VOF、凝固/熔化模型、SIMPLEC、PRESTO!、Compressive、二阶离散和细时间步,这些设置可直接作为你的初始方案。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

以下菜单名称按 Fluent 2025 R2 表述,2023 R2/2024 R2 的位置基本一致。当前 Fluent 官方文档仍支持 VOF 与 Solidification/Melting 联合计算;非凝固相的潜热、固相线和液相线应设为零。(ANSYS Help)


一、先确定仿真层级

建议分两步做。

第一步:基础模型

采用:

  • 二维轴对称;
  • 熔盐为第二相,空气为第一相;
  • 平板直接设为恒温壁面;
  • 固定静态接触角;
  • 忽略蒸发、辐射、化学分解和自然对流;
  • 采用恒定或温度相关物性。

这个模型最容易收敛,适合先跑通液滴铺展和凝固动画。

第二步:论文或工程模型

在基础模型上增加:

  • 有厚度的固体基板;
  • 基板内部瞬态导热;
  • 熔盐—基板热接触热阻;
  • 温度相关黏度、导热系数、密度和表面张力;
  • 必要时采用动态接触角;
  • 实际固相线—液相线区间;
  • 网格及时间步无关性分析。

对于“低温平面”,如果基板在撞击后会明显升温,推荐第二种共轭传热模型;如果平板温度由强制制冷保持不变,可以直接使用恒温壁面。


二、建模前需要准备的参数

熔盐种类不同,物性差别很大。不要直接把文献中的锡物性替换成某一套通用熔盐参数。至少准备下列数据。

参数符号Fluent 中的位置
液态密度ρl(T)\rho_l(T)Materials → Density
固态密度ρs(T)\rho_s(T)可用分段函数表示
液态比热cp,l(T)c_{p,l}(T)Materials → Specific Heat
固态比热cp,s(T)c_{p,s}(T)分段函数
液态导热系数kl(T)k_l(T)Thermal Conductivity
固态导热系数ks(T)k_s(T)分段函数
液态动力黏度μl(T)\mu_l(T)Viscosity
潜热LhL_hLatent Heat
固相线温度TsolT_{\rm sol}Solidus Temperature
液相线温度TliqT_{\rm liq}Liquidus Temperature
表面张力σ(T)\sigma(T)Phase Interaction
接触角θ\thetaWall Adhesion
液滴初温T0T_0初始化 Patch
平板温度TsT_sWall 或基板底面
撞击速度v0v_0初始化 Patch
液滴直径DD几何和初始化
热接触热阻RtcR''_{tc}耦合壁面或等效薄层

对于盐混合物,应优先采用真实的固相线和液相线,而不是人为设置极窄的相变区。文献中的锡接近纯物质,为改善收敛采用了约 2 K 的平滑相变区;这一数值不应直接用于所有熔盐。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

建议先计算几个无量纲数:

We=ρv02DσWe=\frac{\rho v_0^2D}{\sigma} Re=ρv0DμRe=\frac{\rho v_0D}{\mu} Oh=μρσDOh=\frac{\mu}{\sqrt{\rho\sigma D}} Ste=cp(TliqTs)LhSte=\frac{c_p(T_{\rm liq}-T_s)}{L_h}

它们有助于判断是惯性、黏性、表面张力还是凝固效应占主导。


三、推荐的基础算例尺寸

以直径 D=2 mmD=2\ \mathrm{mm} 的球形液滴为例。

项目建议值
计算域高度5D=10 mm5D=10\ \mathrm{mm}
计算域径向宽度5D=10 mm5D=10\ \mathrm{mm}
液滴与平板初始间隙0.02D0.05D0.02D\sim0.05D
基板厚度,若建固体区1D3D1D\sim3D
液滴附近网格D/100D/200D/100\sim D/200
初始时间步1×1062.5×106 s1\times10^{-6}\sim2.5\times10^{-6}\ \mathrm{s}
每时间步迭代次数20~30
总物理时间20~100 ms,直到基本完全凝固

文献中采用 5D×5D5D\times5D 的计算域,液滴附近网格约 10 μm、相当于直径方向约 200 个单元,并采用 2.5μs2.5\,\mu s 时间步。这个分辨率可作为高精度参考,但计算量很大。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

建议先使用 D/100D/100 网格跑通,再增加到 D/150D/150D/200D/200 做无关性分析。


四、几何模型建立

4.1 坐标方向

Fluent 二维轴对称模型通常以:

  • xx:轴向,即竖直方向;
  • yy:径向;
  • y=0y=0:旋转轴。

设置:

  • 平板表面位于 x=0x=0
  • 流体域位于 x>0x>0
  • 重力沿负 xx 方向;
  • 液滴中心位于 y=0y=0

不要让二维轴对称网格进入 y<0y<0 区域。


4.2 方案 A:恒温平板

在 SpaceClaim 或 DesignModeler 中画一个矩形:

0x5D,0y5D0\leq x\leq5D,\qquad 0\leq y\leq5D

矩形边界命名:

  • y=0y=0axis
  • x=0x=0cold-wall
  • x=5Dx=5Dtop-outlet
  • y=5Dy=5Dside-outlet

整个矩形为一个流体面,不需要在几何中直接画出液滴。液滴在 Fluent 中通过 Patch 或 UDF 初始化。


4.3 方案 B:带厚度固体基板

画两个相邻矩形:

流体区:

0x5D,0y5D0\leq x\leq5D,\qquad0\leq y\leq5D

固体基板:

Hsx0,0y5D-H_s\leq x\leq0,\qquad0\leq y\leq5D

执行:

  1. Share Topology;
  2. 保证流体和固体共用 x=0x=0 的界面;
  3. 将上部命名为 fluid-domain
  4. 将下部命名为 substrate-solid
  5. 基板底部命名为 substrate-bottom
  6. 基板外侧命名为 substrate-side

这样可以模拟平板吸热和温升。


五、网格划分

推荐使用 Fluent Meshing、ANSYS Meshing 或 ICEM CFD。

5.1 网格类型

二维轴对称液滴撞击建议:

  • 结构化四边形网格最好;
  • 或 Quad Dominant 四边形主导网格;
  • 尽量减少高扭曲三角形;
  • 界面和壁面附近保持近似等尺寸网格。

文献采用 Quad/Tri 自由面网格和局部细化。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

5.2 局部细化区域

至少细化以下区域:

  1. 液滴初始位置;
  2. 液滴下落路径;
  3. 平板上方高度约 1.5D1.5D 的区域;
  4. 平板径向 02.5D0\sim2.5D 的铺展区域;
  5. 固体基板靠近流固界面的区域。

推荐:

  • 液滴及撞击区:Δ=1020 μm\Delta=10\sim20\ \mu m,以 D=2D=2 mm 为例;
  • 外部空气区:50~100 μm;
  • 增长率:不超过 1.1~1.15;
  • 基板界面附近:与流体侧近似相同尺寸。

5.3 网格质量

检查:

  • Maximum Skewness 尽量小于 0.85;
  • Orthogonal Quality 最小值最好大于 0.15;
  • 平板界面附近不得有突然放大的单元;
  • 轴线附近单元应规则。

六、启动 Fluent

Workbench 中:

  1. 双击 Setup
  2. 选择 2D
  3. 建议使用 Double Precision
  4. 使用 CPU Solver;
  5. 并行核数根据网格量选择。

当前原生 GPU VOF 求解器对多相能量模型仍存在功能限制,因此 VOF、能量方程和凝固联合计算建议使用 CPU Solver。(ANSYS Help)

进入 Fluent 后首先执行:

Mesh → Check

确认:

  • 无负体积;
  • 轴对称轴边界正确;
  • 网格尺度为 m;
  • 几何尺寸与设计一致。

如果导入单位错误,通过:

Mesh → Scale

进行修正。


七、General 设置

进入:

Setup → General

设置:

  • Solver:Pressure-Based
  • Time:Transient
  • 2D Space:Axisymmetric
  • Velocity Formulation:Absolute
  • Gravity:开启
  • gx=9.81 m/s2g_x=-9.81\ \mathrm{m/s^2}
  • gy=0g_y=0

如果你的 xx 方向向下,则重力符号相反。

建议先使用:

  • Gradient:Least Squares Cell Based

后续在 Solution Methods 中开启 Warped-Face Gradient Correction。


八、物理模型设置

8.1 能量方程

路径:

Setup → Models → Energy

选择:

On

没有开启 Energy,就无法计算熔盐冷却和凝固。


8.2 黏性模型

路径:

Setup → Models → Viscous

初始建议选择:

Laminar

对于毫米级、低至中等速度的熔盐液滴,首先采用层流模型通常更稳定。若计算得到的 ReRe 很高,并且液滴出现明显湍动或破碎,再考虑 SST kωk-\omega 等模型。


8.3 VOF 多相模型

路径:

Setup → Models → Multiphase → Edit

设置:

  • Model:Volume of Fluid
  • Number of Eulerian Phases:2
  • Formulation:Implicit
  • Interface Modeling:Sharp
  • 开启 Surface Tension
  • 开启 Wall Adhesion

Implicit VOF 配合 Compressive 离散可使用较稳定的物理时间步;Fluent 官方文档也支持 VOF 表面张力、壁面黏附及固化/熔化模型联合使用。(ANSYS Help)

两相定义

第一相:

  • Primary Phase:air

第二相:

  • Secondary Phase:molten-salt

液滴必须是第二相,后续 Patch 的就是第二相体积分数。


8.4 凝固/熔化模型

路径:

Setup → Models → Solidification & Melting

开启:

Solidification/Melting

推荐初始参数:

  • Mushy Zone Constant:5×1055\times10^5
  • Include Pull Velocities:关闭或设为 0
  • 相变区域由材料的 Solidus 和 Liquidus 决定

文献通过对比实验铺展曲线和最终形状,将糊状区常数选为 5×1055\times10^5。这个值适合作为初始值,但熔盐需要进行 1041055×10510610^4、10^5、5\times10^5、10^6 等灵敏度分析,不能直接认为它对所有盐都准确。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

糊状区常数过大时:

  • 熔盐一进入相变区速度就迅速降为零;
  • 铺展量可能偏小;
  • 容易在很早阶段“冻结”。

糊状区常数过小时:

  • 已凝固区域仍可能发生明显流动;
  • 最终液滴形状不稳定;
  • 凝固层可能被不合理地拖动。

九、材料设置

路径:

Setup → Materials

9.1 空气

空气可采用:

  • Density:1.225 kg/m³,或 Ideal Gas;
  • cpc_p:约 1006 J/(kg·K);
  • kk:约 0.024~0.026 W/(m·K);
  • μ\mu:约 1.8×1051.8\times10^{-5} Pa·s。

对于空气:

  • Latent Heat:0
  • Solidus Temperature:0
  • Liquidus Temperature:0

Fluent 官方明确要求,在 VOF 凝固问题中不发生凝固的相应将这些相变参数设为零。(ANSYS Help)

9.2 熔盐

新建材料 molten-salt,填写实际物性。

最简模型可以先采用常数:

  • Density:ρl\rho_l
  • Specific Heat:cpc_p
  • Thermal Conductivity:klk_l
  • Viscosity:μl\mu_l
  • Latent Heat:LhL_h
  • Solidus:TsolT_{\rm sol}
  • Liquidus:TliqT_{\rm liq}

正式计算建议用:

Piecewise-Linear

输入温度相关物性,尤其是:

  • 黏度;
  • 密度;
  • 导热系数;
  • 比热;
  • 表面张力。

如果固态和液态导热系数不同,可以使用分段函数,例如在固相线以下使用 ksk_s,液相线以上使用 klk_l,中间线性过渡。

9.3 基板材料

若建立了固体基板,在 Materials 中创建:

  • steel
  • aluminum
  • ceramic
  • 或实际基板材料。

至少输入:

  • 密度;
  • 比热;
  • 导热系数。

随后在:

Cell Zone Conditions

substrate-solid 指定为对应固体材料。


十、相间作用设置

路径一般为:

Setup → Phase Interaction

10.1 表面张力

在 air 与 molten-salt 之间输入表面张力系数:

σ=σref\sigma=\sigma_{\rm ref}

基础模型先使用常数,确认模型能稳定运行后,再使用温度相关 UDF。

文献使用温度相关的表面张力:

σ=σ0dσdT(TTref)\sigma=\sigma_0-\frac{d\sigma}{dT}(T-T_{\rm ref})

并采用连续表面力 CSF 模型处理气液界面。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

10.2 壁面黏附和接触角

Wall Adhesion 中开启接触角模型。

初始没有实验值时,可以做三组参数:

  • θ=60\theta=60^\circ:润湿性较好;
  • θ=90\theta=90^\circ:中性润湿;
  • θ=120\theta=120^\circ140140^\circ:润湿性较差。

文献使用静态接触角,并发现接触角越小,液滴通常铺展越充分、接触面积越大、凝固越快。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

实际计算前,建议先做一个无重力静滴算例,检查输入接触角后液滴的平衡形状,以确认相序和接触角测量方向没有设置反。


十一、边界条件

11.1 轴线

axis

类型:

axis

11.2 顶部和外侧边界

top-outletside-outlet

类型:

pressure-outlet

设置:

  • Gauge Pressure:0 Pa
  • Backflow Volume Fraction of molten-salt:0
  • Backflow Temperature:环境空气温度 TT_\infty

边界应离液滴足够远,防止液滴或压力波受边界影响。

11.3 恒温平板方案

cold-wall

类型:

wall

Momentum:

  • Stationary Wall
  • No Slip

Thermal:

  • Temperature
  • 输入 TsT_s

Multiphase 或 Wall Adhesion:

  • 输入静态接触角 θ\theta

11.4 固体基板方案

流体—固体界面:

  • Wall Thermal Condition:Coupled
  • 流体侧设置接触角;
  • 固体侧不需要接触角。

基板底面:

  • 固定温度 TsT_s,或
  • 对流边界 h(TT)h(T-T_\infty)

基板侧面:

  • 一般设为绝热;
  • 若真实平板很宽,外侧边界远离撞击区即可。

热接触热阻

如果已知面积热阻 RtcR''_{tc},优先在耦合壁面的 Thermal 设置中输入 Contact Resistance/Thermal Resistance。

若当前版本或界面中没有直接输入项,可建立一个等效薄层:

δ=keqRtc\delta=k_{\rm eq}R''_{tc}

其中 δ\delta 是等效层厚度,keqk_{\rm eq} 是给定的等效导热系数。

文献锡液滴验证算例采用的热接触热阻约为 3.2×106 m2K/W3.2\times10^{-6}\ \mathrm{m^2K/W},但这个数值取决于材料、粗糙度、氧化层和润湿状态,不能直接用于熔盐。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)


十二、Solution Methods 设置

路径:

Solution → Methods

建议设置为:

项目设置
Pressure-Velocity CouplingSIMPLEC
GradientLeast Squares Cell Based
PressurePRESTO!
MomentumSecond Order Upwind
EnergySecond Order Upwind
Volume FractionCompressive
Transient FormulationBounded Second Order Implicit
Warped-Face Gradient CorrectionOn

这与文献的主要数值离散设置一致。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

如果初始接触阶段发散,可暂时设置:

  • Momentum:First Order Upwind;
  • Transient:First Order Implicit;
  • 时间步减半。

稳定运行 20~100 个时间步后,再切回二阶。


十三、Solution Controls 和残差

建议先保持默认松弛因子。

如果不稳定,可以尝试:

  • Pressure:0.2~0.3
  • Momentum:0.4~0.6
  • Volume Fraction:0.3~0.5
  • Energy:0.9~1.0

残差标准参考文献:

  • Continuity:10410^{-4}
  • Axial Velocity:10410^{-4}
  • Radial Velocity:10410^{-4}
  • Energy:10610^{-6}
  • Volume Fraction:10610^{-6}

每个物理时间步最多迭代 30 次。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

VOF 瞬态计算不能只看残差,还应同步检查:

  • 熔盐总体积是否守恒;
  • 液滴界面是否出现非物理锯齿;
  • 每个时间步内速度和温度是否稳定;
  • 固相区域速度是否接近零。

十四、液滴初始化

有两种方法。


14.1 球形液滴:使用 Cell Register 和 Patch

这是最推荐的基础方法,不需要 UDF。

第一步:全域初始化

路径:

Solution → Initialization

选择:

Hybrid Initialization

设置初始化温度为环境温度 TT_\infty,然后点击:

Initialize

第二步:建立圆形液滴区域

路径通常为:

Domain → Adapt → Region

或:

Results/Domain → Cell Registers → New → Region

选择:

Circle

输入:

  • Center X:x0=D/2+h0x_0=D/2+h_0
  • Center Y:0
  • Radius:D/2D/2

例如:

  • D=2D=2 mm;
  • h0=0.05h_0=0.05 mm;
  • 圆心 x0=1.05x_0=1.05 mm;
  • y0=0y_0=0
  • 半径 1 mm。

二维计算域中只保留 y0y\geq0 的半圆,绕轴旋转后就是完整球形液滴。

Fluent 支持通过圆形 Cell Register 选定区域,再对体积分数、温度和速度进行 Patch。(ANSYS Help)

第三步:Patch 熔盐体积分数

路径:

Solution → Initialization → Patch

设置:

  • Phase:molten-salt
  • Variable:Volume Fraction
  • Value:1
  • Register:刚创建的 droplet-circle
  • 可勾选 Patch Reconstructed Interface

点击:

Patch

Fluent 的 Reconstructed Interface 选项可在几何区域边界单元中设置更合理的初始体积分数,而不是简单地把所有相交单元全部设为 1。(ANSYS Help)

第四步:Patch 液滴温度

  • Variable:Static Temperature
  • Value:T0T_0
  • Register:droplet-circle

点击:

Patch

第五步:Patch 撞击速度

轴向为 xx 方向时:

  • Variable:X Velocity
  • Value:v0-v_0
  • Register:droplet-circle

径向速度:

  • Variable:Y Velocity
  • Value:0

第六步:检查初始化结果

显示:

Results → Graphics → Contours

依次查看:

  • Volume Fraction → molten-salt;
  • Static Temperature;
  • X Velocity。

正确状态应为:

  • 半圆形熔盐液滴;
  • 液滴温度为 T0T_0
  • 液滴内部速度为 v0-v_0
  • 周围空气体积分数为 1;
  • 液滴未与平板重叠。

14.2 椭球液滴:使用 UDF

文献研究的是椭球液滴。对于轴对称椭球,定义:

AR=abAR=\frac{a}{b}

其中:

  • aa:轴向半轴;
  • bb:径向半轴。

保持与直径 DeD_e 的球体体积相同:

a=ARba=AR\,b b=De2AR1/3b=\frac{D_e}{2AR^{1/3}} a=DeAR2/32a=\frac{D_eAR^{2/3}}{2}

下面的 UDF 可同时初始化椭球体积分数、温度和速度。假设:

  • 第一相为空气;
  • 第二相为熔盐;
  • xx 为竖直轴向;
  • yy 为径向;
  • 平板位于 x=0x=0
c
#include "udf.h" #include "math.h" /* ---------- 用户参数,全部使用 SI 单位 ---------- */ #define SALT_PHASE_INDEX 1 /* 主相空气=0,第二相熔盐=1 */ #define DE 0.002 /* 等效直径,m */ #define AR 1.0 /* a/b;1.0 为球形 */ #define GAP 0.00005 /* 液滴底部到平板的间隙,m */ #define V_IMP 0.5 /* 撞击速度绝对值,m/s */ #define T_DROP 800.0 /* 液滴初始温度,K */ DEFINE_INIT(init_molten_salt_drop, domain) { Thread *t; cell_t c; real xc[ND_ND]; real b = 0.5 * DE / pow(AR, 1.0 / 3.0); real a = AR * b; real x0 = GAP + a; thread_loop_c(t, domain) { if (FLUID_THREAD_P(t)) { Thread *salt_thread = THREAD_SUB_THREAD(t, SALT_PHASE_INDEX); begin_c_loop(c, t) { real ellipse_value; C_CENTROID(xc, c, t); ellipse_value = ((xc[0] - x0) * (xc[0] - x0)) / (a * a) + (xc[1] * xc[1]) / (b * b); if (ellipse_value <= 1.0) { C_VOF(c, salt_thread) = 1.0; C_T(c, t) = T_DROP; C_U(c, t) = -V_IMP; C_V(c, t) = 0.0; } else { C_VOF(c, salt_thread) = 0.0; } } end_c_loop(c, t) } } }

编译与挂接

  1. 将代码保存为 init_drop.c

  2. 路径中尽量不要使用中文和空格;

  3. 进入:

    User-Defined → Functions → Compiled

  4. Add 文件;

  5. 点击 Build

  6. 点击 Load

  7. 进入:

    User-Defined → Function Hooks

  8. 在 Initialization 中选择:

    init_molten_salt_drop

  9. 再执行 Hybrid Initialization。

Fluent 官方 UDF 文档提供 DEFINE_INITDEFINE_PROPERTY 等宏用于初始化和自定义物性。(ANSYS Help)


十五、温度相关表面张力 UDF

基础模型稳定后,再使用此 UDF。

c
#include "udf.h" #define SIGMA_REF 0.150 /* 参考温度下的表面张力,N/m */ #define T_REF 800.0 /* 参考温度,K */ #define DSIGMA_DT -1.0e-4 /* d(sigma)/dT,N/(m K) */ #define SIGMA_MIN 1.0e-6 DEFINE_PROPERTY(molten_salt_surface_tension, c, t) { real temperature = C_T(c, t); real sigma; sigma = SIGMA_REF + DSIGMA_DT * (temperature - T_REF); if (sigma < SIGMA_MIN) sigma = SIGMA_MIN; return sigma; }

替换:

  • SIGMA_REF
  • T_REF
  • DSIGMA_DT

为实际熔盐数据。

编译后,在:

Phase Interaction → Surface Tension Coefficient

选择:

molten_salt_surface_tension

不要在缺少可靠数据的情况下随意设置 dσ/dTd\sigma/dT,因为它会改变界面切向应力和液滴内部流动。


十六、时间步设置

对于 D=2D=2 mm、撞击区网格 10~20 μm 的模型,可以从:

Δt=2.5×106 s\Delta t=2.5\times10^{-6}\ \mathrm{s}

开始。

主要检查界面 Courant 数:

CoUΔtΔxCo\sim\frac{U\Delta t}{\Delta x}

建议早期撞击阶段控制在约 0.25~0.5 以下。Fluent 对显式界面处理会根据网格、速度和 Courant 数计算 VOF 子时间步;即使使用隐式 VOF,也不宜把物理时间步设置得过大。(ANSYS Help)

推荐分阶段:

物理时间时间步建议
撞击前至 5 ms12.5 μs1\sim2.5\ \mu s
5~20 ms2.55 μs2.5\sim5\ \mu s
20 ms 以后510 μs5\sim10\ \mu s

只能在液滴界面运动明显减慢后增大时间步。

可选的时间步 UDF:

c
#include "udf.h" DEFINE_DELTAT(adaptive_drop_timestep, domain) { real time = CURRENT_TIME; if (time < 0.005) return 2.5e-6; else if (time < 0.020) return 5.0e-6; else return 1.0e-5; }

挂接到:

User-Defined → Function Hooks → Deltat

或在 Run Calculation 中选择 User-Defined Time Step,具体名称随版本略有差异。


十七、监测量设置

除残差外,至少建立以下监测量。

17.1 熔盐质量或体积守恒

定义熔盐相总体积:

Vsalt=ΩαsaltdVV_{\rm salt}=\int_\Omega \alpha_{\rm salt}\,dV

路径:

Solution → Report Definitions → New → Volume Integral

选择:

  • Field Variable:Volume Fraction of molten-salt
  • Cell Zone:fluid-domain

记录整个计算过程中的变化。

建议:

V(t)V(0)V(0)<1%\frac{|V(t)-V(0)|}{V(0)}<1\%

高精度模型最好小于 0.5%。


17.2 总体固相率

文献定义固相率为凝固熔盐体积占总熔盐体积的比例。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

在 Fluent 中可写为:

Sf=1ΩαsaltϕdVΩαsaltdVS_f= 1-\frac{\int_\Omega \alpha_{\rm salt}\phi\,dV} {\int_\Omega\alpha_{\rm salt}\,dV}

其中:

  • αsalt\alpha_{\rm salt}:熔盐体积分数;
  • ϕ\phi:Liquid Fraction;
  • ϕ=1\phi=1:液态;
  • ϕ=0\phi=0:固态。

建立 Custom Field Function:

text
salt-liquid = VolumeFraction(molten-salt) * Liquid-Fraction

然后分别计算:

  • Volume Integral of salt-liquid
  • Volume Integral of VolumeFraction(molten-salt)

最后计算:

text
solid-fraction = 1 - salt-liquid-volume / salt-total-volume

Sf>0.99S_f>0.99 时,可以认为液滴基本完全凝固。


17.3 铺展因子

平面上的铺展因子:

β=LD=2rwD\beta=\frac{L}{D}=\frac{2r_w}{D}

其中 rwr_w 是最大润湿半径。

实用提取方法:

  1. 创建熔盐体积分数为 0.5 的等值线: Surface → Iso-Surface
  2. Field:Volume Fraction of molten-salt;
  3. Iso-Value:0.5;
  4. 对该等值线计算最大 Y Coordinate;
  5. 得到 rmaxr_{\max}
  6. 计算:
β=2rmaxD\beta=\frac{2r_{\max}}{D}

若发生卫星液滴或飞溅,最大径向坐标可能来自脱离的小液滴,此时应改为在平板上方一个网格高度处建立水平线,只提取与主液滴相连的润湿范围。


17.4 其他监测量

建议同时输出:

  • 液滴最大温度;
  • 液滴最小温度;
  • 熔盐平均温度;
  • 最大速度;
  • 平板热流;
  • 熔盐—基板接触面积;
  • 最大铺展时刻 tβ,maxt_{\beta,\max}
  • 完全凝固时间 tsolidt_{\rm solid}

十八、自动保存

进入:

Solution → Calculation Activities → Autosave

建议:

  • 每 20~50 个时间步保存一次;
  • 保存 Case and Data;
  • 文件名添加 flow-time 或 time-step;
  • 确保磁盘空间充足。

例如总时间 50 ms,动画需要约 200 帧,则每:

0.05/200=0.00025 s0.05/200=0.00025\ \mathrm{s}

保存一帧即可。

若时间步为 2.5 μs,大约每 100 个时间步保存一次。


十九、开始计算

进入:

Solution → Run Calculation

初始建议:

  • Time Step Size:2.5×1062.5\times10^{-6} s
  • Number of Time Steps:先运行 100~200 步测试
  • Max Iterations/Time Step:30
  • Reporting Interval:1

先短跑检查:

  1. 液滴是否向下运动;
  2. 液滴是否正常接触平板;
  3. 是否出现铺展;
  4. 平板接触区域温度是否下降;
  5. Liquid Fraction 是否从底部开始降低;
  6. 已凝固区域速度是否接近零;
  7. 熔盐总体积是否守恒。

短跑没有问题后,再运行完整时间。


二十、结果显示

20.1 液滴界面

路径:

Results → Graphics → Contours

选择:

  • Phases → Volume Fraction → molten-salt
  • Range:0~1
  • Filled:On
  • Draw Mesh:Off

体积分数 0.5 的边界通常视为气液界面。

20.2 温度场

选择:

  • Temperature → Static Temperature
  • 固定颜色范围,例如 TsT_sT0T_0

动画中必须固定 Min/Max,否则每一帧颜色范围变化,会让温度变化产生误导。

20.3 凝固情况

选择:

  • Solidification/Melting → Liquid Fraction
  • Range:0~1

解释:

  • 0:完全固态;
  • 0~1:糊状区;
  • 1:完全液态。

建议同时叠加熔盐体积分数 0.5 的界面线,避免把空气区域误认为熔盐液相。

20.4 速度矢量

选择:

Results → Graphics → Vectors

设置:

  • Vectors of Velocity;
  • Color by Velocity Magnitude;
  • 只在熔盐区域显示,或配合体积分数裁剪;
  • Scale 适当减小;
  • Skip 数量适当提高,防止箭头过密。

二十一、制作仿真动画

Fluent 可以在计算过程中记录动画,也可以利用已保存的瞬态 Data 文件事后制作。官方文档支持对 Contour、Vector、Scene、XY Plot 和 Report Plot 进行记录,并可输出 MP4。(ANSYS Help)

方法一:计算前创建 Solution Animation

第一步:建立温度 Contour 对象

进入:

Results → Graphics → Contours

创建并命名:

temperature-contour

设置:

  • Static Temperature;
  • Fixed Range;
  • Filled;
  • 选中 fluid-domain;
  • 调整视角和缩放;
  • 保存对象。

第二步:建立液相率 Contour 对象

创建:

liquid-fraction-contour

设置:

  • Liquid Fraction;
  • Range 0~1;
  • 固定色标;
  • 保存对象。

第三步:建立体积分数对象

创建:

droplet-interface

设置:

  • Volume Fraction of molten-salt;
  • Range 0~1;
  • 或仅显示 0.5 等值线。

第四步:创建动画定义

路径:

Solution → Activities → Create → Solution Animations

设置:

  • Name:salt_drop_temperature
  • Record After Every:例如 20 Time Steps
  • Storage Type:PNG Image 或 In Memory
  • Storage Directory:指定英文目录
  • Animation Object:temperature-contour
  • View:Use Active
  • Append File Name With:flow-time
  • 点击 OK

分别为:

  • 温度;
  • 液相率;
  • 熔盐体积分数;
  • 速度矢量

创建不同动画。

官方建议计算前选定要记录的图形对象、记录频率和存储类型,可使用 PNG/JPEG/TIFF 等逐帧图像。(ANSYS Help)

第五步:计算结束后播放

进入:

Solution → Calculation Activities → Solution Animations

或在动画对象上右键:

Playback

检查:

  • 帧顺序;
  • 播放速度;
  • 视图是否固定;
  • 色标是否固定;
  • 是否有遗漏帧。

第六步:导出 MP4

在 Playback 窗口:

  1. 选择 Write/Record;
  2. Format:Video File;
  3. Video Name:molten_salt_drop.mp4
  4. Video Options;
  5. Format:MP4;
  6. 设置帧率,例如 20~30 fps;
  7. 设置分辨率,例如 1920×1080;
  8. 点击 Write。

Fluent 官方支持 MP4、AVI、FLV、MOV 和 MPEG 等格式,其中 MP4 最通用。(ANSYS Help)


方法二:计算后利用 Data 文件制作动画

适用于计算前忘记建立动画,但保存了多个瞬态 Data 文件的情况。

  1. 打开最终 Case;
  2. 启用: Transient Postprocessing
  3. Timestep Selector 中选择所有 Data 文件;
  4. 建立温度、液相率或体积分数 Contour;
  5. 点击 Animation...
  6. 选择对应 Graphics Object;
  7. 生成动画;
  8. 导出 MP4。

Transient Postprocessing 可利用已保存的时间步文件事后生成 Contour、Vector、Scene 和 Plot 动画。(ANSYS Help)


二十二、制作类似文献图 6 的动画

文献第 8 页的结果图将:

  • 左半部分显示温度;
  • 右半部分显示液相率和速度矢量;
  • 黑线表示液滴外轮廓;
  • 底部表示冷基板。

(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

在 Fluent 中最稳妥的实现方法是分别导出:

  1. 温度场动画;
  2. 液相率动画;
  3. 速度矢量动画;
  4. 体积分数 0.5 界面动画。

然后在视频编辑软件中左右拼接。

也可以在 Fluent Scene 中叠加:

  • Temperature Contour;
  • Volume Fraction 0.5 Iso-Line;
  • Wall Mesh;
  • Velocity Vectors。

但同一计算域内左右半边显示不同变量操作较复杂,直接输出两个同步动画再拼接更可靠。


二十三、验证流程

正式模拟熔盐前,建议先复现文献中的平板液态锡验证算例:

  • D=2.7D=2.7 mm;
  • 球形液滴;
  • v0=1.0v_0=1.0 m/s;
  • θ=140\theta=140^\circ
  • 液滴温度 513 K;
  • 平板温度 298 K;
  • 平板为水平平面;
  • 糊状区常数从 10510^55×1055\times10^5 比较。

文献给出了液态锡和固态锡的密度、比热、导热系数、黏度、潜热和熔点,并将数值铺展曲线与实验数据进行了比较。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

完成验证后再依次替换:

  1. 锡物性 → 熔盐物性;
  2. 锡熔点 → 熔盐固相线/液相线;
  3. 接触角;
  4. 液滴初温;
  5. 基板温度;
  6. 热接触热阻;
  7. 液滴尺寸和速度。

这样可以区分“软件设置错误”和“熔盐参数导致的物理变化”。


二十四、必须完成的无关性分析

24.1 网格无关性

至少比较:

  • D/50D/50
  • D/100D/100
  • D/200D/200

比较指标:

  • 最大铺展因子 βmax\beta_{\max}
  • 最大铺展时刻;
  • 完全凝固时间;
  • 最终液滴高度;
  • 熔盐质量误差。

相邻两套网格结果差异小于约 2%~5% 后,才能认为基本独立。

24.2 时间步无关性

至少比较:

  • 5 μs5\ \mu s
  • 2.5 μs2.5\ \mu s
  • 1.25 μs1.25\ \mu s

如果液滴界面或铺展曲线差异明显,应继续减小时间步。

24.3 糊状区常数

比较:

Amush=104, 105, 5×105, 106A_{\rm mush}=10^4,\ 10^5,\ 5\times10^5,\ 10^6

观察:

  • 凝固层速度;
  • 最大铺展;
  • 回缩;
  • 最终形貌;
  • 凝固时间。

24.4 接触角和热阻

熔盐接触角和热接触热阻通常是不确定性最大的两个参数。建议至少分别进行三组灵敏度计算。


二十五、常见问题与处理

1. 液滴接触平板时立即发散

处理顺序:

  1. 时间步减半;
  2. Momentum 改一阶;
  3. 瞬态格式改一阶;
  4. 液滴初始间隙减小;
  5. 降低每一步的界面 Courant 数;
  6. 检查平板附近网格质量;
  7. 暂时使用恒定表面张力。

2. 液滴体积不断减少

通常由以下原因导致:

  • 时间步过大;
  • 界面网格过粗;
  • 压力出口太近;
  • 界面离散过度扩散;
  • 网格高扭曲。

处理:

  • 采用 Compressive;
  • 减小时间步;
  • 加密界面;
  • 增大空气计算域;
  • 监控熔盐总体积。

3. 界面附近出现很大非物理速度

这通常是寄生流或表面张力数值误差。

处理:

  • 加密界面;
  • 使用更均匀的网格;
  • 减小时间步;
  • 先用常数表面张力;
  • 平滑物性跳变;
  • 检查接触角;
  • 必要时测试 Coupled Level Set + VOF。

Fluent 官方文档指出,在表面张力作用下可能出现寄生流,并提供了 Coupled Level Set + VOF 和修正表面张力加权方式作为可选处理。(ANSYS Help)

4. 液滴完全不铺展

检查:

  • 接触角是否过大;
  • 糊状区常数是否过大;
  • 基板温度是否过低;
  • 热接触热阻是否过小;
  • 液滴初温是否只略高于液相线;
  • 撞击速度是否 Patch 正确;
  • 重力方向是否正确。

5. 已凝固区域仍高速流动

检查:

  • Solidification/Melting 是否开启;
  • 潜热、固相线和液相线是否输入到熔盐材料;
  • Liquid Fraction 是否确实接近 0;
  • Mushy Zone Constant 是否过小;
  • 是否错误地把相变参数填到了空气材料。

6. 液滴落下前已经明显变形

处理:

  • 减小初始间隙;
  • 直接给定撞击速度,不模拟长距离自由落体;
  • 球形液滴一般影响较小;
  • 椭球液滴可像文献一样采用 UDF 控制飞行阶段表面张力,但要在论文中说明这一人为处理。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

二十六、建议的最终计算顺序

  1. 仅 VOF、关闭能量和凝固,验证液滴能正常撞击和铺展;
  2. 开启 Energy,但暂时关闭凝固,验证温度传热;
  3. 开启 Solidification/Melting;
  4. 使用恒定表面张力和静态接触角跑通;
  5. 增加温度相关物性;
  6. 增加固体基板;
  7. 增加热接触热阻;
  8. 完成网格无关性;
  9. 完成时间步无关性;
  10. 完成糊状区常数、接触角、基板温度和撞击速度的参数分析;
  11. 最后重跑高精度算例并输出动画。

提供具体的熔盐体系、液滴直径、液滴初温、平板温度、撞击速度、接触角和基板材料后,可进一步把上述模板替换成一套完整的物性表、计算域尺寸、时间步以及可直接编译的 UDF 参数。

質問

我现在只想复现这个文献的结果,给出具体步骤

Python

python
import math rho=6970 D=0.002 sigma=0.5263-6.59e-5*(513-505) for We in [0.265,2.386,6.626,12.99,26.51]: v=math.sqrt(We*sigma/(rho*D)) print(We,v)

思考

Computed baseline velocities and initial ellipse dimensions

Worked for 2m 21s

回答

可以。为了避免一开始就同时复现几十组参数,建议按论文的实际验证顺序执行:

  1. 先复现图5的“球形锡液滴撞击低温平板”验证算例;
  2. 再复现图6~8的基准工况;
  3. 最后改变 AR、曲率、We 和接触角,复现图9~17。

以下参数均以文献为依据。论文使用 ANSYS Fluent 2023R2、二维轴对称、VOF 与 Solidification/Melting 耦合,并采用焓—多孔介质法处理凝固。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

需要先说明:论文没有公开三个 UDF 的源代码,也没有给出初始下落间距 h0h_0、局部细化区的精确尺寸、撞击前后自适应时间步的切换规则,以及热接触热阻在 Fluent 面板中的具体实现。因此可以复现其主要曲线、形貌和凝固时间,但仅靠论文不能保证逐像素完全一致。


一、先复现论文图5:平板验证算例

这是最容易检查模型是否正确的算例。

1.1 图5算例参数

参数数值
液滴材料
等效直径 DeD_e2.7 mm
液滴形状球形,AR=1
撞击速度 v0v_01.0 m/s
静态接触角140140^\circ
液滴初温 T0T_0513 K
平板温度 TsT_s298 K
热接触热阻3.2×106 m2K/W3.2\times10^{-6}\ \mathrm{m^2K/W}
糊状区常数5×1055\times10^5
相变温度区间504~506 K
时间步2.5×1062.5\times10^{-6} s
平板类型无厚度恒温壁面

论文用该算例与 Aziz 等人的实验进行了对比,并通过改变糊状区常数,最终选择 Amush=5×105A_{\rm mush}=5\times10^5。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)


二、建立二维轴对称平板模型

2.1 Workbench 操作

  1. 打开 ANSYS Workbench。
  2. 拖入 Fluid Flow (Fluent)
  3. 双击 Geometry
  4. 单位选择 mm
  5. 建立二维矩形流体域。

对于图5验证算例,计算域大小按论文的 5De×5De5D_e\times5D_e

5De=5×2.7=13.5 mm5D_e=5\times2.7=13.5\ {\rm mm}

建立:

  • 轴向长度:13.5 mm;
  • 径向长度:13.5 mm。

Fluent 二维轴对称默认绕 xx 轴旋转,因此建议:

  • xx:竖直撞击方向;
  • yy:径向;
  • y=0y=0:轴线;
  • 液滴沿负 xx 方向撞向平板;
  • 平板位于 x=0x=0

矩形四条边分别命名:

名称类型
y=0y=0axis轴线
x=0x=0cold-wall低温平板
x=13.5x=13.5 mmtop-outlet压力出口
y=13.5y=13.5 mmside-outlet压力出口

生成一个流体面,命名为:

fluid-domain


三、平板算例网格

3.1 网格要求

论文最终采用:

  • 液滴和撞击区:10 μm;
  • 外部区域最大尺寸:100 μm;
  • 时间步:2.5 μs;
  • 液滴直径方向约 200 个单元。

论文指出继续加密后,铺展因子、液滴形貌和凝固前沿变化已经很小。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

图5中的液滴直径是 2.7 mm。严格按 10 μm 划分,相当于直径方向约 270 个单元,计算量很大。

建议分两次:

调试网格

  • 撞击区:20~25 μm;
  • 外部区:100 μm。

最终网格

  • 撞击区:10~13.5 μm;
  • 外部区:100 μm。

3.2 局部细化区域

在平板上方建立一个局部细化矩形:

  • 轴向范围:x=04x=0\sim4 mm;
  • 径向范围:y=05y=0\sim5 mm。

该区域覆盖:

  • 初始液滴;
  • 撞击位置;
  • 液滴最大铺展区;
  • 后续回缩区。

网格设置:

  • Method:Quadrilateral Dominant
  • Face Meshing:Quad/Tri
  • 细化区尺寸:0.01 mm
  • 外部区尺寸:0.1 mm
  • Growth Rate:1.1~1.15

检查:

  • Skewness 尽量小于 0.85;
  • 平板附近单元不能突然变大;
  • 轴线附近单元保持规则。

四、启动 Fluent

在 Workbench 中双击 Setup

  • Dimension:2D
  • Double Precision:开启
  • Solver:CPU
  • 并行核数按电脑配置选择

进入 Fluent 后:

General → Check

确认无负体积和网格错误。

设置:

  • Solver:Pressure-Based
  • Time:Transient
  • 2D Space:Axisymmetric
  • Velocity Formulation:Absolute
  • Gravity:开启
  • gx=9.81 m/s2g_x=-9.81\ {\rm m/s^2}
  • gy=0g_y=0

论文未提到湍流模型,因此选择:

Models → Viscous → Laminar


五、开启物理模型

5.1 Energy

路径:

Models → Energy

设置:

On

5.2 VOF

路径:

Models → Multiphase

选择:

  • Model:Volume of Fluid
  • Number of Phases:2
  • Primary Phase:air
  • Secondary Phase:liquid-tin

推荐:

  • Formulation:Implicit
  • Interface Modeling:Sharp
  • Surface Tension:开启
  • Wall Adhesion:开启

论文明确采用空气为主相、液态锡为第二相,并通过 VOF 捕捉气液界面。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

5.3 Solidification/Melting

路径:

Models → Solidification & Melting

设置:

  • Solidification/Melting:On
  • Mushy Zone Constant:
Amush=5×105A_{\rm mush}=5\times10^5
  • Pull Velocity:0
  • Include Pull Velocities:关闭

相变采用:

Tsolid=504 KT_{\rm solid}=504\ {\rm K} Tliquid=506 KT_{\rm liquid}=506\ {\rm K}

论文以熔点 505 K 为中心,设置 2 K 的相变区间。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)


六、设置空气和锡的物性

6.1 空气

在:

Materials → air

设置:

物性数值
Density1.18 kg/m³
Specific Heat1006 J/(kg·K)
Thermal Conductivity0.0263 W/(m·K)
Viscosity1.85×1051.85\times10^{-5} Pa·s

空气不发生相变:

  • Latent Heat:0
  • Solidus Temperature:0
  • Liquidus Temperature:0

6.2 创建液态锡材料

新建材料:

liquid-tin

设置:

物性数值
Density6970 kg/m³
Specific Heat244 J/(kg·K)
Viscosity1.92×1031.92\times10^{-3} Pa·s
Latent Heat59200 J/kg
Solidus Temperature504 K
Liquidus Temperature506 K

导热系数需要区分固态和液态:

  • 固态锡:66.6 W/(m·K)
  • 液态锡:33.6 W/(m·K)

在 Thermal Conductivity 中选择:

Piecewise-Linear

输入例如:

温度 K导热系数 W/(m·K)
29866.6
50466.6
50633.6
51333.6
60033.6

论文中的锡和空气物性如上。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)


七、表面张力设置

论文采用:

σ=0.52636.59×105(T505)\sigma=0.5263-6.59\times10^{-5}(T-505)

单位为 N/m。

在液滴初温 513 K 时:

σ(513)=0.52636.59×105×80.52577 N/m\sigma(513) =0.5263-6.59\times10^{-5}\times8 \approx0.52577\ {\rm N/m}

调试阶段可先输入常数:

σ=0.52577 N/m\sigma=0.52577\ {\rm N/m}

当常数表面张力模型能够稳定运行后,再使用温度相关 UDF。

c
#include "udf.h" #define SIGMA0 0.5263 #define DSIGMA 6.59e-5 #define TREF 505.0 DEFINE_PROPERTY(tin_surface_tension, c, t) { real temperature; real sigma; temperature = C_T(c, t); sigma = SIGMA0 - DSIGMA * (temperature - TREF); if (sigma < 0.0) sigma = 0.0; return sigma; }

保存为:

tin_surface_tension.c

编译:

  1. User-Defined → Functions → Compiled
  2. Add
  3. Build
  4. Load

然后在空气—锡的 Surface Tension Coefficient 中选择:

tin_surface_tension

论文还在液滴飞行阶段将表面张力设为零,以保持椭球初始形状;图5是球形液滴,初始间隙很小时可以暂时不做这项处理。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)


八、边界条件

8.1 压力出口

top-outletside-outlet

  • Type:Pressure Outlet
  • Gauge Pressure:0 Pa
  • Backflow Temperature:513 K
  • Backflow Volume Fraction of liquid-tin:0

注意:论文不是把空气初温设为 298 K,而是将整个计算域初始温度设为 513 K,低温条件只施加在壁面上。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

8.2 平板

cold-wall

设置:

  • Wall Motion:Stationary
  • Shear Condition:No Slip
  • Contact Angle:140140^\circ

等效实现热接触热阻

论文同时规定:

  • 基板温度:298 K;
  • 热接触热阻:
Rtc=3.2×106 m2K/WR''_{tc}=3.2\times10^{-6}\ {\rm m^2K/W}

Fluent 中可以用等效对流边界实现:

h=1Rtc=13.2×106=312500 W/(m2K)h=\frac{1}{R''_{tc}} =\frac{1}{3.2\times10^{-6}} =312500\ {\rm W/(m^2K)}

在 Thermal 中选择:

  • Thermal Condition:Convection
  • Heat Transfer Coefficient:312500 W/(m²·K)
  • Free Stream Temperature:298 K

这个边界满足:

q=Tfluid298Rtcq=\frac{T_{\rm fluid}-298}{R''_{tc}}

这与“298 K 恒温基板通过接触热阻吸热”等效。论文没有说明它在 Fluent 面板中的具体输入方式,因此这是基于其数学边界条件的等效实现。


九、初始化球形液滴

9.1 初始间隙

论文没有给出 h0h_0。建议先令液滴底部距平板:

h0=0.02 mmh_0=0.02\ {\rm mm}

也就是约两个 10 μm 网格。

对于图5:

  • 液滴半径:1.35 mm;
  • 液滴中心:
x0=1.35+0.02=1.37 mmx_0=1.35+0.02=1.37\ {\rm mm}
  • y0=0y_0=0

9.2 创建圆形 Cell Register

进入:

Adapt → Region

选择 Circle:

  • Center X:0.00137 m
  • Center Y:0
  • Radius:0.00135 m

创建 Register:

droplet-region

9.3 初始化全域

Solution Initialization

设置:

  • X Velocity:0
  • Y Velocity:0
  • Temperature:513 K
  • liquid-tin Volume Fraction:0

点击:

Initialize

9.4 Patch 液滴

进入:

Patch

依次 Patch:

锡体积分数

  • Phase:liquid-tin
  • Variable:Volume Fraction
  • Value:1
  • Register:droplet-region

轴向速度

  • Variable:X Velocity
  • Value:-1.0 m/s

径向速度

  • Variable:Y Velocity
  • Value:0

温度

  • Variable:Static Temperature
  • Value:513 K

显示液态锡体积分数,确认:

  • 液滴为半圆;
  • 绕轴旋转后对应完整球体;
  • 液滴没有与壁面重叠;
  • 液滴向负 xx 方向运动。

十、求解方法

进入:

Solution Methods

设置如下:

项目设置
Pressure-Velocity CouplingSIMPLEC
GradientLeast Squares Cell Based
PressurePRESTO!
MomentumSecond Order Upwind
EnergySecond Order Upwind
Volume FractionCompressive
Transient FormulationBounded Second Order Implicit
Warped-Face Gradient CorrectionOn

残差:

方程收敛标准
Continuity10410^{-4}
X Velocity10410^{-4}
Y Velocity10410^{-4}
Energy10610^{-6}
Volume Fraction10610^{-6}

每个时间步:

  • Max Iterations:30

这些设置与论文完全一致。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)


十一、图5计算设置

进入:

Run Calculation

输入:

  • Time Step Size:
2.5×106 s2.5\times10^{-6}\ {\rm s}
  • Number of Time Steps:8000
  • Max Iterations/Time Step:30

8000 步对应:

8000×2.5 μs=20 ms8000\times2.5\ \mu s=20\ {\rm ms}

论文图5主要比较以下时刻:

  • 0.7 ms
  • 1.7 ms
  • 4.5 ms
  • 8.3 ms
  • 11.3 ms
  • 15.3 ms
  • 17.3 ms

由于你的液滴初始有 0.02 mm 间隙,撞击前约需要:

tc0.02 mm1 m/s=0.02 mst_c\approx\frac{0.02\ {\rm mm}}{1\ {\rm m/s}} =0.02\ {\rm ms}

后处理时应使用:

tpaper=tFluenttct_{\rm paper}=t_{\rm Fluent}-t_c

将液滴首次接触壁面的时刻定义为论文的 t=0t=0


十二、检查图5结果

至少显示:

  1. 液态锡体积分数;
  2. 液相率;
  3. 温度;
  4. 速度;
  5. 铺展因子。

论文图5中的主要特征为:

  • 0.7~4.5 ms:液滴快速铺展;
  • 8.3 ms:中心形成细长锥形结构;
  • 11.3 ms:锥形顶部可能发生脱离;
  • 15.3~17.3 ms:剩余液滴快速回缩并在中心形成凹陷;
  • Amush=5×105A_{\rm mush}=5\times10^5 的最终稳定铺展最接近实验。

如果完全没有锥形和回缩,优先检查:

  • 接触角;
  • 表面张力;
  • 糊状区常数;
  • 时间步;
  • 网格;
  • 热接触热阻。

十三、复现论文基准工况图6~8

图5验证通过后,再建立曲面工况。

13.1 基准工况参数

论文的参考算例为:

参数数值
DeD_e2 mm
D=Ds/DeD^*=D_s/D_e5
球面直径 DsD_s10 mm
球面半径 RsR_s5 mm
球冠高度 hsh_s1 mm
WeWe6.626
v0v_00.5 m/s
接触角140140^\circ
液滴初温513 K
曲面温度298 K
热接触热阻3.2×106 m2K/W3.2\times10^{-6}\ \mathrm{m^2K/W}
网格撞击区 10 μm
时间步2.5 μs

论文以此为参考工况,并分别改变液滴形状、球面曲率、撞击速度和接触角。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)


十四、曲面几何的精确画法

Fluent 坐标仍采用:

  • xx:轴向;
  • yy:径向;
  • y=0y=0:轴线。

基准曲面:

Rs=5 mmR_s=5\ {\rm mm} hs=1 mmh_s=1\ {\rm mm}

令球面最高点位于:

(x,y)=(1,0) mm(x,y)=(1,0)\ {\rm mm}

圆心位于:

(xc,yc)=(hsRs,0)=(4,0) mm(x_c,y_c)=(h_s-R_s,0)=(-4,0)\ {\rm mm}

圆方程:

(x+4)2+y2=52(x+4)^2+y^2=5^2

球冠与 x=0x=0 平面的交点满足:

ye=Rs2(Rshs)2y_e=\sqrt{R_s^2-(R_s-h_s)^2} ye=2516=3 mmy_e=\sqrt{25-16}=3\ {\rm mm}

因此曲面弧线从:

(x,y)=(1,0)(x,y)=(1,0)

延伸至:

(x,y)=(0,3)(x,y)=(0,3)

14.1 建立流体域

对于 De=2D_e=2 mm:

5De=10 mm5D_e=10\ {\rm mm}

流体域边界依次为:

  • 轴线:从 (1,0)(1,0)(10,0)(10,0)
  • 顶部:从 (10,0)(10,0)(10,10)(10,10)
  • 径向出口:从 (10,10)(10,10)(0,10)(0,10)
  • 底部水平壁:从 (0,10)(0,10)(0,3)(0,3)
  • 球冠壁:从 (0,3)(0,3) 沿圆弧到 (1,0)(1,0)

球冠壁命名:

curved-wall

底部远端水平壁可以命名:

bottom-wall

并赋予相同温度。液滴不会铺展到距离轴线 3 mm 以外,因此远端水平壁对结果影响很小。


十五、椭球液滴初始化

论文定义:

AR=a0b0AR=\frac{a_0}{b_0}

其中:

  • a0a_0:轴向半轴;
  • b0b_0:径向半轴。

保持液滴体积与直径 DeD_e 的球体相同:

b0=De2AR1/3b_0=\frac{D_e}{2AR^{1/3}} a0=DeAR2/32a_0=\frac{D_eAR^{2/3}}{2}

对于 De=2D_e=2 mm:

ARa0a_0, mmb0b_0, mm
0.50.63001.2599
0.70.78841.1262
1.01.00001.0000
1.51.31040.8736
2.01.58740.7937

使用以下 UDF 初始化所有椭球工况。

c
#include "udf.h" #include "math.h" /* 第二相为锡,因此索引为 1 */ #define TIN_PHASE_INDEX 1 /* 基准工况参数 */ #define DE 0.002 #define AR 1.0 #define HS 0.001 #define GAP 0.000020 #define V0 0.5 #define T_DROP 513.0 DEFINE_INIT(initialize_tin_droplet, domain) { Thread *mix_thread; Thread *tin_thread; cell_t c; real xyz[ND_ND]; real b0 = DE / (2.0 * pow(AR, 1.0 / 3.0)); real a0 = AR * b0; /* 曲面最高点位于 x=HS */ real x_center = HS + GAP + a0; thread_loop_c(mix_thread, domain) { if (FLUID_THREAD_P(mix_thread)) { tin_thread = THREAD_SUB_THREAD(mix_thread, TIN_PHASE_INDEX); begin_c_loop(c, mix_thread) { real ellipse; C_CENTROID(xyz, c, mix_thread); ellipse = pow((xyz[0] - x_center) / a0, 2.0) + pow(xyz[1] / b0, 2.0); if (ellipse <= 1.0) { C_VOF(c, tin_thread) = 1.0; /* 轴向向下撞击 */ C_U(c, mix_thread) = -V0; C_V(c, mix_thread) = 0.0; C_T(c, mix_thread) = T_DROP; } else { C_VOF(c, tin_thread) = 0.0; } } end_c_loop(c, mix_thread) } } }

修改:

c
#define AR 1.0

即可得到不同液滴形状。

编译并挂接:

  1. User-Defined → Functions → Compiled
  2. Add initialize_tin_droplet.c
  3. Build
  4. Load
  5. User-Defined → Function Hooks
  6. Initialization 选择 initialize_tin_droplet
  7. 执行 Initialization

十六、复现图6~8的第一组算例

依次计算:

算例AR其他参数
10.5D=5, We=6.626,θ=140D^*=5,\ We=6.626,\theta=140^\circ
20.7同上
31.0同上
41.5同上
52.0同上

每个算例至少计算到 32 ms:

N=0.0322.5×106=12800 stepsN=\frac{0.032}{2.5\times10^{-6}} =12800\ {\rm steps}

建议设置:

  • Number of Time Steps:13000
  • Time Step Size:2.5e-6 s
  • Max Iterations:30

十七、基准算例的验收结果

如果模型设置正确,应该接近论文以下数据:

ARβmax\beta_{\max}tβmaxt_{\beta\max}, ms完全凝固时间,ms
0.51.3691.52527.00
0.71.3172.15028.25
1.01.3123.40028.50
1.51.3944.15026.00
2.01.4524.40024.25

对于 AR=1 的基准算例,论文描述为:

  • 1.525 ms:液滴形成球冠形铺展形态;
  • 3.025 ms:中心形成平滑隆起;
  • 3.400 ms:达到最大铺展因子 1.312;
  • 4.525 ms:中心形成凹盆;
  • 6.025 ms:中心流体向上运动;
  • 7.525 ms:中心达到较大高度;
  • 之后发生逐渐衰减的上下振荡;
  • 28.50 ms:完全凝固。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

你的结果满足下面范围时,可认为基本复现:

  • βmax\beta_{\max} 偏差小于 3%~5%;
  • 最大铺展时间偏差小于约 0.3~0.5 ms;
  • 完全凝固时间偏差小于约 1~2 ms;
  • 主要形貌变化顺序相同。

十八、提取铺展因子

论文对曲面定义:

β=LDe\beta=\frac{L}{D_e} L=Dsarcsin(2xwDs)L=D_s\arcsin\left(\frac{2x_w}{D_s}\right)

其中 xwx_w 是液滴接触线相对于轴线的径向坐标。

在当前 Fluent 坐标中:

xw=ywx_w=y_w

操作:

  1. 创建锡体积分数 0.5 的 Iso-Line;
  2. 找到该等值线与曲面壁相交的位置;
  3. 读取最大径向坐标 ywy_w
  4. 计算:
β=DsDearcsin(2ywDs)\beta= \frac{D_s}{D_e} \arcsin\left(\frac{2y_w}{D_s}\right)

基准工况:

  • Ds=0.010D_s=0.010 m
  • De=0.002D_e=0.002 m

所以:

β=5arcsin(2yw0.010)\beta= 5\arcsin\left(\frac{2y_w}{0.010}\right)

注意三角函数使用弧度。


十九、提取总体固相率

论文定义:

Sf=已凝固锡体积锡液滴总体积S_f= \frac{\text{已凝固锡体积}} {\text{锡液滴总体积}}

Fluent 中液相率为 ϕ\phi,所以:

Sf=αtin(1ϕ)dVαtindVS_f= \frac{\int \alpha_{\rm tin}(1-\phi)\,dV} {\int \alpha_{\rm tin}\,dV}

建立 Custom Field Function:

text
tin-solid = Volume-Fraction-tin * (1 - Liquid-Fraction)

再建立两个 Volume Integral:

Vsolid=αtin(1ϕ)dVV_{\rm solid} =\int \alpha_{\rm tin}(1-\phi)dV Vtin=αtindVV_{\rm tin} =\int\alpha_{\rm tin}dV

最后:

Sf=VsolidVtinS_f=\frac{V_{\rm solid}}{V_{\rm tin}}

实际计算中可以用:

Sf0.999S_f\geq0.999

作为“完全凝固”,因为数值上恰好达到 1 有时需要较长时间。


二十、完整参数矩阵

完成基准工况后,按以下矩阵复现整篇论文。

20.1 液滴形状,图6~8

固定:

  • D=5D^*=5
  • We=6.626We=6.626
  • v0=0.5v_0=0.5 m/s
  • θ=140\theta=140^\circ

改变:

AR=0.5, 0.7, 1.0, 1.5, 2.0AR=0.5,\ 0.7,\ 1.0,\ 1.5,\ 2.0

20.2 曲率,图9~11

固定:

  • We=6.626We=6.626
  • v0=0.5v_0=0.5 m/s
  • θ=140\theta=140^\circ

液滴形状:

AR=0.5, 1.0, 2.0AR=0.5,\ 1.0,\ 2.0

曲率:

D=2, 5, 10, 20, +D^*=2,\ 5,\ 10,\ 20,\ +\infty

对应:

DD^*DsD_s, mm半径,mm
242
5105
102010
204020
++\infty平板平板

球冠高度始终保持:

hs=1 mmh_s=1\ {\rm mm}

每个曲面的圆心:

xc=hsDs2x_c=h_s-\frac{D_s}{2}

20.3 Weber 数,图12~14

固定:

  • D=5D^*=5
  • θ=140\theta=140^\circ

液滴形状:

AR=0.5, 1.0, 2.0AR=0.5,\ 1.0,\ 2.0

速度和 Weber 数:

v0v_0, m/sWe
0.10.265
0.32.386
0.56.626
0.712.99
1.026.51

20.4 接触角,图15~17

固定:

  • D=5D^*=5
  • We=6.626We=6.626
  • v0=0.5v_0=0.5 m/s

液滴形状:

AR=0.5, 1.0, 2.0AR=0.5,\ 1.0,\ 2.0

接触角:

θ=40, 70, 90, 140, 160\theta=40^\circ,\ 70^\circ,\ 90^\circ,\ 140^\circ,\ 160^\circ

二十一、建议的实际复现顺序

不要一次建立全部 50 多个算例。按下面顺序最稳妥:

  1. 图5平板验证算例,常数表面张力;
  2. 图5平板验证算例,温度相关表面张力;
  3. 曲面基准算例 AR=1;
  4. 检查是否得到 βmax1.312\beta_{\max}\approx1.312tsolid28.50t_{\rm solid}\approx28.50 ms;
  5. 计算 AR=0.5 和 AR=2;
  6. 确认形貌和凝固时间趋势正确;
  7. 完成全部 AR;
  8. 完成曲率组;
  9. 完成 Weber 数组;
  10. 完成接触角组。

论文报告单个高分辨率算例约消耗 1500 core-hours,因此建议先用 20 μm 网格调试,确认设置后再用 10 μm 网格做最终结果。(关于椭球形金属液滴撞击曲面固体表面时的扩散与凝固过程的数值研究.pdf)

このQ&Aを共有