双调欧几里德旅行商问题

         题目:

       一道经典的动态规划问题,但是路径一个环,  功力太浅,一开始就无从下手,完全找不到问题的最优子结构.最后还是翻看了英文答案,才找到思路.

       思考：

     先把点序列按x坐标排序，花费 (nlgn)时间，P1,P2,P3.......Pn 。把问题看成从pn点出发严格向左到p1，在从p1严格向右到pn，看成一条单线路径。  

       定义p(i,j)表示从i点出发严格向左到p1，再严格向右到pj的一条路径，再定义b[i][j]为路径p(i,j)的最短路径长度，从而可以看到所求问题为 b[n][n]; 定义好这些再来找问题的子结构就容易多了。

       找b[i][j]的最优解结构。 1<= i <= j <=n.

         1.当 i=j 时  ，需要的只有 b[n,n] , 这个问题只有一个子问题 b[n-1][n] ， 所以  b[n][n]=b[n-1][n] +| Pn-1Pn |; 

         2.当i != j时，是不是可以简单的取b[i][j]=b[i][j-1]呢，答案是否定的，因为i可能等于j-1

            a.i =j-1  ， Pj-1 在严格向左的路径上，而在严格向右的路径上必定存在一个点k， Pk为 Pj的直接前驱。 k<j-1 ，这样可以得到多个子问题，   b[j-1][j] =   b[k][j-1] +| PkPj |      ,  1<=k<j-1   其中取 b[k][j-1] +| PkPj |  最小的为最优解       

            b. i< j-1   Pj-1在严格向右的路径上， 可得 b[i][j]=b[i][j-1]+| Pj-1Pj | 

        3.特殊的， i=1 j=2 直接得到 b[1][2]=|P1P2|

 我自己感觉这个问题很抽象，很难想得特别明白，花了一下子问题图，能稍微感觉明白一点

               设一共4个点，  求 b[4][4]=b[3][4]+|P3P4|

                                             b[3][4]= min(     b[1][3]+|P1P4|  , b[2][3]+|P2P4|   )

                                             b[1][3]=b[1][2]+|P2P3|                b[2][3]=b[1][2]+|P1P3|

                     分解成的最小子问题就是 b[1][2] 加上某一段长 推出了的， 其实画在个表格更加明显。在下面程序代码中可以看出来，我是一行一行计算的。

#include<iostream>
#include<cmath>
#include<algorithm>
using namespace std;

class node{
public:
    double x;
    double y;
};

double  computerdis(node a,node b){
double sum=(b.x-a.x)*(b.x-a.x)+(b.y-a.y)*(b.y-a.y);
return   sqrt(sum);
}



void print(int **r,int pre,int j);

void printTour(int **r,int n){  //输出顺序是 Pn Pn-1 Pn-1的左前驱，...一直到 P1 再输出P1的右后继 一直向右到Pn
cout<<'P'<<n  <<"   ";    //输出 Pn
cout<<'P'<<n-1<<"   ";// 输出Pn-1
int k=r[n-1][n];     //Pn 向右路径上的的前驱
print(r,k,n-1);  //先输出向左路径上 Pn-1的前驱
cout<<'P'<<k<<"  ";
}


void print(int **r,int pre,int j)
{                             //无论 pre 和 j怎么取值，都是先输出Pj向左路径上的前驱
    if(pre<j){
        int k=r[pre][j]; // k为Pj的前驱
        if(k!=pre)  cout<<'P'<<k<<"   ";   //k=1的时候发生
        if(k>1) print(r,pre,k);
    }
    else{

        int k=r[j][pre];
        if(k>1){
            print(r,k,j);
            cout<<'P'<<k<<"  ";
        }
    }

    }


void BitonicTours(node *nodes,int n){ //O（n^2）
    double b[n+1][n+1];    //记录P(i,j)的最短双调路径长

    int **route=new int*[n+1];  //用于重构解 r[i][j]等于 P(i,j) 中Pj的直接前驱
    for(int i=0;i!=n+1;++i){
        route[i]=new int[n+1];
    }

 for(int i=1;i!=n+1;++i){   //从第一行开始填表 一直计算到n行 其中 i=j只计算 i=j=n , j>=i       O(n)

    for(int j=i;j!=n+1;++j){                          n+ n-1 +.....+1
        if(i==1&&j==2){ //b[1][2]=|p1p2|   情况3
                b[i][j]=computerdis(nodes[i],nodes[2]);
                continue;
            }
            else if(j==i+1){  //i=j-1  b[i][j]=b[j-1][j]=min{ b[p,i]+|P(p)Pj| }  情况2a
             double mmm=100000;
             for(int p=1;p!=j-1;++p){  //只发生一次
                if(b[p][i]+computerdis(nodes[p],nodes[j])<mmm){
                    mmm=b[p][i]+computerdis(nodes[p],nodes[j]);
                        route[i][j]=p;
                }
             } //for p
b[i][j]=mmm;
            }//2a

            else if(i==j&&i==n){ //i=j=n b[n][n]=b[n-1][n]+|P(n-1)P(n)|   //最后所要求的值
                b[i][j]=b[i-1][j]+computerdis(nodes[i-1],nodes[j]);

            }
            else if(i<j-1){// i<j-1   b[i][j]=b[i][j-1]+|P(j-1)P(j)|     //情况2b
                b[i][j]=b[i][j-1]+computerdis(nodes[j-1],nodes[j]);

            route[i][j]=j-1;
            }

        }//for j

    }//for i


printTour(route,n);  //重构打印路径

cout<<b[n][n]<<endl;
}

bool compa(node a,node b){
return a.x<b.x;
}

int main(){

node shu[8];

double xx[]={0, 5,8,6,0,7,2,1};
double yy[]={0, 4,2,1,6,5,3,0};


for(int i=1;i!=8;++i){
shu[i].x=xx[i];
shu[i].y=yy[i];
}

sort(shu+1,shu+8,compa);  //  nlgn    从小到大排序
BitonicTours(shu,7);

return 0;
}

因为计算表格只计算一般所以 O（n^2）;