遗传算法（GA）解决旅行商问题（TSP）

import numpy as np
import matplotlib.pyplot as plt



def citys_points(m, n=2):  
    # 生成m个城市的坐标，坐标值在[0,10]之间
    citys = np.round(np.random.uniform(0, 10, size=(m, n)), decimals=4)
    np.savetxt('citys.txt', citys, delimiter=',')
    return citys

def init(citys, pop_size=100):  
    [m,n] = np.shape(citys)
    # 初始化种群，每个个体都是一个城市序列
    init_pop = np.zeros((pop_size, m))  
    arr1 = np.arange(0, m)
    
    for i in range(pop_size):
        np.random.shuffle(arr1)
        init_pop[i, :] = arr1
    # 将序列首尾相连，形成一个环路
    init_pop = np.concatenate([init_pop, init_pop[:, 0].reshape(pop_size, 1)], axis=1)
    return init_pop


def get_fitness(pop, citys):
    # Calculate the pairwise Euclidean distance between adjacent rows of the citys variable and sum them up
    # 计算城市之间的距离
    [pop_size, num_genes] = np.shape(pop)  # 100,21
    distance = np.zeros((pop_size))

    for i, row in enumerate(pop):
        for index, item in enumerate(row):
            if index == 20:
                break
            p1 = citys[int(item)]
            # print(f"index:{index}, item:{item}")
            # print(p1, p2)
            p2 = citys[int(row[index + 1])]
            x1, y1 = p1[0], p1[1]
            # print(x1,y1)
            x2, y2 = p2[0], p2[1]
            # print(x2,y2)
            # print(np.sqrt((x2 - x1)**2 + (y2 - y1)**2))
            distance[i] += np.sqrt((x2 - x1)**2 + (y2 - y1)**2)
        #     break
        # break
    
    fitness = 1/distance
    fitness_max_index = np.argmax(fitness)
    fitness_max_item = fitness[fitness_max_index]
    best_route = pop[fitness_max_index]
    shortest_distance = distance[fitness_max_index]
    return fitness,fitness_max_item,shortest_distance,best_route


def roulette_select_pop(pop,citys):
    # 选择种群
    [pop_size, num_genes] = np.shape(pop)     # 100, 21
    res = get_fitness(pop, citys)
    fitness = res[0]
    # print(fitness)
    # Normalize fitness
    # 归一化适应度
    fit_max = fitness.max()
    fit_min = fitness.min()
    fitness_nomal = (fitness - fit_min)/(fit_max - fit_min)
    # print(fitness_nomal)
    fitness_prob = fitness_nomal/(np.sum(fitness_nomal))
    # print(fitness_prob)
    fitness_prob_cumsum = np.cumsum(fitness_prob)
    # print(fitness_prob_cumsum)

    seed = np.random.rand(pop_size)
    choosen_seed = np.sort(seed)
    # # Choose population based on cumulative sum of normalized fitness and random seed
    # # 基于归一化适应度和随机种子选择种群
    choosen_pop = np.zeros((pop_size, 21))
    i = 0
    j = 0

    while i < pop_size:
        # print(choosen_pop[i,:],init_pop[j,:],sep="\n -----------------------")
        if choosen_seed[i] < fitness_prob_cumsum[j]:
            choosen_pop[i, :] = pop[j, :]
            i = i + 1
        else:
            if j == 99:
                break
            j = j + 1
    return choosen_pop


def cross_gene(pop, p):  # pop 为带交换种群， p是交换概率
    # 交叉基因
    
    pop = np.delete(pop, 20, axis=1)
    # print(pop[0:10,:])
    [pop_size,num_genes] = np.shape(pop)
    for i in range(0, pop_size, 2):

        if np.random.rand() < p:

            gene_part = pop[i, 8:12].copy()
            index_item_arr = np.where(np.isin(pop[i + 1, :], gene_part))

            # print(gene_part)
            np.random.shuffle(index_item_arr)
            # print('4-0', np.shape(index_item_arr))
            for j in range(len(gene_part)):
                # print("gene_part",gene_part,"j:",j)
                # print('4:',index_item_arr[0][j],gene_part[j])
                pop[i,8+j] = pop[i+1][index_item_arr[0][j]]
                pop[i + 1, index_item_arr[0][j]] = gene_part[j]

            # print(pop[i],pop[i+1],sep="\n----------------------------\n")

    pop = np.concatenate((pop, pop[:, 0].reshape(100, 1)), axis=1)
    return pop


def mutation(pop, p):
        pop = np.delete(pop, 20, axis=1)
        for index,row in enumerate(pop):
            if np.random.rand(1) < p:
                # print(row)
                points_muta = np.random.randint(0, 20, size=(2))
                gene1 = row[points_muta[0]].copy()
                row[points_muta[0]] = row[points_muta[1]]
                row[points_muta[1]] = gene1
                # print(f"points_muta-----------{points_muta}")
                # index_points_muta = np.where(np.isin(row, points_muta))
                # print(f"here is row:{row}")
                pop[index,:] = row

        return np.concatenate((pop, pop[:, 0].reshape(100, 1)), axis=1)

def next_pop():
    pass

if __name__ == "__main__":
    m = 1000
    p_cross = 0.6
    p_muta  = 0.4
    citys = citys_points(20)
    # citys = np.loadtxt('citys.txt', delimiter=',')
    init_pop = init(citys)
    arr_fitness_max = np.empty(0)
    arr_shortest_distance = np.empty(0)
    
    fig, axs = plt.subplots(3, 1)
    plt.ion()
    for i in range(m):
        # plt.clf()  
        [fitness, fitness_max, shortest_distance, best_route] = get_fitness(init_pop, citys)
        # print(fitness_max,shortest_distance,best_route,"\n-----------------------\n")
        # print(best_route)
        arr_fitness_max = np.append(arr_fitness_max, fitness_max * 100)
        arr_shortest_distance = np.append(arr_shortest_distance, shortest_distance)
        # print(arr_fitness_max,arr_shortest_distance,sep="\n-------------------------\n")
        # ga.draw(arr_fitness_max, arr_shortest_distance, best_route, citys, i,m)
        axs[0].cla()
        axs[0].plot(np.arange(i + 1), arr_fitness_max)
        axs[0].set_xlim((0, m))
        axs[0].set_ylim((0, 4))  # set y-axis range to 0-4
        axs[0].annotate(f'Max Fitness: {fitness_max:.2f}', xy=(i, fitness_max*100), xytext=(i+10, fitness_max*100+0.5), arrowprops=dict(facecolor='black', shrink=0.05))

        # Plot shortest_distance as a line graph in the second subplot
        axs[1].cla()
        axs[1].plot(np.arange(i + 1), arr_shortest_distance)
        axs[1].set_xlim((0, m))
        axs[1].set_ylim((0, 100))  # set y-axis range to 0-100
        axs[1].annotate(f'Shortest Distance: {shortest_distance:.2f}', xy=(i, shortest_distance), xytext=(i+10, shortest_distance+5), arrowprops=dict(facecolor='black', shrink=0.05))

        # Plot best_route as a line graph in the third subplot
        x = [citys[int(i)][0] for i in best_route.astype(int)]
        y = [citys[int(i)][1] for i in best_route.astype(int)]
        # print(x,y,sep="\n--------------\n")
        # axs[2].set_xlim((0, 40))
        # axs[2].set_ylim((0, 40))
        axs[2].cla()

        axs[2].plot(x, y)
        # axs[2].annotate(f'Generation: {i}', xy=(x[0], y[0]), xytext=(x[0]+5, y[0]+5), arrowprops=dict(facecolor='black', shrink=0.05))
        for j, txt in enumerate(best_route.astype(int)):
            axs[2].annotate(txt, (x[j], y[j]))

        plt.pause(0.2)
        roul_pop = roulette_select_pop(init_pop, citys)

        cross_pop = cross_gene(roul_pop, p_cross)
        muta_pop = mutation(cross_pop, p_muta)

        new_pop = np.concatenate((muta_pop[:-1, :], best_route.reshape(1, 21)), axis=0)
        print(new_pop)
        init_pop = new_pop
    plt.ioff()
    plt.show()