TA的每日心情![](source/plugin/dsu_paulsign/img/emot/kx.gif) | 开心 昨天 15:13 |
---|
签到天数: 153 天 [LV.7]常住居民III
管理员
中国膜结构网www.mjgou.com
- 积分
- 2147483647
![QQ](static/image/common/qq_big.gif)
|
- //main.cpp
- #include"Head.h"
- int main() {
- //更新随机种子
- srand((unsigned)time(NULL));
- Graph* graph = new Graph;
- graph->readCity();//启用文件读取,注销掉则为随机构造城市
- graph->innitial_Graph();
- AntColony* ant_colony = new AntColony;
- ant_colony->innitial_AntColony();
- for (int NC_now = 0; NC_now < NC_MAX; ++NC_now) {
- for (int i = 1; i < N; ++i) {
- ant_colony->ChooseNextCity(graph, i);
- }
- ant_colony->Update(graph,NC_now);
- ant_colony->innitial_AntColony();
- }
-
- cout << "The BestLength=" << BestLength << endl;
- ShowByPython(graph);
- int last, next;
- double sum = 0;
- for (int i = 1; i < N; ++i) {
- last = BestPath[i-1];
- next = BestPath[i];
- cout << "Travel " << i << " : " << last << " \t-> " << next << " \tdistante: " << graph->Length[last][next] << endl;
- sum += graph->Length[last][next];
- }
- last = BestPath[N - 1];
- next = BestPath[0];
- cout << "Travel " << N << " : " << last << " \t-> " << next << " \tdistante: " << graph->Length[last][next] << endl;
- sum += graph->Length[last][next];
- cout << "最优解路程和 = " << sum << endl;
- cout << "TSPlib上的已知最优解 :" << graph->KnowBest << endl;
- cout << "与已知最优解的偏差为 :" << ((sum - graph->KnowBest) / graph->KnowBest) * 100 << "%" << endl;//最后两句必须启用文件读取才有效
- delete graph;
- delete ant_colony;
- return 0;
- }
复制代码
- //Head.h
- #include<iostream>
- #include<cstdlib>
- #include<string>
- #include<cmath>
- #include<fstream>
- #include<Python.h>
- using namespace std;
- #define PATH "E:\\TSPlib\\ch150.txt"//TSPlib数据文件的路径
- #define N 150 //城市数量
- #define M 100 //蚂蚁数量
- #define NC_MAX 1200 //迭代次数
- #define α 1 //信息素浓度重要性
- #define β 1 //启发式重要性
- #define Q 100 //蚂蚁所携带信息素
- #define ρ 0.2 //蒸发系数
- //注意:如果从文件中读取City,目前需要手动调整城市数量N,要与文件中保持一致
- //并且修改 PATH 路径与启用main函数中的readCity函数
- //以下数据结构用于最后分析
- int BestPath_PerRound[NC_MAX][N] = { 0 };//记录每轮的最优路径
- double BestLength_PerRound[NC_MAX] = { 0 };//记录每轮最优路径的长度
- int BestPath[N] = { 0 };//全局最短路径(后一轮迭代未必一定比前一轮更优,只是总体是这样,不绝对),不保存每轮的信息
- double BestLength = 1000000;//全局最短路径长度
- double AverageLength_PerRound[NC_MAX] = { 0 };//记录每轮所有蚂蚁所走路径的平均长度
- class City {
- public:
- //无参构造函数,创建城市时默认在[0,999]内随机获得x,y坐标
- City() {
- x = rand() % 1000;
- y = rand() % 1000;
- }
- double x;
- double y;
- };
- class Graph {
- public:
- City citys[N];//无参构造函数创建N个城市
- double Length[N][N] = { 0 };//路径长度表
- double tao[N][N] = { 0 };//信息素浓度表
- double Heuristic[N][N] = { 0 };//启发式表,将城市之间距离的倒数作为启发函数
- public:
- //使用fstream读取文件中的城市坐标,并创建城市
- int KnowBest;//读取文件时存放已知最优解
- void readCity() {
- ifstream CityFile;
- CityFile.open(PATH, ios::in);
- if (!CityFile.is_open()) {
- cout << "Open file failure" << endl;
- exit(0);
- }
- CityFile >> KnowBest;
- cout << "The known best result is :" << KnowBest << endl;
- int a;double b[2];
- while (!CityFile.eof()) {
- CityFile >> a >> b[0] >> b[1];
- citys[a - 1].x = b[0];
- citys[a - 1].y = b[1];
- }
- CityFile.close();
- }
- //初始化城市图,并计算各城市之间的距离
- void innitial_Graph() {
- double x = 0, y = 0;
- for (int i = 0; i < N; ++i) {
- for (int j = 0; j < N; ++j) {
- x = citys[i].x - citys[j].x;
- y = citys[i].y - citys[j].y;
- Length[i][j] = sqrt(pow(x, 2) + pow(y, 2));//就是两点间距离公式
- if (i == j) {
- Length[i][j] = 10000;//默认城市到自己的距离为10000,防止启发式表出现无穷大,并不会计算城市自己到自己的概率
- }
- else {
- if (x == 0 && y == 0) {
- //这种情况是城市坐标完全重合,如果遇到这种情况则不能继续计算了,否则i与j之间的距离为0,启发式就会为无穷大
- //随之p的分子和分母都会变成无穷大,无穷大除无穷大...
- cout << "innitial_Graph ERROR!!! City overlapping, please retry." << endl;
- cout << "City :" << i << " and City :" << j << endl;
- exit(0);
- }
- }
- Heuristic[i][j] = 1 / Length[i][j];//将距离的倒数作为启发式
- tao[i][j] = 1;//信息素tao不能初始化为0
- }
- }
- }
- };
- class AntColony {
- public:
- int tabu[M][N] = { 0 };//禁忌列表,同时记录了蚂蚁的路径,每个元素取值为[0,N-1]
- int allow[M][N] = { 0 };//允许列表,完全依赖于禁忌列表,只是为了算法方便一些而设置,无法记录路径,每个元素取值为[0,1]
- double probability[M][N] = { 0 };//概率表,存放各蚂蚁选择下一作城市的概率
- double cumulative_probability[M][N] = { 0 };//累积概率表,用于轮盘赌算法随机选择城市
- public:
- //每轮迭代前都应该被调用一次
- void innitial_AntColony() {
- //初始化各列表
- for (int i = 0; i < M; ++i) {
- for (int j = 0; j < N; ++j) {
- tabu[i][j] = -1;
- allow[i][j] = 1;
- probability[i][j] = 0;
- cumulative_probability[i][j] = 0;
- }
- }
- //为每一只蚂蚁随机一个初始城市,并将其加入禁忌列表
- for (int i = 0; i < M; i++) {
- tabu[i][0] = rand() % N;
- allow[i][tabu[i][0]] = 0;
- }
- }
- //求取概率表和累积概率表,然后根据轮盘赌算法选择下一个城市,每轮迭代应该被调用N-1次(初始化时已经有一个城市了)
- void ChooseNextCity(Graph* graph, int times) {//times表示此轮迭代已经是第几次调用这个函数了,times应该从1开始,到N-1次结束
- //求概率表
- double sum = 0;
- int city_now = -1;
- for (int i = 0; i < M; ++i) {
- sum = 0;
- city_now = tabu[i][times - 1];//蚂蚁i目前在城市city_now
- for (int j = 0; j < N; ++j) {
- if (allow[i][j] == 1) {
- probability[i][j] = pow(graph->tao[city_now][j], α) * pow(graph->Heuristic[city_now][j], β);
- sum += probability[i][j];
- }
- else {
- probability[i][j] = 0;
- sum += 0;
- }
- }
- //上面求出的probability表并不是真正的概率表,而只是概率的分子,下面求出真正的概率表
- for (int j = 0; j < N; ++j) {
- probability[i][j] = probability[i][j] / sum;
- }
- }
- //用概率表求累积概率表
- for (int i = 0; i < M; ++i) {
- cumulative_probability[i][0] = probability[i][0];
- for (int j = 1; j < N; ++j) {
- cumulative_probability[i][j] = cumulative_probability[i][j - 1] + probability[i][j];
- }
- }
- //依据累积概率表,用轮盘赌算法为每只蚂蚁选择下一个城市
- double p = 0;
- for (int i = 0; i < M; ++i) {
- p = rand() % 1000;
- p /= 1000;
- for (int j = 0; j < N; ++j) {
- if (p <= cumulative_probability[i][j]) {
- //如果满足此式,则让蚂蚁i前往城市j(下面更新tabu表的操作就是从逻辑上把蚂蚁移动到城市j)
- tabu[i][times] = j;
- allow[i][j] = 0;
- break;
- }
- }
- }
- }
- //更新函数每次迭代后都应该被调用一次,用于更新信息素表graph->tabu,并记录本轮迭代的相关信息(最优路径,最优路径的长度,所有蚂蚁所走路径的平均长度)
- //NC_now表示目前的迭代次数,应该从0开始,到NC_MAX-1结束,依次调用此函数
- void Update(Graph* graph, int NC_now) {
- double delta_ktao = 0; //用于保存当前蚂蚁留在所经边的信息素
- double delta_tao[N][N] = { 0 }; //此轮迭代的信息素增量表,表示本轮信息素的增量,配合蒸发系数更新信息素表
- double sum_Length[M] = { 0 }; //保存此轮每只蚂蚁所经过的路径长度
- int last_city, next_city;
- //遍历tabu表,计算每只蚂蚁的路径长度,存放在sum_Length[]中
- for (int i = 0; i < M; ++i) {
- for (int j = 1; j < N; ++j) {
- last_city = tabu[i][j - 1];
- next_city = tabu[i][j];
- sum_Length[i] += graph->Length[last_city][next_city];
- }
- //还要再加上终点到起点的距离
- last_city = tabu[i][N - 1];
- next_city = tabu[i][0];
- sum_Length[i] += graph->Length[last_city][next_city];
- }
- //遍历tabu表,计算信息素增量表delta_tao[][]
- for (int i = 0; i < M; ++i) {
- delta_ktao = Q / sum_Length[i];
- for (int j = 1; j < N; ++j) {
- last_city = tabu[i][j - 1];
- next_city = tabu[i][j];
- delta_tao[last_city][next_city] += delta_ktao;
- }
- }
- //利用信息素增量表和蒸发系数更新信息素表tao[][]
- for (int i = 0; i < N; ++i) {
- for (int j = 0; j < N; ++j) {
- graph->tao[i][j] = graph->tao[i][j] * (1 - ρ) + delta_tao[i][j];
- }
- }
- //计算本轮最优路径,最优路径的长度,所有蚂蚁所走路径的平均长度
- int flag = 0;
- double min = 1000000, sum = 0; //min的初始值必须是一个足够大的值
- for (int i = 0; i < M; ++i) {
- sum += sum_Length[i];
- if (min > sum_Length[i]) {
- min = sum_Length[i];
- flag = i;//标记最好的蚂蚁
- }
- }
- //记录本轮信息至全局变量中
- for (int i = 0; i < N; ++i) {
- BestPath_PerRound[NC_now][i] = tabu[flag][i];
- }
- BestLength_PerRound[NC_now] = min;
- AverageLength_PerRound[NC_now] = sum / M;
- //更新全局最优路径和其长度
- if (BestLength > min) {
- for (int i = 0; i < N; ++i) {
- BestPath[i] = BestPath_PerRound[NC_now][i];
- }
- BestLength = min;
- }
- //用2-opt局部搜索优化本轮最优路径并用信息素强化
- {
- int TempPath[N];
- int flag = false;
- int t;
- for (int i = 0; i < N; ++i) {
- TempPath[i] = BestPath_PerRound[NC_now][i];
- }
- for (int a = 0; a < N-1; ++a) {
- //如果路径被优化为新路径,则再对新路径从头来一次局部搜索
- if (flag == true) {
- a = 0;
- flag = false;
- }
- for (int b = a + 1; b < N; ++b) {
- //逆序a~b的路径
- for (int i = a, j = b; i < j; ++i, --j) {
- t = TempPath[i];
- TempPath[i] = TempPath[j];
- TempPath[j] = t;
- }
- //重新求优化后的路径长度,加上终点到起点的距离
- sum = 0;
- for (int i = 1; i < N; ++i) {
- last_city = TempPath[i - 1];
- next_city = TempPath[i];
- sum += graph->Length[last_city][next_city];
- }
- last_city = TempPath[N - 1];
- next_city = TempPath[0];
- sum += graph->Length[last_city][next_city];
- //如果此解更优则更新本轮最优解
- if (sum < BestLength_PerRound[NC_now]) {
- flag = true;//如果路径被更新,则局部搜索从0开始
- //先更新平均路径再更新最优路径
- AverageLength_PerRound[NC_now] = (AverageLength_PerRound[NC_now] * M - BestLength_PerRound[NC_now] + sum) / M;
- BestLength_PerRound[NC_now] = sum;
- for (int i = 0; i < N; ++i) {
- BestPath_PerRound[NC_now][i] = TempPath[i];
- }
- //如果比全局最优解还优,则更新全局最优解
- if (sum < BestLength) {
- for (int i = 0; i < N; ++i) {
- BestPath[i] = TempPath[i];
- }
- BestLength = sum;
- }
- }
- }
- }
- //使用2opt局部搜索获得更强的解之后,给最优路径追加信息素奖励
- double reward = Q / sum;
- for (int i = 1; i < N; ++i) {
- last_city = BestPath[i-1];
- next_city = BestPath[i];
- graph->tao[last_city][next_city] += reward;
- }
- }
- }
- };
- //此函数的目的是为了用图展示出蚁群算法的结果,用C++调用python的matplotlib库
- bool ShowByPython(Graph* graph) {
- //C++中初始化python环境
- Py_Initialize();
- if (!Py_IsInitialized()) {
- std::cout << "Py_Initialize Error!" << std::endl;
- return false;
- }
- try {
- //执行python语句
- PyRun_SimpleString("import sys");
- PyRun_SimpleString("import matplotlib.pyplot as plt");
- string importPy = "sys.argv = ['suibiantian']";
- PyRun_SimpleString(importPy.c_str());
- double x, y;
- PyRun_SimpleString("plt.figure(num=1)");
- PyRun_SimpleString("plt.title('Ant Colony algorithm for solving TSP')");//蚁群算法解决TSP问题
- PyRun_SimpleString("plt.xlabel('x')");
- PyRun_SimpleString("plt.ylabel('y')");
- importPy= "plt.xlabel('Shortest=" + to_string(BestLength) + "')";
- PyRun_SimpleString(importPy.c_str());
- for (int i = 0; i < N; ++i) {
- x = graph->citys[i].x;//Graph[i].get_x();
- y = graph->citys[i].y;
- importPy = "plt.scatter(" + to_string(x) + "," + to_string(y) + ")";
- PyRun_SimpleString(importPy.c_str());
- }
- int last_city, next_city;
- double X1, X2, Y1, Y2;
- string str;
- //把全局最优路径按顺序画出来,画在第一副图上
- PyRun_SimpleString("plt.figure(num=1)");
- for (int i = 1; i < N; ++i) {
- last_city = BestPath[i - 1];
- next_city = BestPath[i];
- X1 = graph->citys[last_city].x;
- Y1 = graph->citys[last_city].y;
- X2 = graph->citys[next_city].x;
- Y2 = graph->citys[next_city].y;
- str = "plt.plot([" + to_string(X1) + "," + to_string(X2) + "],[" + to_string(Y1) + "," + to_string(Y2) + "])";
- PyRun_SimpleString(str.c_str());
- }
- //把终点到起点的线也画上
- last_city = BestPath[N - 1];
- next_city = BestPath[0];
- X1 = graph->citys[last_city].x;
- Y1 = graph->citys[last_city].y;
- X2 = graph->citys[next_city].x;
- Y2 = graph->citys[next_city].y;
- str = "plt.plot([" + to_string(X1) + "," + to_string(X2) + "],[" + to_string(Y1) + "," + to_string(Y2) + "])";
- PyRun_SimpleString(str.c_str());
- //下面展示随着迭代进行,每轮最优路径长度的变化,这是第二幅图
- double last = 0, next = 0;
- //string str;
- PyRun_SimpleString("plt.figure(num=2)");
- PyRun_SimpleString("plt.title('ShortestLength Convergence graph')");//最优路径收敛图
- PyRun_SimpleString("plt.xlabel('times of iterations')");
- PyRun_SimpleString("plt.ylabel('BestLength')");
- for (int i = 1; i < NC_MAX; ++i) {
- last = BestLength_PerRound[i - 1];
- next = BestLength_PerRound[i];
- str = "plt.plot([" + to_string(i - 1) + "," + to_string(i) + "],[" + to_string(last) + "," + to_string(next) + "])";
- PyRun_SimpleString(str.c_str());
- }
- //下面展示随着迭代进行,每轮平均路径长度的变化,这是第三幅图
- PyRun_SimpleString("plt.figure(num=3)");
- PyRun_SimpleString("plt.title('AverageLength Convergence graph')");//平均路径收敛图
- PyRun_SimpleString("plt.xlabel('times of iterations')");
- PyRun_SimpleString("plt.ylabel('AverageLength')");
- for (int i = 1; i < NC_MAX; ++i) {
- last = AverageLength_PerRound[i - 1];
- next = AverageLength_PerRound[i];
- str = "plt.plot([" + to_string(i - 1) + "," + to_string(i) + "],[" + to_string(last) + "," + to_string(next) + "])";
- PyRun_SimpleString(str.c_str());
- }
- //展示上面描绘的三幅图
- PyRun_SimpleString("plt.show()");
- }
- catch (...) {
- PyErr_Print();
- PyErr_Clear();
- Py_Finalize();
- return false;
- }
- Py_Finalize();
- return true;
- }
复制代码 |
|