船舶横摇运动matlab仿真

根据平稳随机过程理论 分析船舶运动功率谱 模拟船舶横摇运动

clc
clear
syms w we;
w=0.4:0.05:1.6;mu=90;v=6;%kb=1;kt=1;%航速6m/s，航向90度，即横浪，cos（90）=0
T=3;B=10;Sw=0.785;%船宽10米，吃水3米水线面系数假设为0.785
we=w+w.^2.*v/9.8;
for i=1:24
    delta_we（i）=we（i+1）-we（i）;
end
delta_we（25）=delta_we（24）;
z=0.194375;
s_kesi=0.778./（w.^5）.*exp（-z./（w.^4））; %含风速=15m/s即三一波高=4m在内了 %2.47./（w.^6）.*exp（-2.*9.8.^2./24.5.^2./（w.^2））;
s_alpha=w.^4.*s_kesi./（9.8^2）;%波倾角谱
%以下进行波倾角修正
%按照船宽和吃水定kbkt进而求s_alpha_v;
kb=exp（-w.^2.*T./（6.24*pi））;
kt=1-sqrt（Sw）.*B.^2.*w.^4./（1.56^2*16*pi^4）;
s_kesi_v=0.778./（w.^5）.*exp（-0.194375./（w.^4））./（1+2.*w./9.8.*v.*cos（mu*pi/180））;
s_alpha_v=w.^4.*s_kesi_v./（9.8^2）.*kb.^2.*kt.^2;

%figure;hold on;set（gca‘XLim‘[0.22]）;
%plot（ws_kesi‘r‘ws_alpha‘g--‘）;
s_alpha_ap=sqrt（2.*s_alpha_v.*0.05）;%波倾角幅值
s_alpha_v_ap=sqrt（2.*s_alpha_v.*delta_we）;
%产生0~2*pi之间均匀分布的随机变量epsilon
y=zeros（251）;y（1）=1977326743;b=7;M=10^10;
for i=1:25
    y（i+1）=mod（7*y（i）M）;
    x（i）=y（i）/M;
end
epsilon=2*pi*x;

alpha=0;alpha_v=0;
t=0:0.1:30*pi;
for i=1:25
    a=s_alpha_ap（i）.*sin（w（i）*t+epsilon（i））;
    alpha=alpha+a;
    aa=s_alpha_v_ap（i）.*sin（we（i）*t+epsilon（i））;
    alpha_v=alpha_v+aa;
end

miu_roll=0.1614;%横摇阻尼，可以变化，0~2；
omg_roll=0.8;%横摇固有频率给定？根据船型计算？
A=we./omg_roll;%遭遇频率和固有频率之比
W_roll_2=1./（（1-A.^2）.^2+4.*（miu_roll^2）.*（A.^2））;%幅频响应函数
s_roll=W_roll_2.*s_alpha_v;%横摇角能谱用s_alpha_v代替s_alpha
%plot（ws_roll）

s_roll_ap=s_alpha_v_ap.*sqrt（W_roll_2）;
epsilon_roll=atan（2.*miu_roll.*A./（1-A.^2））;
hyj=0;
for i=1:25
    r=s_roll_ap（i）.*sin（we（i）*t-epsilon_roll（i））;
    hyj=hyj+r;
end

hold on;
figure（1）plot（t10.*alpha‘r‘t10*alpha_v‘b‘）;xlabel（‘t （s）‘）;ylabel（‘\alpha‘）;  %绘制波倾角时域图
legend（‘波倾角‘‘有效波倾角‘）;
%plot（t20.*alpha_v）;xlabel（‘t （s）‘）;ylabel（‘\alpha_v‘）;
figure（2）plot（t10.*hyj）;xlabel（‘t （s）‘）;ylabel（‘横摇角‘）;   %绘制横摇角时域图

%subplot（313）;
%plot（t20.*hyj）;xlabel（‘t （s）‘）;ylabel（‘横摇角‘）;   %绘制横摇角时域图

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件       2075  2010-10-14 08:27  roll.m

----------- ---------  ---------- -----  ----

                 2075                    1