JONSWAP 公式 (Joint North Sea Wave Project Spectrum) 是海洋工程和海岸工程中广泛使用的风浪频谱模型。
JONSWAP 数学公式
公式和参数
$$ S(\omega) = \alpha g^2 \omega^{-5} \exp\left[ -\frac{5}{4} \left( \frac{\omega_p}{\omega} \right)^{4} \right] \gamma^{\exp\left[ -\frac{(\omega - \omega_p)^2}{2 \sigma^2 \omega_p^2} \right]} $$| 符号 | 典型值 | 含义 | 备注 |
|---|---|---|---|
| $S(\omega)$ | - | 频谱能量密度$m^2 \cdot s$ | 表示在特定角频率$\omega$下，海浪所包含的能量密度 |
| $\omega$ | - | 角频率 | |
| $\omega_p$ | - | 谱峰频率 | 能量密度最大的频率 |
| $\alpha$ | - | 能量尺度参数 | 决定了波谱的总能量强度 |
| $\gamma$ | $\text{典型值 } 3.3\quad (1.0 \le \gamma \le 7.0)$ | 峰值提升因子 | 控制谱峰尖锐程度 |
| $\sigma$ | $\begin{cases} 0.07, & \text{当 } \omega \le \omega_p \\ 0.09, & \text{当 } \omega > \omega_p \end{cases}$ | 谱宽参数 | 决定了谱峰在 $\omega_p$ 附近的胖瘦程度 |
| $g$ | $9.81$ | 重力加速度 |
相关参数计算
下面用到的参数示意：
| 符号 | 含义 |
|---|---|
| $U_{10}$ | 海平面上方 10 米处的风速 (米/秒) |
| $X$ | 吹程 / 风区长度 (米)，通常在 10000 到 100000 之间 |
谱峰频率 $\omega_p$
经验公式：
$$ \omega_p = 22 \left( \frac{g^2}{U_{10} X} \right)^{1/3} $$能量尺度参数 $\alpha$
$\alpha$ 通常需要实测，如果没有特定海域的实测数据，可通过经验公式计算或直接取参考值。 经验公式：
$$ \alpha = 0.076 \left( \frac{g X}{U_{10}^2} \right)^{-0.22} $$在深海且风力充沛的情况下，α 的典型值通常在 0.0081 到 0.01 之间。 在一些简化计算中，如果不考虑吹程，常沿用 Pierson-Moskowitz 谱的参考值 0.0081。
封装函数
C++ 预算部分
在 C++ 中根据经验公式预先计算 $\omega_p$ 和 $\alpha$
float WindSpeed = 15.0f; // 风速
float Fetch = 50000.0f; // 吹程
float JOmegaP = 22.0f * FMath::Pow((GRAVITY * GRAVITY) / (WindSpeed * Fetch), 1.0f / 3.0f);
float JAlpha = 0.076f * FMath::Pow((WindSpeed * WindSpeed) / (Fetch * GRAVITY), 0.22f);HLSL 部分
float JOmegaP;
float JAlpha;
float JGamma;
float JONSWAP(float omega) {
// 计算 sigma
float sigma = (omega <= JOmegaP) ? 0.07 : 0.09;
// 计算 gamma 的指数
float r = exp(
-(omega - JOmegaP) * (omega - JOmegaP) /
2 * sigma * sigma * JOmegaP * JOmegaP
);
// 优化计算 omega^-5（转为1个倒数指令+3个乘法指令）
float _invOmega = 1.0 / omega;
float _invOmega2 = _invOmega * _invOmega;
float invOmega5 = _invOmega2 * _invOmega2 * _invOmega;
// 优化计算 (omega_p / omega) ^ 4（转为1个倒数复用+3个乘法指令）
float _omegaPOverOmega = JOmegaP * _invOmega; // 复用之前的倒数
float _omegaPOverOmega2 = _omegaPOverOmega * _omegaPOverOmega;
float omegaPOverOmega4 = _omegaPOverOmega2 * _omegaPOverOmega2;
// 最终计算
return JAlpha * GRAVITY * GRAVITY * invOmega5 * exp(-1.25 * omegaPOverOmega4) * pow(JGamma, r);
}在波矢量 $k$ 空间下计算 JONSWAP
在上一章中，我们将像素坐标空间转到了波矢量，然而 JONSWAP 的自变量是角频率，故需要对其进行转换。
波矢量 $k$ 转 角频率 $\omega$ 并计算
这个转换非常简单，在深水的情况下可以直接使用深水色散关系转换，
$$ \omega = \sqrt{g k} $$注意这里的 $k$ 取的是模长，
float omega = sqrt(GRAVITY * length(k));
float jonswapResult = JONSWAP(omega);- TODO 浅水
计算结果转回波矢量 $k$ 空间
这个推导根本看不懂，这里直接贴代码了：
float dOmegaDk = sqrt(g) / (2.0 * pow(kLength, 1.5));
jonswapResult *= dOmegaDk;提示
如果只考虑深水，可将深水波公式代入 $S(\omega)$ 得到 关于波数的 JONSWAP 公式 $S(k)$