Problem:
 Score : 600 points    Problem Statement Takahashi wants to grill N pieces of meat on a grilling net, which can be seen as a two-dimensional plane. The coordinates of the i -th piece of meat are \left(x_i, y_i\right) , and its hardness is c_i .  Takahashi can use one heat source to grill the meat. If he puts the heat source at coordinates \left(X, Y\right) , where X and Y are real numbers, the i -th piece of meat will be ready to eat in c_i \times \sqrt{\left(X - x_i\right)^2 + \left(Y-y_i\right)^2} seconds.  Takahashi wants to eat K pieces of meat. Find the time required to have K or more pieces of meat ready if he put the heat source to minimize this time.      Constraints  All values in input are integers.  1 \leq N \leq 60  1 \leq K \leq N  -1000 \leq x_i , y_i \leq 1000  \left(x_i, y_i\right) \neq \left(x_j, y_j\right) \left(i \neq j \right)  1 \leq c_i \leq 100         Input Input is given from Standard Input in the following format:  N  K  x_1  y_1  c_1  \vdots  x_N  y_N  c_N       Output Print the answer.  It will be considered correct if its absolute or relative error from our answer is at most 10^{-6} .        Sample Input 1 4 3
-1 0 3
0 0 3
1 0 2
1 1 40      Sample Output 1 2.4  If we put the heat source at \left(-0.2, 0\right) , the 1 -st, 2 -nd, and 3 -rd pieces of meat will be ready to eat within 2.4 seconds. This is the optimal place to put the heat source.       Sample Input 2 10 5
-879 981 26
890 -406 81
512 859 97
362 -955 25
128 553 17
-885 763 2
449 310 57
-656 -204 11
-270 76 40
184 170 16      Sample Output 2 7411.2252

Solution (Python):
import math
EPS = 10**-7


class P2():
    def __init__(self, x, y):
        self.x = x
        self.y = y
        self.norm2 = (self.x**2 + self.y**2)**0.5

    def __add__(self, other):
        return P2(self.x + other.x, self.y + other.y)

    def __sub__(self, other):
        return P2(self.x - other.x, self.y - other.y)

    def smul(self, a):
        return P2(self.x * a, self.y * a)

    def dot(self, other):
        return self.x * other.x + self.y * other.y

    def det(self, other):
        return self.x * other.y - self.y * other.x


n, k = map(int, input().split())
p = []
c = []
for _ in range(n):
    x, y, cc = map(float, input().split())
    p.append(P2(x, y))
    c.append(cc)

dist = [[(p[i] - p[j]).norm2 for j in range(n)] for i in range(n)]


def check(r):
    r_c = tuple(r / c[i] for i in range(n))
    
    cand = [pi for pi in p]
    for i in range(n - 1):
        for j in range(i + 1, n):
            d_ij = p[i] - p[j]

            if dist[i][j] < r_c[i] + r_c[j] - EPS:
                
                x = ((r_c[j])**2 -
                     (r_c[i])**2 + d_ij.norm2**2) / (2 * d_ij.norm2)
                h = math.sqrt(max((r_c[j])**2 - x**2, 0.0))
                v = P2(-d_ij.y, d_ij.x).smul(h / d_ij.norm2)
                dx = p[j] + d_ij.smul(x / d_ij.norm2)
                cand.append(dx + v)
                cand.append(dx - v)

    for cand_p in cand:
        cnt = 0
        for i in range(n):
            if (p[i] - cand_p).norm2 < r_c[i] + EPS:
                cnt += 1
                if cnt >= k:
                    return True

    return False


lb = 0.0  
ub = 3000 * 100  
while ub - lb > EPS:
    mid = (ub + lb) / 2
    if check(mid):
        ub = mid
    else:
        lb = mid

print(ub)

Solution (Python):

"""
答えをにブタン
円がKこ以上重なる点が生まれる最小のtを求める
どうやって判定する？
"""

def wantnum(m):

    plis = []

    for i in range(N):
        plis.append([xyc[i][0] , xyc[i][1]])

    for i in range(N):
        for j in range(i):

            xA = xyc[i][0]
            yA = xyc[i][1]
            xB = xyc[j][0]
            yB = xyc[j][1]
            r1 = m / xyc[i][2]
            r2 = m / xyc[j][2]

            X = xB-xA
            Y = yB-yA

            a = (X**2 + Y**2 + r1**2 - r2**2) / 2

            retx1 = ( a*X + Y * (( (X**2+Y**2)*(r1**2) - a**2 )**0.5) ) / (X**2 + Y**2)
            rety1 = ( a*Y - X * (( (X**2+Y**2)*(r1**2) - a**2 )**0.5) ) / (X**2 + Y**2)

            retx1 += xA
            rety1 += yA

            retx2 = ( a*X - Y * (( (X**2+Y**2)*(r1**2) - a**2 )**0.5) ) / (X**2 + Y**2)
            rety2 = ( a*Y + X * (( (X**2+Y**2)*(r1**2) - a**2 )**0.5) ) / (X**2 + Y**2)

            retx2 += xA
            rety2 += yA

            if type(retx1) is not complex and type(rety1) is not complex:
                plis.append([retx1,rety1])
            if type(retx2) is not complex and type(rety2) is not complex:
                plis.append([retx2,rety2])

    ans = 0
    for nx,ny in plis:

        nnum = 0

        for i,j,r in xyc:
            if (i-nx)**2 + (j-ny)**2 < (m/r)**2 + 10**(-7):
                nnum += 1

        ans = max(ans,nnum)

    return ans

            

N,K = map(int,input().split())

xyc = []
l = 0
r = 0

for i in range(N):
    
    x,y,c = map(int,input().split())
    xyc.append([x,y,c])

    r = max(r , c * ((x**2+y**2)**0.5) )

while r - l > 10**(-7) and r-l > l * (10**(-7)):

    m = (r+l) / 2

    ret = wantnum(m)
    

    if ret >= K:
        r = m
    else:
        l = m

print (r)

Solution (Python):
import sys,math
range = xrange
input = raw_input

inp = [int(x) for x in sys.stdin.read().split()]; ii = 0

n = inp[ii]; ii += 1
k = inp[ii]; ii += 1

X = inp[ii + 0: ii + 3 * n: 3]
Y = inp[ii + 1: ii + 3 * n: 3]
C = inp[ii + 2: ii + 3 * n: 3]
ii += 3 * n

def dist(i,j):
    return math.sqrt((X[i] - X[j])**2 + (Y[i] - Y[j])**2)

l = 0.0
r = 100.0 * 1e18
while r - l >= 1e-9:
    t = (l + r) / 2

    R = [t/C[i] for i in range(n)]

    cand = [(X[i], Y[i]) for i in range(n)]

    for i in range(n):
        r0 = R[i]
        x0 = X[i]
        y0 = Y[i]
        for j in range(i):
            r1 = R[j]
            x1 = X[j]
            y1 = Y[j]

            d = dist(i,j)
            if d == 0.0:
                continue
            a = (r0**2 - r1**2 + d**2)/(2*d)
            if 0.0 < d <= r0 + r1 and r0**2 - a**2 >= 0.0:
                h = math.sqrt(r0**2 - a**2)

                x2 = x0 + a * (x1 - x0)/d
                y2 = y0 + a * (y1 - y0)/d
                

                x3 = x2 + h * (y1 - y0)/d
                y3 = y2 - h * (x1 - x0)/d
                cand.append((x3,y3))
                
                x3 = x2 + h * (y1 - y0)/d
                y3 = y2 - h * (x1 - x0)/d
                cand.append((x3,y3))

    for x,y in cand:
        count = 0
        for i in range(n):
            if math.sqrt((X[i] - x)**2 + (Y[i] - y)**2) <= R[i] + 1e-9:
                count += 1

        if count >= k:
            r = t
            break
    else:
        l = t
print (l + r)/2

Solution (Python):
import cmath
import itertools
import math
import os
import random
import sys

INF = float("inf")
PI = cmath.pi
TAU = cmath.pi * 2
EPS = 1e-10

if os.getenv("LOCAL"):
    sys.stdin = open("_in.txt", "r")

sys.setrecursionlimit(10 ** 9)
INF = float("inf")
IINF = 10 ** 18
MOD = 10 ** 9 + 7




class Point:
    """
    2次元空間上の点
    """

    
    CCW_COUNTER_CLOCKWISE = 1
    
    CCW_CLOCKWISE = -1
    
    CCW_ONLINE_BACK = 2
    
    CCW_ONLINE_FRONT = -2
    
    CCW_ON_SEGMENT = 0

    def __init__(self, x: float, y: float):
        self.c = complex(x, y)

    @property
    def x(self):
        return self.c.real

    @property
    def y(self):
        return self.c.imag

    @staticmethod
    def from_complex(c: complex):
        return Point(c.real, c.imag)

    @staticmethod
    def from_polar(r: float, phi: float):
        c = cmath.rect(r, phi)
        return Point(c.real, c.imag)

    def __add__(self, p):
        """
        :param Point p:
        """
        c = self.c + p.c
        return Point(c.real, c.imag)

    def __iadd__(self, p):
        """
        :param Point p:
        """
        self.c += p.c
        return self

    def __sub__(self, p):
        """
        :param Point p:
        """
        c = self.c - p.c
        return Point(c.real, c.imag)

    def __isub__(self, p):
        """
        :param Point p:
        """
        self.c -= p.c
        return self

    def __mul__(self, f: float):
        c = self.c * f
        return Point(c.real, c.imag)

    def __imul__(self, f: float):
        self.c *= f
        return self

    def __truediv__(self, f: float):
        c = self.c / f
        return Point(c.real, c.imag)

    def __itruediv__(self, f: float):
        self.c /= f
        return self

    def __repr__(self):
        return "({}, {})".format(round(self.x, 10), round(self.y, 10))

    def __neg__(self):
        c = -self.c
        return Point(c.real, c.imag)

    def __eq__(self, p):
        return abs(self.c - p.c) < EPS

    def __abs__(self):
        return abs(self.c)

    @staticmethod
    def ccw(a, b, c):
        """
        線分 ab に対する c の位置
        線分上にあるか判定するだけなら on_segment とかのが速い
        Verify: http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_1_C&lang=ja
        :param Point a:
        :param Point b:
        :param Point c:
        """
        b = b - a
        c = c - a
        det = b.det(c)
        if det > EPS:
            return Point.CCW_COUNTER_CLOCKWISE
        if det < -EPS:
            return Point.CCW_CLOCKWISE
        if b.dot(c) < -EPS:
            return Point.CCW_ONLINE_BACK
        if b.dot(b - c) < -EPS:
            return Point.CCW_ONLINE_FRONT
        return Point.CCW_ON_SEGMENT

    def dot(self, p):
        """
        内積
        :param Point p:
        :rtype: float
        """
        return self.x * p.x + self.y * p.y

    def det(self, p):
        """
        外積
        :param Point p:
        :rtype: float
        """
        return self.x * p.y - self.y * p.x

    def dist(self, p):
        """
        距離
        :param Point p:
        :rtype: float
        """
        return abs(self.c - p.c)

    def norm(self):
        """
        原点からの距離
        :rtype: float
        """
        return abs(self.c)

    def phase(self):
        """
        原点からの角度
        :rtype: float
        """
        return cmath.phase(self.c)

    def angle(self, p, q):
        """
        p に向いてる状態から q まで反時計回りに回転するときの角度
        -pi <= ret <= pi
        :param Point p:
        :param Point q:
        :rtype: float
        """
        return (cmath.phase(q.c - self.c) - cmath.phase(p.c - self.c) + PI) % TAU - PI

    def area(self, p, q):
        """
        p, q となす三角形の面積
        :param Point p:
        :param Point q:
        :rtype: float
        """
        return abs((p - self).det(q - self) / 2)

    def projection_point(self, p, q, allow_outer=False):
        """
        線分 pq を通る直線上に垂線をおろしたときの足の座標
        Verify: http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_1_A&lang=ja
        :param Point p:
        :param Point q:
        :param allow_outer: 答えが線分の間になくても OK
        :rtype: Point|None
        """
        diff_q = q - p
        
        r = (self - p).dot(diff_q) / abs(diff_q)
        
        phase = diff_q.phase()

        ret = Point.from_polar(r, phase) + p
        if allow_outer or (p - ret).dot(q - ret) < EPS:
            return ret
        return None

    def reflection_point(self, p, q):
        """
        直線 pq を挟んで反対にある点
        Verify: http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_1_B&lang=ja
        :param Point p:
        :param Point q:
        :rtype: Point
        """
        
        r = abs(self - p)
        
        angle = p.angle(q, self)
        
        angle = (q - p).phase() - angle
        return Point.from_polar(r, angle) + p

    def on_segment(self, p, q, allow_side=True):
        """
        点が線分 pq の上に乗っているか
        Verify: http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_1_C&lang=ja
        :param Point p:
        :param Point q:
        :param allow_side: 端っこでギリギリ触れているのを許容するか
        :rtype: bool
        """
        if not allow_side and (self == p or self == q):
            return False
        
        
        return abs((p - self).det(q - self)) < EPS and (p - self).dot(q - self) < EPS

    @staticmethod
    def circumcenter_of(p1, p2, p3):
        """
        外心
        :param Point p1:
        :param Point p2:
        :param Point p3:
        :rtype: Point|None
        """
        if abs((p2 - p1).det(p3 - p1)) < EPS:
            
            return None
        
        a = (p2.x - p3.x) ** 2 + (p2.y - p3.y) ** 2
        b = (p3.x - p1.x) ** 2 + (p3.y - p1.y) ** 2
        c = (p1.x - p2.x) ** 2 + (p1.y - p2.y) ** 2
        num = p1 * a * (b + c - a) + p2 * b * (c + a - b) + p3 * c * (a + b - c)
        den = a * (b + c - a) + b * (c + a - b) + c * (a + b - c)
        return num / den

    @staticmethod
    def incenter_of(p1, p2, p3):
        """
        内心
        :param Point p1:
        :param Point p2:
        :param Point p3:
        """
        
        d1 = p2.dist(p3)
        d2 = p3.dist(p1)
        d3 = p1.dist(p2)
        return (p1 * d1 + p2 * d2 + p3 * d3) / (d1 + d2 + d3)


class Circle:
    def __init__(self, o, r):
        """
        :param Point o:
        :param float r:
        """
        self.o = o
        self.r = r

    def __eq__(self, other):
        return self.o == other.o and abs(self.r - other.r) < EPS

    def ctc(self, c):
        """
        共通接線 common tangent の数
        4: 離れてる
        3: 外接
        2: 交わってる
        1: 内接
        0: 内包
        inf: 同一
        Verify: http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_7_A&lang=ja
        :param Circle c:
        :rtype: int
        """
        if self.o == c.o:
            return INF if abs(self.r - c.r) < EPS else 0
        
        d = self.o.dist(c.o) - self.r - c.r
        if d > EPS:
            return 4
        elif d > -EPS:
            return 3
        
        elif d + min(self.r, c.r) * 2 > EPS:
            return 2
        elif d + min(self.r, c.r) * 2 > -EPS:
            return 1
        return 0

    def has_point_on_edge(self, p):
        """
        指定した点が円周上にあるか
        :param Point p:
        :rtype: bool
        """
        return abs(self.o.dist(p) - self.r) < EPS

    def contains(self, p, allow_on_edge=True):
        """
        指定した点を含むか
        :param Point p:
        :param bool allow_on_edge: 辺上の点を許容するか
        """
        if allow_on_edge:
            
            return self.o.dist(p) - self.r < EPS
        else:
            
            return self.o.dist(p) - self.r < -EPS

    def area(self):
        """
        面積
        """
        return self.r ** 2 * PI

    def circular_segment_area(self, angle):
        """
        弓形⌓の面積
        :param float angle: 角度ラジアン
        """
        
        sector_area = self.area() * angle / TAU
        
        return sector_area - self.r ** 2 * math.sin(angle) / 2

    def intersection_points(self, other, allow_outer=False):
        """
        :param Segment|Circle other:
        :param bool allow_outer:
        """
        if isinstance(other, Segment):
            return self.intersection_points_with_segment(other, allow_outer=allow_outer)
        if isinstance(other, Circle):
            return self.intersection_points_with_circle(other)
        raise NotImplementedError()

    def intersection_points_with_segment(self, s, allow_outer=False):
        """
        線分と交差する点のリスト
        Verify: http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_7_D&lang=ja
        :param Segment s:
        :param bool allow_outer: 線分の間にない点を含む
        :rtype: list of Point
        """
        
        projection_point = self.o.projection_point(s.p1, s.p2, allow_outer=True)
        
        dist = self.o.dist(projection_point)
        
        if dist - self.r > EPS:
            return []
        if dist - self.r > -EPS:
            if allow_outer or s.has_point(projection_point):
                return [projection_point]
            else:
                return []
        
        diff = Point.from_polar(math.sqrt(self.r ** 2 - dist ** 2), s.phase())
        ret1 = projection_point + diff
        ret2 = projection_point - diff
        ret = []
        if allow_outer or s.has_point(ret1):
            ret.append(ret1)
        if allow_outer or s.has_point(ret2):
            ret.append(ret2)
        return ret

    def intersection_points_with_circle(self, other):
        """
        円と交差する点のリスト
        Verify: http://judge.u-aizu.ac.jp/onlinejudge/description.jsp?id=CGL_7_E&langja
        :param circle other:
        :rtype: list of Point
        """
        ctc = self.ctc(other)
        if not 1 <= ctc <= 3:
            return []
        if ctc == 3:
            
            return [Point.from_polar(self.r, (other.o - self.o).phase()) + self.o]
        if ctc == 1:
            
            if self.r > other.r:
                return [Point.from_polar(self.r, (other.o - self.o).phase()) + self.o]
            else:
                return [Point.from_polar(self.r, (self.o - other.o).phase()) + self.o]
        
        assert ctc == 2

        a = other.r
        b = self.r
        c = self.o.dist(other.o)
        
        cos_a = (b ** 2 + c ** 2 - a ** 2) / (2 * b * c)
        angle = math.acos(cos_a)
        phi = (other.o - self.o).phase()
        return [
            self.o + Point.from_polar(self.r, phi + angle),
            self.o + Point.from_polar(self.r, phi - angle),
        ]

    def tangent_points_with_point(self, p):
        """
        p を通る接線との接点
        :param Point p:
        :rtype: list of Point
        """
        dist = self.o.dist(p)
        
        if dist - self.r < -EPS:
            
            return []
        if dist - self.r < EPS:
            
            return [Point(p.x, p.y)]

        a = math.sqrt(dist ** 2 - self.r ** 2)
        b = self.r
        c = dist
        
        cos_a = (b ** 2 + c ** 2 - a ** 2) / (2 * b * c)
        angle = math.acos(cos_a)
        phi = (p - self.o).phase()
        return [
            self.o + Point.from_polar(self.r, phi + angle),
            self.o + Point.from_polar(self.r, phi - angle),
        ]

    def tangent_points_with_circle(self, other):
        """
        other との共通接線との接点
        :param Circle other:
        :rtype: list of Point
        """
        ctc = self.ctc(other)
        if ctc > 4:
            raise ValueError('2つの円が同一です')
        if ctc == 0:
            return []
        if ctc == 1:
            return self.intersection_points_with_circle(other)

        assert ctc in (2, 3, 4)
        ret = []
        
        
        if abs(self.r - other.r) < EPS:
            
            phi = (other.o - self.o).phase()
            ret.append(self.o + Point.from_polar(self.r, phi + PI / 2))
            ret.append(self.o + Point.from_polar(self.r, phi - PI / 2))
        else:
            
            intersection = self.o + (other.o - self.o) / (self.r - other.r) * self.r
            ret += self.tangent_points_with_point(intersection)

        
        
        intersection = self.o + (other.o - self.o) / (self.r + other.r) * self.r
        ret += self.tangent_points_with_point(intersection)
        return ret

    @staticmethod
    def circumscribed_of(p1, p2, p3):
        """
        p1・p2・p3 のなす三角形の外接円
        Verify:
        :param Point p1:
        :param Point p2:
        :param Point p3:
        """
        if p1.on_segment(p2, p3):
            return Circle((p2 + p3) / 2, p2.dist(p3) / 2)
        if p2.on_segment(p1, p3):
            return Circle((p1 + p3) / 2, p1.dist(p3) / 2)
        if p3.on_segment(p1, p2):
            return Circle((p1 + p2) / 2, p1.dist(p2) / 2)
        o = Point.circumcenter_of(p1, p2, p3)
        return Circle(o, o.dist(p1))

    @staticmethod
    def min_enclosing(points):
        """
        points をすべて含む最小の円
        計算量の期待値は O(N)
        https://www.jaist.ac.jp/~uehara/course/2014/i481f/pdf/ppt-7.pdf
        Verify: https://www.spoj.com/problems/QCJ4/
        :param list of Point points:
        :rtype: Circle
        """
        points = points[:]
        random.shuffle(points)

        
        ret = Circle(points[0], 0)
        for i, p in enumerate(points):
            if ret.contains(p):
                continue
            ret = Circle(p, 0)
            for j, q in enumerate(points[:i]):
                if ret.contains(q):
                    continue
                
                ret = Circle((p + q) / 2, abs(p - q) / 2)
                for k, r in enumerate(points[:j]):
                    if ret.contains(r):
                        continue
                    
                    ret = Circle.circumscribed_of(p, q, r)
        return ret


def argsort(li, key=None, reverse=False):
    return [i for _, i in sorted(
        [(a, i) for i, a in enumerate(li)], key=(lambda t: key(t[0])) if key else None, reverse=reverse
    )]


N, K = list(map(int, sys.stdin.buffer.readline().split()))
XYC = [list(map(int, sys.stdin.buffer.readline().split())) for _ in range(N)]


def test(t: float):
    
    
    circles = []
    for x, y, c in XYC:
        circles.append(Circle(Point(x, y), t /c))
    
    points = []
    for c1, c2 in itertools.combinations(circles, r=2):
        points += c1.intersection_points_with_circle(c2)
    for c in circles:
        points.append(c.o)

    for p in points:
        cnt = 0
        for c in circles:
            cnt += c.contains(p)
        if cnt >= K:
            return True
    return False


























ok = 1 << 100
ng = 0
while abs(ok - ng) > EPS:
    mid = (ok + ng) / 2
    if test(mid):
        ok = mid
    else:
        ng = mid
print(ok)