Abstract
ural 1913?Titan Ruins: Old Generators Are Fine Too
implementation 幾何 圓交
Source
http://acm.timus.ru/problem.aspx?space=1&num=1913
Solution
主要思路就是判斷距離到某兩個點是r，另一個點是2r的點。
比賽時候用三分亂搞，大概是被卡了精度還是什么的，沒有通過。
賽后用圓交重寫了一遍，一開始還是沒通過。
然后發(fā)現(xiàn)自己的圓交模板是針對面積問題的，單點的情況沒有考慮清楚，亂改的時候改錯了。大致就是 處理圓覆蓋的情況時，在預(yù)處理統(tǒng)計時算被覆蓋了一次，在算兩圓交點的時候又被統(tǒng)計了一次 。還是沒弄熟悉圓交的原理所致。
Code
#include <iostream>
#include
<cassert>
#include
<cmath>
#include
<cstdio>
#include
<cstring>
#include
<algorithm>
using
namespace
std;
const
double
EPS = 1e-
8
;
const
double
PI = acos(-
1.0
);
inline
int
sgn(
double
x) {
if
(fabs(x)<EPS)
return
0
;
return
x>
0
?
1
:-
1
;
}
struct
sp {
double
x, y;
double
pa;
int
cnt;
sp() {}
sp(
double
a,
double
b): x(a), y(b) {}
sp(
double
a,
double
b,
double
c,
int
d): x(a), y(b), pa(c), cnt(d) {}
bool
operator
<(
const
sp &rhs)
const
{
if
(sgn(pa-rhs.pa)==
0
)
return
cnt>
rhs.cnt;
return
pa<
rhs.pa;
}
void
read() {scanf(
"
%lf%lf
"
, &x, &
y);}
void
write() {printf(
"
%.10f %.10f\n
"
, x, y);}
}t, s1, s2, p[
5
], e[
20
];
double
r, rad[
5
];
int
cover[
5
];
sp
operator
*(
double
d,
const
sp &
v) {
return
sp(d*v.x, d*
v.y);
}
sp
operator
-(
const
sp &u,
const
sp &
v) {
return
sp(u.x-v.x, u.y-
v.y);
}
sp
operator
+(
const
sp &u,
const
sp &
v) {
return
sp(u.x+v.x, u.y+
v.y);
}
double
dot(
const
sp &u,
const
sp &
v) {
return
u.x*v.x+u.y*
v.y;
}
double
det(
const
sp &u,
const
sp &
v) {
return
u.x*v.y-u.y*
v.x;
}
double
dis(
const
sp &u,
const
sp &
v) {
double
dx = u.x-
v.x;
double
dy = u.y-
v.y;
return
sqrt(dx*dx+dy*
dy);
}
double
dissqr(
const
sp &u,
const
sp &
v) {
double
dx = u.x-
v.x;
double
dy = u.y-
v.y;
return
dx*dx+dy*
dy;
}
void
write(sp u, sp v) {
puts(
"
Now we have enough power
"
);
u.write(); v.write();
}
bool
cirint(
const
sp &u,
double
ru,
const
sp &v,
double
rv, sp &a, sp &
b) {
double
d =
dis(u, v);
if
(sgn(d-(ru+rv))>
0
|| sgn(d-fabs(ru-rv))<=
0
)
return
0
;
sp c
= u-
v;
double
ca, sa, cb, sb, csum, ssum;
ca
= c.x/d, sa = c.y/
d;
cb
= (rv*rv+d*d-ru*ru)/(
2
*rv*d), sb = sqrt(
1
-cb*
cb);
csum
= ca*cb-sa*
sb;
ssum
= sa*cb+sb*
ca;
a
= sp(rv*csum, rv*
ssum);
a
= a+
v;
sb
= -
sb;
csum
= ca*cb-sa*
sb;
ssum
= sa*cb+sb*
ca;
b
= sp(rv*csum, rv*
ssum);
b
= b+
v;
return
1
;
}
bool
cu(
int
N, sp p[],
double
r[], sp &
res) {
int
i, j, k;
memset(cover,
0
,
sizeof
cover);
for
(i =
0
; i < N; ++
i)
for
(j =
0
; j < N; ++
j) {
double
rd = r[i]-
r[j];
if
(i!=j && sgn(rd)>
0
&& sgn(dis(p[i], p[j])-rd)<=
0
)
cover[j]
++
;
}
sp a, b;
for
(i =
0
; i < N; ++
i) {
int
ecnt =
0
;
e[ecnt
++] = sp(p[i].x-r[i], p[i].y, -PI,
1
);
e[ecnt
++] = sp(p[i].x-r[i], p[i].y, PI, -
1
);
for
(j =
0
; j < N; ++
j) {
if
(j==i)
continue
;
if
(cirint(p[i], r[i], p[j], r[j], a, b)) {
e[ecnt
++] = sp(a.x, a.y, atan2(a.y-p[i].y, a.x-p[i].x),
1
);
e[ecnt
++] = sp(b.x, b.y, atan2(b.y-p[i].y, b.x-p[i].x), -
1
);
if
(sgn(e[ecnt-
2
].pa-e[ecnt-
1
].pa)>
0
) {
e[
0
].cnt++
;
e[
1
].cnt--
;
}
}
}
sort(e, e
+
ecnt);
int
cnt = e[
0
].cnt;
for
(j =
1
; j < ecnt; ++
j) {
if
(cover[i]+cnt==
N) {
double
a = (e[j].pa+e[j-
1
].pa)/
2
;
res
= p[i]+sp(r[i]*cos(a), r[i]*
sin(a));
return
1
;
}
cnt
+=
e[j].cnt;
}
}
return
0
;
}
bool
solve(
const
sp &o,
const
sp &u,
const
sp &v,
double
r) {
p[
0
] = o, rad[
0
] = r+
r;
p[
1
] = u; p[
2
] = v; rad[
1
] = rad[
2
] =
r;
sp a, b;
if
(cu(
3
, p, rad, a)) {
b
=
0.5
*(a-
o);
b
= o+
b;
write(a, b);
return
1
;
}
else
return
0
;
}
int
main() {
cin
>>
r;
t.read(); s1.read(); s2.read();
sp a, b, c, d;
if
(cirint(s1, r, s2, r, a, b)) {
if
(solve(t, s1, s2, r))
return
0
;
}
else
{
if
(cirint(s1, r, t, r, a, b) &&
cirint(s2, r, t, r, c, d)) {
write(a, c);
return
0
;
}
else
{
if
(solve(s1, t, s2, r))
return
0
;
if
(solve(s2, t, s1, r))
return
0
;
}
}
puts(
"
Death
"
);
return
0
;
}
順便貼圓交模板，這回的應(yīng)該沒問題了。
#include <cmath>
#include
<algorithm>
using
namespace
std;
const
double
PI = acos(-
1.0
);
const
double
EPS = 1e-
8
;
int
sgn(
double
d) {
if
(fabs(d)<EPS)
return
0
;
return
d>
0
?
1
:-
1
;
}
struct
sp {
double
x, y;
double
pa;
int
cnt;
sp() {}
sp(
double
a,
double
b): x(a), y(b) {}
sp(
double
a,
double
b,
double
c,
int
d): x(a), y(b), pa(c), cnt(d) {}
bool
operator
<(
const
sp &rhs)
const
{
/*
需要處理單點的情況時
if (sgn(pa-rhs.pa)==0) return cnt>rhs.cnt;
*/
return
pa<
rhs.pa;
}
void
read() {scanf(
"
%lf%lf
"
, &x, &
y);}
void
write() {printf(
"
%lf %lf\n
"
, x, y);}
};
sp
operator
+(
const
sp &u,
const
sp &
v) {
return
sp(u.x+v.x, u.y+
v.y);
}
sp
operator
-(
const
sp &u,
const
sp &
v) {
return
sp(u.x-v.x, u.y-
v.y);
}
double
det(
const
sp &u,
const
sp &
v) {
return
u.x*v.y-v.x*
u.y;
}
double
dir(
const
sp &p,
const
sp &u,
const
sp &
v) {
return
det(u-p, v-
p);
}
double
dis(
const
sp &u,
const
sp &
v) {
double
dx = u.x-
v.x;
double
dy = u.y-
v.y;
return
sqrt(dx*dx+dy*
dy);
}
//
計算兩圓交點
//
輸入圓心坐標(biāo)和半徑
//
返回兩圓是否相交(相切不算)
//
如果相交，交點返回在a和b(從a到b為u的交弧)
bool
cirint(
const
sp &u,
double
ru,
const
sp &v,
double
rv, sp &a, sp &
b) {
double
d =
dis(u, v);
/*
需要處理單點的情況時 覆蓋的情況在這里不統(tǒng)計
if (sgn(d-(ru+rv))>0 || sgn(d-fabs(ru-rv))<=0) return 0;
*/
if
(sgn(d-(ru+rv))>=
0
|| sgn(d-fabs(ru-rv))<=
0
)
return
0
;
sp c
= u-
v;
double
ca, sa, cb, sb, csum, ssum;
ca
= c.x/d, sa = c.y/
d;
cb
= (rv*rv+d*d-ru*ru)/(
2
*rv*d), sb = sqrt(
1
-cb*
cb);
csum
= ca*cb-sa*
sb;
ssum
= sa*cb+sb*
ca;
a
= sp(rv*csum, rv*
ssum);
a
= a+
v;
sb
= -
sb;
csum
= ca*cb-sa*
sb;
ssum
= sa*cb+sb*
ca;
b
= sp(rv*csum, rv*
ssum);
b
= b+
v;
return
1
;
}
sp e[
222
];
int
cover[
111
];
//
求圓并
//
輸入點數(shù)N，圓心數(shù)組p，半徑數(shù)組r，答案數(shù)組s
//
s[i]表示被至少i個圓覆蓋的面積(最普通的圓并就是s[1])
void
circle_union(
int
N, sp p[],
double
r[],
double
s[]) {
int
i, j, k;
memset(cover,
0
,
sizeof
cover);
for
(i =
0
; i < N; ++
i)
for
(j =
0
; j < N; ++
j) {
double
rd = r[i]-
r[j];
//
覆蓋的情況在這里統(tǒng)計
if
(i!=j && sgn(rd)>
0
&& sgn(dis(p[i], p[j])-rd)<=
0
)
cover[j]
++
;
}
for
(i =
0
; i < N; ++
i) {
int
ecnt =
0
;
e[ecnt
++] = sp(p[i].x-r[i], p[i].y, -PI,
1
);
e[ecnt
++] = sp(p[i].x-r[i], p[i].y, PI, -
1
);
for
(j =
0
; j < N; ++
j) {
if
(j==i)
continue
;
if
(cirint(p[i], r[i], p[j], r[j], a, b)) {
e[ecnt
++] = sp(a.x, a.y, atan2(a.y-p[i].y, a.x-p[i].x),
1
);
e[ecnt
++] = sp(b.x, b.y, atan2(b.y-p[i].y, b.x-p[i].x), -
1
);
if
(sgn(e[ecnt-
2
].pa-e[ecnt-
1
].pa)>
0
) {
e[
0
].cnt++
;
e[
1
].cnt--
;
}
}
}
sort(e, e
+
ecnt);
int
cnt = e[
0
].cnt;
for
(j =
1
; j < ecnt; ++
j) {
double
pad = e[j].pa-e[j-
1
].pa;
s[cover[i]
+cnt] += (det(e[j-
1
], e[j])+r[i]*r[i]*(pad-sin(pad)))/
2
;
cnt
+=
e[j].cnt;
}
}
}
更多文章、技術(shù)交流、商務(wù)合作、聯(lián)系博主
微信掃碼或搜索：z360901061
微信掃一掃加我為好友
QQ號聯(lián)系： 360901061
您的支持是博主寫作最大的動力，如果您喜歡我的文章，感覺我的文章對您有幫助，請用微信掃描下面二維碼支持博主2元、5元、10元、20元等您想捐的金額吧，狠狠點擊下面給點支持吧，站長非常感激您！手機微信長按不能支付解決辦法：請將微信支付二維碼保存到相冊，切換到微信，然后點擊微信右上角掃一掃功能，選擇支付二維碼完成支付。
【本文對您有幫助就好】元
