What are your chances of hitting a fly with a tennis racquet?

To start with, ignore the racquet's handle. Assume the racquet is a perfect
ring, of outer radius **R** and thickness **t** (so the inner radius of
the ring is **R**−**t**).

The ring is covered with horizontal and vertical strings. Each string is a
cylinder of radius **r**. Each string is a chord of the ring (a straight
line connecting two points of the circle). There is a gap of length **g**
between neighbouring strings. The strings are symmetric with respect to the
center of the racquet i.e. there is a pair of strings whose centers meet at
the center of the ring.

The fly is a sphere of radius **f**. Assume that the racquet is moving in
a straight line perpendicular to the plane of the ring. Assume also that the
fly's center is inside the outer radius of the racquet and is equally likely
to be anywhere within that radius. Any overlap between the fly and the
racquet (the ring or a string) counts as a hit.

One line containing an integer **N**, the number of test cases in the
input file.

The next **N** lines will each contain the numbers **f**, **R**,
**t**, **r** and **g** separated by exactly one space. Also the
numbers will have exactly 6 digits after the decimal point.

**N** lines, each of the form "Case #**k**: **P**", where **k**
is the number of the test case and **P** is the probability of hitting the
fly with a piece of the racquet.

Answers with a relative or absolute error of at most 10^{-6} will be
considered correct.

Time limit: 60 seconds per test set.

Memory limit: 1GB.

**f**, **R**, **t**, **r** and **g** will be positive and smaller or equal to 10000.

**t** < **R**

**f** < **R**

**r** < **R**

1 <= **N** <= 30

The total number of strings will be at most 60 (so at most 30 in each direction).

1 <= **N** <= 100

The total number of strings will be at most 2000 (so at most 1000 in each direction).

5 0.250000 1.000000 0.100000 0.010000 0.500000 0.250000 1.000000 0.100000 0.010000 0.900000 0.000010 10000.000000 0.000010 0.000010 1000.000000 0.400000 10000.000000 0.000010 0.000010 700.000000 1.000000 100.000000 1.000000 1.000000 10.000000

Case #1: 1.000000 Case #2: 0.910015 Case #3: 0.000000 Case #4: 0.002371 Case #5: 0.573972

This was the hardest problem of the set, as we can see from the submission
statistics, but 652 people succeeded in solving it correctly.

The first step in solving this problem involves a standard trick. Instead of
looking for the region that will hit the fly's circular body, we look for the
region where the center of the fly must avoid. The problem is thus
transformed into an equivalent one where we thicken the size of the racquet
to **t** + **f** and the radius of the strings from **r** to
**r** + **f**. And we consider the fly to be simply a single point,
only focusing on its center.

Now the main problem is to find the area of the holes (the areas between the
strings, or between the strings and the racquet), and divide it by the area
of the disc corresponding to the original racquet. It is easy to see this
gives the probability that the racquet *does not* hit the fly.

A simplifying observation is that we can just look at the first quadrant of
the circle, because the racquet is symmetrical.

Now we can go through each square gap in the racquet and check a few cases.
Each case consists in computing the area of intersection between the square
and the circle representing the inside of the racquet. The case where the
square is totally inside or totally outside the circle is easy. The cases
left when one corner, two or three corners of the square are inside the
circle and the other corners are outside. There are no other cases because we
just look at the intersection of a quarter of the circle with small squares.

Each of these cases can be solved by splitting the intersection area into
triangles and a circular segment.
Here's some code that does this:

double circle_segment(double rad, double th) {
return rad*rad*(th - sin(th))/2;
}
double rad = R-t-f;
double ar = 0.0;
for (double x1 = r+f; x1 < R-t-f; x1 += g+2*r)
for (double y1 = r+f; y1 < R-t-f; y1 += g+2*r) {
double x2 = x1 + g - 2*f;
double y2 = y1 + g - 2*f;
if (x2 <= x1 || y2 <= y1) continue;
if (x1*x1 + y1*y1 >= rad*rad) continue;
if (x2*x2 + y2*y2 <= rad*rad) {
// All points are inside circle.
ar += (x2-x1)*(y2-y1);
} else if (x1*x1 + y2*y2 >= rad*rad &&
x2*x2 + y1*y1 >= rad*rad) {
// Only (x1,y1) inside circle.
ar += circle_segment(rad, acos(x1/rad) - asin(y1/rad)) +
(sqrt(rad*rad - x1*x1)-y1) *
(sqrt(rad*rad - y1*y1)-x1) / 2;
} else if (x1*x1 + y2*y2 >= rad*rad) {
// (x1,y1) and (x2,y1) inside circle.
ar += circle_segment(rad, acos(x1/rad) - acos(x2/rad)) +
(x2-x1) * (sqrt(rad*rad - x1*x1)-y1 +
sqrt(rad*rad - x2*x2)-y1) / 2;
} else if (x2*x2 + y1*y1 >= rad*rad) {
// (x1,y1) and (x1,y2) inside circle.
ar += circle_segment(rad, asin(y2/rad) - asin(y1/rad)) +
(y2-y1) * (sqrt(rad*rad - y1*y1)-x1 +
sqrt(rad*rad - y2*y2)-x1) / 2;
} else {
// All except (x2,y2) inside circle.
ar += circle_segment(rad, asin(y2/rad) - acos(x2/rad)) +
(x2-x1)*(y2-y1) -
(y2-sqrt(rad*rad - x2*x2)) *
(x2-sqrt(rad*rad - y2*y2)) / 2;
}
}
printf("Case #%d: %.6lf\n", prob++, 1.0 - ar / (PI*R*R/4));

This solution takes O(**S**^{2}) time, where **S** is the
number of vertical strings of the racquet. It's not hard to come up with an
O(**S**) solution because there are at most 4**S** border squares which
can be found efficiently, but the previous solution was fast enough.

Instead of solving the problem exactly, an iterative solution which
approximates the area to the needed precision would have also worked. One
such solution uses divide and conquer by splitting the square into four
smaller squares and then checking the simple cases where the squares are
totally inside or totally outside the square. In the cases where the circle
and square intersect just recurse if the current square is larger than some
chosen precision. An observation is that we can divide every length by the
radius of the racquet because it gets canceled in the division between the
area of the gaps in the racquet and the disc area. This observation helps
the iterative solution since we can make the number of iterations smaller.
Here's some sample code:

double intersection(double x1, double y1,
double x2, double y2) {
// the normalized radius is 1
if (x1*x1 + y1*y1 > 1) {
return 0;
}
if (x2*x2 + y2*y2 < 1) {
return (x2-x1) * (y2-y1);
}
// EPS2 = 1e-6 * 1e-6
if ((x2-x1)*(y2-y1) < EPS2) {
// this trick helps in doing 10 or 100 times less
// iterations than we would need to get the same
// precision if we just return 0;
return (x2-x1) * (y2-y1) / 2;
}
double mx = (x1 + x2) / 2;
double my = (y1 + y2) / 2;
return intersection(x1, y1, mx, my) +
intersection(mx, y1, x2, my) +
intersection(x1, my, mx, y2) +
intersection(mx, my, x2, y2);
}

```
use rustgym_util::*;
use std::fmt::Write;
use std::io::*;
fn solve(case_no: usize, reader: &mut impl BufRead, writer: &mut impl Write) {
let args: Vec<f64> = reader.parse_vec();
let f = args[0];
let r0 = args[1];
let t = args[2];
let r = args[3];
let g = args[4];
let res = probability(f, r0, t, r, g);
writeln!(writer, "Case #{}: {:.6}", case_no, res).unwrap();
}
fn probability(f: f64, r0: f64, t: f64, r: f64, g: f64) -> f64 {
let mut area = 0.0;
let r1 = r0 - t - f;
let r2 = r1 * r1;
let mut x1 = r + f;
while x1 < r1 {
let mut y1 = r + f;
while y1 < r1 {
let x2 = x1 + g - 2.0 * f;
let y2 = y1 + g - 2.0 * f;
area += square(x1, x2, y1, y2, r1, r2);
y1 += g + 2.0 * r;
}
x1 += g + 2.0 * r;
}
1.0 - area / (r0 * r0 * std::f64::consts::PI / 4.0)
}
fn square(x1: f64, x2: f64, y1: f64, y2: f64, r1: f64, r2: f64) -> f64 {
if x2 <= x1 || y2 <= y1 || x1 * x1 + y1 * y1 >= r2 {
0.0
} else {
let c = (x2 - x1) * (y2 - y1);
if x2 * x2 + y2 * y2 <= r2 {
c
} else if x1 * x1 + y2 * y2 >= r2 && x2 * x2 + y1 * y1 >= r2 {
let a = circle_segment(r1, (x1 / r1).acos() - (y1 / r1).asin());
let b = ((r2 - x1 * x1).sqrt() - y1) * ((r2 - y1 * y1).sqrt() - x1) / 2.0;
a + b
} else if x1 * x1 + y2 * y2 >= r2 {
let a = circle_segment(r1, (x1 / r1).acos() - (x2 / r1).acos());
let b = ((r2 - x1 * x1).sqrt() - y1 + (r2 - x2 * x2).sqrt() - y1) * (x2 - x1) / 2.0;
a + b
} else if x2 * x2 + y1 * y1 >= r2 {
let a = circle_segment(r1, (y2 / r1).asin() - (y1 / r1).asin());
let b = ((r2 - y1 * y1).sqrt() - x1 + (r2 - y2 * y2).sqrt() - x1) * (y2 - y1) / 2.0;
a + b
} else {
let a = circle_segment(r1, (y2 / r1).asin() - (x2 / r1).acos());
let b = (y2 - (r2 - x2 * x2).sqrt()) * (x2 - (r2 - y2 * y2).sqrt()) / 2.0;
a - b + c
}
}
}
fn circle_segment(radius: f64, theta: f64) -> f64 {
radius * radius * (theta - theta.sin()) / 2.0
}
google_test_gen!(test, "input.txt", "output.txt");
```