【计算几何】【bitset】Gym

三维空间中有一些(<=2000)气球，一些光源(<=15)，给定一个目标点，问你在移除不超过K个气球的前提下，目标点所能接受到的最大光照。

枚举每个光源，预处理其若要照射到光源，需要移走哪些气球，构建成一个bitset。

然后再2^15枚举光源集合，看看要让集合中所有光源照到目标点所要移走的气球是否在K以内，尝试更新答案。

需要注意的一点是，三维叉积叉出来的向量的长度的绝对值，就是原来两个向量所张成的平行四边形面积的大小。

#include<cstdio>
#include<cmath>
#include<bitset>
#include<iostream>
using namespace std;
bitset<2001>S[16],Ss;
const double EPS=0.0000001;
struct Point{
	int x,y,z,t;
	Point(const int &x,const int &y,const int &z,const int &t){
		this->x=x;
		this->y=y;
		this->z=z;
		this->t=t;
	}
	Point(const int &x,const int &y,const int &z){
		this->x=x;
		this->y=y;
		this->z=z;
	}
	Point(){}
	void read(){
		scanf("%d%d%d%d",&x,&y,&z,&t);
	}
	double length(){
		return sqrt((double)x*(double)x+(double)y*(double)y+(double)z*(double)z);
	}
	int length2(){
		return x*x+y*y+z*z;
	}
}ba[2010],lig[17],aim;
typedef Point Vector;
Vector operator - (const Point &a,const Point &b){
	return Vector(a.x-b.x,a.y-b.y,a.z-b.z);
}
bool in(const Point &BA,const Point &p){
	return (BA-p).length2()<BA.t*BA.t;
}
Vector Cross(const Vector &a,const Vector &b){
	return Vector(a.y*b.z-a.z*b.y,a.z*b.x-a.x*b.z,a.x*b.y-a.y*b.x);
}
int dot(const Vector &a,const Vector &b){
	return a.x*b.x+a.y*b.y+a.z*b.z;
}
double DisToSegment(Point P,Point A,Point B)
{
    Vector v1=B-A,v2=P-A,v3=P-B;
    if(dot(v1,v2)<0) return (double)v2.length();
    else if(dot(v1,v3)>0) return (double)v3.length();
    else return fabs((double)Cross(v1,v2).length())/(double)v1.length();
}
int n,m,K;
double ans;
int main(){
//	freopen("g.in","r",stdin);
	while(1){
		scanf("%d%d%d",&n,&m,&K);
		if(!n && !m && !K){
			return 0;
		}
		ans=0;
		for(int i=1;i<=m;++i){
			S[i].reset();
		}
		for(int i=1;i<=n;++i){
			ba[i].read();
		}
		for(int i=1;i<=m;++i){
			lig[i].read();
		}
		scanf("%d%d%d",&aim.x,&aim.y,&aim.z);
		for(int i=1;i<=m;++i){
			for(int j=1;j<=n;++j){
				bool ina=in(ba[j],aim),inl=in(ba[j],lig[i]);
				if((ina^inl) || ((!ina && !inl) && DisToSegment(ba[j],aim,lig[i])-(double)ba[j].t<EPS)){
					S[i].set(j);
				}
			}
//			for(int j=1;j<=n;++j){
//				cout<<S[i][j];
//			}
//			puts("");
		}
		for(int i=0;i<(1<<m);++i){
			double tmp=0;
			Ss.reset();
			for(int j=1;j<=m;++j){
				if(i&(1<<(j-1))){
					Ss|=S[j];
					if(Ss.count()>K){
						goto OUT;
					}
					tmp+=(double)lig[j].t/(double)(aim-lig[j]).length2();
				}
			}
			ans=max(ans,tmp);
			OUT:;
		}
		printf("%.10lf
",ans);
	}
	return 0;
}