大家好呀！我们你们好久不见的。。。咳咳，初次见面的小编！
之前重新整理了ILS的代码，有人留言问能不能提供java版。
正好最近在学启发式算法和java，为了造福人类小编打算提供模拟退火法和迭代局部搜索求解TSP的java版本，方便一些不喜欢C++的同鞋~~
代码是基于我自己写的版本，但我是学习了公众号推文之后写的，同时有参照原文代码，因此虽然与C++代码有些许区别，但总体是类似的。
本文不再赘述TSP或者SA，ILS了，有需要请点击下方链接学习（一看就懂的那种哦！）：
干货 | 用模拟退火(SA, Simulated Annealing)算法解决旅行商问题
干货|迭代局部搜索算法(Iterated local search)探幽（附C++代码及注释）
不多说了，开始看代码吧~！
SA求解TSP的JAVA代码
SA分为四个类：MainRun，Data，Path，SimulatedAnnealing。
MainRun是程序的入口。
package SA;
/** * 主函数运行类 */
public class MainRun {
public static void main (String args[]) { long begintime = System.nanoTime();
Data.PrintData(); SimulatedAnnealing.SA(); SimulatedAnnealing.PrintPath();
long endtime = System.nanoTime(); double usedTime= (endtime - begintime)/(1e9); System.out.println(); System.out.println("程序耗时："+usedTime+"s"); }}
Data是数据类，存放SA和TSP的数据。
package SA;
/** * 数据类，包括： * TSP城市坐标，采用柏林52城 * SA系数。 */
public class Data { public static final double T0=50000.0, // 初始温度 T_min=1e-8, // 结束温度 q=0.98, // 退火系数 K=1; //公式中的常数K public static final int L=1000, // 每个温度时的迭代次数，即链长 N=52; // 城市数量
public static double[][] city_pos= //柏林52城城市坐标，最优解7542 { { 565,575 },{ 25,185 },{ 345,750 },{ 945,685 },{ 845,655 }, { 880,660 },{ 25,230 },{ 525,1000 },{ 580,1175 },{ 650,1130 }, { 1605,620 },{ 1220,580 },{ 1465,200 },{ 1530,5 }, { 845,680 },{ 725,370 },{ 145,665 }, { 415,635 },{ 510,875 },{ 560,365 },{ 300,465 },{ 520,585 },{ 480,415 }, { 835,625 },{ 975,580 },{ 1215,245 },{ 1320,315 },{ 1250,400 },{ 660,180 }, { 410,250 },{ 420,555 },{ 575,665 },{ 1150,1160 },{ 700,580 },{ 685,595 }, { 685,610 },{ 770,610 },{ 795,645 },{ 720,635 },{ 760,650 },{ 475,960 }, { 95,260 },{ 875,920 },{ 700,500 },{ 555,815 },{ 830,485 },{ 1170,65 }, { 830,610 },{ 605,625 },{ 595,360 },{ 1340,725 },{ 1740,245 } };
public static void PrintData(){ System.out.println("模拟退火算法，初始温度T0="+T0+ ",降温系数q="+q+",每个温度迭代"+L+"次"); }
}
Path是路径类，打包处理路径的静态方法。
package SA;import static SA.Data.*;import static java.lang.Math.*;
/** * 路径类，打包处理路径的静态方法： * 计算两点间距离 distance * 计算路径长度 path_len * 产生新解（新路径）create_new */
public class Path { public static double distance(double[] p1,double[] p2) //计算两点间距离 { double dis=0; dis=sqrt(pow(p1[0]-p2[0],2)+pow(p1[1]-p2[1],2)); return dis; }
public static double path_len(int[] city_list) //计算路径长度 { double path=0; for (int i=0;i<(N-1);i++) path+=distance(city_pos[city_list[i]],city_pos[city_list[i+1]]); path+=distance(city_pos[city_list[0]],city_pos[city_list[N-1]]) ; return path; }
public static void create_new(int[] city_list) //产生新解{ int i=(int) (random()*N); int j=(int) (random()*N); while(j==i) j=(int) (random()*N);
int temp; temp=city_list[i]; city_list[i]=city_list[j]; city_list[j]=temp; }}
SimulatedAnnealing开始模拟退火。
package SA;import static SA.Data.*; import static SA.Path.*;import static java.lang.Math.*;import java.util.*;
/** * 模拟退火过程 */
public class SimulatedAnnealing {
private static int[] bestpath=new int[N]; private static int count=0; private static double f2;
public static void SA() //模拟退火{ for(int i=0;i<N;i++) // 初始化总最优路径 bestpath[i]=i;
int[] newpath=Arrays.copyOf(bestpath,bestpath.length); // 新解
double f1; f2=path_len(bestpath); double T=T0;
while(T>T_min) { for (int i=0;i<L;i++) { create_new(newpath); f1=path_len(newpath); double de=f1-f2; if (de<0) { System.arraycopy(newpath, 0, bestpath, 0, bestpath.length); f2=f1; } else { double r = random(); if(exp(de/(K*T))<r) { System.arraycopy(newpath, 0, bestpath, 0,bestpath.length); f2=f1; } else System.arraycopy(bestpath, 0, newpath, 0, bestpath.length); } } T*=q; count++; } }
public static void PrintPath() //打印最优解{ System.out.println("共降温"+count+"次，得到的TSP最优距离为"+f2+"路径为"); for(int j=0;j<N;j++) { if (j==0) System.out.print(bestpath[j]); else System.out.print("--->"+bestpath[j]); } }}
ILS求解TSP的JAVA代码
ILS分为MainRun，City，Delta，Perturbation，LocalSearch，Solution。
主函数运行类，包括了ILS总方法。
package ILSTSP;import static ILSTSP.City.*; import static ILSTSP.Solution.*; import static ILSTSP.Perturbation.*;import static ILSTSP.LocalSearch.*;
/** * 主函数运行类，包括了ILS总方法以及计时功能。 */
public class MainRun { public static void main(String args[]) { int max_iterations = 600; int max_no_improve = 50; long begintime = System.nanoTime();
Solution best_solution=new Solution(); iterated_local_search(best_solution, max_iterations, max_no_improve);
System.out.println(); System.out.println("搜索完成！最优路线总长度 = "+best_solution.cost); System.out.println("最优访问城市序列如下：" ); for (int i = 0; i < CITY_SIZE;i++) System.out.print(best_solution.permutation[i]+"\t");
long endtime = System.nanoTime(); double usedTime= (endtime - begintime)/(1e9); System.out.println(); System.out.println("程序耗时："+usedTime+"s"); }
static void iterated_local_search(Solution best_solution, int max_iterations, int max_no_improve) { Solution current_solution = new Solution();
//获得初始随机解 random_permutation(best_solution.permutation);
best_solution.cost = cost_total(best_solution.permutation); local_search(best_solution, max_no_improve); //初始搜索
for (int i = 0; i < max_iterations; i++) { perturbation(best_solution,current_solution); //扰动+判断是否接受新解 local_search(current_solution, max_no_improve);//继续局部搜索
//找到更优解 if (current_solution.cost < best_solution.cost) { for (int j = 0; j < CITY_SIZE; j++) { best_solution.permutation[j] = current_solution.permutation[j]; } best_solution.cost = current_solution.cost; }
System.out.println("迭代搜索 "+i+ " 次\t" + "最优解 = "+ best_solution.cost+ "当前解 = "+ current_solution.cost ); } }}
City类存放城市数据，柏林52城坐标。
package ILSTSP;
import static java.lang.Math.*;
/** * TSP数据类， * 存放柏林52城坐标， * 两点间距离计算distance_2city。 */
public class City {
public static int CITY_SIZE=52; // 城市数量
public static int[][] cities= { { 565,575 },{ 25,185 },{ 345,750 },{ 945,685 },{ 845,655 }, { 880,660 },{ 25,230 },{ 525,1000 },{ 580,1175 },{ 650,1130 }, { 1605,620 },{ 1220,580 },{ 1465,200 },{ 1530,5 }, { 845,680 },{ 725,370 },{ 145,665 }, { 415,635 },{ 510,875 },{ 560,365 },{ 300,465 },{ 520,585 },{ 480,415 }, { 835,625 },{ 975,580 },{ 1215,245 },{ 1320,315 },{ 1250,400 },{ 660,180 }, { 410,250 },{ 420,555 },{ 575,665 },{ 1150,1160 },{ 700,580 },{ 685,595 }, { 685,610 },{ 770,610 },{ 795,645 },{ 720,635 },{ 760,650 },{ 475,960 }, { 95,260 },{ 875,920 },{ 700,500 },{ 555,815 },{ 830,485 },{ 1170,65 }, { 830,610 },{ 605,625 },{ 595,360 },{ 1340,725 },{ 1740,245 } };
public static double distance_2city(int[] c1,int[] c2) //计算两点间距离 { double dis=0; dis=sqrt(pow(c1[0]-c2[0],2)+pow(c1[1]-c2[1],2)); return dis; }}
Solution类，存放解，即路径。
package ILSTSP;
import static java.lang.Math.*; import static ILSTSP.City.*;
public class Solution { public int[] permutation=new int[CITY_SIZE]; //城市排列 public double cost;
//获取随机城市排列 public static void random_permutation(int[] cities_permutation){ int n=CITY_SIZE; int[] temp=new int[CITY_SIZE]; for(int i=0;i<CITY_SIZE;i++) temp[i]=i; for(int i=0;i<CITY_SIZE-1;i++) { int r=(int) random()*n; cities_permutation[i]=temp[r]; temp[r]=temp[n-1]; n--; } cities_permutation[CITY_SIZE-1]=temp[0]; }
public static double cost_total(int[] cities_permutation){ double total_distance = 0; for (int i = 0; i < CITY_SIZE-1; i++) total_distance += distance_2city(cities[cities_permutation[i]], cities[cities_permutation[i+1]]); total_distance += distance_2city(cities[cities_permutation[0]], cities[cities_permutation[CITY_SIZE-1]]); //最后一个城市和第一个城市计算距离 return total_distance; }}
扰动类。
package ILSTSP;import static ILSTSP.Solution.*; import static ILSTSP.City.*;import static java.lang.Math.*;
/** * 扰动类，跳出局部。 */
public class Perturbation { //扰动 public static void perturbation(Solution best_solution, Solution current_solution){ double_bridge_move(best_solution.permutation, current_solution.permutation); current_solution.cost = cost_total(current_solution.permutation); }
//将城市序列分成4块，然后按块重新打乱顺序。 //用于扰动函数 private static void double_bridge_move(int[] cities_permutation, int[] new_cities_permutation){ int[] pos=new int[5]; pos[0]= 0; pos[1]= (int) random()*(CITY_SIZE/3)+1; pos[2]= (int) (random()*(CITY_SIZE/3)+CITY_SIZE/3); pos[3]= (int) (random()*(CITY_SIZE/3)+(CITY_SIZE/3)*2); pos[4]= CITY_SIZE;
int n=4; int[] random_order=new int[4], temp=new int[4]; for(int i=0;i<4;i++) temp[i]=i; for(int i=0;i<3;i++) { int r=(int) (random()*n); random_order[i]=temp[r]; temp[r]=temp[n-1]; n--; } random_order[3]=temp[0];
int deadprotect=0; do { int i=0; for (int j=pos[random_order[0]];j<pos[random_order[0]+1];j++) { new_cities_permutation[i]=cities_permutation[j]; i++; }
for (int j=pos[random_order[1]];j<pos[random_order[1]+1];j++) { new_cities_permutation[i]=cities_permutation[j]; i++; }
for (int j=pos[random_order[2]];j<pos[random_order[2]+1];j++) { new_cities_permutation[i]=cities_permutation[j]; i++; }
for (int j=pos[random_order[3]];j<pos[random_order[3]+1];j++) { new_cities_permutation[i]=cities_permutation[j]; i++; }
deadprotect++; if (deadprotect==5) break;
}while(AcceptanceCriterion(cities_permutation, new_cities_permutation)); }
//判断接受准则 private static boolean AcceptanceCriterion(int[] cities_permutation, int[] new_cities_permutation){ int AcceptLimite=500; double c1=cost_total(cities_permutation); double c2=cost_total(new_cities_permutation)-AcceptLimite; if (c1<c2) return false; else return true; }}
局部搜索类
package ILSTSP;
import static ILSTSP.City.*; import static ILSTSP.Delta.*;
/** * 局部搜索类， * 采用反转i到j之间城市序列的邻域动作。 */
public class LocalSearch { //本地局部搜索，边界条件 max_no_improve //best_solution最优解 //current_solution当前解 public static void local_search(Solution best_solution, int max_no_improve){ int count = 0; int i, k; double inital_cost = best_solution.cost; //初始花费 double now_cost = 0;
Solution current_solution = new Solution(); //为了防止爆栈……直接new了，你懂的
for (i = 0; i < CITY_SIZE - 1; i++) { for (k = i + 1; k < CITY_SIZE; k++) { Delta[i][k] = calc_delta(i, k, best_solution.permutation); } }
do { //枚举排列 for (i = 0; i < CITY_SIZE - 1; i++) { for (k = i + 1; k < CITY_SIZE; k++) { //邻域动作 two_opt_swap(best_solution.permutation, current_solution.permutation, i, k); now_cost = inital_cost + Delta[i][k]; current_solution.cost = now_cost; if (current_solution.cost < best_solution.cost) { count = 0; //better cost found, so reset for (int j = 0; j < CITY_SIZE; j++) { best_solution.permutation[j] = current_solution.permutation[j]; } best_solution.cost = current_solution.cost; inital_cost = best_solution.cost; Update(i, k, best_solution.permutation); } } } count++; } while (count <= max_no_improve); }
//邻域动作 反转index_i private static void two_opt_swap(int[] cities_permutation, int[] new_cities_permutation, int index_i, int index_j){ for (int i = 0; i < CITY_SIZE; i++) { new_cities_permutation[i] = cities_permutation[i]; }
swap_element(new_cities_permutation, index_i, index_j); }
//颠倒数组中下标begin到end的元素位置 private static void swap_element(int[] p, int begin, int end){ int temp; while (begin < end) { temp = p[begin]; p[begin] = p[end]; p[end] = temp; begin++; end--; } }}
Delta类，在原推文中未讲解，这里略微讲解一下。
Delta实际上是对局部搜索的过程进行去重优化。
Delta[i][j]中存放的数字表示反转i，j间路径后对总距离的改变量。（我们之前没有定义计算总距离的方法，因为这次改为了由Delta计算总距离）
对calc_delta部分，每次翻转以后没必要再次重新计算Delta值，只需要在翻转的头尾做个小小处理。
比如：
有城市序列 1-2-3-4-5 总距离 = d_12 + d_23 + d_34 + d_45 + d_51 = A
翻转后的序列 1-4-3-2-5 总距离 = d_14 + d_43 + d_32 + d_25 + d_51 = B
由于 d_ij 与 d_ji是一样的，所以B也可以表示成 B = A – d_12 – d_45 + d_14 + d_25。
对Update部分，每次翻转以后不需要依次更新所有Delta值，有一部分是可以忽略的。
比如：
对于城市序列1-2-3-4-5-6-7-8-9-10，如果对3-5应用了邻域操作2-opt （即翻转）， 事实上对于Delta 1-2、6-10是不需要重复计算的。
package ILSTSP;
import static ILSTSP.City.*;
/** * Delta类， * 对局部搜索的过程进行去重优化。 * Delta[i][j]数组中存放的数字表示反转i，j间路径后对总距离的改变量。 */
public class Delta { public static double[][] Delta=new double[CITY_SIZE][CITY_SIZE];
public static double calc_delta(int i, int k, int[] tmp) { /* 以下计算说明： 对于每个方案，翻转以后没必要再次重新计算总距离 只需要在翻转的头尾做个小小处理
比如： 有城市序列 1-2-3-4-5 总距离 = d12 + d23 + d34 + d45 + d51 = A 翻转后的序列 1-4-3-2-5 总距离 = d14 + d43 + d32 + d25 + d51 = B 由于 dij 与 dji是一样的，所以B也可以表示成 B = A - d12 - d45 + d14 + d25 下面的优化就是基于这种原理 */ double delta=0; if((i==0)&&(k==CITY_SIZE-1)) delta=0; else { int i2=(i-1+CITY_SIZE)%CITY_SIZE; int k2=(k+1)%CITY_SIZE; delta = 0 - distance_2city(cities[tmp[i2]], cities[tmp[i]]) + distance_2city(cities[tmp[i2]], cities[tmp[k]]) - distance_2city(cities[tmp[k]], cities[tmp[k2]]) + distance_2city(cities[tmp[i]], cities[tmp[k2]]); } return delta; }
public static void Update(int i, int k, int[] tmp) { /* 去重处理，对于Delta数组来说，对于城市序列1-2-3-4-5-6-7-8-9-10，如果对3-5应用了邻域操作2-opt ， 事实上对于 7-10之间的翻转是不需要重复计算的。所以用Delta提前预处理一下。
当然由于这里的计算本身是O（1） 的，事实上并没有带来时间复杂度的减少（更新操作反而增加了复杂度） 如果delta计算 是O（n）的，这种去重操作效果是明显的。 */ if (i!=0 && k != CITY_SIZE - 1){ i --; k ++; for (int j = i; j <= k; j ++){ for (int l = j + 1; l < CITY_SIZE; l ++){ Delta[j][l] = calc_delta(j, l, tmp); } }
for (int j = 0; j < k; j ++){ for (int l = i; l <= k; l ++){ if (j >= l) continue; Delta[j][l] = calc_delta(j, l, tmp); } } }// 如果不是边界，更新(i-1, k + 1)之间的 else{ for (i = 0; i < CITY_SIZE - 1; i++) { for (k = i + 1; k < CITY_SIZE; k++) { Delta[i][k] = calc_delta(i, k, tmp); } } }// 边界要特殊更新 }}