competitive_programming_rs/geometry/
minimum_bounding_circle.rs1pub mod minimum_bounding_circle {
2 pub fn make_circle(ps: &[Point<f64>]) -> (Point<f64>, f64) {
3 let n = ps.len();
4 assert!(n >= 2);
5
6 let mut c = make_circle2(ps[0], ps[1]);
7 for i in 2..n {
8 if is_included(ps[i], c) {
9 continue;
10 }
11 c = make_circle2(ps[0], ps[i]);
12 for j in 1..i {
13 if is_included(ps[j], c) {
14 continue;
15 }
16 c = make_circle2(ps[i], ps[j]);
17 for k in 0..j {
18 if is_included(ps[k], c) {
19 continue;
20 }
21 c = make_circle3(ps[i], ps[j], ps[k]);
22 }
23 }
24 }
25 c
26 }
27
28 fn make_circle3(a: Point<f64>, b: Point<f64>, c: Point<f64>) -> (Point<f64>, f64) {
29 let ea = (b - c).norm();
30 let eb = (c - a).norm();
31 let ec = (a - b).norm();
32 let s = (b - a).det(&(c - a));
33
34 let center = (a * ea * (eb + ec - ea) + b * eb * (ec + ea - eb) + c * ec * (ea + eb - ec))
35 / (s * s * 4.0);
36 let r2 = (center - a).norm();
37 (center, r2)
38 }
39
40 fn make_circle2(a: Point<f64>, b: Point<f64>) -> (Point<f64>, f64) {
41 let c = (a + b) / 2.0;
42 let r2 = (a - c).norm();
43 (c, r2)
44 }
45
46 fn is_included(a: Point<f64>, circle: (Point<f64>, f64)) -> bool {
47 let (center, r2) = circle;
48 (a - center).norm() <= r2
49 }
50
51 pub struct Point<T> {
52 pub x: T,
53 pub y: T,
54 }
55
56 impl<T> Copy for Point<T> where T: Copy {}
57
58 impl<T> Clone for Point<T>
59 where
60 T: Clone,
61 {
62 fn clone(&self) -> Point<T> {
63 Point {
64 x: self.x.clone(),
65 y: self.y.clone(),
66 }
67 }
68 }
69
70 impl<T> std::ops::Sub for Point<T>
71 where
72 T: std::ops::Sub<Output = T>,
73 {
74 type Output = Point<T>;
75 fn sub(self, other: Point<T>) -> Point<T> {
76 Point {
77 x: self.x - other.x,
78 y: self.y - other.y,
79 }
80 }
81 }
82
83 impl<T> std::ops::Mul<T> for Point<T>
84 where
85 T: Copy + std::ops::Mul<Output = T>,
86 {
87 type Output = Point<T>;
88 fn mul(self, rhs: T) -> Point<T> {
89 Point {
90 x: self.x * rhs,
91 y: self.y * rhs,
92 }
93 }
94 }
95
96 impl<T> std::ops::Div<T> for Point<T>
97 where
98 T: Copy + std::ops::Div<Output = T>,
99 {
100 type Output = Point<T>;
101 fn div(self, rhs: T) -> Point<T> {
102 Point {
103 x: self.x / rhs,
104 y: self.y / rhs,
105 }
106 }
107 }
108
109 impl<T> std::ops::Add<Point<T>> for Point<T>
110 where
111 T: std::ops::Add<Output = T>,
112 {
113 type Output = Point<T>;
114
115 fn add(self, rhs: Point<T>) -> Point<T> {
116 Point {
117 x: self.x + rhs.x,
118 y: self.y + rhs.y,
119 }
120 }
121 }
122
123 impl<T> Point<T>
124 where
125 T: Copy + std::ops::Mul<Output = T> + std::ops::Sub<Output = T> + std::ops::Add<Output = T>,
126 {
127 pub fn det(&self, other: &Point<T>) -> T {
128 self.x * other.y - self.y * other.x
129 }
130 pub fn norm(&self) -> T {
131 self.x * self.x + self.y * self.y
132 }
133 }
134}
135
136#[cfg(test)]
137mod tests {
138 use super::minimum_bounding_circle::*;
139 use crate::utils::test_helper::Tester;
140
141 #[test]
142 fn solve_aoj2423() {
143 let tester = Tester::new("./assets/AOJ2423/in/", "./assets/AOJ2423/out/");
144 tester.test_solution(|sc| {
145 let circles: usize = sc.read();
146 let people: usize = sc.read();
147 let rs: Vec<f64> = sc.vec(circles);
148 let mut people_r2 = vec![0.0; people];
149 let mut people_list = vec![];
150 for i in 0..people {
151 let n: usize = sc.read();
152 let mut ps = vec![];
153 for _ in 0..n {
154 let x: f64 = sc.read();
155 let y: f64 = sc.read();
156 ps.push(Point { x, y });
157 }
158 let (_, r2) = make_circle(&ps);
159 people_list.push((r2, i));
160 people_r2[i] = r2;
161 }
162
163 let mut circle_list = rs
164 .iter()
165 .enumerate()
166 .map(|(i, &r)| (r, i))
167 .collect::<Vec<_>>();
168 circle_list.sort_by(|a, b| a.0.partial_cmp(&b.0).unwrap());
169 people_list.sort_by(|a, b| a.0.partial_cmp(&b.0).unwrap());
170
171 let mut ans = vec![0; people];
172
173 for _ in 0..people {
174 let mut circle_buf = 0;
175 let mut people_buf = people_list.len();
176 for p_pos in (0..people_list.len()).rev() {
177 for c_pos in (0..(circle_list.len() - circle_buf)).rev() {
178 if circle_list[c_pos].0 * circle_list[c_pos].0 >= people_list[p_pos].0 {
179 circle_buf += 1;
180 } else {
181 break;
182 }
183 }
184 if people_list.len() - p_pos > circle_buf {
185 sc.write("NG\n");
186 return;
187 } else if people_list.len() - p_pos == circle_buf {
188 people_buf = people_list.len() - p_pos;
189 break;
190 }
191 }
192
193 let youngest_people = people_list
194 .iter()
195 .rev()
196 .take(people_buf)
197 .map(|p| p.1)
198 .min()
199 .unwrap();
200 let youngest_circle = circle_list
201 .iter()
202 .rev()
203 .take(circle_buf)
204 .map(|c| c.1)
205 .filter(|&i| rs[i] * rs[i] >= people_r2[youngest_people])
206 .min()
207 .unwrap();
208 ans[youngest_people] = youngest_circle;
209 people_list = people_list
210 .into_iter()
211 .filter(|p| p.1 != youngest_people)
212 .collect::<Vec<_>>();
213 circle_list = circle_list
214 .into_iter()
215 .filter(|c| c.1 != youngest_circle)
216 .collect::<Vec<_>>();
217 }
218
219 for i in ans.into_iter() {
220 sc.write(format!("{}\n", i + 1));
221 }
222 });
223 }
224}