中国邮递员问题

2018-04-19

中国邮递员问题

1. 问题描述

中国邮递员问题是邮递员在某一地区的信件投递路程问题。邮递员每天从邮局出发，走遍该地区所有街道再返回邮局，问题是他应如何安排送信的路线可以使所走的总路程最短。

2. 问题解决

图本身就是一个欧拉回路图，那么直接走一个欧拉回路就访问了所有路径并回到了起点

（欧拉回路图的前提条件是每个点的度数都是偶数）
图不是一个欧拉回路图，需要对某些边重复走多次回到起点，等价于添加了一些已经存在的重边，构建了欧拉回路。

不是欧拉回路需要通过增加一些边使得图变成一个欧拉回路，并能保证问题的最优解，添加的边权值和一定是尽可能短的。
- 计算度数为奇数的点
- 算出所有点之间的最短路径
- 奇度数的点一定为偶数，点与点之间构建二分图，权值为两点之间的最短路，找到一个最小权值匹配集
- 对于匹配集，找到每一个匹配集两点之间形成的最短路径，那么这条最短路径就需要加入到额外的边中，且一定能保证这条路径上除了这两个匹配的奇度数点度数加了1变成了偶数，其他所有中间点都是加了2，不影奇偶性。
- 所有边都添加好后，就是一个欧拉回路了，计算从起点开始的欧拉回路路径

3.相关算法

Fleury算法：求解欧拉回路路径，可以写成非常简洁的形式

dfs(source , euler_path);
    
void dfs(int u , vector<pii> &euler_path){
    for(int i=0 ; i<edges[u].size() ; i++){
        int v = edges[u][i];
        if(0 < euler_edge_degree[u][v]){
            //                cout<<"here debug: "<<u<<" -- "<<v<<endl;
            euler_edge_degree[u][v] -- ;
            euler_edge_degree[v][u] -- ;
            dfs(v , euler_path);
            euler_path.push_back(make_pair(u, v));
        }
    }
}

Floyd计算最短路算法

void floyd(){
        for(int i=0 ; i<n ; i++)
            for(int j=0 ; j<n ; j++)
                for (int k=0 ; k<n ; k++)
                    min_dis[i][j] = min(min_dis[i][j] , min_dis[i][k]+min_dis[k][j]);
    }

最大权值匹配KM算法(算最小权匹配，只要将权值设置为负数就可以了)

// 这里主要通过mp[i][j]存储连接的权值，下标从0开始，这部分代码没有考虑左右两点存在无权值的情况
int weight[N][N],match[N],lx[N],ly[N],visx[N],visy[N],slack[N];
int n , m;  //n表示左侧二分图的节点个数，m表示右侧二分图节点个数。
bool dfs(int x){
    visx[x]=1;
    for(int y=0 ; y<n ; y++){
        if(visy[y]) continue;
        int t=lx[x]+ly[y]-weight[x][y];
        if(t==0){
            visy[y]=1;
            if(match[y]==-1||dfs(match[y])){
                match[y]=x;
                return true;
            }
        }
        else if(slack[y]>t) slack[y]=t;
    }
    return false;
}
// KM() 返回的是最大权匹配，匹配结果保存在match[]中，可通过get_match_pair得到，要计算最小权匹配，初始化的时候权值加个负号就可以了
int KM(){
    memset(match,-1,sizeof(match));
    memset(lx,-INF,sizeof(lx));
    memset(ly,0,sizeof(ly));
    for(int i=0;i<n;i++){
        for(int j=0;j<n;j++){
            if(weight[i][j]>lx[i]) lx[i]=weight[i][j];
        }
    }
    for(int i=0;i<n;i++){
        for(int y=0;y<n;y++)
            slack[y]=INF;
        while(true){
            memset(visx,0,sizeof(visx));
            memset(visy,0,sizeof(visy));
            if(dfs(i)) break;
            int d=INF;
            for(int y=0;y<n;y++){
                if(!visy[y]&&d>slack[y]) d=slack[y];
            }
            for(int x=0;x<n;x++){
                if(visx[x]) lx[x]-=d;
            }
            for(int y=0;y<n;y++){
                if(visy[y]) ly[y]+=d;
                else slack[y]-=d;
            }
        }
    }
    int res=0;
    for(int i=0;i<n;i++){
        if(match[i]>-1) res+=weight[match[i]][i];
    }
    return res;
}

找到两点之间的最短路径Dijkstra算法

#define X first
#define Y second
int dis[N][N] , n;
// 找到x和y之间的最短路径，用 <<x,y>,w>的方式存下来
vector<piii> find_shortest_path(int x , int y){
    cout<<"debug find shortest path: "<<x<<" -- "<<y<<endl;
    vector<piii> shortest_path = vector<piii>();
    int lastpos[N];
    int optimal[N];
    memset(lastpos , -1 , sizeof(lastpos));
    memset(optimal , 0x3f , sizeof(optimal));
    set<pii> st;
    st.insert(make_pair(0 , x));
    optimal[x] = 0;
    while(!st.empty()){
        pii u = *(st.begin());
        st.erase(st.begin());
        //            cout<<"debug dijkstra: "<<u.Y<<" : "<<u.X<<endl;
        if(u.Y == y) break;
        if(u.X > optimal[u.Y]) continue;
        for(int v=0 ; v<n ; v++){
            if(dis[u.Y][v] != INF){
                if(u.X + dis[u.Y][v] < optimal[v]){
                    optimal[v] = u.X + dis[u.Y][v];
                    lastpos[v] = u.Y;
                    st.insert(make_pair(optimal[v] , v));
                }
            }
        }
    }
    int cur = y;
    while(lastpos[cur] != -1){
        // cout<<"debug last pos: "<<cur<<endl;
        shortest_path.push_back(make_pair(make_pair(lastpos[cur], cur), dis[lastpos[cur]][cur]));
        cur = lastpos[cur];
    }
    //        cout<<"debug 3 section finish"<<endl;
    return shortest_path;
}

4.整体代码实现

//
//  main.cpp
//  ChinesePostmanProblem
//
//  Created by ganjun on 2018/4/10.
//  Copyright © 2018年 ganjun. All rights reserved.
//
#include <iostream>
#include <vector>
#include <math.h>
#include <algorithm>
#include <map>
#include <set>
#include "stdio.h"
using namespace std;
#define N 105
#define INF 0x3f3f3f3f
#define pii pair<int,int>
#define piii pair<pii,int>
class Match{
public:
    // 这里主要通过mp[i][j]存储连接的权值，下标从0开始，这部分代码没有考虑左右两点存在无权值的情况
    int weight[N][N],match[N],lx[N],ly[N],visx[N],visy[N],slack[N];
    int n , m;  //n表示左侧二分图的节点个数，m表示右侧二分图节点个数。
    Match(){
        
    }
    Match(int weight[][N] , int n , int m){
        for(int i=0 ; i<n ; i++){
            for(int j=0 ; j<m ; j++){
                this->weight[i][j] = weight[i][j];
            }
        }
        this->n = n;
        this->m = m;
    }
    bool dfs(int x){
        visx[x]=1;
        for(int y=0 ; y<n ; y++){
            if(visy[y]) continue;
            int t=lx[x]+ly[y]-weight[x][y];
            if(t==0){
                visy[y]=1;
                if(match[y]==-1||dfs(match[y])){
                    match[y]=x;
                    return true;
                }
            }
            else if(slack[y]>t) slack[y]=t;
        }
        return false;
    }
    // KM() 返回的是最大权匹配，匹配结果保存在match[]中，可通过get_match_pair得到，要计算最小权匹配，初始化的时候权值加个负号就可以了
    int KM(){
        memset(match,-1,sizeof(match));
        memset(lx,-INF,sizeof(lx));
        memset(ly,0,sizeof(ly));
        for(int i=0;i<n;i++){
            for(int j=0;j<n;j++){
                if(weight[i][j]>lx[i]) lx[i]=weight[i][j];
            }
        }
        for(int i=0;i<n;i++){
            for(int y=0;y<n;y++)
                slack[y]=INF;
            while(true){
                memset(visx,0,sizeof(visx));
                memset(visy,0,sizeof(visy));
                if(dfs(i)) break;
                int d=INF;
                for(int y=0;y<n;y++){
                    if(!visy[y]&&d>slack[y]) d=slack[y];
                }
                for(int x=0;x<n;x++){
                    if(visx[x]) lx[x]-=d;
                }
                for(int y=0;y<n;y++){
                    if(visy[y]) ly[y]+=d;
                    else slack[y]-=d;
                }
            }
        }
        int res=0;
        for(int i=0;i<n;i++){
            if(match[i]>-1) res+=weight[match[i]][i];
        }
        return res;
    }
    
    vector<pii> get_match_pair(){
        KM();
        vector<pii> match_pairs = vector<pii>();
        for(int i=0 ; i<n ; i++){
            if(match[i] >= 0){
                cout<<i<<" -- " <<match[i]<<endl;
                match_pairs.push_back(make_pair(i , match[i]));
            }
        }
        return match_pairs;
    }
};
class ShortestPath{
public:
#define X first
#define Y second
    int dis[N][N] , n;
    int min_dis[N][N];
    ShortestPath(){
        
    }
    ShortestPath(int dis[][N] , int n){
        for(int i=0 ; i<n ; i++)
            for(int j=0 ; j<n ; j++)
                this->dis[i][j] = this->min_dis[i][j] = dis[i][j];
        
        this->n = n;
    }
    
    void floyd(){
        for(int i=0 ; i<n ; i++)
            for(int j=0 ; j<n ; j++)
                for (int k=0 ; k<n ; k++)
                    min_dis[i][j] = min(min_dis[i][j] , min_dis[i][k]+min_dis[k][j]);
    }
    // 找到x和y之间的最短路径，用 <<x,y>,w>的方式存下来
    vector<piii> find_shortest_path(int x , int y){
        cout<<"debug find shortest path: "<<x<<" -- "<<y<<endl;
        vector<piii> shortest_path = vector<piii>();
        int lastpos[N];
        int optimal[N];
        memset(lastpos , -1 , sizeof(lastpos));
        memset(optimal , 0x3f , sizeof(optimal));
        set<pii> st;
        st.insert(make_pair(0 , x));
        optimal[x] = 0;
        while(!st.empty()){
            pii u = *(st.begin());
            st.erase(st.begin());
//            cout<<"debug dijkstra: "<<u.Y<<" : "<<u.X<<endl;
            if(u.Y == y) break;
            if(u.X > optimal[u.Y]) continue;
            for(int v=0 ; v<n ; v++){
                if(dis[u.Y][v] != INF){
                    if(u.X + dis[u.Y][v] < optimal[v]){
                        optimal[v] = u.X + dis[u.Y][v];
                        lastpos[v] = u.Y;
                        st.insert(make_pair(optimal[v] , v));
                    }
                }
            }
        }
        int cur = y;
        while(lastpos[cur] != -1){
           // cout<<"debug last pos: "<<cur<<endl;
            shortest_path.push_back(make_pair(make_pair(lastpos[cur], cur), dis[lastpos[cur]][cur]));
            cur = lastpos[cur];
        }
//        cout<<"debug 3 section finish"<<endl;
        return shortest_path;
    }
};
class ChinesePostmanProblem{
public:
    // 请勿加入重边
    int dis[N][N] , n;
    int degree[N];
    ShortestPath path;
    map<int , int> vertex;
    
    // 后面部分用于处理欧拉回路的
    vector<int> edges[N];  // edge[x][i] = (y , w) 表示与x相连的第i条边的点是y，这里不需要权值了，权值可以直接去dis中获取
    int euler_edge_degree[N][N]; // euler_edge_degree[i][j]表示 i和j之间连了多少条边
    
    // 最终的路径以及总长度
    vector<pii> postman_path;
    int postman_path_len;
    
    ChinesePostmanProblem(int n){
        vertex.clear();
        memset(degree , 0 , sizeof(int)*n);
        memset(dis , 0x3f , sizeof(dis));
        postman_path.clear();
        this->n = n;
    }
    void add_edge(int x , int y , int w){
        dis[x][y] = w;
        dis[y][x] = w;
        degree[x] ++;
        degree[y] ++;
    }
    
    // 第一部分计算，找到度数为奇数的点
    vector<int> odd_degree_nodes(){
        vector<int> nodes = vector<int>();
        for(int i=0 ; i<n ; i++){
            if (degree[i] & 1)
                nodes.push_back(i);
        }
        return nodes;
    }
    
    // 第二部分计算，所有点之间的最短距离，并构建距离矩阵用于权值匹配算法，因为要映射到左右两边的点，所以用map构建映射关系
    // 注意一点 match是有一倍的重复的，因为1—n在左右两边都有，匹配的时候答案多了一倍
    Match optimal_match(vector<int> odd_nodes){
        int weight[N][N];
        int cnt = odd_nodes.size();
        path = ShortestPath(dis , n);
        path.floyd();
        
        for(int i=0 ; i<cnt ; i++){
            vertex[i] = odd_nodes[i];
        }
        for(int i=0 ; i<cnt ; i++){
            for(int j=0 ; j<cnt ; j++){
                if(i == j) weight[i][j] = -INF;  // 因为求的是最小权匹配，所以用的是负数
                else weight[i][j] = path.min_dis[odd_nodes[i]][odd_nodes[j]];
            }
        }
        Match match = Match(weight , cnt , cnt);
        return match;
    }
    
    // 第三部分，通过匹配的点映射回原来的点后找到需要额外添加的边 piii 表示<<x,y>,w>
    vector<piii> more_edges(vector<int> odd_nodes , Match match){
        vector<piii> added_edges = vector<piii>();
        vector<pii> match_pairs = match.get_match_pair();
        for(int i=0 ; i<match_pairs.size() ; i++){
            if(match_pairs[i].X > match_pairs[i].Y) continue; // 匹配有两倍的重复，只取一半
            vector<piii> optimal_path = path.find_shortest_path(odd_nodes[match_pairs[i].X], odd_nodes[match_pairs[i].Y]);
            for(int j=0 ; j<optimal_path.size() ; j++)
                added_edges.push_back(optimal_path[j]);
        }
        return added_edges;
    }
    
    // 第四部分，根据现有的无重边的dis[][]邻接矩阵表示，无边dis[i][j] = INF，和添加的边vector<piii> added_edges，找到欧拉回路 ， source表示起点
    vector<pii> euler_path(vector<piii> added_edges , int source){
        vector<pii> euler_path = vector<pii>();
        // 所有边的邻接表表示
        for(int i=0 ; i<n ; i++){
            edges[i].clear();
            for(int j=0 ; j<n ; j++){
                if(dis[i][j] != INF){
                    edges[i].push_back(j);
                }
            }
        }
        for(int i=0 ; i<added_edges.size() ; i++){
            int u = added_edges[i].X.X;
            int v = added_edges[i].X.Y;
            edges[u].push_back(v);
            edges[v].push_back(u);
        }
        // 边集处理完毕，进行euler回路操作(fleury 算法)
        memset(euler_edge_degree , 0 , sizeof(euler_edge_degree));
        for(int i=0 ; i<n ; i++)
            for(int j=0 ; j<edges[i].size() ; j++){
//                cout<<"debug 2: "<<i<<" -- "<<edges[i][j]<<endl;
                euler_edge_degree[i][edges[i][j]] += 1;
            }
        
        dfs(source , euler_path);
        return euler_path;
    }
    
    void dfs(int u , vector<pii> &euler_path){
        for(int i=0 ; i<edges[u].size() ; i++){
            int v = edges[u][i];
            if(0 < euler_edge_degree[u][v]){
//                cout<<"here debug: "<<u<<" -- "<<v<<endl;
                euler_edge_degree[u][v] -- ;
                euler_edge_degree[v][u] -- ;
                dfs(v , euler_path);
                euler_path.push_back(make_pair(u, v));
            }
        }
    }
    
    void postman_solve(int source){
        // 1
        vector<int> odd_degree_nodes = this->odd_degree_nodes();
        cout<<"first section : "<<endl;
        for(int i=0 ; i<odd_degree_nodes.size() ; i++){
            cout<<"odd_degree_nodes: "<<odd_degree_nodes[i]<<endl;
        }
        // 2
        Match _match = this->optimal_match(odd_degree_nodes);
        // 3
        cout<<"third section : "<<endl;
        vector<piii> more_edges = this->more_edges(odd_degree_nodes, _match);
        for(int i=0 ; i<more_edges.size() ; i++){
            cout<<"more edges: "<<more_edges[i].X.X<<" "<<more_edges[i].X.Y<<endl;
        }
        // 4
        postman_path = this->euler_path(more_edges, source);
        postman_path_len = 0;
        cout<<"fourth section :"<<endl;
        for(int i=0 ; i<postman_path.size() ; i++){
            postman_path_len += dis[postman_path[i].first][postman_path[i].second];
            cout<<"final path of index "<<i<<" : "<<postman_path[i].first<<" -- "<<postman_path[i].second<<endl;
        }
    }
    
};
int main(int argc, const char * argv[]) {
    // insert code here...
//    int weight[3][N] = {
//        {5, 3, 4},
//        {4, 2, 3},
//        {3, 5, 2}
//    };
//    Match match = Match(weight, 3, 3);
//    printf("%d\n" , match.KM());
//    match.get_match_pair();
//
//    int dis[5][N] = {
//        {INF , 3 , 4 , INF , INF},
//        {3 , INF , 2 , INF , 2},
//        {4 , 2 , INF , 4 , INF},
//        {INF , INF , 4 , INF , 1},
//        {INF , 2 , INF , 1 , INF}
//    };
//    ShortestPath path = ShortestPath(dis , 5);
//    vector<piii> find_optimal_path = path.find_shortest_path(0 , 3);
//    for(int i= 0 ; i < find_optimal_path.size() ; i++){
//        cout<<find_optimal_path[i].X.X<<" -- "<<find_optimal_path[i].X.Y<<" : "<<find_optimal_path[i].Y<<endl;
//    }
    /*
     1 5
     1 9
     5 3
     5 4
     5 8
     2 3
     2 4
     4 6
     6 7
     6 8
     7 8
     8 9
     */
    //test ChinesePostmanProblem
    ChinesePostmanProblem pro = ChinesePostmanProblem(9);
    pro.add_edge(0 , 4 , 1);
    pro.add_edge(0 , 8 , 1);
    pro.add_edge(4 , 2 , 1);
    pro.add_edge(4 , 3 , 1);
    pro.add_edge(4 , 7 , 1);
    pro.add_edge(1 , 2 , 1);
    pro.add_edge(1 , 3 , 1);
    pro.add_edge(3 , 5 , 1);
    pro.add_edge(5 , 6 , 1);
    pro.add_edge(5 , 7 , 1);
    pro.add_edge(6 , 7 , 1);
    pro.add_edge(7 , 8 , 1);
//    vector<pii> euler_path = pro.euler_path(vector<pair<pair<int, int>, int> > (), 3);
//    for(int i=0 ; i<euler_path.size() ; i++){
//        cout<<euler_path[i].first<<" -- "<<euler_path[i].second<<endl;
//    }
    pro.postman_solve(3); // 这里3表示邮递员从3号位置出发，最后回到3号位置
    cout<<"sum len: "<<pro.postman_path_len<<endl;
    return 0;
}