题解 P9141 【[THUPC 2023 初赛] 乱西星上的空战】
liuzhangfeiabc · · 题解
大家好我又来了，这次带给大家的依然是 THUPC 的传统艺能大模拟~
个人点评：这题是我个人写完并通过的算法竞赛中的大模拟题中最长的一道，但可能并非最难写的。甚至在我看来，在你有一个相对完整的三维计算几何模板的情况下，码长只有这题三分之一的 P7147 [THUPC2021 初赛] 麻将模拟器 的实现难度都不低于这题。
究其原因，这题的模块化特别强，各个模块之间基本上相对独立，特别适合进行抽象及使用各种 oop 技巧，尤其是在题目已经天然地给你划分好阶段的前提下。这样下来的整个编码逻辑是十分清晰的，尤其不会出现像我写 P7610 [THUPC2021] 群星连结 时各种逻辑线搅乱在一起根本不知道下一步该开哪一条的现象。
甚至于我在赛前预测，说不定会有队伍专门来写这题。如果队友们能有效地并行工作，每人写一个模块的话，应该是有希望在赛时写完并调过的，但是为什么没有队这么干呢 qwq。
我们一点点来看吧：
首先，正如麻将模拟器那题的真正难点在于 dp，这题的真正难点大概在于计算几何部分吧。除了常规的空间向量运算外，还有点到线段的距离、点在平面上的投影等。如果第一次写可能会绕不过弯来，但是熟练了或者有现成模板的话应该问题不大。
这里一个致命的小细节是我调试时才发现的：一开始误以为全都是整点的情况下应该有办法避免精度误差，因此自信地直接写实数比较，结果分分钟被教做人了。后来发现主要原因出自求角度时的 acos 函数，即使只有机器精度级别的误差，在后面判相等之类的地方都是致命的。所以不该偷懒的地方千万别偷懒啊。
namespace geometry{
#define EPS 1e-12l
inline ldb my_acosl(ldb x){return acosl(max(-1.0l,min(1.0l,x)));}
inline bool chkeq(ldb x,ldb y){return fabsl(x - y) < EPS;}
struct vec2{
ldb x,y;
ldb dis()const {return sqrtl(x * x + y * y);} //长度
ldb dis2()const {return x * x + y * y;} //长度的平方
ldb dot(const vec2 &p)const {return x * p.x + y * p.y;}
ldb cross(const vec2 &p)const {return x * p.y - y * p.x;}
ldb dis_to_rectangle(ldb lx,ldb hy)const {return min(fabsl(x - lx),fabsl(x + lx)) + min(fabsl(y - hy),fabsl(y + hy));}
};
struct vec3{
ldb x,y,z;
vec3 to_int(){x = round(x);y = round(y);z = round(z);return *this;}
ldb dis()const {return sqrtl(x * x + y * y + z * z);}
ldb dis2()const {return x * x + y * y + z * z;}
vec3 get_norm()const {ldb d = dis();return {x / d,y / d,z / d};}//获取向量单位化后的结果，但向量本身不单位化
vec3 norm(){ldb d = dis();x /= d;y /= d;z /= d;return *this;}//获取向量单位化后的结果并将向量本身单位化
ldb dot(const vec3 &p)const {return x * p.x + y * p.y + z * p.z;}
vec3 cross(const vec3 &p)const {return {y * p.z - z * p.y,z * p.x - x * p.z,x * p.y - y * p.x};}
ldb get_angle(const vec3 &p)const {return my_acosl(dot(p) / dis() / p.dis());}
void input(){x = read_ldb();y = read_ldb();z = read_ldb();}
};
vec3 operator + (const vec3 &a,const vec3 &b){return {a.x + b.x,a.y + b.y,a.z + b.z};}
vec3 operator - (const vec3 &a,const vec3 &b){return {a.x - b.x,a.y - b.y,a.z - b.z};}
vec3 operator * (ldb k,const vec3 &a){return {k * a.x,k * a.y,k * a.z};}
bool operator == (const vec3 &a,const vec3 &b){return chkeq(a.x,b.x) && chkeq(a.y,b.y) && chkeq(a.z,b.z);}
struct line{
vec3 p,v;
vec3 projection(const vec3 &x)const {return p + v.dot(x - p) * v;}
ldb get_min_dis(const vec3 &x)const {return (x - projection(x)).dis();}
};
struct segment{
vec3 p,q;
ldb len()const {return (p - q).dis();} //长度
ldb len2()const {return (p - q).dis2();} //长度的平方
ldb get_min_dis(const vec3 &x)const {//求点到线段的最近距离
if(p == q) return (x - p).dis();
ldb tmp = (x - p).dot(q - p);
if(tmp <= 0) return (x - p).dis();
if(tmp >= len2()) return (x - q).dis();
line y = {p,(q - p).norm()};
return y.get_min_dis(x);
}
};
struct plain{
vec3 p,n;
vec3 projection(const vec3 &x)const {return x - n.dot(x - p) * n;}
};
};接下来是无人机的定义。因为要预判无人机的走位，这里采用的方法是在一个对象里分别定义了无人机当前的运动状态和预判的运动状态。比较麻烦的部分大概就是求敌机在自己雷达平面上的投影，以及判断一个位置是否能飞到。前者在已有的计算几何板子下也并不难，后者主要就是把各种是否要滚转、向哪个方向滚转，以及到底是正杆还是负杆搞清楚就好。
struct plane{
int id;
int status;
bool team;
vec3 p,d,u,l; //当前的运动状态向量
vec3 nxtp,nxtd,nxtu,nxtl; //这一回合即将移动到的运动状态向量
ldb tu,td,r,vm,lx,hy;
int target;
bool tar_in_radar;
vector<int> exploded;
void getl(){l = u.cross(d);}
void input(int _id,bool _team){
id = _id;
team = _team;
p.input();d.input();u.input();getl();
tu = read_ldb();td = read_ldb();r = read_ldb();
vm = read_ldb();lx = read_ldb();hy = read_ldb();
status = ALIVE;
target = 0;
exploded.clear();
}
bool in_horizon(const plane &x)const {return d.dot(x.p - p) > EPS;}
bool in_nxt_horizon(const plane &x)const {return nxtd.dot(x.p - nxtp) > EPS;}
vec2 get_radar_r(const plane &x)const {//求x在自己的雷达平面上的投影
plain pl = {p,d};
vec3 nd = pl.projection(x.p);
return {l.dot(nd - p),u.dot(nd - p)};
}
vec2 get_nxt_radar_r(const plane &x)const {//求移动后x在自己的雷达平面上的投影
plain pl = {nxtp,nxtd};
vec3 nd = pl.projection(x.p);
return {nxtl.dot(nd - p),nxtu.dot(nd - p)};
}
bool in_radar_r(const vec2 &r)const {return fabsl(r.x) <= lx + EPS && fabsl(r.y) <= hy + EPS;}
bool can_reach(const vec3 &x){//飞机能否飞行x距离（飞到p+x位置），能的话更新nxtp,nxtu,nxtd
nxtp = p + x;
nxtd = x.get_norm();
if(d == nxtd) nxtl = l;//不需要滚转
else if(d == -1 * nxtd) return 0;
else{//需要滚转
nxtl = d.cross(nxtd).norm();//要通过滚转把l转到与d和nxtd所在的平面垂直的方向
if(l.dot(nxtl) < 0) nxtl = -1 * nxtl;//滚转只能在90度以内
}
nxtu = nxtd.cross(nxtl);
return l.get_angle(nxtl) / r //滚转
+ d.get_angle(nxtd) / (u.dot(nxtd) >= 0 ? tu : td) //俯仰
+ x.dis() / vm <= 1 + EPS; //直线飞行
}
}a[210]; 导弹的定义与无人机类似，鉴于同一时刻存在的来自于同一架无人机的导弹最多只会存在一枚，我这里的写法是一个导弹对象就是固定由某架无人机发射的，那么需要注意多次初始化的问题；另外导弹会比无人机多一个功能，即预判是否锁定目标及相应的锁定角，不过这不难写。
struct missile{
int id;
int status;
bool team;
vec3 p,d; //当前的运动状态向量
vec3 nxtp,nxtd; //这一回合即将移动到的运动状态向量
ldb tr,vm,ds,dp,bs;
int tz,timer,target;
void input(int _id,bool _team){
id = _id;
team = _team;
tr = read_ldb();vm = read_ldb();ds = read_ldb();
dp = read_ldb();bs = read_ldb();tz = read();
status = DEAD;
timer = target = 0;
}
void init(const vec3 &_p,const vec3 &_d,int _target){
p = _p;d = _d;
timer = 0;
target = _target;
status = INACTIVE;
}
bool nxt_can_lock(const plane &x,ldb &angle){//导弹即将飞到的位置能否锁定目标即将飞到的位置，能的话锁定角是多少（传给angle）
if(x.nxtp == nxtp){//两者目标位置相同
angle = 0;
return 1;
}
ldb dott = nxtd.dot(x.nxtp - nxtp);
angle = nxtd.get_angle(x.nxtp - nxtp);//锁定角
return dott > 0 && angle <= bs + EPS;//在前方且锁定角不超过最大锁定角
}
bool can_reach(const vec3 &x){//导弹能否飞行x距离（飞到p+x位置），能的话更新nxtp,nxtd
nxtp = p + x;
nxtd = x.get_norm();
return d.get_angle(nxtd) / tr //偏航
+ x.dis() / vm <= 1 + EPS; //直线飞行
}
}b[210];接下来，抽象出几个比较麻烦的步骤：
无人机目标选择：这一部分需要注意优先级为“先前选定的目标——在雷达范围内的目标（按到自己的距离排序）——在视野范围内但不在雷达范围内的目标（按到雷达范围的曼哈顿距离排序）”，然后分类讨论即可。为了后续实现方便，可以在无人机对象里记录这一步求出的“目标敌机是否在雷达范围内”信息。
void find_target(plane &x){//无人机选定目标
if(x.status != ALIVE) return;
int last_target = x.target;x.target = 0;x.tar_in_radar = 0;
ldb min_dis = INF,min_dis_r = INF;
for(int i = 1;i <= m;++i){
plane &y = a[i];
if(y.status != ALIVE || y.team == x.team || !x.in_horizon(y)) continue;
vec2 r = x.get_radar_r(y);
if(last_target == i){//先前的锁定目标，现在还能锁定
x.target = i;
x.tar_in_radar = x.in_radar_r(r);
return;
}
if(x.in_radar_r(r)){//在雷达范围内，优先选
ldb nw_dis = (x.p - y.p).dis2();
if(nw_dis < min_dis - EPS){//取距离自己最近的
x.target = i;
min_dis = nw_dis;
x.tar_in_radar = 1;
}
}
else{//不在雷达范围内
if(min_dis < INF) continue;
ldb nw_dis = r.dis_to_rectangle(x.lx,x.hy);
if(nw_dis < min_dis_r - EPS){//取距离雷达范围最近的
x.target = i;
min_dis_r = nw_dis;
}
}
}
}无人机和导弹的下一步飞行策略：这里出题人期待的做法就是大力枚举所有附近的整点，因为本身情况特别多，想给出一个比较方便的排除掉大量点的方案还是很困难的。本题时限之所以开得很大也正是希望确保所有在这里即使是最暴力枚举而不加任何优化的实现也可以轻松通过而不卡常。
在枚举终点后，是否可达由先前定义中已经给出的函数来判断；然后是优先级，这里千万要看清楚，比如无人机的移动在“目标敌机都仍落在视野内”的前提下，应当先比较“与目标敌机的距离”而非“目标敌机是否在雷达范围内”。
void get_plane_fly_info(plane &x){//求出无人机的飞行策略
if(x.status != ALIVE) return;
if(x.target){//如果无人机有选定目标
int vi = floor(x.vm);
vec3 nxtp = {0,0,0},nxtd = {0,0,0},nxtu = {0,0,0};
ldb min_tar_dis = INF,min_len_rq = INF,min_dis_radar = INF,min_straight_fly = INF;
for(int i = -vi;i <= vi;++i){
for(int j = -vi;j <= vi;++j){
for(int k = -vi;k <= vi;++k){
if(!i && !j && !k) continue;//不能原地不动
if(i * i + j * j + k * k > x.vm * x.vm) continue;
vec3 np = {(ldb)i,(ldb)j,(ldb)k};
if(!x.can_reach(np)) continue;//必须要合法到达
if(x.in_nxt_horizon(a[x.target])){//目标在视野内
ldb ds = (x.p + np - a[x.target].p).dis2();//到目标的距离
if(ds - EPS > min_tar_dis) continue;
vec2 r = x.get_nxt_radar_r(a[x.target]);
ldb tmp1 = INF,tmp2 = INF;
if(x.in_radar_r(r)) tmp1 = r.dis();//在雷达范围内
else tmp2 = r.dis_to_rectangle(x.lx,x.hy);//不在雷达范围内
if((ds < min_tar_dis - EPS) //优先：距离最小
|| (chkeq(ds,min_tar_dis) && tmp1 < min_len_rq - EPS) //其次：在雷达范围内且距视野中心尽可能近
|| (chkeq(ds,min_tar_dis) && min_len_rq == INF && tmp2 < min_dis_radar - EPS)){ //再次：不在雷达范围内，距雷达范围尽可能近
min_tar_dis = ds;
min_len_rq = tmp1;
min_dis_radar = tmp2;
nxtp = x.nxtp;
nxtd = x.nxtd;
nxtu = x.nxtu;
}
}
else{//目标不在视野内
if(min_tar_dis < INF) continue;
ldb ds = (np - x.vm * x.d).dis();
if(ds < min_straight_fly - EPS){//找距离直飞最近的整点
min_straight_fly = ds;
nxtp = x.nxtp;
nxtd = x.nxtd;
nxtu = x.nxtu;
}
}
}
}
}
x.nxtp = nxtp;
x.nxtd = nxtd;
x.nxtu = nxtu;
}
else{//眼镜蛇机动
x.nxtp = x.p;
x.nxtd = x.u;
x.nxtu = -1 * x.d;
}
}
void get_missile_fly_info(missile &x){//求出导弹的飞行策略
if(x.status == DEAD) return;
int vi = floor(x.vm);
vec3 nxtp = {0,0,0},nxtd = {0,0,0};
ldb min_tar_dis = INF,min_lock_angle = INF,min_straight_fly = INF;
for(int i = -vi;i <= vi;++i){
for(int j = -vi;j <= vi;++j){
for(int k = -vi;k <= vi;++k){
if(!i && !j && !k) continue;//不能原地不动
if(i * i + j * j + k * k > x.vm * x.vm) continue;
vec3 np = {(ldb)i,(ldb)j,(ldb)k};
if(!x.can_reach(np)) continue;//必须要合法到达
ldb tmp_angle = INF;
if(x.target && x.nxt_can_lock(a[x.target],tmp_angle)){//有目标，且位移之后仍能锁定
ldb ds = (x.p + np - a[x.target].nxtp).dis2();//到目标即将飞到的位置的距离
if((ds < min_tar_dis - EPS) //优先：距离最小
|| (chkeq(ds,min_tar_dis) && tmp_angle < min_lock_angle - EPS)){//其次：锁定角最小
min_tar_dis = ds;
min_lock_angle = tmp_angle;
nxtp = x.nxtp;
nxtd = x.nxtd;
}
}
else{//没有目标，或位移之后脱锁
if(min_tar_dis < INF) continue;
ldb ds = (np - x.vm * x.d).dis();
if(ds < min_straight_fly - EPS){//找距离直飞最近的整点
min_straight_fly = ds;
nxtp = x.nxtp;
nxtd = x.nxtd;
}
}
}
}
}
if(min_tar_dis == INF) x.target = 0;//脱锁
x.nxtp = nxtp;
x.nxtd = nxtd;
}
接下来是对各个阶段的拆解：
第一阶段，无人机选定目标并确定飞行策略：已经搞定了。
void stage1(){//所有无人机选定目标，并确定当前时刻内的飞行策略
for(int i = 1;i <= m;++i) find_target(a[i]);
for(int i = 1;i <= m;++i) get_plane_fly_info(a[i]);
}第二阶段，所有能发射导弹的无人机发射导弹：简单判断即可。
void launch_missile(plane &x){//发射导弹
if(x.status != ALIVE || !x.tar_in_radar) return;
if(b[x.id].status != DEAD) return;
b[x.id].init(x.p,(a[x.target].p - x.p).norm(),x.target);
}
void stage2(){//所有能发射导弹的无人机发射导弹
for(int i = 1;i <= m;++i) launch_missile(a[i]);
}第三阶段，所有导弹确定飞行策略并位移：注意判断是否有导弹会爆炸，这里要用到点到线段的最短距离，还要注意导弹是否激活，以及直接撞上去的特殊情况。建议把需要判断的内容都放在这里进行，这样接下来处理导弹爆炸就会很容易。
bool chk_missile_explode(const vec3 &st,const vec3 &ed,const vec3 &tar,ldb dp){
//从st飞到ed，距离tar的最近距离是否不超过dp
if(st == ed) return (tar - st).dis2() <= dp * dp;
segment seg = {st,ed};
return seg.get_min_dis(tar) <= dp;
}
void missile_fly(missile &x){//导弹飞行
if(x.status == DEAD) return;
for(int i = 1;i <= m;++i){
plane &y = a[i];
if(y.status == DEAD) continue;
if(x.nxtp == y.p //导弹直接撞到飞机上，无论是否激活都摧毁
|| (x.status == ACTIVE && chk_missile_explode(x.p,x.nxtp,y.p,x.dp))){//判断激活的导弹x能否摧毁飞机y
y.status = HANDLING;
y.exploded.pb(x.id);
}
}
x.p = x.nxtp;
x.d = x.nxtd;
}
void stage3(){//所有导弹确定飞行策略并位移，该过程中部分无人机可能被摧毁；
for(int i = 1;i <= m;++i) get_missile_fly_info(b[i]);
for(int i = 1;i <= m;++i) missile_fly(b[i]);
}第四阶段，所有可空爆的导弹爆炸并消失：这里可以与第六阶段一起写成一个函数，难度不大，主要在于准备输出信息的部分。
void crash(plane &x,bool flag){//飞机坠毁
if(x.status != HANDLING) return;
if(!flag){
++out_en.p1;
out_en.q1.pb(mp(x.id,x.exploded));
}
else{
++out_en.p2;
out_en.q2.pb(mp(x.id,x.exploded));
}
x.status = DEAD;
for(int i = 0;i < x.exploded.size();++i) b[x.exploded[i]].status = DEAD;
for(int i = 1;i <= m;++i){//瞄准它的导弹立刻脱锁
if(b[i].target == x.id) b[i].target = 0;
}
}
void stage4(){//所有可空爆的导弹爆炸并消失
for(int i = 1;i <= m;++i) crash(a[i],0);
}第五阶段，所有无人机按第一阶段中确定的飞行策略位移：与第三阶段类似，要处理导弹引爆的问题。
void plane_fly(plane &x){//无人机飞行
if(x.status == DEAD) return;
for(int i = 1;i <= m;++i){//检查飞机能否被别的导弹摧毁
missile &y = b[i];
if(y.status == DEAD) continue;
if(x.nxtp == y.p //导弹直接撞到飞机上，无论是否激活都摧毁
|| (y.status == ACTIVE && chk_missile_explode(x.p,x.nxtp,y.p,y.dp))){//判断激活的导弹y能否摧毁飞机x
x.status = HANDLING;
x.exploded.pb(i);
}
}
x.p = x.nxtp;
x.d = x.nxtd;
x.u = x.nxtu;
x.getl();
}
void stage5(){//所有无人机按 1. 中确定的飞行策略位移，该过程中部分无人机可能被摧毁；
for(int i = 1;i <= m;++i) plane_fly(a[i]);
}第六阶段，所有可空爆的导弹爆炸并消失：已经搞定了。
void stage6(){//所有可空爆的导弹爆炸并消失
for(int i = 1;i <= m;++i) crash(a[i],1);
}第七阶段，所有位置相同的无人机发生碰撞并坠毁：这里主要是一个碰撞检测，难度并不大，主要是准备输出内容，注意不要重复统计就好。
void chk_collide(plane &x){//碰撞检测
if(x.status != ALIVE) return;
bool flag = 0;
vector<int> vt;vt.clear();
for(int j = x.id + 1;j <= m;++j){
plane &y = a[j];
if(y.status != ALIVE) continue;
if(y.p == x.p){
x.status = y.status = DEAD;
if(!flag){
flag = 1;
vt.pb(x.id);
for(int k = 1;k <= m;++k){//瞄准它的导弹立刻脱锁
if(b[k].target == x.id) b[k].target = 0;
}
}
vt.pb(j);
for(int k = 1;k <= m;++k){//瞄准它的导弹立刻脱锁
if(b[k].target == j) b[k].target = 0;
}
}
}
if(flag){
++out_en.p3;
out_en.q3.pb(mp(x.id,vt));
}
}
void stage7(){//所有位置相同的无人机发生碰撞并坠毁
for(int i = 1;i <= m;++i) chk_collide(a[i]);
}第八阶段，导弹消失：因为导弹消失不会带走别的飞机，因此简单判断就好。
void self_explosion(missile &x){//导弹自爆
if(x.status == DEAD) return;
if(x.timer++ == x.tz //超时
|| (!x.target && x.status == ACTIVE)) //脱锁并且已经激活
x.status = DEAD;
}
void stage8(){//所有超过制导时长和脱锁且已激活的导弹消失。
for(int i = 1;i <= m;++i) self_explosion(b[i]);
}第九阶段，导弹激活：要注意导弹激活的前提除了远离发射它的无人机外，还有可能是发射它的无人机已经坠毁了。也都不难判断。
void activate(missile &x){//导弹激活
if(x.status != INACTIVE) return;
if(a[x.id].status == DEAD || (x.p - a[x.id].p).dis2() > x.ds * x.ds) x.status = ACTIVE;
}
void stage9(){//所有可激活的导弹被激活
for(int i = 1;i <= m;++i) activate(b[i]);
}大功告成！我们来看一下完整代码。
#include<bits/stdc++.h>
#define gc getchar()
#define pc putchar
#define li long long
#define uli unsigned li
#define fi first
#define se second
#define pb push_back
#define mp make_pair
#define pii pair<int,int>
#define md int mid = l + r >> 1
#define ls q << 1
#define rs q << 1 | 1
#define ln ls,l,mid
#define rn rs,mid + 1,r
#define ldb long double
#define INF 1e18l
using namespace std;
inline li read(){
li x = 0;
int y = 0,c = gc;
while(c < '0' || c > '9') y = c,c = gc;
while(c >= '0' && c <= '9') x = x * 10 + c - '0',c = gc;
return y == '-' ? -x : x;
}
inline ldb read_ldb(){
double x;scanf("%lf",&x);return (ldb)x;
}
inline void prt(li x){
if(x >= 10) prt(x / 10);
pc(x % 10 + '0');
}
inline void print(li x){
if(x < 0) pc('-'),x = -x;
prt(x);
}
int n,m,t;
/**************************GEOMETRY BEGIN********************************/
namespace geometry{
#define EPS 1e-12l
inline ldb my_acosl(ldb x){return acosl(max(-1.0l,min(1.0l,x)));}
inline bool chkeq(ldb x,ldb y){return fabsl(x - y) < EPS;}
struct vec2{
ldb x,y;
ldb dis()const {return sqrtl(x * x + y * y);} //长度
ldb dis2()const {return x * x + y * y;} //长度的平方
ldb dot(const vec2 &p)const {return x * p.x + y * p.y;}
ldb cross(const vec2 &p)const {return x * p.y - y * p.x;}
ldb dis_to_rectangle(ldb lx,ldb hy)const {return min(fabsl(x - lx),fabsl(x + lx)) + min(fabsl(y - hy),fabsl(y + hy));}
};
struct vec3{
ldb x,y,z;
vec3 to_int(){x = round(x);y = round(y);z = round(z);return *this;}
ldb dis()const {return sqrtl(x * x + y * y + z * z);}
ldb dis2()const {return x * x + y * y + z * z;}
vec3 get_norm()const {ldb d = dis();return {x / d,y / d,z / d};}//获取向量单位化后的结果，但向量本身不单位化
vec3 norm(){ldb d = dis();x /= d;y /= d;z /= d;return *this;}//获取向量单位化后的结果并将向量本身单位化
ldb dot(const vec3 &p)const {return x * p.x + y * p.y + z * p.z;}
vec3 cross(const vec3 &p)const {return {y * p.z - z * p.y,z * p.x - x * p.z,x * p.y - y * p.x};}
ldb get_angle(const vec3 &p)const {return my_acosl(dot(p) / dis() / p.dis());}
void input(){x = read_ldb();y = read_ldb();z = read_ldb();}
};
vec3 operator + (const vec3 &a,const vec3 &b){return {a.x + b.x,a.y + b.y,a.z + b.z};}
vec3 operator - (const vec3 &a,const vec3 &b){return {a.x - b.x,a.y - b.y,a.z - b.z};}
vec3 operator * (ldb k,const vec3 &a){return {k * a.x,k * a.y,k * a.z};}
bool operator == (const vec3 &a,const vec3 &b){return chkeq(a.x,b.x) && chkeq(a.y,b.y) && chkeq(a.z,b.z);}
struct line{
vec3 p,v;
vec3 projection(const vec3 &x)const {return p + v.dot(x - p) * v;}
ldb get_min_dis(const vec3 &x)const {return (x - projection(x)).dis();}
};
struct segment{
vec3 p,q;
ldb len()const {return (p - q).dis();} //长度
ldb len2()const {return (p - q).dis2();} //长度的平方
ldb get_min_dis(const vec3 &x)const {//求点到线段的最近距离
if(p == q) return (x - p).dis();
ldb tmp = (x - p).dot(q - p);
if(tmp <= 0) return (x - p).dis();
if(tmp >= len2()) return (x - q).dis();
line y = {p,(q - p).norm()};
return y.get_min_dis(x);
}
};
struct plain{
vec3 p,n;
vec3 projection(const vec3 &x)const {return x - n.dot(x - p) * n;}
};
};
using namespace geometry;
/**************************GEOMETRY END**********************************/
/**************************DEFINITION BEGIN******************************/
struct out_entry{
int p1,p2,p3;
vector<pair<int,vector<int> > >q1,q2,q3;
void init(){
q1.clear();q2.clear();q3.clear();
p1 = p2 = p3 = 0;
}
void out_sort(){
for(int i = 0;i < p1;++i) sort(q1[i].second.begin(),q1[i].second.end());
sort(q1.begin(),q1.end());
for(int i = 0;i < p2;++i) sort(q2[i].second.begin(),q2[i].second.end());
sort(q2.begin(),q2.end());
for(int i = 0;i < p3;++i){
sort(q3[i].second.begin(),q3[i].second.end());
q3[i].first = q3[i].second[0];
}
sort(q3.begin(),q3.end());
}
void output(){
out_sort();
print(p1);pc(' ');print(p2);pc(' ');print(p3);pc('\n');
for(int i = 0;i < p1;++i){
print(q1[i].first);pc(' ');print(q1[i].second.size());
for(int j = 0;j < q1[i].second.size();++j) pc(' '),print(q1[i].second[j]);
pc('\n');
}
for(int i = 0;i < p2;++i){
print(q2[i].first);pc(' ');print(q2[i].second.size());
for(int j = 0;j < q2[i].second.size();++j) pc(' '),print(q2[i].second[j]);
pc('\n');
}
for(int i = 0;i < p3;++i){
print(q3[i].second.size());
for(int j = 0;j < q3[i].second.size();++j) pc(' '),print(q3[i].second[j]);
pc('\n');
}
init();
}
}out_en;
#define DEAD 0
#define ALIVE 1
#define INACTIVE 2
#define ACTIVE 3
#define HANDLING 4
struct plane{
int id;
int status;
bool team;
vec3 p,d,u,l;
vec3 nxtp,nxtd,nxtu,nxtl;
ldb tu,td,r,vm,lx,hy;
int target;
bool tar_in_radar;
vector<int> exploded;
void getl(){l = u.cross(d);}
void input(int _id,bool _team){
id = _id;
team = _team;
p.input();d.input();u.input();getl();
tu = read_ldb();td = read_ldb();r = read_ldb();
vm = read_ldb();lx = read_ldb();hy = read_ldb();
status = ALIVE;
target = 0;
exploded.clear();
}
bool in_horizon(const plane &x)const {return d.dot(x.p - p) > EPS;}
bool in_nxt_horizon(const plane &x)const {return nxtd.dot(x.p - nxtp) > EPS;}
vec2 get_radar_r(const plane &x)const {//求x在自己的雷达平面上的投影
plain pl = {p,d};
vec3 nd = pl.projection(x.p);
return {l.dot(nd - p),u.dot(nd - p)};
}
vec2 get_nxt_radar_r(const plane &x)const {//求移动后x在自己的雷达平面上的投影
plain pl = {nxtp,nxtd};
vec3 nd = pl.projection(x.p);
return {nxtl.dot(nd - p),nxtu.dot(nd - p)};
}
bool in_radar_r(const vec2 &r)const {return fabsl(r.x) <= lx + EPS && fabsl(r.y) <= hy + EPS;}
bool can_reach(const vec3 &x){//飞机能否飞行x距离（飞到p+x位置），能的话更新nxtp,nxtu,nxtd
nxtp = p + x;
nxtd = x.get_norm();
if(d == nxtd) nxtl = l;//不需要滚转
else if(d == -1 * nxtd) return 0;
else{//需要滚转
nxtl = d.cross(nxtd).norm();//要通过滚转把l转到与d和nxtd所在的平面垂直的方向
if(l.dot(nxtl) < 0) nxtl = -1 * nxtl;//滚转只能在90度以内
}
nxtu = nxtd.cross(nxtl);
return l.get_angle(nxtl) / r //滚转
+ d.get_angle(nxtd) / (u.dot(nxtd) >= 0 ? tu : td) //俯仰
+ x.dis() / vm <= 1 + EPS; //直线飞行
}
}a[210];
struct missile{
int id;
int status;
bool team;
vec3 p,d;
vec3 nxtp,nxtd;
ldb tr,vm,ds,dp,bs;
int tz,timer,target;
void input(int _id,bool _team){
id = _id;
team = _team;
tr = read_ldb();vm = read_ldb();ds = read_ldb();
dp = read_ldb();bs = read_ldb();tz = read();
status = DEAD;
timer = target = 0;
}
void init(const vec3 &_p,const vec3 &_d,int _target){
p = _p;d = _d;
timer = 0;
target = _target;
status = INACTIVE;
}
bool nxt_can_lock(const plane &x,ldb &angle){//导弹即将飞到的位置能否锁定目标即将飞到的位置，能的话锁定角是多少（传给angle）
if(x.nxtp == nxtp){//两者目标位置相同
angle = 0;
return 1;
}
ldb dott = nxtd.dot(x.nxtp - nxtp);
angle = nxtd.get_angle(x.nxtp - nxtp);//锁定角
return dott > 0 && angle <= bs + EPS;//在前方且锁定角不超过最大锁定角
}
bool can_reach(const vec3 &x){//导弹能否飞行x距离（飞到p+x位置），能的话更新nxtp,nxtd
nxtp = p + x;
nxtd = x.get_norm();
return d.get_angle(nxtd) / tr //偏航
+ x.dis() / vm <= 1 + EPS; //直线飞行
}
}b[210];
/**************************DEFINITION END*********************************/
/**************************STAGE 1 BEGIN**********************************/
void find_target(plane &x){//无人机选定目标
if(x.status != ALIVE) return;
int last_target = x.target;x.target = 0;x.tar_in_radar = 0;
ldb min_dis = INF,min_dis_r = INF;
for(int i = 1;i <= m;++i){
plane &y = a[i];
if(y.status != ALIVE || y.team == x.team || !x.in_horizon(y)) continue;
vec2 r = x.get_radar_r(y);
if(last_target == i){//先前的锁定目标，现在还能锁定
x.target = i;
x.tar_in_radar = x.in_radar_r(r);
return;
}
if(x.in_radar_r(r)){//在雷达范围内，优先选
ldb nw_dis = (x.p - y.p).dis2();
if(nw_dis < min_dis - EPS){//取距离自己最近的
x.target = i;
min_dis = nw_dis;
x.tar_in_radar = 1;
}
}
else{//不在雷达范围内
if(min_dis < INF) continue;
ldb nw_dis = r.dis_to_rectangle(x.lx,x.hy);
if(nw_dis < min_dis_r - EPS){//取距离雷达范围最近的
x.target = i;
min_dis_r = nw_dis;
}
}
}
}
void get_plane_fly_info(plane &x){//求出无人机的飞行策略
if(x.status != ALIVE) return;
if(x.target){//如果无人机有选定目标
int vi = floor(x.vm);
vec3 nxtp = {0,0,0},nxtd = {0,0,0},nxtu = {0,0,0};
ldb min_tar_dis = INF,min_len_rq = INF,min_dis_radar = INF,min_straight_fly = INF;
for(int i = -vi;i <= vi;++i){
for(int j = -vi;j <= vi;++j){
for(int k = -vi;k <= vi;++k){
if(!i && !j && !k) continue;//不能原地不动
if(i * i + j * j + k * k > x.vm * x.vm) continue;
vec3 np = {(ldb)i,(ldb)j,(ldb)k};
if(!x.can_reach(np)) continue;//必须要合法到达
if(x.in_nxt_horizon(a[x.target])){//目标在视野内
ldb ds = (x.p + np - a[x.target].p).dis2();//到目标的距离
if(ds - EPS > min_tar_dis) continue;
vec2 r = x.get_nxt_radar_r(a[x.target]);
ldb tmp1 = INF,tmp2 = INF;
if(x.in_radar_r(r)) tmp1 = r.dis();//在雷达范围内
else tmp2 = r.dis_to_rectangle(x.lx,x.hy);//不在雷达范围内
if((ds < min_tar_dis - EPS) //优先：距离最小
|| (chkeq(ds,min_tar_dis) && tmp1 < min_len_rq - EPS) //其次：在雷达范围内且距视野中心尽可能近
|| (chkeq(ds,min_tar_dis) && min_len_rq == INF && tmp2 < min_dis_radar - EPS)){ //再次：不在雷达范围内，距雷达范围尽可能近
min_tar_dis = ds;
min_len_rq = tmp1;
min_dis_radar = tmp2;
nxtp = x.nxtp;
nxtd = x.nxtd;
nxtu = x.nxtu;
}
}
else{//目标不在视野内
if(min_tar_dis < INF) continue;
ldb ds = (np - x.vm * x.d).dis();
if(ds < min_straight_fly - EPS){//找距离直飞最近的整点
min_straight_fly = ds;
nxtp = x.nxtp;
nxtd = x.nxtd;
nxtu = x.nxtu;
}
}
}
}
}
x.nxtp = nxtp;
x.nxtd = nxtd;
x.nxtu = nxtu;
}
else{//眼镜蛇机动
x.nxtp = x.p;
x.nxtd = x.u;
x.nxtu = -1 * x.d;
}
}
void stage1(){//所有无人机选定目标，并确定当前时刻内的飞行策略
for(int i = 1;i <= m;++i) find_target(a[i]);
for(int i = 1;i <= m;++i) get_plane_fly_info(a[i]);
}
/**************************STAGE 1 END************************************/
/**************************STAGE 2 BEGIN**********************************/
void launch_missile(plane &x){//发射导弹
if(x.status != ALIVE || !x.tar_in_radar) return;
if(b[x.id].status != DEAD) return;
b[x.id].init(x.p,(a[x.target].p - x.p).norm(),x.target);
}
void stage2(){//所有能发射导弹的无人机发射导弹
for(int i = 1;i <= m;++i) launch_missile(a[i]);
}
/**************************STAGE 2 END************************************/
/**************************STAGE 3 BEGIN**********************************/
bool chk_missile_explode(const vec3 &st,const vec3 &ed,const vec3 &tar,ldb dp){
//从st飞到ed，距离tar的最近距离是否不超过dp
if(st == ed) return (tar - st).dis2() <= dp * dp;
segment seg = {st,ed};
return seg.get_min_dis(tar) <= dp;
}
void get_missile_fly_info(missile &x){//求出导弹的飞行策略
if(x.status == DEAD) return;
int vi = floor(x.vm);
vec3 nxtp = {0,0,0},nxtd = {0,0,0};
ldb min_tar_dis = INF,min_lock_angle = INF,min_straight_fly = INF;
for(int i = -vi;i <= vi;++i){
for(int j = -vi;j <= vi;++j){
for(int k = -vi;k <= vi;++k){
if(!i && !j && !k) continue;//不能原地不动
if(i * i + j * j + k * k > x.vm * x.vm) continue;
vec3 np = {(ldb)i,(ldb)j,(ldb)k};
if(!x.can_reach(np)) continue;//必须要合法到达
ldb tmp_angle = INF;
if(x.target && x.nxt_can_lock(a[x.target],tmp_angle)){//有目标，且位移之后仍能锁定
ldb ds = (x.p + np - a[x.target].nxtp).dis2();//到目标即将飞到的位置的距离
if((ds < min_tar_dis - EPS) //优先：距离最小
|| (chkeq(ds,min_tar_dis) && tmp_angle < min_lock_angle - EPS)){//其次：锁定角最小
min_tar_dis = ds;
min_lock_angle = tmp_angle;
nxtp = x.nxtp;
nxtd = x.nxtd;
}
}
else{//没有目标，或位移之后脱锁
if(min_tar_dis < INF) continue;
ldb ds = (np - x.vm * x.d).dis();
if(ds < min_straight_fly - EPS){//找距离直飞最近的整点
min_straight_fly = ds;
nxtp = x.nxtp;
nxtd = x.nxtd;
}
}
}
}
}
if(min_tar_dis == INF) x.target = 0;//脱锁
x.nxtp = nxtp;
x.nxtd = nxtd;
}
void missile_fly(missile &x){//导弹飞行
if(x.status == DEAD) return;
for(int i = 1;i <= m;++i){
plane &y = a[i];
if(y.status == DEAD) continue;
if(x.nxtp == y.p //导弹直接撞到飞机上，无论是否激活都摧毁
|| (x.status == ACTIVE && chk_missile_explode(x.p,x.nxtp,y.p,x.dp))){//判断激活的导弹x能否摧毁飞机y
y.status = HANDLING;
y.exploded.pb(x.id);
}
}
x.p = x.nxtp;
x.d = x.nxtd;
}
void stage3(){//所有导弹确定飞行策略并位移，该过程中部分无人机可能被摧毁；
for(int i = 1;i <= m;++i) get_missile_fly_info(b[i]);
for(int i = 1;i <= m;++i) missile_fly(b[i]);
}
/**************************STAGE 3 END************************************/
/**************************STAGE 4 BEGIN**********************************/
void crash(plane &x,bool flag){//飞机坠毁
if(x.status != HANDLING) return;
if(!flag){
++out_en.p1;
out_en.q1.pb(mp(x.id,x.exploded));
}
else{
++out_en.p2;
out_en.q2.pb(mp(x.id,x.exploded));
}
x.status = DEAD;
for(int i = 0;i < x.exploded.size();++i) b[x.exploded[i]].status = DEAD;
for(int i = 1;i <= m;++i){//瞄准它的导弹立刻脱锁
if(b[i].target == x.id) b[i].target = 0;
}
}
void stage4(){//所有可空爆的导弹爆炸并消失
for(int i = 1;i <= m;++i) crash(a[i],0);
}
/**************************STAGE 4 END************************************/
/**************************STAGE 5 BEGIN**********************************/
void plane_fly(plane &x){//无人机飞行
if(x.status == DEAD) return;
for(int i = 1;i <= m;++i){//检查飞机能否被别的导弹摧毁
missile &y = b[i];
if(y.status == DEAD) continue;
if(x.nxtp == y.p //导弹直接撞到飞机上，无论是否激活都摧毁
|| (y.status == ACTIVE && chk_missile_explode(x.p,x.nxtp,y.p,y.dp))){//判断激活的导弹y能否摧毁飞机x
x.status = HANDLING;
x.exploded.pb(i);
}
}
x.p = x.nxtp;
x.d = x.nxtd;
x.u = x.nxtu;
x.getl();
}
void stage5(){//所有无人机按 1. 中确定的飞行策略位移，该过程中部分无人机可能被摧毁；
for(int i = 1;i <= m;++i) plane_fly(a[i]);
}
/**************************STAGE 5 END************************************/
/**************************STAGE 6 BEGIN**********************************/
void stage6(){//所有可空爆的导弹爆炸并消失
for(int i = 1;i <= m;++i) crash(a[i],1);
}
/**************************STAGE 6 END************************************/
/**************************STAGE 7 BEGIN**********************************/
void chk_collide(plane &x){//碰撞检测
if(x.status != ALIVE) return;
bool flag = 0;
vector<int> vt;vt.clear();
for(int j = x.id + 1;j <= m;++j){
plane &y = a[j];
if(y.status != ALIVE) continue;
if(y.p == x.p){
x.status = y.status = DEAD;
if(!flag){
flag = 1;
vt.pb(x.id);
for(int k = 1;k <= m;++k){//瞄准它的导弹立刻脱锁
if(b[k].target == x.id) b[k].target = 0;
}
}
vt.pb(j);
for(int k = 1;k <= m;++k){//瞄准它的导弹立刻脱锁
if(b[k].target == j) b[k].target = 0;
}
}
}
if(flag){
++out_en.p3;
out_en.q3.pb(mp(x.id,vt));
}
}
void stage7(){//所有位置相同的无人机发生碰撞并坠毁
for(int i = 1;i <= m;++i) chk_collide(a[i]);
}
/**************************STAGE 7 END************************************/
/**************************STAGE 8 BEGIN**********************************/
void self_explosion(missile &x){//导弹自爆
if(x.status == DEAD) return;
if(x.timer++ == x.tz //超时
|| (!x.target && x.status == ACTIVE)) //脱锁并且已经激活
x.status = DEAD;
}
void stage8(){//所有超过制导时长和脱锁且已激活的导弹消失。
for(int i = 1;i <= m;++i) self_explosion(b[i]);
}
/**************************STAGE 8 END************************************/
/**************************STAGE 9 BEGIN**********************************/
void activate(missile &x){//导弹激活
if(x.status != INACTIVE) return;
if(a[x.id].status == DEAD || (x.p - a[x.id].p).dis2() > x.ds * x.ds) x.status = ACTIVE;
}
void stage9(){//所有可激活的导弹被激活
for(int i = 1;i <= m;++i) activate(b[i]);
}
/**************************STAGE 9 END************************************/
int main(){
int i;
n = read();t = read();m = n + n;
for(i = 1;i <= m;++i){
a[i].input(i,i <= n);b[i].input(i,i <= n);
}
out_en.init();
for(i = 1;i <= t;++i){
stage1();
stage2();
stage3();
stage4();
stage5();
stage6();
stage7();
stage8();
stage9();
out_en.output();
}
return 0;
}后记：
从造题者的角度讲，最大的难点大概不是写出完整代码，而是确保代码的正确性与题面的准确和无歧义性。因为这种题目即使是出题人也很难保证代码一开始就是对的而且完美契合题意的，许多 bug，无论是程序上的还是题目叙述上的，都是出题人和验题人反复对拍、调试、磨合后才能发现的。作为验题人，写完整个代码大概只有零零碎碎的几个小时，但是“与出题人拍过+检验题面正确性”这事花了我们整整一个晚上……
算了，心累了，以后不搞这玩意了（大概率下次会真香吧）。