三种海浪谱仿真

海浪谱：PM谱、JONSWAP谱、布氏谱

三种海浪仿真，并设置坐标轴字体和字号

%%JONSWAP谱

[x,S_Jon] = JONSWAP1(5 ,1 * 10^6); %5m/s
[x1,S_Jon1] = JONSWAP1(10,1 * 10^6); %10m/s
[x2,S_Jon2] = JONSWAP1(15 ,1 * 10^6); %15m/s

figure;
hold on;
plot(x,S_Jon,'LineWidth',1.5);
hold on;
plot(x1,S_Jon1,'LineWidth',1.5);
plot(x2,S_Jon2,'LineWidth',1.5);
box on
hold off
set(gca,'fontsize',14,'FontName','Times New Roman');%设置刻度字体大小
legend({'5m/s', '10m/s','15m/s'},'Location','northeast','FontSize',14);
xlabel('\omega/Hz','fontsize',14,'FontName','Times New Roman');
ylabel('S(\omega)/m^2/Hz','fontsize',14,'FontName','Times New Roman');

%% 布氏谱
[w,S_B1] = Bretsch(20);
[w,S_B2] = Bretsch(25);
[w,S_B3] = Bretsch(30);
figure;

hold on;
plot(w,S_B1,'LineWidth',1.5);
hold on;
plot(w,S_B2,'LineWidth',1.5);
plot(w,S_B3,'LineWidth',1.5);
box on

hold off
set(gca,'fontsize',14,'FontName','Times New Roman');
legend({'20m/s', '25m/s','30m/s'},'Location','northeast','FontSize',14);
xlabel('\omega/Hz','fontsize',14,'FontName','Times New Roman');
ylabel('S(\omega)/m^2/Hz','fontsize',14,'FontName','Times New Roman');
% saveas(gca,'步式谱1.png');
%% PM谱
[w,S_PM1] = PM1(15);
[w,S_PM2] = PM1(18);
[w,S_PM3] = PM1(20);
figure;

hold on;
plot(w,S_PM1,'LineWidth',1.5);
hold on;
plot(w,S_PM2,'LineWidth',1.5);
plot(w,S_PM3,'LineWidth',1.5);
box on
hold off
set(gca,'fontsize',14,'FontName','Times New Roman');
legend({'15m/s', '18m/s','20m/s'},'Location','northeast','FontSize',14);
xlabel('\omega/Hz','fontsize',14,'FontName','Times New Roman');
ylabel('S(\omega)/m^2/Hz','fontsize',14,'FontName','Times New Roman');
% saveas(gca,'PM谱1.png');

%% 几个谱函数

%% PM谱
function [w,s] = PM1(U)
i = 0;
for w = 0:0.01:1
    i = i+ 1;
    s(i) = 0.0081 * 9.81^2 / w^5 * exp((-0.74) * (9.81 / ( U * w))^4);
end
w = 0:0.01:1;

end

%% JONSWAP谱
function [w,S] = JONSWAP1(U,X)
i = 0;
 gamma = 3.3;
 for omega = 0.01:0.01:1
     
     i = i+1;
     g = 9.81;
     alpha = 0.076 * (g * X/ U^2)^-0.22;
     omega_p = 22 * (g/U) * (g * X/U^2)^-0.33;
     if omega > omega_p
         sigma = 0.09;
     else
         sigma = 0.07;
     end
     K = exp(-(omega - omega_p)^2 /(2 * sigma^2 * omega_p^2));
     
     S(i) = alpha * g^2 / omega^2 * exp(-5/4 * (omega_p / omega)^4) * gamma^K;
 end
    w = 0.01:0.01:1;
end

function [w,S] = Bretsch(U)
%% 布氏谱
i = 0;
W_m = 9.81 / U;
H = 0.22 * U^2 / 9.81;

for w = 0:0.01:1
    i = i+1;
    S(i) = 1.25 / 4 * W_m^4 / w^5 * H * exp(-1.25*(W_m/w)^4);
end
    
w = 0:0.01:1;

end

如有错误请指正！