fly swatter

Problem

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 Rt).

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.

Input

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.

Output

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.

Limits

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

Small dataset (Test set 1 - Visible)

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

Large dataset (Test set 2 - Hidden)

1 <= N <= 100
The total number of strings will be at most 2000 (so at most 1000 in each direction).

Sample Input & Output

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

Rust Solution

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");

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(S2) 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 4S 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);
}

Having problems with this solution? Click here to submit an issue on github.