五阶stokes波的MATLAB代码

斯托克斯波是一种用于描述有限振幅波动的数学模型，它考虑了波的非线性和色散特性。斯托克斯波理论是由英国科学家乔治·斯托克斯在19世纪提出的，用于解释和预测海洋波浪的行为。

五阶斯托克斯波是斯托克斯波理论的一个扩展，它包括了直到第五阶的非线性项。这个模型可以更精确地描述波浪的形态，包括波峰、波谷以及波的陡峭程度，该理论在海洋工程、海岸工程和船舶设计等领域有着重要的应用，它可以帮助工程师更准确地预测波浪对结构物的作用力，以及波浪对船舶运动的影响......(省略公式，自己查书找）

话不多说，直接上代码：

% 五阶stokes波高示意图
% write by YANG Liu 2024/10/22
clc;
clear;
close all;
% 定义参数
k = 10.04;
T = 0.6;
d = 0.5;
% [x, t] = meshgrid(-60:0.1:60, -60:0.2:60);
t= -60:0.2:60;
lamda = 1;

% 计算 hyperbolic cosine 和 sine
c = cosh(k * d);
s = sinh(k * d);

% 系数
A11 = 1 / s;
A13 = -c^2 * (5 * c^2 + 1) / (8 * s^5);
A15 = -(1184 * c^10 - 1440 * c^8 - 1992 * c^6 + 2641 * c^4 - 249 * c^2 + 18) / ...
      (1536 * s^11);
A22 = 3 / (8 * s^4);
A24 = (192 * c^8 - 424 * c^6 - 312 * c^4 + 480 * c^2 - 17) / (768 * s^10);
A33 = (13 - 4 * c^2) / (64 * s^7);
A35 = (512 * c^12 + 4224 * c^10 - 6800 * c^8 - 12808 * c^6 + 16704 * c^4 - 3154 * ...
       c^2 + 107) / (4096 * s^13) / (6 * c^2 - 1);
A44 = (80 * c^6 - 816 * c^4 + 1338 * c^2 - 197) / (1536 * s^10) / (6 * c^2 - 1);
A55 = -(2880 * c^10 - 72480 * c^8 + 324000 * c^6 - 432000 * c^4 + 163470 * ...
        c^2 - 16245) / (61440 * s^11) / (6 * c^2 - 1) / (8 * c^4 - 11 * c^2 + 3);

B22 = (2 * c^2 + 1) * c / (4 * s^3);
B24 = c * (272 * c^8 - 504 * c^6 - 192 * c^4 + 322 * c^2 + 21) / (384 * s^9);
B44 = c * (768 * c^10 - 488 * c^8 - 48 * c^6 + 48 * c^4 + 106 * c^2 - 21) / ...
      (384 * s^9) / (6 * c^2 - 1);
B33 = 3 * (8 * c^6 + 1) / (64 * s^6);
B35 = (88128 * c^14 - 208224 * c^12 + 70848 * c^10 + 54000 * c^8 - ...
       21816 * c^6 + 6264 * c^4 - 54 * c^2 - 81) / (12288 * s^12) / (6 * c^2 - 1);
B55 = (192000 * c^16 - 262720 * c^14 + 83680 * c^12 + 20160 * c^10 - 7280 * c^8 + ...
       7160 * c^6 - 1800 * c^4 - 1050 * c^2 + 225) / (12288 * s^10) / (6 * c^2 - 1) / ...
       (8 * c^4 - 11 * c^2 + 3);

% 计算波速和角频率
C0s = 9.8 * tanh(k * d);
wc = sqrt(C0s * (1 + lamda^2 * (8 * c^4 - 8 * c^2 + 9) / 8 / s * 4 + ...
                lamda^4 * (3840 * c^12 - 4096 * c^10 + 2592 * c^8 - 1008 * c^6 + ...
                           5944 * c^4 - 1830 * c^2 + 147) / 512 / s * 10 / (6 * c^2 - 1)));
omega = 2 * pi / T;

% 计算波高
% theta = (k * x - omega * t);
theta = (k * 0 - omega * t);
E = zeros(1, 5);
E(1) = lamda / k;
E(2) = (lamda^2 * B22 + lamda^4 * B24) / k;
E(3) = (lamda^3 * B33 + lamda^5 * B35) / k;
E(4) = lamda^4 * B44 / k;
E(5) = lamda^5 * B55 / k;

% eta = zeros(size(x));
eta = zeros(size(t));
for n = 1:5
    eta = eta + E(n) * cos(n * theta);
end

% 绘制波高图
figure;
% surf(x, t, eta, 'EdgeColor', 'none');
plot(t,eta)
% xlabel('x (m)');ylabel('t(s)');zlabel('\eta (m)');
 xlabel('t (s)');ylabel('\eta (m)');
 axis tight
title('5th Order Stokes Wave Elevation');