1 importnumpy as np2 importmatplotlib.pyplot as plt3
4
5 #建立“蚂蚁”类
6 classAnt(object):7 def __init__(self, path):8 self.path = path #蚂蚁当前迭代整体路径
9 self.length = self.calc_length(path) #蚂蚁当前迭代整体路径长度
10
11 def calc_length(self, path_): #path=[A, B, C, D, A]注意路径闭环
12 length_ =013 for i in range(len(path_)-1):14 delta = (path_[i].x - path_[i+1].x, path_[i].y - path_[i+1].y)15 length_ +=np.linalg.norm(delta)16 returnlength_17
18 @staticmethod19 def calc_len(A, B): #静态方法,计算城市A与城市B之间的距离
20 return np.linalg.norm((A.x - B.x, A.y -B.y))21
22
23 #建立“城市”类
24 classCity(object):25 def __init__(self, x, y):26 self.x =x27 self.y =y28
29
30 #建立“路径”类
31 classPath(object):32 def __init__(self, A): #A为起始城市
33 self.path =[A, A]34
35 def add_path(self, B): #追加路径信息,方便计算整体路径长度
36 self.path.append(B)37 self.path[-1], self.path[-2] = self.path[-2], self.path[-1]38
39
40 #构建“蚁群算法”的主体
41 classACO(object):42 def __init__(self, ant_num=50, maxIter=300, alpha=1, beta=5, rho=0.1, Q=1):43 self.ants_num = ant_num #蚂蚁个数
44 self.maxIter = maxIter #蚁群最大迭代次数
45 self.alpha = alpha #信息启发式因子
46 self.beta = beta #期望启发式因子
47 self.rho = rho #信息素挥发速度
48 self.Q = Q #信息素强度
49 ###########################
50 self.deal_data('coordinates.dat') #提取所有城市的坐标信息
51 ###########################
52 self.path_seed = np.zeros(self.ants_num).astype(int) #记录一次迭代过程中每个蚂蚁的初始城市下标
53 self.ants_info = np.zeros((self.maxIter, self.ants_num)) #记录每次迭代后所有蚂蚁的路径长度信息
54 self.best_path = np.zeros(self.maxIter) #记录每次迭代后整个蚁群的“历史”最短路径长度
55 ###########################
56 self.solve() #完成算法的迭代更新
57 self.display() #数据可视化展示
58
59 defdeal_data(self, filename):60 with open(filename, 'rt') as f:61 temp_list = list(line.split() for line in f) #临时存储提取出来的坐标信息
62 self.cities_num = len(temp_list) #1. 获取城市个数
63 self.cities = list(City(float(item[0]), float(item[1])) for item in temp_list) #2. 构建城市列表
64 self.city_dist_mat = np.zeros((self.cities_num, self.cities_num)) #3. 构建城市距离矩阵
65 for i inrange(self.cities_num):66 A =self.cities[i]67 for j inrange(i, self.cities_num):68 B =self.cities[j]69 self.city_dist_mat[i][j] = self.city_dist_mat[j][i] =Ant.calc_len(A, B)70 self.phero_mat = np.ones((self.cities_num, self.cities_num)) #4. 初始化信息素矩阵
71 #self.phero_upper_bound = self.phero_mat.max() * 1.2 ###信息素浓度上限
72 self.eta_mat = 1/(self.city_dist_mat + np.diag([np.inf]*self.cities_num)) #5. 初始化启发函数矩阵
73
74 defsolve(self):75 iterNum = 0 #当前迭代次数
76 while iterNum <77 self.random_seed>
78 delta_phero_mat = np.zeros((self.cities_num, self.cities_num)) #初始化每次迭代后信息素矩阵的增量
79 ##########################################################################
80 for i inrange(self.ants_num):81 city_index1 = self.path_seed[i] #每只蚂蚁访问的第一个城市下标
82 ant_path = Path(self.cities[city_index1]) #记录每只蚂蚁访问过的城市
83 tabu = [city_index1] #记录每只蚂蚁访问过的城市下标,禁忌城市下标列表
84 non_tabu = list(set(range(self.cities_num)) -set(tabu))85 for j in range(self.cities_num-1): #对余下的城市进行访问
86 up_proba = np.zeros(self.cities_num-len(tabu)) #初始化状态迁移概率的分子
87 for k in range(self.cities_num-len(tabu)):88 up_proba[k] = np.power(self.phero_mat[city_index1][non_tabu[k]], self.alpha) *\89 np.power(self.eta_mat[city_index1][non_tabu[k]], self.beta)90 proba = up_proba/sum(up_proba) #每条可能子路径上的状态迁移概率
91 while True: #提取出下一个城市的下标
92 random_num =np.random.rand()93 index_need = np.where(proba >random_num)[0]94 if len(index_need) >0:95 city_index2 =non_tabu[index_need[0]]96 break
97 ant_path.add_path(self.cities[city_index2])98 tabu.append(city_index2)99 non_tabu = list(set(range(self.cities_num)) -set(tabu))100 city_index1 =city_index2101 self.ants_info[iterNum][i] =Ant(ant_path.path).length102 if iterNum == 0 and i == 0: #完成对最佳路径城市的记录
103 self.best_cities =ant_path.path104 else:105 if self.ants_info[iterNum][i] < Ant(self.best_cities).length: self.best_cities =ant_path.path106 tabu.append(tabu[0]) #每次迭代完成后,使禁忌城市下标列表形成完整闭环
107 for l inrange(self.cities_num):108 delta_phero_mat[tabu[l]][tabu[l+1]] += self.Q/self.ants_info[iterNum][i]109
110 self.best_path[iterNum] =Ant(self.best_cities).length111
112 self.update_phero_mat(delta_phero_mat) #更新信息素矩阵
113 iterNum += 1
114
115 defupdate_phero_mat(self, delta):116 self.phero_mat = (1 - self.rho) * self.phero_mat +delta117 #self.phero_mat = np.where(self.phero_mat > self.phero_upper_bound, self.phero_upper_bound, self.phero_mat) # 判断是否超过浓度上限
118
119 def random_seed(self): #产生随机的起始点下表,尽量保证所有蚂蚁的起始点不同
120 if self.ants_num <= self.cities_num: #蚂蚁数 <= 城市数
121 self.path_seed[:] =np.random.permutation(range(self.cities_num))[:self.ants_num]122 else: #蚂蚁数 > 城市数
123 self.path_seed[:self.cities_num] =np.random.permutation(range(self.cities_num))124 temp_index =self.cities_num125 while temp_index + self.cities_num <=self.ants_num:126 self.path_seed[temp_index:temp_index + self.cities_num] =np.random.permutation(range(self.cities_num))127 temp_index +=self.cities_num128 temp_left = self.ants_num %self.cities_num129 if temp_left !=0:130 self.path_seed[temp_index:] =np.random.permutation(range(self.cities_num))[:temp_left]131
132 def display(self): #数据可视化展示
133 plt.figure(figsize=(6, 10))134 plt.subplot(211)135 plt.plot(self.ants_info, 'g.')136 plt.plot(self.best_path, 'r-', label='history_best')137 plt.xlabel('Iteration')138 plt.ylabel('length')139 plt.legend()140 plt.subplot(212)141 plt.plot(list(city.x for city in self.best_cities), list(city.y for city in self.best_cities), 'g-')142 plt.plot(list(city.x for city in self.best_cities), list(city.y for city in self.best_cities), 'r.')143 plt.xlabel('x')144 plt.ylabel('y')145 plt.savefig('ACO.png', dpi=500)146 plt.show()147 plt.close()148
149
150 ACO()
77>