PSO算法求解（上）-基本PSO

粒子群算法是智能算法之一，应用还是比较广泛的，尤其是在毕设当中。

前一次关于“PSO求解整数规划问题”的举例说明比较突兀，并没有完整的依据，且求解方法是笨拙的，根本不适应一般求解。

这里，将从基本粒子群算法出发，然后过渡到一般问题的求解。首先，使用基本粒子群求解无约束问题。
图片中文字来源：粒子群算法及应用_纪震著_北京：科学出版社_2009.01

基本粒子群算法求解无约束问题：

#include<iostream>
#include<stdlib.h>   
#include<stdio.h>
#include<fstream>
#include <string>
#include<ctime>
#include <algorithm>
#include <vector>
#include<cmath>

using namespace std;

/*测试用例：
min Z
	Z=(1.5-x_1+x_1*x_2)^2+(2.25-x_1+x_1*(x_2)^2)^2+(2.625-x_1+x_1*(x_2)^3)^2
	x->[-4.5,4.5]
最优解：Z(3,0.5)=0
*/

//PSO参数
const int person_number = 10;//种群个体数目
const int dimension = 2;//搜索空间维度->目标变量数
double weight;//惯性权重
double LDW_ini = 0.9;//线性递减权值策略（LDW）
double LDW_end = 0.4;
double c1 = 2, c2 = 2;//学习因子
const int maxgen = 250;//迭代次数 

//PSO变量
double X_pos[person_number][dimension];//粒子位置变量
double V_speed[person_number][dimension];//粒子速度变量
double fitness[person_number];//粒子适应值
double fitness_tem[person_number];
double pbest[person_number][dimension];//个体极值->位置
double pbest_tem;//临时个体极值
int pbest_index;//粒子编号
double gbest[dimension];//粒子最优值
double gbest_tem;//临时最优
double gbest_fit;

//数据记录
double iter_arr[maxgen];//记录迭代次数
double gbest_arr[maxgen];//记录迭代过程中的最优值
double gbest_X[dimension];//记录变量取值

//适应度函数
double Fitness(double x_pos[person_number][dimension],int pers) {//输入粒子位置矩阵、粒子编号
	double fit = 0;
	fit = pow(1.5 - x_pos[pers][0] + x_pos[pers][0] * x_pos[pers][1],2) + pow(2.25 - x_pos[pers][0] + x_pos[pers][0] * pow(x_pos[pers][1],2),2) + pow(2.625 - x_pos[pers][0] + x_pos[pers][0] * pow(x_pos[pers][1],3),2);
	return fit;
}

//权重设置
double Weight(int ITER) {//输入迭代次数
	int pso_weight = 0;
	pso_weight = ((LDW_ini - LDW_end)*(maxgen - ITER) / maxgen) + LDW_end;
	return pso_weight;
}

//初始化
void initial() {
	//初始化位置和速度
	srand((unsigned int)time(NULL));
	for (int i = 0; i < person_number;i++) {
		for (int j = 0; j < dimension;j++) {
			X_pos[i][j] = (rand() % (90 + 1)-45) / 10.0;//-4.5-4.5
			V_speed[i][j] = (rand() % (90 + 1)-45) / 10.0;
		}
	}

	//初始化个体极值
	for (int i = 0; i < person_number; i++) {
		for (int j = 0; j < dimension; j++) {
			pbest[i][j] = X_pos[i][j];
		}
	}

	//初始化全局极值
	for (int i = 0; i < person_number; i++) {//计算适应值，定位最小适应值的粒子位置和编号
		fitness[i] = Fitness(X_pos,i);
	}
	pbest_tem = 10000;
	pbest_index = -1;
	for (int i = 0; i < person_number; i++) {//定位个体极值和粒子编号
		if (fitness[i] < pbest_tem) {
			pbest_tem = fitness[i];
			pbest_index = i;
		}
		else {
			pbest_tem = pbest_tem;
			pbest_index = pbest_index;
		}
	}

	for (int j = 0; j < dimension;j++) {//初始化全局极值
		gbest[j] = pbest[pbest_index][j];
	}
	gbest_fit = Fitness(X_pos,pbest_index);
}

//PSO
void pso() {
	int iter = 0;//迭代计数
	srand((unsigned int)time(NULL));
	while (iter < maxgen) {
		//更新权重
		weight = Weight(iter);

		//更新粒子位置和速度
		for (int i = 0; i < person_number;i++) {
			for (int j = 0; j < dimension;j++) {
				//速度更新
				V_speed[i][j] = weight * V_speed[i][j] + c1 * rand() * (pbest[i][j] - X_pos[i][j]) + c2 * rand() * (gbest[j] - X_pos[i][j]);
				if(V_speed[i][j]<-4.5|| V_speed[i][j]>4.5)
					V_speed[i][j] = (rand() % (90 + 1) - 45) / 10.0;
				//位置更新
				X_pos[i][j] = X_pos[i][j] + V_speed[i][j];
				if (X_pos[i][j]<-4.5 || X_pos[i][j]>4.5)
					X_pos[i][j] = (rand() % (90 + 1) - 45) / 10.0;
			}
		}

		//计算适应值
		for (int i = 0; i < person_number; i++) {
			fitness_tem[i] = Fitness(X_pos, i);
		}
		//更新历史个体极值
		for (int i = 0; i < person_number; i++) {
			if (fitness_tem[i] < fitness[i]) {
				for (int j = 0; j < dimension;j++) {
					pbest[i][j] = X_pos[i][j];
				}
				fitness[i] = fitness_tem[i];
			}
		}

		//定位个体极值和粒子编号
		pbest_tem = 10000;
		pbest_index = -1;
		for (int i = 0; i < person_number; i++) {
			if (fitness[i] < pbest_tem) {
				pbest_tem = fitness[i];
				pbest_index = i;
			}
			else {
				pbest_tem = pbest_tem;
				pbest_index = pbest_index;
			}
		}
		gbest_tem = pbest_tem;//全局极值

		//更新全局极值
		if (gbest_tem < gbest_fit) {
			gbest_fit = gbest_tem;
			for (int j = 0; j < dimension;j++) {
				gbest[j] = pbest[pbest_index][j];
			}
		}

		//记录
		iter_arr[iter]=iter+1;
		gbest_arr[iter]=gbest_fit;
		iter++;
	}
	


}

int main() {
	initial();
	pso();
	for (int i = 0; i < maxgen; i++) {
		cout << iter_arr[i] << " ";
	}
	cout << endl;
	for (int i = 0; i < maxgen; i++) {
		cout << gbest_arr[i] << " ";
	}
	system("pause");
	return 0;
}

迭代过程如图：

接下来，引入罚函数概念，将约束问题转化为无约束问题，运用PSO解决一般的数学规划问题。