三维空间中有一些(<=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\n",ans);
}
return 0;
}