精华内容
下载资源
问答
  • 数学建模源码集锦-粒子群算法的寻优算法应用实例
  • 数学建模源码集锦-基于混合粒子群算法的TSP搜索算法应用实例
  • 数学建模源码集锦-基于粒子群算法的多目标搜索算法应用实例
  • 数学建模源码集锦-基于动态粒子群算法的动态环境寻优算法应用实例
  • 数学建模源码集锦-基于粒子群算法的PID控制器优化方法应用实例
  • 数学建模:现代优化算法之粒子群算法 《数学建模算法与应用》在现代优化算法中没有提及,依照组内进度要求,所以需要在此节补充有关“粒子群算法”的内容。 数学建模数学建模:现代优化算法之粒子群算法前言一、...

    数学建模:现代优化算法之粒子群算法

    《数学建模算法与应用》在现代优化算法中没有提及,依照组内进度要求,所以需要在此节补充有关“粒子群算法”的内容。



    前言

    先来看一段视频:Matlab 粒子群算法PSO实例学习(有代码和详细注释)

    来自 “百科”: 粒子群优化算法(Particle Swarm optimization,PSO)又翻译为粒子群算法、微粒群算法、或微粒群优化算法。是通过模拟鸟群觅食行为而发展起来的一种基于群体协作的随机搜索算法。通常认为它是群集智能 (Swarm intelligence, SI) 的一种。它可以被纳入多主体优化系统(Multiagent Optimization System, MAOS)。粒子群优化算法是由Eberhart博士和kennedy博士发明。

    一、算法简介

    参考CSDN博主「lx青萍之末」的原创文章
    最优化算法之粒子群算法(PSO)

    二、算法步骤

    1.初始化

    PSO初始化为一群随机粒子(随机解),然后通过迭代找到最优解,在每一次迭代中,粒子通过跟踪两个“极值”来更新自己。第一个就是粒子本身所找到的最优解,这个解叫做个体极值pBest,另一个极值是整个种群找到的最优解,这个极值是全局极值gBest。另外也可以不用整个种群而只是用其中一部分最优粒子的邻居,那么在所有邻居中的极值就是局部极值。

    2.更新规则

    在这里插入图片描述
    公式(1)的第一部分称为【记忆项】,表示上次速度大小和方向的影响;公式(1)的第二部分称为【自身认知项】,是从当前点指向粒子自身最好点的一个矢量,表示粒子的动作来源于自己经验的部分;公式(1)的第三部分称为【群体认知项】,是一个从当前点指向种群最好点的矢量,反映了粒子间的协同合作和知识共享。粒子就是通过自己的经验和同伴中最好的经验来决定下一步的运动。以上面两个公式为基础,形成了PSO的标准形式。
    在这里插入图片描述

    三、PSO算法的流程图和伪代码

    在这里插入图片描述

    四、PSO算法举例

    【例一】
    在这里插入图片描述
    【解】
    在这里插入图片描述

    【例二】 本部分参考CSDN博主「nightmare_dimple」的原创文章
    原文链接:https://blog.csdn.net/nightmare_dimple/article/details/74331679

    对于函数 f=xsin(x)cos(2x)-2xsin(3x) ,求其在区间[0,20]上该函数的最大值。
    【解】

    clc;clear;close all;
    %% 初始化种群
    f= @(x)x .* sin(x) .* cos(2 * x) - 2 * x .* sin(3 * x); % 函数表达式
    figure(1);ezplot(f,[0,0.01,20]);
    N = 50;                         % 初始种群个数
    d = 1;                          % 空间维数
    ger = 100;                      % 最大迭代次数     
    limit = [0, 20];                % 设置位置参数限制
    vlimit = [-1, 1];               % 设置速度限制
    w = 0.8;                        % 惯性权重
    c1 = 0.5;                       % 自我学习因子
    c2 = 0.5;                       % 群体学习因子 
    for i = 1:d
        x = limit(i, 1) + (limit(i, 2) - limit(i, 1)) * rand(N, d);%初始种群的位置
    end
    v = rand(N, d);                  % 初始种群的速度
    xm = x;                          % 每个个体的历史最佳位置
    ym = zeros(1, d);                % 种群的历史最佳位置
    fxm = zeros(N, 1);               % 每个个体的历史最佳适应度
    fym = -inf;                      % 种群历史最佳适应度
    hold on
    plot(xm, f(xm), 'ro');title('初始状态图');
    figure(2)
    %% 群体更新
    iter = 1;
    record = zeros(ger, 1);          % 记录器
    while iter <= ger
         fx = f(x) ; % 个体当前适应度   
         for i = 1:N      
            if fxm(i) < fx(i)
                fxm(i) = fx(i);     % 更新个体历史最佳适应度
                xm(i,:) = x(i,:);   % 更新个体历史最佳位置
            end 
         end
    if fym < max(fxm)
            [fym, nmax] = max(fxm);   % 更新群体历史最佳适应度
            ym = xm(nmax, :);      % 更新群体历史最佳位置
     end
        v = v * w + c1 * rand * (xm - x) + c2 * rand * (repmat(ym, N, 1) - x);% 速度更新
        % 边界速度处理
        v(v > vlimit(2)) = vlimit(2);
        v(v < vlimit(1)) = vlimit(1);
        x = x + v;% 位置更新
        % 边界位置处理
        x(x > limit(2)) = limit(2);
        x(x < limit(1)) = limit(1);
        record(iter) = fym;%最大值记录
    %     x0 = 0 : 0.01 : 20;
    %     plot(x0, f(x0), 'b-', x, f(x), 'ro');title('状态位置变化')
    %     pause(0.1)
        iter = iter+1;
    end
    figure(3);plot(record);title('收敛过程')
    x0 = 0 : 0.01 : 20;
    figure(4);plot(x0, f(x0), 'b-', x, f(x), 'ro');title('最终状态位置')
    disp(['最大值:',num2str(fym)]);
    disp(['变量取值:',num2str(ym)]);
    

    结果如下:
    在这里插入图片描述
    由上图可以看出算法已成功找出了最优解,其最优解为18.3014,而其最大值为32.1462。
    如果想看粒子群算法中粒子的搜索过程的可以将代码中注释掉的三行代码放上去。效果是这样的:

    在这里插入图片描述
    【例三】
    版权声明:本部分参考CSDN博主「箫韵」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
    原文链接:https://blog.csdn.net/qq_41053605/article/details/88924287

    (1)构建目标函数。这里使用只有两个参数的函数,这样便于画出三维图型。

    function y=A11_01(x)
    y=x(1)^2+x(2)^2-x(1)*x(2)-10*x(1)-4*x(2)+60;
    

    首先从直观上看看这个函数:

    x1=-15:1:15;
    x2=-15:1:15;
    [x1,x2]=meshgrid(x1,x2);
    y=x1.^2+x2.^2-x1.*x2-10.*x1-4.*x2+60;
    mesh(x1,x2,y);
    

    在这里插入图片描述
    图 1 目标函数的三维网格图

    (2)整体实现代码

    clear 
    clc
    
    %绘制原图        图1目标函数的三维网格图
    x1=-15:1:15;
    x2=-15:1:15;
    [x1,x2]=meshgrid(x1,x2);
    y=x1.^2+x2.^2-x1.*x2-10.*x1-4.*x2+60;
    mesh(x1,x2,y);    
    hold on;
    
    %%预设参数
    n=100; %粒子群的规模
    d=2; %变量个数
    c1=2;
    c2=2;
    w=0.9;%权重一般设为0.9
    K=50; %迭代次数
    
    %%分布粒子的取值范围及速度的取值范围
    x=-10+20*rand(n,d);  %x在[-10,10]中取值
    v=-5+10*rand(n,d); %v在[-5,5]中取值
    
    %%计算适应度
    fit=zeros(n,1);
    for j=1:n
        fit(j)=A11_01(x(j,:));
    end
    pbest=x;
    ind=find(min(fit)==fit);
    gbest=x(ind,:);
    h=scatter3(x(:,1),x(:,2),fit,'o');  %2 粒子的初始分布图
    
    %%更新速度与位置
    for i=1:K
        for m=1:n
           v(m,:)=w*v(m,:) + c1*rand*(pbest(m,:)-x(m,:)) + c2*rand*(gbest-x(m,:));%rand是[0,1]随机数
           
           v(m,find(v(m,:)<-5))=-5;%这里发现速度小于-5时取-5
           v(m,find(v(m,:)>5))=5;%这里发现速度大于5时取5
           
           x(m,:)=x(m,:)+0.5*v(m,:);
           x(m,find(x(m,:)<-10))=-10;%这里发现位置小于-10时取-10
           x(m,find(x(m,:)>10))=10;%这里发现位置大于10时取10
           
           %重新计算适应度
           fit(m)=A11_01(x(m,:));
           if x(m,:)<A11_01(pbest(m,:))
               pbest(m,:)=x(m,:);
           end
           if A11_01(pbest(m,:))<A11_01(gbest)
               gbest=pbest(m,:);
           end
        end
        fitnessbest(i)=A11_01(gbest);
        pause(0.01);    %为了直观,每更新一次,暂停0.01秒
        h.XData=x(:,1);
        h.YData=x(:,2);
        h.ZData=fit;
    end
    % hold off;
    % plot(fitnessbest);
    % xlabel('迭代次数');
    

    (3)粒子的初始分布图

    在这里插入图片描述
    图2 粒子的初始分布图

    (4)粒子群移动图
    在这里插入图片描述
    图 3 粒子群移动图

    (5)迭代过程
    在这里插入图片描述
    图 4 迭代过程图
    从图中可以看出,迭代次数基本上在13,14的时候就达到最优,即目标函数取得最小值为8.

    五、PSO算法的demo

    #include <iostream>
    #include <vector>
    #include <cmath>
    #include <map>
    #include <algorithm>
    #include <random>
    #include <ctime>
    #include <Eigen/Dense>
    using namespace Eigen;
    using namespace std;
    
    const int dim = 1;//维数
    const int p_num = 10;//粒子数量
    const int iters = 100;//迭代次数
    const int inf = 999999;//极大值
    const double pi = 3.1415;
    //定义粒子的位置和速度的范围
    const double v_max = 4;
    const double v_min = -2;
    const double pos_max = 2;
    const double pos_min = -1;
    //定义位置向量和速度向量
    vector<double> pos;
    vector<double> spd;
    //定义粒子的历史最优位置和全局最优位置
    vector<double> p_best;
    double g_best;
    //使用eigen库定义函数值矩阵和位置矩阵
    Matrix<double, iters, p_num> f_test;
    Matrix<double, iters, p_num> pos_mat;
    
    //定义适应度函数
    double fun_test(double x)
    {
        double res = x * x + 1;
        return res;
    }
    
    //初始化粒子群的位置和速度
    void init()
    {
        //矩阵中所有元素初始化为极大值
        f_test.fill(inf);
        pos_mat.fill(inf);
        //生成范围随机数
        static std::mt19937 rng;
        static std::uniform_real_distribution<double> distribution1(-1, 2);
        static std::uniform_real_distribution<double> distribution2(-2, 4);
        for (int i = 0; i < p_num; ++i)
        {
            pos.push_back(distribution1(rng));
            spd.push_back(distribution2(rng));
        }
        vector<double> vec;
        for (int i = 0; i < p_num; ++i)
        {
            auto temp = fun_test(pos[i]);//计算函数值
            //初始化函数值矩阵和位置矩阵
            f_test(0, i) = temp;
            pos_mat(0, i) = pos[i];
            p_best.push_back(pos[i]);//初始化粒子的历史最优位置
        }
        std::ptrdiff_t minRow, minCol;
        f_test.row(0).minCoeff(&minRow, &minCol);//返回函数值矩阵第一行中极小值对应的位置
        g_best = pos_mat(minRow, minCol);//初始化全局最优位置
    }
    
    void PSO()
    {
        static std::mt19937 rng;
        static std::uniform_real_distribution<double> distribution(0, 1);
        for (int step = 1; step < iters; ++step)
        {
            for (int i = 0; i < p_num; ++i)
            {
                //更新速度向量和位置向量
                spd[i] = 0.5 * spd[i] + 2 * distribution(rng) * (p_best[i] - pos[i]) +
                    2 * distribution(rng) * (g_best - pos[i]);
                pos[i] = pos[i] + spd[i];
                //如果越界则取边界值
                if (spd[i] < -2 || spd[i] > 4)
                    spd[i] = 4;
                if (pos[i] < -1 || pos[i] > 2)
                    pos[i] = -1;
                //更新位置矩阵
                pos_mat(step, i) = pos[i];
            }
            //更新函数值矩阵
            for (int i = 0; i < p_num; ++i)
            {
                auto temp = fun_test(pos[i]);
                f_test(step, i) = temp;
            }
            for (int i = 0; i < p_num; ++i)
            {
                MatrixXd temp_test;
                temp_test = f_test.col(i);//取函数值矩阵的每一列
                std::ptrdiff_t minRow, minCol;
                temp_test.minCoeff(&minRow, &minCol);//获取每一列的极小值对应的位置
                p_best[i] = pos_mat(minRow, i);//获取每一列的极小值,即每个粒子的历史最优位置
            }
            g_best = *min_element(p_best.begin(), p_best.end());//获取全局最优位置
        }
        cout << fun_test(g_best);
    }
    
    int main()
    {
        init();
        PSO();
        system("pause");
        return 0;
    }
    

    总结

    在这里插入图片描述

    参考内容:

    版权声明:本文部分为CSDN博主「lx青萍之末」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
    原文链接:https://blog.csdn.net/daaikuaichuan/article/details/81382794

    ————————————————
    版权声明:本文部分为CSDN博主「nightmare_dimple」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
    原文链接:https://blog.csdn.net/nightmare_dimple/article/details/74331679

    展开全文
  • 1.1 粒子群算法建模粒子群算法首先在给定的解空间中随机初始化粒子群,待优化问题的变量数决定了解空间的维数。每个粒子有了初始位置与初始速度,然后迭代寻优。每一次迭代中,每个粒子通过跟踪两个...

    8f4917f943304f605dd2ed802c9649bf.png

    来喽!

    一、粒子群算法理论

    粒子群算法来源于鸟类集体活动的规律性,进而利用群体智能建立简化模型。它模拟的是鸟类的觅食行为,将求解问题的空间比作鸟类飞行的时间,每只鸟抽象成没有体积和质量的粒子,来表征一个问题的可行解。

    1.1 粒子群算法建模

    粒子群算法首先在给定的解空间中随机初始化粒子群,待优化问题的变量数决定了解空间的维数。每个粒子有了初始位置与初始速度,然后迭代寻优。每一次迭代中,每个粒子通过跟踪两个极值来更新自己的解空间中的位置和速度,一个是单个粒子本身在迭代中找到的最优粒子(个体极值),一个是所有粒子在迭代过程中的最优解粒子(全局极值)。

    1.2 粒子群算法特点

    (1)基于群体智能理论的优化算法,高效的并行算法。

    (2)粒子群算法随机初始化种群,使用适应值来评价个体的优劣程度,进行一系列的随机搜索。但粒子群算法根自己的速度来决定搜索,不需要进行交叉变异等操作,避免了复杂的操作。

    (3)每个粒子在算法结束时仍保持个体极值,所以除了得到问题的最优解以外,还能得到若干次优解,给出更多方案。

    (4)粒子群算法特有的记忆使其可以动态的跟踪当前搜索情况并调整搜索方向。

    二、粒子群算法种类

    2.1 基本粒子群算法

    假设在一个D维的目标搜索空间中,有N个粒子组成群落,其中第i个粒子表示为一个D维向量:

    第i个粒子飞行速度也是一个D维向量:

    第i个粒子迄今为止搜索到最优位置称为个体极值:

    整个粒子群迄今为止搜索到的最优位置为全局极值:

    在找到这两个最优值后,粒子群更新自己的速度和位置:

    其中,

    是学习因子(加速常数);
    为[0,1]范围内的均匀随机数;
    是粒子速度;
    增加了粒子飞行的随机性。更新速度的式子由三部分组成,第一部分为‘惯性’,反应了粒子群运动的‘习惯’;第二部分是‘认知’,反映了对自身历史经验的记忆;第三部分是‘社会’,反映了粒子间协同合作知识共享的群体历史经验。

    2.2 标准粒子算法

    引入了两个新的概念:

    探索:指粒子在一定程度上离开原先的轨迹,向新的方向搜索,体现了开拓未知区域的能力,探索能力是全局搜索能力。

    开发:指粒子在一定程度上继续在原先的轨迹上进行进一步搜索,是局部搜索能力。

    1998年,Shi Yuhui等人提出了带有惯性权重的改进粒子群算法:

    第一个式子有三个部分,第一部分保证算法的全局收敛性,第二、三部分保证局部收敛能力。w较大,全局收敛能力强局部收敛能力减弱;w较小,全局收敛能力减弱,局部收敛能力较强。当w=1时,与基本粒子算法相同,实验结果表明:w在0.8~1.2之间有更快的收敛速度,w>1.2时,容易陷入局部极值。

    在搜索时还可以对w进行动态调整,采用较多的是Shi提出的线性递减权值策略:

    表示最大进化代数;
    表示惯性最大权重;
    表示惯性最小权重。

    2.3 压缩粒子群算法

    Clerc等人利用约束因子来控制最终的收敛,可以有效搜索不同区域:

    为压缩因子:

    实验结果表明:与使用惯性权重的粒子群算法比较,使用具有压缩因子的粒子群算法有更快的收敛速度。

    2.4 离散粒子群算法

    Kennefy和Eberhart提出了离散二进制版的粒子群算法。将离散的问题映射到了连续的粒子运动空间,计算上扔保留速度-位置更新的运算规则。粒子在空间状态的取值只限于0,1两个值,而速度的每一维

    代表每一位
    取值为1的可能性,因此
    的公式保持不变,但是
    只能在[0,1]之间取值:

    686d14c77d38f36471bcf6a936e3f811.png

    r是从

    分布产生的随机数。

    三、粒子群算法流程

    4813d2d11f4e7b486c2bbdabd61d3e3a.png

    四、MATLAB实例及仿真

    参数说明:

    粒子种群规模N:一般为100~200。

    惯性权重w:控制算法的开发和探索能力,一般取值为[0.8,1.2]。

    加速常数

    :调节向
    方向飞行的最大步长,分别决定了粒子个体经验和群体经验对粒子运动轨迹的影响。
    粒子缺乏认知能力,
    个体之间没有信息交互,所以一般设置
    ,通常取
    。这样个体经验和群体经验有了同样的影响力。

    粒子最大速度

    :粒子速度在每一维都有一个最大速度限制值,用来对粒子速度进行钳制,该值太大,粒子也许会飞过优秀区域,太小粒子们可能会陷入局部最优,无法移动足够远跳出局部最优。研究发现 设定
    和调整惯性权重的作用是等效的,所以
    一般用于初始化进行设定,二不再对最大速度进行细致的选择和调节。

    邻域结构设定:全局版的粒子群算法将整个群体作为粒子的邻域,具有搜索快的优点但是容易陷入局部最优。局部版本粒子群算法将位置相近的个体作为粒子的邻域,收敛速度慢,不易陷入局部最优。实际中可先采用全局粒子群算法寻找最优解方向,再采用局部粒子群算法细致搜索。

    边界条件处理:边界吸收的方法。

    例1:计算函数

    的最小值,其中个体x的维数n=10,这是一个简单的平方和函数,只有一个极小点x=(0,0,...,0),理论上最小值分f(0,0,...,0)=0。

    解:

    clear all;
    close all;
    clc;
    N=100;                         %群体粒子个数
    D=10;                          %粒子维数
    T=200;                         %最大迭代次数
    c1=1.5;                        %学习因子1
    c2=1.5;                        %学习因子2
    w=0.8;                         %惯性权重
    Xmax=20;                       %位置最大值
    Xmin=-20;                      %位置最小值
    Vmax=10;                       %速度最大值
    Vmin=-10;                      %速度最小值
    %初始化个体
    x=rand(N,D)*(Xmax-Xmin)+Xmin;
    v=rand(N,D)*(Vmax-Vmin)+Vmin;
    %初始化个体最优位置和最优值
    p=x;
    pbest=ones(N,1);
    for i=1:N
        pbest(i)=func1(x(i,:));
    end
    %初始化全局最优位置和最优值
    g=ones(1,D);
    gbest=inf;
    for i=1:N
        if (pbest(i)<gbest)
            g=p(i,:);
            gbest=pbest(i);
        end
    end
    gb=ones(1,T);
    %按照公式依次迭代直到满足精度或者迭代次数
    for i=1:T
        for j=1:N
            if (func1(x(j,:))<pbest(j))
                p(j,:)=x(j,:);
                pbest(j)=func1(x(j,:));
            end
            if (pbest(j)<gbest)
                g=p(j,:);
                gbest=pbest(j);
            end
            v(j,:)=w*v(j,:)+c1*rand*(p(j,:)-x(j,:))+c2*rand*(g-x(j,:));
            x(j,:)=x(j,:)+v(j,:);
            %边界条件处理
            for ii=1:D
                if (v(j,ii)<Vmin)||(v(j,ii)>Vmax)
                    v(j,ii)=rand*(Vmax-Vmin)+Vmin;
                end
                if (x(j,ii)<Xmin)|(x(j,ii)>Xmax)
                    x(j,ii)=rand*(Xmax-Xmin)+Xmin;
                end
            end
        end
        %记录全局最优值
        gb(i)=gbest;
    end
    g;                         %最优个体
    gb(end);                   %最优值
    figure
    plot(gb)
    xlabel('迭代次数')
    ylabel('适应度值')
    title('适应度进化曲线')
    %适应度函数
    function result=func1(x)
    summ=sum(x.^2);
    result=summ;
    end

    eabed64558542dcba992c8a20e7bec2f.png

    例2:求函数

    的最小值,其中x的取值范围是[-4,4],y的取值范围是[-4,4]。这是一个有多个局部极值的函数。函数图形如图:
    clear all;
    close all;
    clc;
    x=-4:0.02:4;
    y=-4:0.02:4;
    N=size(x,2);
    for i=1:N
        for j=1:N
            z(i,j)=3*cos(x(i)*y(j))+x(i)+y(j)*y(j);
        end
    end
    mesh(x,y,z)
    xlabel('x')
    ylabel('y')

    f5a7e45103cc84328cfbdb2e5d224088.png

    解:

    clear all;
    close all;
    clc;
    N=100;                         %群体粒子个数
    D=2;                           %粒子维数
    T=200;                         %最大迭代次数
    c1=1.5;                        %学习因子1
    c2=1.5;                        %学习因子2
    Wmax=0.8;                      %惯性权重最大值
    Wmin=0.4;                      %惯性权重最小值
    Xmax=4;                        %位置最大值
    Xmin=-4;                       %位置最小值
    Vmax=1;                        %速度最大值
    Vmin=-1;                       %速度最小值
    %初始化个体
    x=rand(N,D)*(Xmax-Xmin)+Xmin;
    v=rand(N,D)*(Vmax-Vmin)+Vmin;
    %初始化个体最优位置和最优值
    p=x;
    pbest=ones(1,D);
    for i=1:N
        pbest(i)=func2(x(i,:));
    end
    %初始化全局最优位置和最优值
    g=ones(1,D);
    gbest=inf;
    for i=1:N
        if (pbest(i)<gbest)
            g=p(i,:);
            gbest=pbest(i);
        end
    end
    gb=ones(1,T);
    %按照公式依次迭代直到满足精度或者迭代次数
    for i=1:T
        for j=1:N
            %更新个体最优位置和最优值
            if (func2(x(j,:))<pbest(j))
                p(j,:)=x(j,:);
                pbest(j)=func2(x(j,:));
            end
            %更新全局最优位置和最优值
            if (pbest(j)<gbest)
                g=p(j,:);
                gbest=pbest(j);
            end
            %计算动态惯性权重值
            w=Wmax-(Wmax-Wmin)*i/T;
            %更新位置和速度
            v(j,:)=w*v(j,:)+c1*rand*(p(j,:)-x(j,:))+c2*rand*(g-x(j,:));
            x(j,:)=x(j,:)+v(j,:);
            %边界条件处理
            for ii=1:D
                if (v(j,ii)<Vmin)||(v(j,ii)>Vmax)
                   v(j,ii)=rand*(Vmax-Vmin)+Vmin;
                end
                if (x(j,ii)<Xmin)||(x(j,ii)>Xmax)
                    x(j,ii)=rand*(Xmax-Xmin)+Xmin;
                end
            end
        end
        gb(i)=gbest;
    end
    g;%最优个体
    gb(end);%最优值
    figure
    plot(gb)
    xlabel('迭代次数')
    ylabel('适应度值')
    title('适应度进化曲线')
    %适应度函数
    function value=func2(x)
    value=3*cos(x(1)*x(2))+x(1)+x(2)^2;
    end
            

    32150e88a0a8e765431b0aaa37c06bf6.png

    例3:用离散粒子群算法求函数

    的最小值,其中x的取值范围是[0,9],这是一个有局部多个极值的函数,图形如图所示:
    clear all;
    close all;
    clc;
    x=0:0.01:9;
    y=x+6*sin(4*x)+9*cos(5*x);
    plot(x,y)
    xlabel('x')
    ylabel('f(x)')

    e7b71574f3ece27f7be56e4a6eced536.png

    解:

    clear all;
    close all;
    clc;
    N=100;                         %群体粒子个数
    D=20;                           %粒子维数
    T=200;                         %最大迭代次数
    c1=1.5;                        %学习因子1
    c2=1.5;                        %学习因子2
    Wmax=0.8;                      %惯性权重最大值
    Wmin=0.4;                      %惯性权重最小值
    Xs=9;                        %位置最大值
    Xx=0;                       %位置最小值
    Vmax=10;                        %速度最大值
    Vmin=-10;                       %速度最小值
    %初始化个体
    x=randi([0,1],N,D);
    v=rand(N,D)*(Vmax-Vmin)+Vmin;
    %初始化个体最优位置和最优值
    p=x;
    pbest=ones(N,1);
    for i=1:N
        pbest(i)=func3(x(i,:),Xs,Xx);
    end
    %初始化全局最优位置和最优值
    g=ones(1,D);
    gbest=inf;
    for i=1:N
        if (pbest(i)<gbest)
            g=p(i,:);
            gbest=pbest(i);
        end
    end
    gb=ones(1,T);
    %按照公式依次迭代直到满足精度或者迭代次数
    for i=1:T
        for j=1:N
            %更新个体最优位置和最优值
            if (func3(x(j,:),Xs,Xx)<pbest(j))
                p(j,:)=x(j,:);
                pbest(j)=func3(x(j,:),Xs,Xx);
            end
               %更新全局最优位置和最优值
            if (pbest(j)<gbest)
                g=p(j,:);
                gbest=pbest(j);
            end
              %计算动态惯性权重值
              w=Wmax-(Wmax-Wmin)*i/T;
              %更新位置和速度
            v(j,:)=w*v(j,:)+c1*rand*(p(j,:)-x(j,:))+c2*rand*(g-x(j,:));
           %边界条件处理
            for ii=1:D
                if (v(j,ii)<Vmin)||(v(j,ii)>Vmax)
                   v(j,ii)=rand*(Vmax-Vmin)+Vmin;
                end
            end
            vx(j,:)=1./(1+exp(-v(j,:)));
            for jj=1:D
               if vx(j,jj)>rand
                   x(j,jj)=1;
               else
                   x(j,jj)=0;
               end
            end
        end
        gb(i)=gbest;
    end
    g;%最优个体
    m=0;
    for j=1:D
        m=g(j)*2^(j-1)+m;
    end
    f1=Xx+m*(Xs-Xx)/(2^D-1);%最优值
    figure
    plot(gb)
    xlabel('迭代次数')
    ylabel('适应度值')
    title('适应度进化曲线')
    %适应度函数
    function result=func3(x,Xs,Xx)
    m=0;
    D=length(x);
    for j=1:D
        m=x(j)*2^(j-1)+m;
    end
    f=Xx+m*(Xs-Xx)/(2^D-1);%译码成十进制数
    fit=f+6*sin(4*f)+9*cos(5*f);
    result=fit;
    end
      

    3e4179b64fd905aa4ce474a9f5ba9d37.png

    例4:0-1背包问题,有N件物品和容积为V的包,第i件物品的容积是c(i),价值是w(i),求将这些物品放入包中,使物体的总容积不超过背包容积,且总价值和最大。假设物品数量为10,背包容量为300.每件物品的体积为[95,75,23,73,50,22,6,57,89,98];价值为[89,59,19,43,100,72,44,16,7,64];

    解:

    clear all;
    close all;
    clc;
    N=100;                          %群体粒子个数
    D=10;                           %粒子维数
    T=200;                          %最大迭代次数
    c1=1.5;                         %学习因子1
    c2=1.5;                         %学习因子2
    Wmax=0.8;                       %惯性权重最大值
    Wmin=0.4;                       %惯性权重最小值
    Vmax=10;                        %速度最大值
    Vmin=-10;                       %速度最小值
    V=300;                          %背包容量
    C=[95,75,23,73,50,22,6,57,89,98];            %物品体积
    W=[89,59,19,43,100,72,44,16,7,64];           %物品价值
    afa=2;                                       %惩罚函数系数
    %初始化个体
    x=randi([0,1],N,D);
    v=rand(N,D)*(Vmax-Vmin)+Vmin;
    %初始化个体最优位置和最优值
    p=x;
    pbest=ones(N,1);
    for i=1:N
        pbest(i)=func4(x(i,:),C,W,V,afa);
    end
    %初始化全局最优位置和最优值
    g=ones(1,D);
    gbest=eps;
    for i=1:N
        if (pbest(i)>gbest)
            g=p(i,:);
            gbest=pbest(i);
        end
    end
    gb=ones(1,T);
    %按照公式依次迭代直到满足精度或者迭代次数
    for i=1:T
        for j=1:N
             %更新个体最优位置和最优值
            if (func4(x(j,:),C,W,V,afa)>pbest(j))
                p(j,:)=x(j,:);
                pbest(j)=func4(x(j,:),C,W,V,afa);
            end
            %更新全局最优位置和最优值
            if (pbest(j)>gbest)
                g=p(j,:);
                gbest=pbest(j);
            end
             %计算动态惯性权重值
              w=Wmax-(Wmax-Wmin)*i/T;
              %更新位置和速度
              v(j,:)=w*v(j,:)+c1*rand*(p(j,:)-x(j,:))+c2*rand*(g-x(j,:));
         %边界条件处理
            for ii=1:D
                 if (v(j,ii)<Vmin)||(v(j,ii)>Vmax)
                 v(j,ii)=rand*(Vmax-Vmin)+Vmin;
                 end
            end
            vx(j,:)=1./(1+exp(-v(j,:)));
            for jj=1:D
                 if vx(j,jj)>rand
                 x(j,jj)=1;
                 else
                 x(j,jj)=0;
                 end
            end
        end
        gb(i)=gbest;
    end
    g;%最优个体
    figure
    plot(gb)
    xlabel('迭代次数')
    ylabel('适应度值')
    title('适应度变化曲线')
    %适应度函数
    function result=func4(f,C,W,V,afa)
    fit=sum(f.*W);
    TotalSize=sum(f.*C);
    if TotalSize<=V
        fit=fit;
    else
        fit=fit-afa*(TotalSize-V);
    end
    result=fit;
    end
            

    f81bbfbb038657ea56255f7c16d77998.png

    希望我明天继续爱学习!!!

    展开全文
  • 粒子群算法-讲解+实例

    千次阅读 2019-07-11 21:44:56
    今天给大家讲解的时粒子群算法,首先先牢记以下的基本公式: 1.简单的来讲,粒子群算法是这个样子,当你在一个全解的范围内,想要去找最优解,可以先派出像四面八方而去的小兵去搜索,他们向四面八方去探索时在坚守...

    今天给大家讲解的时粒子群算法,首先先牢记以下的基本公式:在这里插入图片描述
    在这里插入图片描述
    1.简单的来讲,粒子群算法是这个样子,当你在一个全解的范围内,想要去找最优解,可以先派出像四面八方而去的小兵去搜索,他们向四面八方去探索时在坚守自己本身的航线的同时,会通过自身之前所找到的一些局部最优解轨迹和全局最优解即其他所有小兵找到的最优解的最最优解来修正自身轨迹,最终所有小兵都会聚集在全局最优解周围。
    –因此我们就可以看到在这个式子里面分为这么几个部分
    Ⅰ:他本身的航线
    Ⅱ:他历史所走过的找到的局部最优解
    叁:其它包括它本身找到的最优的解
    其中Ⅰ主要代表他对于整个解的空间的探索能力,防止他不会因为一些局部最优的情况下就停止探寻,因此也常引入状态变量W来,一般W取0.8~1.2之间,当W>1.2时容易陷入局部最优值。因此可以对其进行动态调整。让w随着时间的变化而产生变化,在不同的时候给他不同的收敛值。
    在这里插入图片描述
    Tmax最大迭代代数,wmax=0.9 wmin=0.4 t当前迭代次数

    展开全文
  • 配电网粒子群算法实例

    千次阅读 2015-12-18 16:30:23
    粒子群算法求解配电网故障定位评价函数的极小值即可。 Matlab源代码:  clc;  clear all;  global Ij;  global Vmax; %允许误差  MaxNum=120; %粒子最大迭代次数  D=12; %...

    算例:

            



    一、故障信息的数学表示

    在上图中K表示断路器,每一个断路器上均有一个FTU装置,可以反馈断路器开关是否过流,用表示上传的故障信息,反映的是各分段开关处是否流过故障电流有故障电流为1,否则为0)。即:


    因为FTU上传的信息可分为有故障信息及无故障信息两类,对于分段区间来讲也只能是有故障及无故障两种情况,所以我们可以用二进制编码规则对配电网故障定位问题进行数学建模。以上图所示辐射状配电网为例,系统拥有12个分段开关,我们可以用一串12位的二进制代码表示FTU的上传信息,作为程序的输入,1代表对应的开关有过流信息,0代表对应的开关无过流信息。同时用另一串12位的二进制代码作为程序的输出,代表对应馈线区间发生故障,代表无故障。

    二、二进制PSO算法的配电网故障定位原理

    应用二进制PSO算法求解配电网故障定位问题,粒子的位置代表配电网中馈线区段的状态,粒子的维数代表配电网的馈线区段总数。每一馈线区段存在0和1两种状态,0表示正常状态,1表示故障状态,馈线区段的状态为待求量。因此,N段馈线区段的状态求解就转化成N维粒子群优化求解,每个粒子的N维位置都表示为配电网N段馈线区段的潜在状态。每次迭代过程中,通过评价函数评价各粒子位置优劣,更新粒子的当前最优位置和全体粒子的最优位置,进而更新粒子的速度和位置,直到满足程序终止条件为止。最终得出的粒子群的全局最优位置就是所求的各馈线区段的实际状态。

    三、构造配电网故障定位评价函数

        基于待求各馈线区段实际状态下所对应的信息应与实际上传故障信息偏差最小的原则,构造如下评价函数:


    表达式的值为每个潜在解对应的适应度值,值越小表示解越优良,因此评价函数应取极小值。

    式中:为第j个开关FTU上传的故障信息,取值为1认为该开关流过了故障电流,为0则未流过;为各开关节点的期望状态,若该开关流过了故障电流,则期望状态为1,相反,期望状态为0,表达式为各区段状态的函数。N为配电网中馈线区段的总数,为配电网中各设备状态,为1表明设备故障,取0则设备正常。表示一个权系数乘以故障设备数,。是根据故障诊断理论中“最小集”概念设置的权重系数,取值范围介于与1之间,表明故障区间数越少解越优,避免出现误诊断,本文中权系数取0.5。

    四、求期望函数

    针对不同的网络,期望函数的表达式不同,无法用统一的格式表示出来。仅以上述算例为例来说明如何求取期望函数。

    在辐射型配电网中,任何区间发生故障,必将导致其上游开关流过故障过电流。这是求取期望函数的原则。假设图一的馈线区段6处发生故障,那么分段开关K1,K2,K3 ,K6,都将有过电流若馈线区段11处发生故障,则K1,K2,K3 ,K6,K10都将有过电流。由此类推我们可得到各个分段开关的期望函数:




    用粒子群算法求解配电网故障定位评价函数的极小值即可。


    Matlab源代码:

        clc;

        clear all;

        global Ij;

        global Vmax;                       %允许误差

        MaxNum=120;                    %粒子最大迭代次数

        D=12; %粒子的维数                    

        Particlesize=30;                    %粒子群规模

        c1=2.1;                            %每个粒子的个体学习因子,也称为加速常数

        c2=2.1;                            %每个粒子的社会学习因子,也称为加速常数

        w=1;                           %惯性因子

        Vmax=4;                        %粒子的最大飞翔速度

        S_b=randint(Particlesize,D);     %粒子所在的位置

        V=Vmax*rand(Particlesize,D);         %粒子的飞翔速度

         Ij=[1 1 1 1 0 1 1 1 0 0 1 1];%FTU 反馈的值这里就是输入,请修改这里,看结果

        fitness=F(S_b');%这里注意下,要转置

        %inline定义的适应度函数会使程序运行速度大大降低

        for i=1:Particlesize

                f(i)=F(S_b(i,:));

        end

        personalbest_x=S_b;%这里是种群里的每个粒子都认为自己是最好的

        personalbest_faval=f;

        [globalbest_favali]=min(personalbest_faval);%全局最优

        globalbest_x=personalbest_x(i,:);

        k=1;

        Start=0;

        while k<=MaxNum

            Remember=globalbest_faval;

            for i=1:Particlesize

                    f(i)=F(S_b(i,:));%计算适应度函数值

                if f(i)<personalbest_faval(i) %判断当前位置是否是历史上最佳位置

                    personalbest_faval(i)=f(i);

                    personalbest_x(i,:)=S_b(i,:);

                end

            end%这里就是把整个种群中的所有粒子都检查了一遍,保证当前最优

            [globalbest_favali]=min(personalbest_faval);%更新种群最优

            globalbest_x=personalbest_x(i,:);

            for i=1:Particlesize %更新粒子群里每个个体的最新位置

               V(i,:)=w*V(i,:)+c1*rand*(personalbest_x(i,:)-S_b(i,:))...

                   +c2*rand*(globalbest_x-S_b(i,:));%更新了速度

               for j=1:length(Ij)

                   if rand(1)<Sigmoid(V(i,j))%这里体现了速度决定位置

                       S_b(i,j)=1;

                   else

                        S_b(i,j)=0;

                   end

               end

            end

        if globalbest_faval-Remember==0

            Start=Start+1;

            if Start==20

                disp(globalbest_x);

                break;

            end

        end

            k=k+1;

    end

    F.m

    function fitness=F(SB)

    global Ij;

    IE=zeros(1,length(Ij));

    IE(1)=SB(1)|SB(2)|SB(3)|SB(4)|SB(5)|SB(6)|SB(7)|SB(8)|SB(9)|SB(10)|SB(11)|SB(12);

    IE(2)=SB(2)|SB(3)|SB(4)|SB(5)|SB(6)|SB(7)|SB(8)|SB(9)|SB(10)|SB(11)|SB(12);

    IE(3)=SB(3)|SB(4)|SB(5)|SB(6)|SB(7)|SB(8)|SB(9)|SB(10)|SB(11)|SB(12);

    IE(4)=SB(4)|SB(5);

    IE(5)=SB(5);

    IE(6)=SB(6)|SB(7)|SB(8)|SB(9)|SB(10)|SB(11)|SB(12);

    IE(7)=SB(7)|SB(8)|SB(9)|SB(10);

    IE(8)=SB(8)|SB(9)|SB(10);

    IE(9)=SB(9);

    IE(10)=SB(10);

    IE(11)=SB(11)|SB(12);

    IE(12)=SB(12);

    fitness=sum(abs(Ij-IE))+0.5*sum(abs(SB));

    end

    Sigmoid.m

    function AN=Sigmoid(V)

    global Vmax;

        if V>Vmax

            AN=0.98;

        elseif V>-Vmax

           AN=1/(1+exp(-V));

        else

            AN=-0.98;

        end

    end


    展开全文
  • 粒子群算法的matlab实现(一)

    万次阅读 多人点赞 2017-07-04 15:45:39
    粒子群算法(Particle Swarm Optimization,PSO)是20世纪90年代兴起的一门学科,因其概念简明、实现方便、收敛速度快而为人所知。粒子群算法的基本思想是模拟鸟群随机搜寻食物的捕食行为,鸟群通过自身经验和种群...
  • 数学建模粒子群算法详解

    千次阅读 2019-08-31 19:24:03
    粒子群(鸟群算法算法是典型的寻优算法。分为全局最优和局部最优。 基本思想:主要模拟自然界生物捕食的策略,群体迭代,粒子在解空间追随最优的例子进行搜索。是智能算法的一种。 算法特点: 简单易行。 收敛...
  • 文章目录粒子群算法1 粒子群算法概述2 相关变量3 固定的参数4 粒子群算法求解优化问题5 实例6 python实现7 特点 粒子群算法 1 粒子群算法概述 粒子群算法,也称粒子群优化(Particle SwarmOptimization,PSO)算法,...
  • 粒子群算法 车间调度 甘特图 mt06!!!!!!!!!!
  • 1995年,美国社会心理学家James Kennedy和电气工程师Russell Eberhart共同提出了粒子群算法(Particle Swarm Optimization,PSO),该算法的提出是受对鸟类群体行为进行建模与仿真的研究结果的启发。 粒子群优化...
  • C语言实现粒子群算法(PSO)一

    千次阅读 2016-12-03 16:38:31
    粒子群算法的C语言代码实现。
  • 粒子群算法 自然界中的鸟群和鱼群的群体行为一直是科学家的研究兴趣所在。生物学家Craig Reynolds在1987年提出了一个非常有影响的鸟群聚集模型,在他的仿真中,每一个个体都遵循:避免与邻域个体相冲撞;匹配邻域...
  • 本文记录了有约束情况下的粒子群算法函数寻优的解决方法
  • 基于粒子群算法优化的ELMAN动态递归神经网络预测及其MATAB实现 文章目录基于粒子群算法优化的ELMAN动态递归神经网络预测及其MATAB实现1. 模型与算法描述1.1 ELMAN神经网络预测模型介绍1.2 粒子群算法PSO介绍2. 粒子...
  • http://www.cnblogs.com/lyrichu/p/6151272.htmlC语言实现粒子群算法(PSO)一最近在温习C语言,看的书是《C primer Plus》,忽然想起来以前在参加数学建模的时候,用过的一些智能算法,比如遗传算法、粒子群算法、蚁...
  • 粒子群算法求解物流配送路线问题(python)

    千次阅读 多人点赞 2020-02-15 17:05:26
    1.Matlab实现粒子群算法的程序代码:https://www.cnblogs.com/kexinxin/p/9858664.html matlab代码求解函数最优值:https://blog.csdn.net/zyqblog/article/details/80829043 讲解通俗易懂,有数学实例的博文:...
  • 粒子群算法matlab以求解函数最优解为例 clear; clc; close all; N=100; %粒子个数 D=2; %粒子维数 MaxIter=500; %最大迭代次数 C1max=1.8; %权重参数,自适应 C2max=1.8; C1=1.2; C2=1.2; w=0.79; Wmax=0.8; %对...
  • 在灰色Verhulst模型的基础上对等间隔和非等间隔GM(1,1)幂模型...最后实例表明,基于粒子群算法参数辨识的GM(1,1)幂模型建模精度高于灰色Verhulst模型,同时也表明了该方法的有效性和可行性,具有重要的理论意义。
  • 数学建模方法-蚁群算法

    千次阅读 2019-09-26 03:51:40
    它的诞生比粒子群算法还要早3年,在1992年的某一天,一位叫Marco Dorigo的在他的博士论文中提出了蚁群算法,并称其灵感来源于蚂蚁在寻找食物过程中发现路径的行为。。。好了历史背景介绍到止,接下来就认真讲一下...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 465
精华内容 186
关键字:

粒子群算法建模实例