【数学建模】基于matlab三维海浪模型仿真【含Matlab源码 1159期】
【摘要】
一、获取代码方式
获取代码方式1: 完整代码已上传我的资源:【数学建模】基于matlab三维海浪模型仿真【含Matlab源码 1159期】
获取代码方式2: 通过订阅紫极神光博客付费专栏,凭支付凭证,...
一、获取代码方式
获取代码方式1:
完整代码已上传我的资源:【数学建模】基于matlab三维海浪模型仿真【含Matlab源码 1159期】
获取代码方式2:
通过订阅紫极神光博客付费专栏,凭支付凭证,私信博主,可获得此代码。
备注:
订阅紫极神光博客付费专栏,可免费获得1份代码(有效期为订阅日起,三天内有效);
二、部分源代码
clear all;
close all;
nhFig = 0; % Figure Number;
SeaRegLx = 40e+3; % Sea region length, unit: m
SeaRegLy = 40e+3; % Sea region length, unit: m
% 根据频域采样定理,在确定采样周期即波束的足迹宽度之后,
% 可以利用时域采样的采样点数目(通过空间频谱的宽度即根据时域采样定理计算出采样点数),实现频域采样.
% 对于海洋波谱而言,波谱宽度是非常宽的,因此其波谱宽度一般是指包括了大部分能量(或者说主要的能量部分)的不完全谱宽度
% 所以这里事先根据波谱函数计算主要能量部分对应的谱宽度
xSampleStep = 50; % unit: m
ySampleStep = 50;
xNs = round( SeaRegLx / xSampleStep );
yNs = round( SeaRegLy / ySampleStep );
x = linspace( -SeaRegLx / 2, SeaRegLx / 2, xNs );
y = linspace( -SeaRegLy/ 2, SeaRegLy / 2, yNs );
% 海浪谱仿真
g = 9.8; % gravity acceleration
% swell wave spectrum parameter
SwellWaveLength = 1000; % 涌浪波长
KswellWavePeak = 2 * pi / SwellWaveLength;
SwellAngled=0;
SwellAngle = SwellAngled/ 180 * pi; % 涌浪与观测向夹角
KxswellWavePeak = KswellWavePeak * cos( SwellAngle );
KyswellWavePeak = KswellWavePeak * sin( SwellAngle );
SigmaKx = 2.5e-3; % 涌浪谱宽度
SigmaKy = 2.5e-3;
SigmaHSwell = 2; %涌浪波高
% wind wave spectrum parameter
WindAngle = 45 / 180 * pi; % 风向与观测向夹角
U10 = 12; % 10米高处海面风速 % 5m/s, 10m/s, 15m/s 的风浪谱宽度分别为 1.5, 0.4, 0.15.对应的最大空间采样间隔为2m, 10m, 20m
KwindPeak = g / ( 1.2 * U10 )^2;
% sea wave spectrum parameter
NxSeaWave = xNs; % 频域采样点数,方便傅立叶变换
NySeaWave = yNs;
KxSeaWave = 2 * pi / SeaRegLx * ( -xNs / 2 : 1 : xNs / 2 - 1 );
KySeaWave = 2 * pi / SeaRegLy * ( -yNs / 2 : 1 : yNs / 2 - 1 );
KxSeaWaveTicks = KxSeaWave;
KySeaWaveTicks = KySeaWave;
KxSeaWave = ( KxSeaWave == 0 ) * ( max( KxSeaWave ) * 1e-16 ) + KxSeaWave; % avoid divided by 0
KySeaWave = ( KySeaWave == 0 ) * ( max( KySeaWave ) * 1e-16 ) + KySeaWave; % avoid divided by 0
% swell wave spectrum
SpectrumSwell = zeros( NySeaWave, NxSeaWave );
Temp1 = ones( NySeaWave, 1 ) * ( KxSeaWave - KxswellWavePeak ) / SigmaKx;
Temp2 = ( KySeaWave - KyswellWavePeak )' / SigmaKy * ones( 1, NxSeaWave );
SpectrumSwell = SigmaHSwell^2 / 2 / pi / SigmaKx / SigmaKy * exp( -0.5 * ( Temp1.^2 + Temp2.^2 ) );
clear Temp1;
clear Temp2;
figure;
colormap(gray(256));
image( KxSeaWave, KySeaWave, 256 - 255 / ( max( max( abs( SpectrumSwell ) ) ) - min( min ( abs( SpectrumSwell ) ) ) ) * ( abs( SpectrumSwell ) - min( min ( abs( SpectrumSwell ) ) ) ) );
axis('xy');
xlabel( 'kx:X方向波数');
ylabel( 'ky:Y方向波数');
title( '涌浪谱');
KxSeaWaveMatrix = ones( NySeaWave, 1 ) * KxSeaWave;
KySeaWaveMatrix = KySeaWave' * ones( 1, NxSeaWave );
KSeaWaveTemp1 = ( sqrt( KxSeaWaveMatrix.^2 + KySeaWaveMatrix.^2 ) );
Fwk1 = exp( - 1.22 * ( sqrt( KSeaWaveTemp1 ./ KwindPeak ) -1 ).^2 );
HKKpx1 = 1.24 * ( ( KSeaWaveTemp1 / KwindPeak < 0.31 ) & ( KSeaWaveTemp1 / KwindPeak >= 0 ) );
HKKpx2 = 2.61 * ( ( KSeaWaveTemp1 / KwindPeak ).^0.65 ) .* ( ( KSeaWaveTemp1 / KwindPeak < 0.9 ) & ( KSeaWaveTemp1 / KwindPeak >= 0.31 ) );
HKKpx3 = 2.28 * ( ( KSeaWaveTemp1 / KwindPeak ).^( -0.65 ) ) .* ( KSeaWaveTemp1 / KwindPeak >= 0.9 );
HKKp1 = HKKpx1 + HKKpx2 + HKKpx3;
Temp1x = 1.62 * 1e-3 * U10 / ( g^0.5 ) ./ ( KSeaWaveTemp1 ).^3.5;
Temp2x = exp( -( KwindPeak ./ KSeaWaveTemp1 ).^2 ) .* ( 1.7 .^ Fwk1 );
Temp3x = ( HKKp1 .* ( sech( ( HKKp1 .* ( atan( ( KySeaWaveMatrix ./ KxSeaWaveMatrix ) ) - WindAngle ) ) ) ).^2 );
% SpectrumWind1 = ( Temp1x .* Temp2x .* Temp3x )';
fid=fopen('sea_top_wl1000_wh50m.dat','w')
fprintf(fid,'%10.5f\n',real(wh));
status=fclose(fid);
fid=fopen('sea_top_wl1000_wh50m.dat','r');
sea_dem=fscanf(fid,'%g',[800,800]);
fclose(fid);
%image( x, y, 256 - 255 / ( max( max( real( sea_dem ) ) ) - min( min ( real( sea_dem ) ) ) ) * ( real( sea_dem ) - min( min ( real( sea_dem ) ) ) ) );
%
mesh(sea_dem)
colorbar
% view(0,90)
%colormap(gray)
% This function computes the non-coherent integration improvment
% factor using the empirical formula defind in Eq. (2.48)
fact1 = 1.0 + log10( 1.0 / pfa) / 46.6;
fact2 = 6.79 * (1.0 + 0.235 * pd);
fact3 = 1.0 - 0.14 * log10(np) + 0.0183 * (log10(np))^2;
impr_of_np = fact1 * fact2 * fact3 * log10(np);
return
% SpectrumWind = SpectrumWind1';
SpectrumWind = Temp1x .* Temp2x .* Temp3x;
clear KSeaWaveTemp1 Fwk1 HKKpx1 HKKpx2 HKKpx3 HKKp1 Temp1x Temp2x Temp3x;
- 1
- 2
- 3
- 4
- 5
- 6
- 7
- 8
- 9
- 10
- 11
- 12
- 13
- 14
- 15
- 16
- 17
- 18
- 19
- 20
- 21
- 22
- 23
- 24
- 25
- 26
- 27
- 28
- 29
- 30
- 31
- 32
- 33
- 34
- 35
- 36
- 37
- 38
- 39
- 40
- 41
- 42
- 43
- 44
- 45
- 46
- 47
- 48
- 49
- 50
- 51
- 52
- 53
- 54
- 55
- 56
- 57
- 58
- 59
- 60
- 61
- 62
- 63
- 64
- 65
- 66
- 67
- 68
- 69
- 70
- 71
- 72
- 73
- 74
- 75
- 76
- 77
- 78
- 79
- 80
- 81
- 82
- 83
- 84
- 85
- 86
- 87
- 88
- 89
- 90
- 91
- 92
- 93
- 94
- 95
- 96
三、运行结果
四、matlab版本及参考文献
1 matlab版本
2014a
2 参考文献
[1]李昕.MATLAB数学建模[M].清华大学出版社.2017
[2]王健,赵国生.MATLAB数学建模与仿真[M].清华大学出版社.2016
[3]余胜威.MATLAB数学建模经典案例实战[M].清华大学出版社.2015
文章来源: qq912100926.blog.csdn.net,作者:海神之光,版权归原作者所有,如需转载,请联系作者。
原文链接:qq912100926.blog.csdn.net/article/details/119220599
【版权声明】本文为华为云社区用户转载文章,如果您发现本社区中有涉嫌抄袭的内容,欢迎发送邮件进行举报,并提供相关证据,一经查实,本社区将立刻删除涉嫌侵权内容,举报邮箱:
cloudbbs@huaweicloud.com
- 点赞
- 收藏
- 关注作者
评论(0)