Gerstner Waves 简介
Gerstner Waves 最初是作为流体力学的一个精确解析解（严格推导）被提出，假设每一个水质点都在做圆周运动，那么在忽略粘性和假设流体不可压缩的情况下，这组方程可以完美满足流体力学的基本方程。
数学公式
单个波的位移公式
单个波的位移：
$$ P' = \left\{ \begin{aligned} x' &= Q A D_x \cos(\theta) \\ y' &= A \sin(\theta) \\ z' &= Q A D_z \cos(\theta) \end{aligned} \right. $$其中相位 $\theta$ ：
$$ \begin{aligned} \theta &= k(D \cdot P_{xz}) - \omega t \\ &= k(D_x x + D_z z) - \omega t \end{aligned} $$其中 $\omega$ 由色散公式得：
$$ \omega^2 = g k \tanh(kh) \implies \left\{ \begin{aligned} &\text{深水: } &\lim_{h \to 0} w^2 &= gk &\implies w &= \sqrt{gk} \\ &\text{浅水: } &\lim_{h \to \infty} w^2 &= gk^2h &\implies w &= k\sqrt{gh} \end{aligned} \right. $$| 参数 | 含义 |
|---|---|
| $Q$ | 锐利度（波陡度） |
| $A$ | 波高（振幅） |
| $D=(D_x,D_z)$ | 波传播方向 |
| $k$ | 波数（由 $k = \frac{2\pi}{\lambda}$ 得，$\lambda$ 为波长） |
| $h$ | 水深 |
| $\omega$ | 角频率（可通过色散公式计算） |
| $t$ | 时间 |
| $g$ | 重力加速度（ $g \approx 9.81$ ） |
写成绝对位置公式：
$$ P(x,z) = \begin{pmatrix} x &+& Q A D_x \cos(\theta) \\ y &+& A \sin(\theta) \\ z &+& Q A D_z \cos(\theta) \end{pmatrix} $$多波叠加
把多个波的位移累加即可， 写成绝对位置公式：
$$ P(x,z) = \begin{pmatrix} x &+& \sum_{i=1}^{N} Q_i A_i D_{x_i} \cos(\theta_i) \\ y &+& \sum_{i=1}^{N} A_i \sin(\theta_i) \\ z &+& \sum_{i=1}^{N} Q_i A_i D_{z_i} \cos(\theta_i) \end{pmatrix} $$在实际工程中，更多的是多波叠加。
重新计算法线
单波
顶点移动后，为了正确光照，需要重新计算法线。原理：
计算偏导、构建切线：
$$ \begin{aligned} T_x &= \frac{\partial P'}{\partial x} = \begin{pmatrix} 1 - k Q A D_x^2 sin(\theta) \\ k A D_x \cos(\theta) \\ - k Q A D_x D_z \sin(\theta) \end{pmatrix} \\ T_z &= \frac{\partial P'}{\partial z} = \begin{pmatrix} - k Q A D_x D_z \sin(\theta) \\ k A D_z \cos(\theta) \\ 1 - k Q A D_z^2 \sin(\theta) \end{pmatrix} \end{aligned} $$叉乘并归一化得到法线：
$$ N = \text{normalize}(T_x \times T_z) $$根据向量叉乘公式：
$$ \begin{aligned} N_x &= (T_z)_y (T_x)_z - (T_z)_z (T_x)_y \\ N_y &= (T_z)_z (T_x)_x - (T_z)_x (T_x)_z \\ N_z &= (T_z)_x (T_x)_y - (T_z)_y (T_x)_x \end{aligned} $$带入计算 X 分量 和 Z 分量，复杂项都会被抵消：
$$ \begin{aligned} N_x &= - k A D_x \cos(\theta) \\ N_z &= - k A D_z \cos(\theta) \end{aligned} $$带入计算 Y 分量：
$$ N_y = 1 - k Q A \sin(\theta) \cdot (D_x^2 + D_z^2) $$由于 $D$ 是波的传播方向单位向量，所以 $D_x^2 + D_z^2 = 1$，代入后得到：
$$ N_y = 1 - k Q A \sin(\theta) $$对于单个波，它的未归一化法线：
$$ N = \begin{pmatrix} - k D_x A \cos(\theta) \\ 1 - k Q A \sin(\theta) \\ - k D_z A \cos(\theta) \end{pmatrix} $$提示
在实际的工程中，有些项目会仅计算 X 分量和 Z 分量，Y 分量直接取1。虽然数学上不对，但在波浪不极端陡峭的情况下，人眼看不出区别。
多波
对于多波叠加的法线，由于微积分的线性性质，只需要将每个波求出的偏导数项（即偏移量）累加，然后再执行一次最终的计算即可。
根据和的导数等于导数的和，计算偏导、构建切线：
$$ \begin{aligned} T_x &= \frac{\partial P'}{\partial x} = \begin{pmatrix} 1 - \sum_{i=1}^{N} k_i Q_i A_i D_{x_i}^2 sin(\theta_i) \\ \sum_{i=1}^{N} k_i A_i D_{x_i} \cos(\theta_i) \\ - \sum_{i=1}^{N} k_i Q_i A_i D_{x_i} D_{z_i} \sin(\theta_i) \end{pmatrix} \\ T_z &= \frac{\partial P'}{\partial z} = \begin{pmatrix} - \sum_{i=1}^{N} k_i Q_i A_i D_{x_i} D_{z_i} \sin(\theta_i) \\ \sum_{i=1}^{N} k_i A_i D_{z_i} \cos(\theta_i) \\ 1 - \sum_{i=1}^{N} k_i Q_i A_i D_{z_i}^2 \sin(\theta_i) \end{pmatrix} \end{aligned} $$叉乘并归一化得到法线：
$$ N = \text{normalize}(T_x \times T_z) $$提示
在实际的工程中，有些项目会直接将每个波在 X 和 Z 轴的法线偏转量加起来，Y 分量直接取1。虽然数学上不对，但在波浪不极端陡峭的情况下，人眼看不出区别。
代码示例
注意
以下两个示例均为数学精确解，实际工程中请看情况根据上面的 tip 优化法线计算。
波数据结构
// 波属性
struct Wave
{
float amplitude; // 振幅 A
float wavenumber; // 波数 k
float2 direction; // 方向 D
float steepness; // 锐利度 Q
float speed; // 角频率 ω
};
// 计算的结果
struct WaveResult
{
float3 position; // 偏移后的位置
float3 normal; // 偏移后的法线
};单波函数
WaveResult ComputeGerstnerWave(Wave wave, float3 gridPosition, float time)
{
WaveResult result;
// 1. 预计算基础参数
/// 相位 Phase/θ
float phase = wave.wavenumber * dot(wave.direction, gridPosition.xz) - wave.speed * time;
float cosP = cos(phase);
float sinP = sin(phase);
/// Q * A
float QA = wave.steepness * wave.amplitude;
// 2. 计算顶点偏移
result.position.x = gridPosition.x + QA * wave.direction.x * cosP;
result.position.z = gridPosition.z + QA * wave.direction.y * cosP;
result.position.y = wave.amplitude * sinP;
// 3. 计算法线
result.normal.x = -wave.wavenumber * wave.direction.x * wave.amplitude * cosP;
result.normal.z = -wave.wavenumber * wave.direction.y * wave.amplitude * cosP;
result.normal.y = 1.0 - wave.wavenumber * wave.steepness * wave.amplitude * sinP;
result.normal = normalize(result.normal);
return result;
}多波叠加
#define N 8;
WaveResult ComputeMultGerstnerWaves(Wave waves[N], float3 gridPosition, float time)
{
// 1. 初始化累加用于累加的变量
/// 累积偏导（注意有初始值）
float3 dPdx = float3(1.0, 0.0, 0.0);
float3 dPdz = float3(0.0, 0.0, 1.0);
/// 累积位置偏移
float3 totalOffset = float3(0.0, 0.0, 0.0);
// 2. 遍历波
for (int i = 0; i< N; i++)
{
Wave wave = waves[i];
// 1. 预计算基础参数
/// 相位 Phase/θ
float phase = wave.wavenumber * dot(wave.direction, gridPosition.xz) - wave.speed * time;
float cosP = cos(phase);
float sinP = sin(phase);
/// Q * A
float QA = wave.steepness * wave.amplitude;
// 2. 计算顶点偏移
totalOffset.x += QA * wave.direction.x * cosP;
totalOffset.z += QA * wave.direction.y * cosP;
totalOffset.y += wave.amplitude * sinP;
// 3. 计算法线，累积偏导数
/// 预计算部分数据
float a_k_cos = wave.amplitude * wave.wavenumber * cosP;
float qa_k_sin = QA * wave.wavenumber * sinP;
/// 累积 X 方向的偏导
dPdx.x -= wave.direction.x * wave.direction.x * qa_k_sin;
dPdx.y += wave.direction.x * a_k_cos;
dPdx.z -= wave.direction.x * wave.direction.y * qa_k_sin;
/// 累积 Z 方向的偏导
dPdz.x -= wave.direction.x * wave.direction.y * qa_k_sin;
dPdz.y += wave.direction.y * a_k_cos;
dPdz.z -= wave.direction.y * wave.direction.y * qa_k_sin;
}
// 3. 计算最终结果
WaveResult result;
result.normal = normalize(cross(dPdz, dPdx));
result.position = gridPosition + totalOffset;
return result;
}提示
在实际工程中，建议使用 cbuffer 代替 waves 数组。
精简优化法线计算
直接按单波的方法计算并累积，最后归一化。
float3 normal = float3(0, 1.0, 0);
for (...) {
...
normal.x -= wave.wavenumber * wave.direction.x * wave.amplitude * cosP;
normal.z -= wave.wavenumber * wave.direction.y * wave.amplitude * cosP;
}
normal = normalize(normal);程序化构建波信息
经典 Gerstner Wave 实现里，每个 wave 的参数都需要手动定义。为避免麻烦，多数项目会随机生成 wave 参数，但遵循海洋统计规律。例如： 方向：随机 ±30° （模拟风向） 波长：指数分布 振幅：与波长相关 速度：由深水波公式
for (int i = 0; i < N; i++)
{
wave.direction = normalize(windDir + random2(-0.3,0.3));
wave.wavelength = lerp(minWave, maxWave, random());
wave.amplitude = wave.wavelength * 0.03;
wave.speed = sqrt(9.81 * (2 * PI / wave.wavelength));
wave.steepness = 0.7 / (wave.amplitude * N);
waves[i] = wave;
}这样得到的海浪：
- 方向相似（风向）
- 尺度不同
- 不会整齐排列 视觉上就自然很多。
提示
可以使用 fBm 调制各个参数，这样会更符合现实，甚至可以大力出奇迹做出海面的粗糙质感。
波破碎
真实海浪在陡度过高时会破。
$$ Q k A > \text{某个阈值} $$可以据此判断是否发生了波破碎，可以加一些 平化处理 或 产生白沫。