VOF 波公式

一阶 VOF 波

使用斯托克斯波理论的一阶近似对一阶波建模。此近似生成具有规则周期性正弦分布的波。

水平速度的方程为：

υ_{h} = a ω \cos (K \cdot x - ω t) e^{K z}

(260)

垂直速度的方程为：

υ_{v} = a ω \sin (K \cdot x - ω t) e^{K z}

(261)

表面高度的方程为：

η = a \cos (K \cdot x - ω t)

(262)

其中：

$a$ 为波幅值
$ω$ 为波频率
$K$ 为波矢量
$K$ 为波矢量的幅值
$z$ 为与平均水位的垂直距离。

波周期 $T$ 定义如下：

T = \frac{2 π}{ω}

(263)

波长 $λ$ 定义如下：

λ = \frac{2 π}{K}

(264)

一阶波在有限水深 $d$ 中的耗散关系（波周期 $T$ 与波长 $λ$ 之间）为：

T = {[\frac{g}{2 π λ} \tanh (\frac{2 π d}{λ})]}^{- 1 / 2}

(265)

而对于无限水深，耗散关系为：

λ = \frac{g T^{2}}{2 π}

(266)

波形与深度无关。

五阶 VOF 波

使用斯托克斯波理论的五阶近似对五阶波建模。与通过一阶方法生成的波相比，此波更接近于真正的波。波形和波相速度取决于水深、波高和水流。

五阶 VOF 波基于 Fenton 的作品 [113]。

Ursell 数 $U_{R}$ 定义为 [110]：

U_{R} = \frac{H λ^{2}}{d^{3}}

(267)

其中， $H$ 为波高， $λ$ 为波长， $d$ 为水深。

此波理论只对小于 40 的 Ursell 数有效。

椭圆余弦波

椭圆余弦波是 Korteweg–de Vries 方程的非线性精确周期波求解，该方程描述了平层上的波传播。Korteweg & de Vries [117] 获得了他们称之为“椭圆余弦”的周期求解，原因在于表面高度与雅可比椭圆函数 $c n ()$ 的平方成比例。椭圆余弦波用于描述与流体深度相比具有较长波长的表面重力波。椭圆余弦求解显示在自然浅水波中观察到的长平波谷和窄波峰。随着波长趋于无穷大，椭圆余弦求解描述了一种孤波。

椭圆余弦求解采用波高扩大的形式呈现结果。所使用的参数是相对于流体谷深 $H / h$ 的波高，以 $ϵ$ 表示。

级数以 $ϵ / m$ 的幂级数表示，而不是以 $ϵ$ 的幂级数表示 [112]。由于 $m$ 可小于 1，因此，监视 $ϵ / m$ 的幅值要优于在带有系数（为 $1 / m$ 中的多项式）的 $ϵ$ 中拥有一个幂级数。

存在完整的三阶求解。此求解更适用于较短和幅值更小的波，其中的参数 $m$ 小于 0.96。

对于较小的 Ursell 数， $m$ 的求解方案可能无法收敛并导致错误，指示特定椭圆余弦波具有无效的模型参数 $m$ 。文献 [118] 中提供了有关 $m$ 求解方法的更多详细信息。

符号 $c n$ 用于表示 $c n (α X / h | m) = c n (α (x - c t) / h | m)$ 。

表面高度计算

\frac{η}{h} = 1 + (\frac{ϵ}{m}) m {c n}^{2} + {(\frac{ϵ}{m})}^{2} (- \frac{3}{4} m^{2} {c n}^{2} + \frac{3}{4} m^{2} {c n}^{4}) + {(\frac{ϵ}{m})}^{3} ((- \frac{61}{80} m^{2} + \frac{111}{80} m^{3}) {c n}^{2} + (\frac{61}{80} m^{2} - \frac{53}{20} m^{3}) {c n}^{4} + \frac{101}{80} m^{3} {c n}^{6})

(268)

波计算框架中的水平流体速度

假设 $x$ 方向上的变化相对较小，且可采用缩放的无量纲变量 $α x / h$ 表示。 $α$ 为一个较小的物理量，表示 $x$ 方向上的相对较慢变化， $h$ 为谷深。

级数表示为 $α^{2}$ （具体为 $δ = 4 α^{2} / 3$ ）[114]。即使对于高波，结果也是精确的。

\frac{U}{\sqrt{g h}} = - 1 + δ (\frac{1}{2} - m + m {c n}^{2}) + δ^{2} (- \frac{19}{40} + \frac{79}{40} m - \frac{79}{40} m^{2} + {c n}^{2} (- \frac{3}{2} m + 3 m^{2}) - m^{2} {c n}^{4} + {(\frac{Y}{h})}^{2} (- \frac{3}{4} m + \frac{3}{4} m^{2} + {c n}^{2} (\frac{3}{2} m - 3 m^{2}) + \frac{9}{4} m^{2} {c n}^{4})) + δ^{3} (\frac{55}{112} - \frac{3471}{1120} m + \frac{7113}{1120} m^{2} - \frac{2371}{560} m^{3} + {c n}^{2} (\frac{71}{40} m - \frac{339}{40} m^{2} + \frac{339}{40} m^{3}) + {c n}^{4} (\frac{27}{10} m^{2} - \frac{27}{5} m^{3}) + \frac{6}{5} m^{4} {c n}^{6} + {(\frac{Y}{h})}^{2} (\frac{9}{8} m - \frac{27}{8} m^{2} + \frac{9}{4} m^{3} + {c n}^{2} (- \frac{9}{4} m + \frac{27}{2} m^{2} - \frac{27}{2} m^{3}) + {c n}^{4} (- \frac{75}{8} m^{2} + \frac{75}{4} m^{3}) - \frac{15}{2} m^{3} {c n}^{6}) + {(\frac{Y}{h})}^{4} (- \frac{3}{16} m + \frac{9}{16} m^{2} - \frac{3}{8} m^{3} + {c n}^{2} (\frac{3}{8} m - \frac{51}{16} m^{2} + \frac{51}{16} m^{3}) + {c n}^{4} (\frac{45}{16} m^{2} - \frac{45}{8} m^{3}) + \frac{45}{16} m^{3} {c n}^{6}))

(269)

产生前导项 $- 1$ ，原因在于将此波视为沿 $x$ 正向行进。流体沿相对于波的 $x$ 负向以波速阶的速度在波下流动。

垂直流体速度计算

垂直流体速度可通过质量守恒方程 $\partial U / \partial X + \partial V / \partial Y = 0$ 和 $d (c n (θ | m)) / d θ = - s n (θ | m) d n (θ | m)$ 从 Eqn. (269) 算得。每个含有 ${(Y / h)}^{i} {c n}^{j} (α X / h | m)$ （其中 $j > 0$ ）的项都替换为 $α s n () d n () (\frac{j}{i + 1}) \times {(Y / h)}^{i + 1} {c n}^{j - 1} ()$ 。因此，如果 Eqn. (269) 写为：

\frac{U}{\sqrt{g h}} = - 1 + \sum_{i = 1}^{5} δ^{i} \sum_{j = 0}^{i - 1} {(\frac{Y}{h})}^{2 j} \sum_{k = 0}^{i} {c n}^{2 k} () Φ_{i j k}

(270)

其中，每个系数 $Φ_{i j k}$ 均为参数 $i$ 中度数 $m$ 的多项式，则垂直速度分量为：

\frac{V}{\sqrt{g h}} = 2 α c n () s n () d n () \sum_{i = 1}^{5} δ^{i} \sum_{j = 0}^{i - 1} {(\frac{Y}{h})}^{2 j + 1} \sum_{k = 1}^{i} {c n}^{2 (k - 1)} () \frac{k}{2 j + 1} Φ_{i j k}

(271)

椭圆余弦波理论适用于大于 40 的 Ursell 数。

不规则 VOF 波

不规则波使用波谱（即海面垂直位移的功率谱密度函数）对风浪（短期海况）进行建模。不规则波可使用 Pierson-Moskowitz 频谱或 JONSWAP 频谱 [111]。这两种频谱均描述大多数恶劣海况下通常会发生的风浪状况。

Simcenter STAR-CCM+ 提供的 Pierson-Moskowitz 频谱（在文献中也称为双参数 Pierson-Moskowitz 或修正 Pierson-Moskowitz 频谱）是 Bretschneider 频谱 [108] 的改编表示：

S (ω) = A ω^{- 5} exp (- B ω^{- 5})

(272)

Pierson-Moskowitz 频谱 $S_{P M} (ω)$ 随后通过将 Eqn. (272) 中的参数 A 和 B 替换为 $A = \frac{5}{16} (H_{S}^{2} ω_{p}^{4})$ 和 $B = \frac{5}{4} (ω_{p}^{4})$ 获得。

其中：

$ω_{p} = (2 π) / (T_{p})$ 为角谱峰值频率
$H_{S}$ 为有效波高。

JONSWAP 频谱 $S_{J} (ω)$ 表示为：

S_{J} (ω) = A_{γ} {S_{P M} (ω) γ}^{exp (- 0.5 {(\frac{ω - ω_{p}}{{σ ω}_{p}})}^{2})}

(273)

其中：

$S_{P M} (ω)$ 为 Pierson-Moskowitz 频谱
$γ$ 为无量纲峰形参数
$A_{γ} = 1 - 0.287 ln (γ)$ 为归一化因子
$σ$ 为频谱宽度参数：
$σ_{a}$ 对于 $ω \leq ω_{p}$

$σ_{b}$ 对于 $ω > ω_{p}$ 。

波阻尼

可通过对垂直运动施加阻力来产生波阻尼。Simcenter STAR-CCM+ 实施 Choi 和 Yoon [109] 的方法，为 $w$ 速度方程添加了一个阻力项：

S_{z}^{d} = - ρ (f_{1} + f_{2} | w |) \frac{e^{κ} - 1}{e^{1} - 1} w

(274)

且：

κ = {(\frac{x - x_{sd}}{x_{ed} - x_{sd}})}^{n_{d}}

(275)

其中：

$x_{sd}$ 为波阻尼的起点（在 x 方向上传播）
$x_{ed}$ 为波阻尼的终点（边界）
$f_{1}$ 、 $f_{2}$ 和 $n_{d}$ 为阻尼模型的参数
$w$ 为垂直速度分量。

从 Eqn. (274) 中，阻尼强度指示器 $r^{d}$ 由以下公式给出：

r^{d} = \frac{e^{κ} - 1}{e^{1} - 1}

(276)

波力

根据 Kim 等人的方法 ([116])，Simcenter STAR-CCM+ 可将 3D 流体模拟与由 VOF 波指定的理论求解结合起来。

在指定距离内，迫使离散化纳维-斯托克斯方程的解趋向于另一种解（例如理论解或简化数值解），从而使用更小的求解域来减少计算工作量。这种力作用还能消除因逐步力作用的阻尼特性而导致的边界处表面波反射相关问题。

波力仅适用于动量。没有添加相源或湍流源。通过向以下形式的传输（动量）方程添加源项来获得波力：

q_{ϕ} = - γ ρ (ϕ - ϕ^{*})

(277)

其中：

$γ$ 为力系数
$ρ$ 为流体密度
$ϕ$ 为输运方程的当前求解
$ϕ^{*}$ 为力求解所接近的值。

当需要将求解固定为特定的 $ϕ^{*}$ 值时，将较大数值用于 $γ$ 。随后，可以忽略离散方程的其余部分。

来自 Eqn. (277) 的源项与一个可变力系数一同应用到指定的力区域。例如，下图是两个重叠求解域的示意图。3D 纳维-斯托克斯方程在阴影区域内的一个域中求解。内部区域（3D 纳维-斯托克斯）中没有施加任何力，但在外部区域（力区域）中，沿求解域边界激活了力源项。力区域的宽度因边界而异。力区域的最佳宽度取决于正在建模的问题。

力系数平稳地改变，从力区域内缘处的零到边界（力区域外缘）处的最大值 $γ$ 。

γ = - γ_{0} \cos^{2} (π x^{*} / 2)

(278)

由于力区域通常不重合，因此体积耦合可以为单向或双向。

从 Eqn. (278) 中，力强度指示器 $r^{f}$ 由以下公式给出：

r^{f} = \cos^{2} (\frac{π x^{*}}{2})

(279)