简单二维平面的微粒群算法实现效果如下图所示,上升到多维空间的核心思想类似。
准备工作
新建一个MFC单文档应用程序,具体步骤如下图所示:
新建一个类,如下图所示:
我们用微粒群来求解最优值问题,求解的是一个函数,比如y=x^2。
x给几个初始值,代表初始的鸟群(比如10只鸟),
每只鸟位置是随机的,x=0的位置是最优位置。
通过一定的迭代,把最优位置找到。
首先定义一些需要的变量,以及后面需要用到的成员函数,如下图所示:
1、添加所要求解函数的成员函数,比如这里以f(x) = x0² + x1²为例,如下图所示。
代码如下:
float CWeiLiQun::f(float x[])
{
return (x[0])*(x[0]) + (x[1])*(x[1]);//这里以函数f(x) = x0² + x1²为例
}
2、添加运动的成员函数,如下图所示。
3、添加画鸟群的成员函数,如下图所示。
最终变量以及成员函数如下:
// WeiLiQun.h: interface for the CWeiLiQun class.
//
//
#if !defined(AFX_WEILIQUN_H__C52936AD_87DC_4442_A966_93B5BB20BCC8__INCLUDED_)
#define AFX_WEILIQUN_H__C52936AD_87DC_4442_A966_93B5BB20BCC8__INCLUDED_
#if _MSC_VER > 1000
#pragma once
#endif // _MSC_VER > 1000
#define N 10 //定义N只鸟(这里取10)
typedef struct
{
float x[2];//微粒位置
float v[2];//微粒速度
float p[2];//曾经经历过的最好位置
}SWeiLi;//用结构体表示微粒
class CWeiLiQun
{
public:
void Draw(CDC *pDC);
void Move();
float f(float x[]);
CWeiLiQun();
virtual ~CWeiLiQun();
SWeiLi m_wl[N];//定义微粒数组
float m_pg[2];//全局的最优(不是一个微粒,而是一个位置)
float c1,c2;//c1是调节微粒飞向自身最好位置的步长,c2是调节微粒飞向全局最好位置的步长。
float r1,r2;//r1,r2~u(0,1)
CPoint m_YD;//原点
float m_kx,m_ky;//比例尺
int m_times;//迭代次数
};
#endif // !defined(AFX_WEILIQUN_H__C52936AD_87DC_4442_A966_93B5BB20BCC8__INCLUDED_)
紧接着初始化,代码如下:
CWeiLiQun::CWeiLiQun()
{
int i,j;
for( i = 0; i < N; i++)
{
for (j = 0; j < 2; j++)
{
m_wl[i].x[j] = rand()%20 -10;//N只鸟的初始位置,-10到10之间的范围
m_wl[i].v[j] = rand()%6 -3; //N只鸟的初始速度,-3到3之间的范围
m_wl[i].p[j] = m_wl[i].x[j]; //N只鸟当前的最好位置
}
}
m_pg[0] = m_wl[0].p[0];
m_pg[1] = m_wl[0].p[1];//全局的最好位置
for( i = 1; i < N; i++)
{
if(f(m_wl[i].p) < f(m_pg))
{
m_pg[0] = m_wl[i].p[0];
m_pg[1] = m_wl[i].p[1];
}
}
c1 = 1.2;
c2 = 1.2;
r1 = 0.6;
r2 = 0.6;//按照经验给定的初始值
m_YD.x = 400;
m_YD.y = 250;//原点初始值位置
m_kx = 10;
m_ky = -10;//比例尺放大10倍
m_times = 0;
}
最关键的运动代码如下:
void CWeiLiQun::Move()
{
int i,j;
int vMax = 10;
for(i = 0; i < N; i++)
{
for (j = 0; j < 2; j++)
{
m_wl[i].v[j] += c1*r1*(m_wl[i].p[j] - m_wl[i].x[j]) + c2*r2*(m_pg[j] - m_wl[i].x[j]);//给定的公式
if( m_wl[i].v[j] > vMax )
m_wl[i].v[j] = vMax;
if( m_wl[i].v[j] < -vMax)
m_wl[i].v[j] = -vMax;//如果不写这个限制条件,若X范围很大,则粒子群震荡得非常厉害
m_wl[i].x[j] += m_wl[i].v[j];//给定的公式
if(f(m_wl[i].x) < f(m_wl[i].p))//当前位置是否比曾经经历过的最优位置要好
{
m_wl[i].p[0] = m_wl[i].x[0];
m_wl[i].p[1] = m_wl[i].x[1];
}
}
}
for( i = 0; i < N; i++)
{
if(f(m_wl[i].p) < f(m_pg))//和全局最优位置再做比较
{
m_pg[0] = m_wl[i].p[0];
m_pg[1] = m_wl[i].p[1];
}
}
m_times++;//迭代次数
}
再将鸟群画出来,代码如下:
void CWeiLiQun::Draw(CDC *pDC)
{
int i,j;
int x ,y,r;
CString str1,str2;
r = 20;//坐标原点半径
x = m_YD.x;
y = m_YD.y;
pDC->Ellipse(x-r,y-r,x+r,y+r);
for(i = 0; i < N; i++)
{
x = m_YD.x + m_wl[i].x[0] * m_kx;
y = m_YD.y + m_wl[i].x[1] * m_ky;
r = 10;
pDC->Ellipse(x-r,y-r,x+r,y+r);
}
str1.Format("迭代次数:%d",m_times);
str2.Format("全局最优值位置:(%.2f,%.2f)",m_pg[0],m_pg[1]);
pDC->TextOut(600,50,str1);
pDC->TextOut(600,100,str2);
for(i = 0; i < N; i++)
{
str1.Format("第%d只鸟的当前位置:(%.2f,%.2f);曾经最优值位置:(%.2f,%.2f)",i+1,m_wl[i].x[0],m_wl[i].x[1],m_wl[i].p[0],m_wl[i].p[1]);
pDC->TextOut(600,150+50*i,str1);
}
}
新建一个菜单,如下图所示:
再添加一个“连续动画”,如下图所示:
添加窗口消息响应句柄,如下图所示:
在CWeiLiQunSuanFaView里引入
代码如下:
void CWeiLiQunSuanFaView::OnMNext()
{
m_WeiLiQun.Move();
Invalidate(TRUE);
}
void CWeiLiQunSuanFaView::OnMDongHua()
{
SetTimer(1,100,NULL);
}
void CWeiLiQunSuanFaView::OnTimer(UINT nIDEvent)
{
m_WeiLiQun.Move();
Invalidate(TRUE);
CView::OnTimer(nIDEvent);
}
代码如下:
void CWeiLiQunSuanFaView::OnDraw(CDC* pDC)
{
CWeiLiQunSuanFaDoc* pDoc = GetDocument();
ASSERT_VALID(pDoc);
// TODO: add draw code for native data here
m_WeiLiQun.Draw(pDC);
}
实现效果图:
上述案例只是实现了简单的二维平面的微粒群算法,三维等多维空间的核心思想类似,可自行尝试。