BZOJ 1094 粒子运动

发表于 | 更新于

Link

Solution

模拟即可 写了一上午还好意思这么说 在外层循环枚举点对，转化成求点a,ba,b在整个过程中的最近距离。 可以根据粒子撞壁把整个过程划分成O(K)O(K)个阶段。 假设一个时刻，aa撞了壁，转向，bb正在去撞壁的路上。 用它们的速度算出来谁会先撞壁，求出来时间tt。在tt时间内，两个点的轨迹都是一个直线，所以可以用二次函数求出在这一端tt时间内两点距离的最小值。 如果用直线的一般式方程，化出的式子比较麻烦。而用直线的点法式方程，可以很简单地求出最值。 =(va−vb)2t2+2(va−vb)(Pa−Pb)t+(Pa−Pb)2\sqrt{(v_a-v_b)^2t^2+2(v_a-v_b)(P_a-P_b)t+(P_a-P_b)^2} 撞壁的时候需要计算转向，画个图感受感受就会了。。

Tips

没错，OI题也可以用二次函数求最值。 细节很多。。难写。。

Code

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
//Code by Lucida
#include<bits/stdc++.h>
#define red(x) scanf("%d",&x)
#define fred(x) scanf("%lf",&x)
template <class T> inline bool chkmx(T &a,const T &b){return a<b?a=b,1:0;}
template <class T> inline bool chkmn(T &a,const T &b){return a>b?a=b,1:0;}
typedef double ld;
const int MAXN=100+10;
typedef double ld;
const ld eps=1e-5,INF=1e100,pi=acos(-1.0);
using std::min;
using std::pair;
using std::make_pair;
int fcmp(ld x)
{
	if(-eps<x && x<eps) return 0;
	return x<0?-1:1;
}
template <class T> T sqr(T x){return x*x;}
template <class T> T abs(T x){return x>0?x:-x;}
struct vec
{
	ld x,y;
	vec(){}
	vec(ld _x,ld _y):x(_x),y(_y){}
	ld norm(){return sqrt(x*x+y*y);}//sqrt(y*y)???
	vec& operator +=(vec a){x+=a.x,y+=a.y;return *this;}
	vec& operator -=(vec a){x-=a.x,y-=a.y;return *this;}
	vec& operator *=(ld t){x*=t,y*=t;return *this;}
	vec& operator /=(ld t){x/=t,y/=t;return *this;}
};
typedef vec point;
vec operator +(vec a,vec b){return vec(a.x+b.x,a.y+b.y);}
vec operator -(vec a,vec b){return vec(a.x-b.x,a.y-b.y);}
vec operator *(vec a,ld t){return vec(a.x*t,a.y*t);}
vec operator /(vec a,ld t){return vec(a.x/t,a.y/t);}
ld inner(vec a,vec b){return a.x*b.x+a.y*b.y;}
ld outer(vec a,vec b){return a.x*b.y-a.y*b.x;}
struct circle
{
	point o;ld r;
	circle(){}
	circle(point _o,ld _r):o(_o),r(_r){}
};
struct line
{
	point a;vec v;
	line();
	line(point _a,vec _v):a(_a),v(_v){}
};
ld dist(point a,line l){return abs(outer(a-l.a,l.v))/l.v.norm();}
ld dist(point a,point b){return sqrt(sqr(a.x-b.x)+sqr(a.y-b.y));}
vec rotate(vec a,ld rad){return vec(a.x*cos(rad)-a.y*sin(rad),a.x*sin(rad)+a.y*cos(rad));}
point cross(line a,line b){return a.a+a.v*outer(a.a-b.a,b.v)/outer(b.v,a.v);}
pair<point,point> cross(line a,circle c)
{
	point p=cross(a,line(c.o,rotate(a.v,pi/2)));
	ld t=sqrt(sqr(c.r)-sqr(dist(c.o,a)))/a.v.norm();//t只是个长度 需要比norm才是比值。。naive
	return make_pair(p+a.v*t,p-a.v*t);
}
point inflex(line a,circle c)//有向直线的交点
{
	pair<point,point> res=cross(a,c);
	return fcmp(inner(a.v,res.first-a.a))<0?res.second:res.first;//a.v......
}
//ld sd(line a,line b){return min(dist(a.a,b.a),dist(a.a+a.v,b.a+b.v));} naive
ld sd(line a,line b,ld t)
{
	ld A=inner(a.v-b.v,a.v-b.v),B=2*inner(a.v-b.v,a.a-b.a),C=inner(a.a-b.a,a.a-b.a);
	ld axis=0;
	if(fcmp(A))
	{
		axis=-B/(2*A);
		if(fcmp(axis)<0) axis=0;
		else if(fcmp(axis-t)>0) axis=t;
	}
	else if(fcmp(B))
		axis=fcmp(B)<0?0:t;
	return sqrt(A*axis*axis+B*axis+C);
}
int K;
bool coincid(point a,point b){return fcmp(a.x-b.x)==0 && fcmp(a.y-b.y)==0;}
void turn(line &l,circle c)
{
	line axis=line(l.a,rotate(l.a-c.o,pi/2));
	line bot=line(l.a+l.v,c.o-l.a);
	point p=cross(axis,bot);
	l.v=l.a+l.v+(p-(l.a+l.v))*2-l.a;
}
ld calc(line li,line lj,circle c)
{
	ld res=INF;
	for(int i=0,j=0;i<K && j<K;)//已经碰撞了i j次
	{
		point pi=inflex(li,c),pj=inflex(lj,c);
		ld ti=dist(li.a,pi)/li.v.norm(),tj=dist(lj.a,pj)/lj.v.norm(),t; 
		chkmn(res,sd(li,lj,t=min(ti,tj)));
		if(coincid(li.a+=li.v*t,pi)) turn(li,c),i++;
		if(coincid(lj.a+=lj.v*t,pj)) turn(lj,c),j++;
	}
	return res;
}
int main()
{
	//freopen("input","r",stdin);
	static vec v[MAXN];
	static circle c;
	static point p[MAXN];
	fred(c.o.x),fred(c.o.y),fred(c.r);
	ld ANS=INF;int n;red(n),red(K);
	for(int i=1;i<=n;i++)
		fred(p[i].x),fred(p[i].y),fred(v[i].x),fred(v[i].y);
	for(int i=1;i<=n-1;i++)
		for(int j=i+1;j<=n;j++)
			chkmn(ANS,calc(line(p[i],v[i]),line(p[j],v[j]),c));
	printf("%.3lf\n",ANS);
	return 0;
}