Problem p02764 (python)

Fixed submission:
from math import atan2, cos, sin

EPS = 1e-7


def hypot(a):
    return sum(x ** 2 for x in a) ** 0.5


def dot(a, b):
    return a[0] * b[0] + a[1] * b[1] + a[2] * b[2]


def add(a, b):
    return [a[0] + b[0], a[1] + b[1], a[2] + b[2]]


def sub(a, b):
    return [a[0] - b[0], a[1] - b[1], a[2] - b[2]]


def matmul(mat, v):
    return [dot(row, v) for row in mat]


def rotZ(theta):
    return [[cos(theta), -sin(theta), 0], [sin(theta), cos(theta), 0], [0, 0, 1]]


def _circleCircleIntersection(R, d, r):
    # For circles centered at (0, 0) and (d, 0) with radius R and r respectively
    # Equations are:
    #   x^2       + y^2 = R^2
    #   (x - d)^2 + y^2 = r^2
    # Subtracting and expanding:
    #   x^2 - (x^2 - 2 * d * x + d^2) = R^2 - r^2
    #                 2 * d * x - d^2 = R^2 - r^2
    #                               x = (R^2 - r^2 + d^2) / (2*d)
    if d == 0:
        return [(0, 0)]
    x = (R ** 2 - r ** 2 + d ** 2) / (2 * d)
    ySq = R ** 2 - x ** 2
    if ySq == 0:
        return [(x, 0)]
    elif ySq > 0:
        y = ySq ** 0.5
        # assert abs((x ** 2 + y ** 2) - R ** 2) < EPS
        # assert abs(((x - d) ** 2 + y ** 2) - r ** 2) < EPS
        return [(x, y), (x, -y)]
    else:
        raise AssertionError()


def circleCircleIntersection(x1, y1, r1, x2, y2, r2):
    # Either 1 point, 2 points
    # 0 points because one circle contains the other
    # 0 points because the circles don't intersect
    v1 = [x1, y1, 0]
    v2 = [x2, y2, 0]
    d = hypot(sub(v2, v1))
    if r1 + r2 < d:
        # Circles don't intersect
        return []
    if min(r1, r2) + d < max(r1, r2):
        # One circle contains the other
        return []

    # Translate
    offset = list(v1)
    v1 = sub(v1, offset)
    v2 = sub(v2, offset)
    assert v1[0] == 0 and v1[1] == 0

    # Rotate
    angle = -atan2(v2[1], v2[0])
    v2 = matmul(rotZ(angle), v2)
    assert abs(v2[1]) < EPS

    ret = []
    for pt in _circleCircleIntersection(r1, v2[0], r2):
        pt = [pt[0], pt[1], 0]
        # Un-rotate
        pt = matmul(rotZ(-angle), pt)
        # Un-translate
        pt = add(pt, offset)
        ret.append((pt[0], pt[1]))

    return ret


N, K = list(map(int, input().split()))
XYC = []
for i in range(N):
    XYC.append(list(map(int, input().split())))


def numMeat(t):
    # At time t, each x, y can be heated by anything within a circle of radius radius t/c
    # Gather all intersecting points, some of which must be part of the region with most overlaps
    intersect = []
    for i in range(len(XYC)):
        x1, y1, c1 = XYC[i]
        r1 = t / c1
        for j in range(i + 1, len(XYC)):
            x2, y2, c2 = XYC[j]
            r2 = t / c2
            pts = circleCircleIntersection(x1, y1, r1, x2, y2, r2)
            intersect.extend(pts)
        intersect.append((x1, y1))

    best = 0
    for pt in intersect:
        count = 0
        for x, y, c in XYC:
            if c * ((x - pt[0]) ** 2 + (y - pt[1]) ** 2) ** 0.5 <= t + EPS:
                count += 1
        best = max(best, count)

    return best


lo = 0
# Get upperbound (by arbitrarily placing heat source at origin)
hi = max(c * (x ** 2 + y ** 2) ** 0.5 for x, y, c in XYC)
while hi - lo > EPS:
    mid = (hi + lo) / 2
    if numMeat(mid) >= K:
        hi = mid
    else:
        lo = mid
print(hi)

Buggy submission (verdict: Runtime Error):
from math import atan2, cos, sin

EPS = 1e-7


def hypot(a):
    return sum(x ** 2 for x in a) ** 0.5


def dot(a, b):
    return a[0] * b[0] + a[1] * b[1] + a[2] * b[2]


def add(a, b):
    return [a[0] + b[0], a[1] + b[1], a[2] + b[2]]


def sub(a, b):
    return [a[0] - b[0], a[1] - b[1], a[2] - b[2]]


def matmul(mat, v):
    return [dot(row, v) for row in mat]


def rotZ(theta):
    return [[cos(theta), -sin(theta), 0], [sin(theta), cos(theta), 0], [0, 0, 1]]


def _circleCircleIntersection(R, d, r):
    # For circles centered at (0, 0) and (d, 0) with radius R and r respectively
    # Equations are:
    #   x^2       + y^2 = R^2
    #   (x - d)^2 + y^2 = r^2
    # Subtracting and expanding:
    #   x^2 - (x^2 - 2 * d * x + d^2) = R^2 - r^2
    #                 2 * d * x - d^2 = R^2 - r^2
    #                               x = (R^2 - r^2 + d^2) / (2*d)
    x = (R ** 2 - r ** 2 + d ** 2) / (2 * d)
    ySq = R ** 2 - x ** 2
    if ySq == 0:
        return [(x, 0)]
    elif ySq > 0:
        y = ySq ** 0.5
        assert abs((x ** 2 + y ** 2) - R ** 2) < EPS
        assert abs(((x - d) ** 2 + y ** 2) - r ** 2) < EPS
        return [(x, y), (x, -y)]
    else:
        raise AssertionError()


def circleCircleIntersection(x1, y1, r1, x2, y2, r2):
    # Either 1 point, 2 points
    # 0 points because one circle contains the other
    # 0 points because the circles don't intersect
    v1 = [x1, y1, 0]
    v2 = [x2, y2, 0]
    d = hypot(sub(v2, v1))
    if r1 + r2 < d:
        # Circles don't intersect
        return []
    if min(r1, r2) + d < max(r1, r2):
        # One circle contains the other
        return []

    # Translate
    offset = list(v1)
    v1 = sub(v1, offset)
    v2 = sub(v2, offset)
    assert v1[0] == 0 and v1[1] == 0

    # Rotate
    angle = -atan2(v2[1], v2[0])
    v2 = matmul(rotZ(angle), v2)
    assert abs(v2[1]) < EPS

    ret = []
    for pt in _circleCircleIntersection(r1, v2[0], r2):
        pt = [pt[0], pt[1], 0]
        # Un-rotate
        pt = matmul(rotZ(-angle), pt)
        # Un-translate
        pt = add(pt, offset)
        ret.append((pt[0], pt[1]))

    return ret


N, K = list(map(int, input().split()))
XYC = []
for i in range(N):
    XYC.append(list(map(int, input().split())))


def numMeat(t):
    # At time t, each x, y can be heated by anything within a circle of radius radius t/c
    # Gather all intersecting points, some of which must be part of the region with most overlaps
    intersect = []
    for i in range(len(XYC)):
        x1, y1, c1 = XYC[i]
        r1 = t / c1
        for j in range(i + 1, len(XYC)):
            x2, y2, c2 = XYC[j]
            r2 = t / c2
            pts = circleCircleIntersection(x1, y1, r1, x2, y2, r2)
            intersect.extend(pts)
        intersect.append((x1, y1))

    if False:
        import matplotlib.pyplot as plt

        plt.figure()
        ax = plt.gca()
        ax.cla()
        bound = t // 5
        ax.set_xlim((-bound, bound))
        ax.set_ylim((-bound, bound))

        for i in range(len(XYC)):
            x1, y1, c1 = XYC[i]
            r1 = t / c1
            ax.add_artist(plt.Circle((x1, y1), r1, color="blue", fill=False))
        for x, y in intersect:
            ax.plot([x], [y], ".", color="red")
        plt.show()

    best = 0
    for pt in intersect:
        count = 0
        for x, y, c in XYC:
            if c * ((x - pt[0]) ** 2 + (y - pt[1]) ** 2) ** 0.5 <= t + EPS:
                count += 1
        best = max(best, count)

    return best


lo = 0
# Get upperbound (by arbitrarily placing heat source at origin)
hi = max(c * (x ** 2 + y ** 2) ** 0.5 for x, y, c in XYC)
while hi - lo > EPS:
    mid = (hi + lo) / 2
    if numMeat(mid) >= K:
        hi = mid
    else:
        lo = mid
print(hi)