EDIT2:

• I was wrong, and the test cases were correct.

• Added a section discussing my errors.

EDIT1:

• Changed the threshold of Taylor series from to (thanks to Swistakk)

• Fixed use of Fraction in Python: Fraction(str) instead of Fraction(float)

• No change in the outcome.

TL;DR

600D - Area of Two Circles' Intersection requires careful treatment of loss of significance, more than the geometry itself. The test is plagued by the numerical errors, and rejects "more accurate" solution. EDIT2: Solving this problem with double precision turned out very hard!

Geometry

Let's review the geometry, first. We need three numbers (r₁, r₂, d); where d is the distance between two centers,

When two circles do not intersect, which is the case d ≥ r₁ + r₂, the area of intersection is zero.

When one circle is contained in the other, which is the case max(r₁, r₂) ≥ d + min(r₁, r₂), the intersection equals to the smaller circle.

Otherwise, intersection of two circles is determined from two circular sectors subtracted by two triangles, yielding

where angles θ₁ and θ₂ are determined from the cosine formula

So far, so good. I wrote Python code (with double representation), and failed. Rewrote in C++ using long double, and succeeded. But..., what was wrong?

Loss of significance

The culprit is loss of significance (wikipedia). It occurs when you take difference (y - x) of two close numbers x ≈ y. For example,

The difference has only one significant digit although you started with numbers of five significant digits!

This type of numerical error can happen in two parts. One in inside of , the other in the function (θ - sinθ cos θ ) when θ is small .

Inside of the is just arithmetic operations so we can represent numbers as Fractions to avoid the error. (d is a floating point number but I still convert it into a fraction as well.)

(θ - sinθ cos θ ) involves non-arithmetic operation, so I use Taylor expansion. Wolfram|Alpha tells me like this.

The first three terms are good enough. I replaced the function f(x) with this form when ( EDIT1: changed the threshold to 0.01.)

Am I wrong, or test is wrong?

Here is my Python implementation.

"""
Codeforces Educational Round #2
Problem 600 D. Area of Two Circles' Intersection

@author yamaton
@date 2015-11-30
      2015-12-02
"""

import math
import fractions


def f(x):
    """
    compute  x - sin(x) cos(x) 
    without loss of significance
    """
    cos = math.cos
    sin = math.sin

    if abs(x) < 0.01:
        return 2 * x**3 / 3 - 2 * x**5 / 15 + 4 * x**7 / 315
    return x - sin(x) * cos(x)


def solve(r1, r2, d):
    acos = math.acos

    r1, r2 = min(r1, r2), max(r1, r2)
    if d >= r1 + r2:  # circles are far apart
        return 0.0     
    if r2 >= d + r1:  # whole circle is contained in the other
        return math.pi * r1 ** 2

    # Convert numbers into fractions before calculation.
    # Numbers are fractions until `acos()` is applied.
    df = fractions.Fraction("%.20f" % d)
    r1f, r2f = map(fractions.Fraction, [r1, r2])
    r1sq, r2sq, dsq = map(lambda i: i * i, [r1f, r2f, df])
    theta1 = acos((r1sq + dsq - r2sq) / (2 * r1f * df))
    theta2 = acos((r2sq + dsq - r1sq) / (2 * r2f * df))

    return r1 * r1 * f(theta1) + r2 * r2 * f(theta2)


def main():
    x1, y1, r1 = map(int, input().split())
    x2, y2, r2 = map(int, input().split())
    d = math.hypot(x1-x2, y1-y2)
    result = solve(r1, r2, d)
    print("%.20f" % result)


if __name__ == '__main__':
    main()

But this STILL fails on test 36. Here is a screenshot.

Am I wrong? Or, the test is wrong?

 →  ... It turned out I was wrong.

My Errors

I had two errors.

1. d was computed first as hypotenuse, but the double value was insufficient in significance.

2. acos returns inaccurate value when the argument is close to 1. According to Wikipedia,

For angles near 0 and π, arccosine is ill-conditioned and will thus calculate the angle with reduced accuracy in a computer implementation (due to the limited number of digits).

To fix them, I computed θ by keeping exact representation with fractions from scratch,

where σ₁ is the sign (1 or -1) of r₁² + d² - r₂². Also acos(sigma * sqrt(x)) part is computed with series when abs(x) is close to 1.

This is how I finally got AC with double precision: 14663921.

Lesson Learned:

• Just use C++/Java so we can use long double or BigDecimal