2016国赛A题——系泊系统问题粒子群算法求解

这篇具有很好参考价值的文章主要介绍了2016国赛A题——系泊系统问题粒子群算法求解。希望对大家有所帮助。如果存在错误或未考虑完全的地方,请大家不吝赐教,您也可以点击"举报违法"按钮提交疑问。

前言

        本文参考了一篇有关系泊系统的知乎博文,链接如下:2016年全国大学生数学建模竞赛A题系泊系统的设计详解(含建模和程序详解) - 知乎 (zhihu.com)

        因为是练习一下写代码,建模的话基本是按照该博文的方法建的,也推荐大家对数模机理问题有兴趣的同学关注一下这个博主(b站有号)。本文只在未知数的求解过程上做一个创新。(原文使用的是随机搜索法。)

        本文就第一问写一下求解的过程,重点分析粒子群求解吃水深度的部分(其他部分与知乎博主的大致相同)。写作本文的目的是为了记录下自己用粒子群算法求解一些复杂问题的过程,也给大家一些启发,让大家感受到启发式算法除了求解高维函数外的一些其他应用。

题目描述

粒子群算法试题,matlab,启发式算法,算法,线性代数

 粒子群算法试题,matlab,启发式算法,算法,线性代数

 系统的受力分析方程建立

浮标的受力方程和锚链钢管的受力方程如下图所示:

粒子群算法试题,matlab,启发式算法,算法,线性代数

            通过Matlab的实时脚本可以求解该参数方程组。在假定和都是以知量的情况下,该方程组是有解的。(但其实是未知量,也是后文粒子群算法的优化目标)。

 

粒子群算法试题,matlab,启发式算法,算法,线性代数

        对于方程2而言,只要知道了前一个受力对象的拉力和受力角度和倾斜角度,就可以解出下一个对象的拉力和受力角度和倾斜角度。是一个迭代方程组。具体的公式求解可以通过Matlab实时脚本求解。

优化目标

其中是每节钢管或锚链的长度,为水深,都是固定值。为受力方程计算结果。 

程序求解过程

浮力函数和风力函数

        风力的公式题中给出了,浮力的公式大家应该也都会,不做赘述。

function F = windy(h,v_wind)
    D = 2; 
    H = 2;
    S = D * (H-h);
    F = 0.625 * S * v_wind^2;
end
function f_float = float(h)
    g = 9.8;
    rho_sea = 1025;
    d = 2;
    f_float = rho_sea*pi*d^2/4*h*g;
end

 初始方程和迭代方程函数

        求解方程组中的量,并写出各个变量的迭代关系式。

function [theta1,T1] = equation_initial(h,v_wind)
    m = 1000; g = 9.8;
    F_float = float(h);
    F_wind = windy(h,v_wind);
    theta1 = atand(F_wind/(F_float-m*g));
    T1 = F_wind/sind(theta1);
end
function [theta2,T2,phi1] = equation_iteration(T1,theta1,m,V)
    g = 9.8;
    rho_sea = 1025;
    F_float = rho_sea*V*g;
    theta2 = atand(T1*sind(theta1)/(F_float+T1*cosd(theta1)-m*g));
    T2 = T1*sind(theta1)/sind(theta2);
    phi1 = atand(T1*sind(theta1)/(T1*cosd(theta1)+0.5*(F_float-m*g)));
end

 

锚链选型函数

function [dL,dm] = mooring_choose(i)
    mooring = [0.078 3.2;
        0.105 7;
        0.12 12.5;
        0.15 19.5;
        0.18 28.12];
    dL = mooring(i,1);
    dm = mooring(i,2);
end

方程求解函数 

function [deep,r,x,y,phi] = mooring_system(h,v_wind,M,i,L_all)
    rho_Fe = 7800; % 锚链密度
    [dL,dm] = mooring_choose(i); % 锚链选型
    N = round(L_all/dL); % 锚链节数
    m = [10*ones(1,4),100,dL*dm*ones(1,N)]; % 各节质量
    m(6) = m(6) + M; % 球的质量
    D = [0.05*ones(1,4),0.3];  % 钢管和钢桶
    L = ones(1,5); % 钢管和钢桶的长度
    V = zeros(size(m)); V(1:5) = pi*0.5*D.^2.*L; V(6:end) = m(6:end)/rho_Fe; % 体积
    L = [L,dL*ones(1,N)]; % 钢管和锚链长度整合
    N = N + 5; 
    T = zeros(1,N+1); % 受力
    theta = zeros(1,N+1); % 角度
    phi = zeros(1,N); % 倾斜角度
    [theta(1),T(1)] = equation_initial(h,v_wind); % 初始方程计算
    for i = 1:N-1
        [theta(i+1),T(i+1),phi(i)] = equation_iteration(T(i),theta(i),m(i),V(i)); % 迭代方程计算
        if (phi(i) < 0||phi(i)>90) % 角度超出范围跳出循环
            break;
        end
    end
    phi(i:end) = 90; % 超出范围后的倾斜角度视为90度
    deep = sum(L.*cosd(phi)); % 深度
    r = sum(L.*sind(phi)); % 浮游半径
    L = fliplr(L); % L倒置
    phi = fliplr(phi); % phi倒置
    x = cumsum(L.*sind(phi)); % x坐标计算
    y = cumsum(L.*cosd(phi)); % y坐标计算
end

* 粒子群算法

优化函数

        根据模型应该是求最小值,但由于该算法用于求最大值,因此添加一个负号。

function y = Obj_fun1(x,v_wind,M,i,L_all,H)
  deep = mooring_system(x,v_wind,M,i,L_all); % 计算总深度
  y = -abs(deep - H); % 计算误差 
end    

算法局限性

       粒子群算法作为一种群体智能算法,与遗传算法相比有一些不同之处,主要在于该方法只能求解连续性问题。因此整数解的优化问题就不适合用粒子群算法求解。

参数设置

        粒子数量、变量个数、函数上下界根据优化函数进行调整;最大速度可以根据自变量的定义域进行调整。其他参数可以基本不改。

    %% 粒子群算法中的预设参数
    n = 100; % 粒子数量
    narvs = 1; % 变量个数
    c1 = 2;  % 每个粒子的个体学习因子,也称为个体加速常数
    c2 = 2;  % 每个粒子的社会学习因子,也称为社会加速常数
    w = 0.9;  % 惯性权重
    K = 50;  % 迭代的次数
    vmax = 1.2; % 粒子的最大速度
    x_lb = 0; % x的下界
    x_ub = 2; % x的上界

粒子的速度和位置初始化

        位置初始化其实就是在定义域内随机产生初始解。解的个数为定义的粒子的个数。粒子群法通过定义多个粒子进行优化迭代,可以很好地避免陷入局部最优解。

    %% 初始化粒子的位置和速度
    x = zeros(n,narvs);
    for m = 1: narvs
        x(:,m) = x_lb(m) + (x_ub(m)-x_lb(m))*rand(n,1);   
    end
    v = -vmax + 2*vmax .* rand(n,narvs);  % 随机初始化粒子的速度

    %% 计算适应度
    fit = zeros(n,1); 
    for m = 1:n  
        fit(m) = Obj_fun1(x(m,:),v_wind,M,i,L_all,H);  
    end
    pbest = x;  
    ind = find(fit == max(fit), 1); 
    gbest = x(ind,:); 

粒子速度和位置更新

        粒子速度的变化是粒子群的核心点。粒子速度和位置的求解公式如下:

粒子群算法试题,matlab,启发式算法,算法,线性代数
摘自清风老师的课件

        其实,粒子的速度就是函数中自变量在一次位置更新时候移动的范围,与牛顿迭代法中新解的产生方式是相同的。 但不同之处在于,牛顿迭代是基于倒数确定移动的步长,是一种纯数学的方法;而粒子群算法是基于上述核心公式,模拟了鸟的认知过程,是一种社会仿生学的方法。

        然后就是新解与旧解适应度的比较,保留适应度高的那个解作为最新的最佳位置。

 %% 迭代K次来更新速度与位置
    fitnessbest = ones(K,1); 
    for d = 1:K  % 开始迭代,一共迭代K次
        for m = 1:n   % 依次更新第i个粒子的速度与位置
            v(m,:) = w*v(m,:) + c1*rand(1)*(pbest(m,:) - x(m,:)) + c2*rand(1)*(gbest - x(m,:)); 
            for j = 1: narvs
                if v(m,j) < -vmax(j)
                    v(m,j) = -vmax(j);
                elseif v(m,j) > vmax(j)
                    v(m,j) = vmax(j);
                end
            end
            x(m,:) = x(m,:) + v(m,:); % 更新第i个粒子的位置
            % 如果粒子的位置超出了定义域,就对其进行调整
            for j = 1: narvs
                if x(m,j) < x_lb(j)
                    x(m,j) = x_lb(j);
                elseif x(m,j) > x_ub(j)
                    x(m,j) = x_ub(j);
                end
            end
            fit(m) = Obj_fun1(x(m,:),v_wind,M,i,L_all,H);  
            if fit(m) > Obj_fun1(pbest(m,:),v_wind,M,i,L_all,H)  
                pbest(m,:) = x(m,:);  
            end
            if  fit(m) > Obj_fun1(gbest,v_wind,M,i,L_all,H)                
                gbest = pbest(m,:);   
            end
        end
        fitnessbest(d) = Obj_fun1(gbest,v_wind,M,i,L_all,H); 
    end
    h = gbest;

总代码

        粒子群优化函数总代码如下:

function h = scope(H,v_wind,M,i,L_all)
    %% 粒子群算法中的预设参数
    n = 100; % 粒子数量
    narvs = 1; % 变量个数
    c1 = 2;  % 每个粒子的个体学习因子,也称为个体加速常数
    c2 = 2;  % 每个粒子的社会学习因子,也称为社会加速常数
    w = 0.9;  % 惯性权重
    K = 50;  % 迭代的次数
    vmax = 1.2; % 粒子的最大速度
    x_lb = 0; % x的下界
    x_ub = 2; % x的上界

    %% 初始化粒子的位置和速度
    x = zeros(n,narvs);
    for m = 1: narvs
        x(:,m) = x_lb(m) + (x_ub(m)-x_lb(m))*rand(n,1);    
    end
    v = -vmax + 2*vmax .* rand(n,narvs);  % 随机初始化粒子的速度

    %% 计算适应度
    fit = zeros(n,1);  
    for m = 1:n  
        fit(m) = Obj_fun1(x(m,:),v_wind,M,i,L_all,H);   
    end
    pbest = x;   % 初始化这n个粒子迄今为止找到的最佳位置
    ind = find(fit == max(fit), 1);  % 找到适应度最大的那个粒子的下标
    gbest = x(ind,:);  % 定义所有粒子迄今为止找到的最佳位置

    %% 迭代K次来更新速度与位置
    fitnessbest = ones(K,1);  % 初始化每次迭代得到的最佳的适应度
    for d = 1:K  
        for m = 1:n   % 依次更新第i个粒子的速度与位置
            v(m,:) = w*v(m,:) + c1*rand(1)*(pbest(m,:) - x(m,:)) + c2*rand(1)*(gbest - x(m,:));  % 更新第i个粒子的速度
            % 如果粒子的速度超过了最大速度限制,就对其进行调整
            for j = 1: narvs
                if v(m,j) < -vmax(j)
                    v(m,j) = -vmax(j);
                elseif v(m,j) > vmax(j)
                    v(m,j) = vmax(j);
                end
            end
            x(m,:) = x(m,:) + v(m,:); % 更新第i个粒子的位置
            % 如果粒子的位置超出了定义域,就对其进行调整
            for j = 1: narvs
                if x(m,j) < x_lb(j)
                    x(m,j) = x_lb(j);
                elseif x(m,j) > x_ub(j)
                    x(m,j) = x_ub(j);
                end
            end
            fit(m) = Obj_fun1(x(m,:),v_wind,M,i,L_all,H);  % 重新计算第i个粒子的适应度
            if fit(m) > Obj_fun1(pbest(m,:),v_wind,M,i,L_all,H)  
                pbest(m,:) = x(m,:); 
            end
            if  fit(m) > Obj_fun1(gbest,v_wind,M,i,L_all,H)  
                gbest = pbest(m,:);   
            end
        end
        fitnessbest(d) = Obj_fun1(gbest,v_wind,M,i,L_all,H); 
    end
    h = gbest;
end

系统绘图函数

        得到的图像如下:

        以风速36m/s的情况为例:

粒子群算法试题,matlab,启发式算法,算法,线性代数

 

function [] = mooring_plot(deep,phi,long,x,y,v_wind,h0,H)
    figure
    %% 锚链绘制
    plot(x,y,'k') % 锚链
    
    hold on 
    plot(linspace(-5,30,100),y(end)*ones(1,100),'b:') % 水面
    hold on 
    plot(linspace(-5,30,100),H*ones(1,100),'b-') % 理想水面
    axis([0,x(end)+5,0,y(end)+5]) % 标定图像尺寸
    title(['v_{wind}=',num2str(v_wind),'m/s']) % 以风速标题
    hold on 
    
 
    %% 钢球绘制
    y0=linspace(y(end-5),y(end-5)-2,100);
    x0=x(end-5)*ones(size(y0));
    plot(x0,y0,'k');
    t=0:0.01:2*pi;
    x2=x0(end)+0.5*sin(t);y2=y0(end)+0.5*cos(t);
    fill(x2,y2,'k')
    
    %% 标注拖地长度和锚链与海底夹角
    f=find(y>0);
    x0=x(f(1));
    if(phi(3)==90&&x0>1)
        text(-1,2,['拖地长度=',num2str(x0),'m']);% 计算拖地长度
    else
        phi_end=90-phi(3);
        f1=find(y>1);
        x0=x(f1(1));
        text(x0+3,1,['锚链海底夹角=',num2str(phi_end),'°'])% 不拖地计算锚链与海底的夹角
    end
  
    %% 标注水深,游动半径和吃水深度
    phi_bucket=phi(end-4);
    text(x(end-30),y0(end)-2,['倾斜角度=',num2str(phi_bucket),'°'],'fontsize',10)
    text(-1,y(end)-1,['水面深度误差=',num2str(H-deep),'m'],'fontsize',10)
    text(x(end-140),y(end)+2,['浮游半径=',num2str(long),'m'],'fontsize',10)
    text(x(end)+1,y(end)-0.5,['吃水深度=',num2str(h0),'m'],'fontsize',10)
    
    %% 浮标绘制
    x1=[x(end)-0.5,x(end)-0.5,x(end)+0.5,x(end)+0.5];
    y1=[y(end)-h0,y(end)+2-h0,y(end)+2-h0,y(end)-h0];
    patch(x1,y1,'k')
    
    %% 锚绘制
    x3=[-0.5,-0.5,0.5,0.5];
    y3=[0,1,1,0];
    patch(x3,y3,'k')
    axis([-5,30,0,25])
end

主函数

clear;close all;clc;
v_wind = [12,24,36]; M = 1200; H = 18;L_all = 22.05;
h = scope(H,v_wind(3),M,2,L_all); 
[deep,r,x,y,phi] = mooring_system(h,v_wind(3),M,2,L_all);
mooring_plot(deep,phi,r,x,y,v_wind(3),h,H);

总结

       第一问的求解就是这样。对于后面2问,第二问中,只需在粒子群算法中增加一个变量M即可。对于第三问,是个多目标规划问题,用粒子群算法显然不合适,可以采用NSGA2算法,但那个与粒子群算法无关了,此处不探讨了。希望大家看了这篇文章可以受到启发。文章来源地址https://www.toymoban.com/news/detail-627703.html

到了这里,关于2016国赛A题——系泊系统问题粒子群算法求解的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处: 如若内容造成侵权/违法违规/事实不符,请点击违法举报进行投诉反馈,一经查实,立即删除!

领支付宝红包 赞助服务器费用

相关文章

  • 单目标应用:粒子群优化算法(PSO)求解微电网优化MATLAB

    微电网优化模型介绍: 微电网多目标优化调度模型简介_IT猿手的博客-CSDN博客 (1)部分代码 (2)部分结果

    2024年02月06日
    浏览(64)
  • 多目标粒子群算法求解帕累托前沿Pareto,Pareto的原理,测试函数100种求解之21

    目录 背影 parte前沿的定义 注意事项 基于多目标粒子群的帕累托前沿求解 主要参数 MATLAB代码 效果图 结果分析 展望 在目标优化过程种,很多时候都两个或者多个目标,并且目标函数不能同时达到最优,鱼与熊掌不可兼得,这个时候可以通过求解帕累托前沿,通过帕累托前沿

    2024年02月11日
    浏览(42)
  • 【路径规划】粒子群算法求解机器人障碍物环境的Voronoi图路径规划【含GUI Matlab源码 3748期】

    粒子群算法(Particle Swarm Optimization, PSO)可以用于栅格地图上机器人的最短路径规划。在这种问题中,栅格地图被划分为离散的单元格,每个单元格可以是阻挡或可通过的区域。机器人需要从起始位置移动到目标位置,避免碰到阻挡。 PSO算法中,通过使用一群粒子来搜索最优

    2024年01月17日
    浏览(63)
  • 路径之谜 2016年国赛 深度优先搜索

    目录 解题思路 AC代码: 题目描述 小明冒充 XX 星球的骑士,进入了一个奇怪的城堡。 城堡里边什么都没有,只有方形石头铺成的地面。 假设城堡地面是 n×n 个方格。如下图所示。 按习俗,骑士要从西北角走到东南角。可以横向或纵向移动,但不能斜着走,也不能跳跃。每

    2024年02月08日
    浏览(34)
  • 第十二届蓝桥杯国赛试题及解析

    第一题 *选择题严禁使用程序验证设s =HiLanQiao\\\',运行以下哪个选项代码可以输出“LanQiao”子串( A )。 A、print(s[-7:]) B、print(s/-6:-11) C、print(s1-7:01) D、print(s[-7:-1]) 第二题 *选择题严禁使用程序验证已知a=2021.0529,运行以下哪个选项代码可以输出“2021.05” ( B )A、print( .2f1\\\'.format(a

    2024年02月08日
    浏览(46)
  • 基于TSP(旅行商)问题的混合粒子群算法 附直接运行代码

    如果对粒子群一点都不知道的可以看看上文标准粒子群算法, 想看代码的直接去下面1.4标题 即可 链接:(105条消息) 自己对粒子群算法的理解(附matlab直接运行代码)(二维)_吕浩轩的博客-CSDN博客_二维粒子群算法​​​​​​h 好现在开始正文: 标准粒子群通过追随个体极值和

    2023年04月16日
    浏览(76)
  • 关于旅行商问题的多种算法求解

    一、问题描述 旅行商问题是指旅行家要旅行n个城市,要求每个城市经历一次且仅经历一次然后回到出发城市,并要求所走路程最短。 首先通过所给出的一个无向图,即n个顶点,m个无向边,每条边有一个权值代表两个点之间的距离,要求把每一个点都走一遍并回到原点,求

    2024年02月04日
    浏览(50)
  • 遗传算法及基于该算法的典型问题的求解实践

        遗传算法是一个很有用的工具,它可以帮我们解决生活和科研中的诸多问题。最近在看波束形成相关内容时了解到可以用这个算法来优化阵元激励以压低旁瓣,于是特地了解和学习了一下这个算法,觉得蛮有意思的,于是把这两天关于该算法的学习和实践的内容总结成了

    2024年03月21日
    浏览(48)
  • Matlab蚁群算法求解旅行商问题

    目录 问题展现 解决代码 代码1 输出结果 代码2 输出结果 代码3 输出结果 假设有一个旅行商人要拜访全国 31 个省会城市,他需要选择所要走的路径,路径的限制是每个城市只能拜访一次,而且最后要回到原来出发的城市。路径的选择要求是:所选路径的路程为所有路径之中的

    2023年04月15日
    浏览(61)
  • 面试项目算法 - 数桥问题python求解

    本项目基于一个流行的日本谜题--\\\"Hashiwokakero\\\"、\\\"Hashi \\\"或 \\\"Bridges\\\"。你需要编写一个程序来解决这个谜题,并简要说明你所使用的算法和数据结构。 程序的输入将是一个由数字和点组成的矩形数组,例如: 每个数字代表一个 \\\"岛屿\\\",而点代表岛屿之间的空隙(水域)。大于 9

    2024年03月20日
    浏览(42)

觉得文章有用就打赏一下文章作者

支付宝扫一扫打赏

博客赞助

微信扫一扫打赏

请作者喝杯咖啡吧~博客赞助

支付宝扫一扫领取红包,优惠每天领

二维码1

领取红包

二维码2

领红包