[HDOJ 4892] Defence of the Trees [计算几何+最短路+状态压缩]

平面上有一些点，不超过300个，每个点有一个种类，种类个数不超过6。另外还有一些特殊点，不超过40个，问以这些特殊点为顶点(可以不全用)的多边形，内部包含所有种类的点，多边形的周长最小是多少。

首先我们可以知道满足条件的多边形一定是一个凸包，因为如果不是凸包的话，把它补成一个凸包不会让里边包含的点的种类减少，而且会让它的周长减小。

然后我们知道可以把任意的一个凸包转化成若干个三角形组成的一个链式结构(链式三角剖分)，即从凸包的任意一顶点开始，向所有其他顶点连线段。

我们从40个特殊点里边任取3个点，组成一个三角形，把这个三角形作为图G中的一个节点，记录每个三角形内部的点的种类(状态压缩，种类个数不超过6，所以压缩到64以内即可)。

我们记录相邻两个三角形的情况，对图G添加边。如果三角形ABC和三角形BCD相邻，即A和D在BC的异侧，那么从节点ABC到节点BCD的距离为BD+CD-BC。

这样我们给每个节点赋初值为其代表的三角形的周长。任意的一个多边形的周长都可由图G的一条链表示。

对这个图求最短路，找到最短的包含所有种类点的链，链的长度即为结果。如果没有这样的链，则输出Impossible

#include <cstdio>
#include <cmath>
#include <cstring>
#include <queue>

using namespace std;

const double eps=1e-9;

struct Point {
	double x,y;
	Point () {}
	Point (double xx,double yy) {
		x=xx;y=yy;
	}
	void read() {
		scanf("%lf%lf",&x,&y);
	}
	friend Point operator - (const Point &a,const Point &b) {
		return Point(a.x-b.x,a.y-b.y);
	}
};
double crossProduct(const Point &a,const Point &b) {
	return -a.x*b.y+b.x*a.y;
}
double dis(const Point &a,const Point &b) {
	double x=a.x-b.x,y=a.y-b.y;
	return sqrt(x*x+y*y);
}
bool inTriangle(Point &x,Point &a,Point &b,Point &c) {
	double sa=crossProduct(a-x,b-a),sb=crossProduct(b-x,c-b),sc=crossProduct(c-x,a-c);
	if (sa>-eps&&sb>-eps&&sc>-eps) return true;
	if (sa<eps&&sb<eps&&sc<eps) return true;
	return false;
}

int n,m,ki,kk,nn;
Point tree[300];
int treeKind[300];
Point stump[40];

struct Node {
	int fe,k;
	double v[64];
};
struct Edge {
	int t,ne;
	double v;
};
struct QueNode {
	int i,k;
	double v;
	QueNode() {}
	QueNode(int ii,int kk,double vv) {
		i=ii;k=kk;v=vv;
	}
	friend bool operator < (const QueNode &a,const QueNode &b) {
		return a.v>b.v;
	}
};

Node a[9880];
Edge b[548340];
priority_queue<QueNode> c;
int bp;

int getIndex(int i,int j,int k) {//get node number,0<=i<j<k<n
	return (n-i-1)*(n-i-2)*(n-i-3)/6+(n-j-1)*(n-j-2)/2+n-k-1;
}
void putEdge(int x,int y,double z) {
	//printf("%d %d %lf\n",x,y,z);
	b[bp].t=y;
	b[bp].ne=a[x].fe;
	b[bp].v=z;
	a[x].fe=bp++;
}
void add(int idx1,int idx2,Point &x,Point &y,Point &a,Point &b) {
	double sa=crossProduct(x-a,a-b),sb=crossProduct(y-a,a-b);
	if ((sa<eps&&sb<eps)||(sa>-eps&&sb>-eps)) return;
	putEdge(idx1,idx2,dis(a,x)+dis(b,x)-dis(a,b));
}
double dij() {
	while (!c.empty()) {
		int i=c.top().i,k=c.top().k;
		double v=c.top().v;
		//printf("top: %d %d %lf\n",i,k,v);
		if (k==kk-1) return v;
		c.pop();
		if (a[i].v[k]==v) {
			for (int j=a[i].fe;j!=-1;j=b[j].ne) {
				int ii=b[j].t,kk=k|a[ii].k;
				double vv=v+b[j].v;
				if (a[ii].v[kk]<-eps||a[ii].v[kk]>vv) {
					a[ii].v[kk]=vv;
					//printf("--in: %d %d %lf\n",ii,kk,vv);
					c.push(QueNode(ii,kk,vv));
				}
			}
		}
	}
	return -1;
}

int main() {
	int i,j,k,l;
	while (scanf("%d%d%d",&m,&n,&ki)!=EOF) {
		kk=1<<ki;nn=getIndex(0,0,0)+1;bp=0;
		while (!c.empty()) c.pop();
		for (i=0;i<nn;i++) {
			a[i].fe=-1;
			for (j=0;j<kk;j++) a[i].v[j]=-1;
			a[i].k=0;
		}
		for (i=0;i<m;i++) tree[i].read();
		for (i=0;i<m;i++) {
			scanf("%d",&treeKind[i]);
			treeKind[i]=1<<treeKind[i]-1;
		}
		for (i=0;i<n;i++) stump[i].read();
		bp=0;
		for (i=0;i<n;i++) 
			for (j=i+1;j<n;j++) 
				for (k=j+1;k<n;k++) {
					int idx=getIndex(i,j,k);
					for (l=0;l<m;l++) 
						if (inTriangle(tree[l],stump[i],stump[j],stump[k])) 
							a[idx].k|=treeKind[l];
					int tmpk=a[idx].k;
					a[idx].v[tmpk]=dis(stump[i],stump[j])+dis(stump[i],stump[k])+dis(stump[j],stump[k]);
					//printf("--%d %lf\n",idx,a[idx].v[tmpk]);
					c.push(QueNode(idx,tmpk,a[idx].v[tmpk]));
				}
		for (i=0;i<n;i++) 
			for (j=i+1;j<n;j++)
				for (k=j+1;k<n;k++) {
					int idx=getIndex(i,j,k);
					//printf("%d %d %d %d %d\n",i,j,k,idx,a[idx].k);
					for (l=0;l<n;l++) {
						if (l!=i&&l!=j&&l!=k) {
							int tmpIdx1,tmpIdx2,tmpIdx3;
							if (l<i) {
								tmpIdx1=getIndex(l,i,j);
								tmpIdx2=getIndex(l,j,k);
								tmpIdx3=getIndex(l,i,k);
							} else if (l<j) {
								tmpIdx1=getIndex(i,l,j);
								tmpIdx2=getIndex(l,j,k);
								tmpIdx3=getIndex(i,l,k);
							} else if (l<k) {
								tmpIdx1=getIndex(i,j,l);
								tmpIdx2=getIndex(j,l,k);
								tmpIdx3=getIndex(i,l,k);
							} else {
								tmpIdx1=getIndex(i,j,l);
								tmpIdx2=getIndex(j,k,l);
								tmpIdx3=getIndex(i,k,l);
							}
							add(idx,tmpIdx1,stump[l],stump[k],stump[i],stump[j]);
							add(idx,tmpIdx2,stump[l],stump[i],stump[j],stump[k]);
							add(idx,tmpIdx3,stump[l],stump[j],stump[k],stump[i]);
						}
					}
				}
		double ans=dij();
		if (ans>-eps) printf("%.12lf\n",ans);
		else printf("Impossible\n");
	}
	return 0;
}