备战数学建模十图论之二边逐次修正法

找最小圈，然后类似破圈法求最小生成树解决最小旅行商为题的方案

import numpy as np
import networkx as nx
from numpy import inf
from matplotlib import pyplot as plt

def ebzcxz(D, H):
    # 计算初始H圈权值
    W = 0
    for i, j in zip(H, H[1:] + [H[0]]):
        W += D[i, j]
    print('初始H圈为:', [i + 1 for i in H], '权值和为;', W)

    # 计算顶点数
    N = D.shape[0]

    # 初始化迭代次数
    k_iter = 1

    while 1:
        # 取出H圈中一对相邻点
        for i1, j1 in zip(H, H[1:] + [H[0]]):
            # 取出H圈中另一对相邻点
            for i2, j2 in zip(H, H[1:] + [H[0]]):
                if (H.index(i2) - H.index(j1)) >= 1:
                    indi1, indj1, indi2, indj2 = H.index(i1), H.index(j1), H.index(i2), H.index(j2)

                    # 计算原H圈和改H圈不同路径位置的权和
                    ## 原H圈
                    origin = 0
                    for k in range(indi1, indi2):
                        origin += D[H[k], H[k + 1]]
                    origin += D[i2, j2]

                    ## 改H圈
                    new = 0
                    for k in range(indi2, indj1, -1):
                        new += D[H[k], H[k - 1]]
                    new += (D[i1, i2] + D[j1, j2])

                    # 改后的权和小于原权和则更新H圈和权
                    if new < origin:
                        W = W - origin + new
                        H = [H[indi1], H[indi2]] + list(reversed(H[indj1:indi2])) + (
                            H[indj2:] if indj2 != 0 else []) + H[:indi1]

                        print('第', k_iter, '次迭代产生的新H圈:', [i + 1 for i in H], '权值和为:', W)
                        k_iter += 1
                        break
            else:
                continue
            break
        else:
            break
def drawMat(A):
    '''
    画无向赋权图
    created by yzy
    2022.8.3
    '''
    # A=np.mat([[0,3,0,2,0],
    #    [3,0,2,1,0],
    #    [0,2,0,0,5],
    #    [2,1,0,0,4],
    #    [0,0,5,4,0]]) #邻接矩阵
    G=nx.Graph(A) #无向赋权图
    pos=nx.shell_layout(G) #布局设置
    w=nx.get_edge_attributes(G,"weight") #获取权重
    nx.draw_networkx(G=G,pos=pos,node_size=260) #画网络图
    nx.draw_networkx_edge_labels(G=G,pos=pos,font_size=12,edge_labels=w) #标注权重
    plt.show()

# 输出距离矩阵D，这里是无向图
D = np.array([
    [0, 56, 35, 2, 51, 60],
    [56, 0, 21, 57, 78, 70],
    [35, 21, 0, 36, 68, 68],
    [2, 57, 36, 0, 51, 61],
    [51, 78, 68, 51, 0, 13],
    [60, 70, 68, 61, 13, 0]
])
H = list(range(D.shape[0]))
ebzcxz(D, H)
drawMat(D)
# 输出距离矩阵D，此处是有向图即U12=25，U21=5
C = np.array([
    [0, 25, 40, 31, 27],
    [5, 0, 17, 30, 25],
    [19, 15, 0, 6, 1],
    [9, 50, 24, 0, 6],
    [22, 8, 7, 10, 0]
])
Q = list(range(C.shape[0]))
ebzcxz(C, Q)
# drawMat(D)画图函数画不了有向图

运行结果为

 

初始H圈为: [1, 2, 3, 4, 5, 6] 权值和为; 237
第 1 次迭代产生的新H圈: [1, 4, 3, 2, 5, 6] 权值和为: 210
第 2 次迭代产生的新H圈: [2, 6, 5, 1, 4, 3] 权值和为: 193
第 3 次迭代产生的新H圈: [5, 4, 1, 3, 2, 6] 权值和为: 192
初始H圈为: [1, 2, 3, 4, 5] 权值和为; 76
第 1 次迭代产生的新H圈: [3, 5, 4, 1, 2] 权值和为: 62

有个地方不完善，求初始圈的方法太固定了，应该多找几种初始圈然后分别二边逐次修正取最优方案。