Problem:
 Score : 600 points    Problem Statement Given are N points (x_i, y_i) in a two-dimensional plane.  Find the minimum radius of a circle such that all the points are inside or on it.      Constraints  2 \leq N \leq 50  0 \leq x_i \leq 1000  0 \leq y_i \leq 1000  The given N points are all different.  The values in input are all integers.         Input Input is given from Standard Input in the following format:  N  x_1  y_1  :  x_N  y_N       Output Print the minimum radius of a circle such that all the N points are inside or on it.  Your output will be considered correct if the absolute or relative error from our answer is at most 10^{-6} .        Sample Input 1 2
0 0
1 0      Sample Output 1 0.500000000000000000  Both points are contained in the circle centered at (0.5,0) with a radius of 0.5 .       Sample Input 2 3
0 0
0 1
1 0      Sample Output 2 0.707106781186497524       Sample Input 3 10
10 9
5 9
2 0
0 0
2 7
3 3
2 5
10 0
3 7
1 9      Sample Output 3 6.726812023536805158  If the absolute or relative error from our answer is at most 10^{-6} , the output will be considered correct.

Solution (Python):
import math, random










def make_circle(points):
	
	shuffled = [(float(x), float(y)) for (x, y) in points]
	random.shuffle(shuffled)
	
	
	c = None
	for (i, p) in enumerate(shuffled):
		if c is None or not is_in_circle(c, p):
			c = _make_circle_one_point(shuffled[ : i + 1], p)
	return c



def _make_circle_one_point(points, p):
	c = (p[0], p[1], 0.0)
	for (i, q) in enumerate(points):
		if not is_in_circle(c, q):
			if c[2] == 0.0:
				c = make_diameter(p, q)
			else:
				c = _make_circle_two_points(points[ : i + 1], p, q)
	return c



def _make_circle_two_points(points, p, q):
	circ = make_diameter(p, q)
	left  = None
	right = None
	px, py = p
	qx, qy = q
	
	
	for r in points:
		if is_in_circle(circ, r):
			continue
		
		
		cross = _cross_product(px, py, qx, qy, r[0], r[1])
		c = make_circumcircle(p, q, r)
		if c is None:
			continue
		elif cross > 0.0 and (left is None or _cross_product(px, py, qx, qy, c[0], c[1]) > _cross_product(px, py, qx, qy, left[0], left[1])):
			left = c
		elif cross < 0.0 and (right is None or _cross_product(px, py, qx, qy, c[0], c[1]) < _cross_product(px, py, qx, qy, right[0], right[1])):
			right = c
	
	
	if left is None and right is None:
		return circ
	elif left is None:
		return right
	elif right is None:
		return left
	else:
		return left if (left[2] <= right[2]) else right


def make_diameter(a, b):
	cx = (a[0] + b[0]) / 2.0
	cy = (a[1] + b[1]) / 2.0
	r0 = math.hypot(cx - a[0], cy - a[1])
	r1 = math.hypot(cx - b[0], cy - b[1])
	return (cx, cy, max(r0, r1))


def make_circumcircle(a, b, c):
	
	ox = (min(a[0], b[0], c[0]) + max(a[0], b[0], c[0])) / 2.0
	oy = (min(a[1], b[1], c[1]) + max(a[1], b[1], c[1])) / 2.0
	ax = a[0] - ox;  ay = a[1] - oy
	bx = b[0] - ox;  by = b[1] - oy
	cx = c[0] - ox;  cy = c[1] - oy
	d = (ax * (by - cy) + bx * (cy - ay) + cx * (ay - by)) * 2.0
	if d == 0.0:
		return None
	x = ox + ((ax*ax + ay*ay) * (by - cy) + (bx*bx + by*by) * (cy - ay) + (cx*cx + cy*cy) * (ay - by)) / d
	y = oy + ((ax*ax + ay*ay) * (cx - bx) + (bx*bx + by*by) * (ax - cx) + (cx*cx + cy*cy) * (bx - ax)) / d
	ra = math.hypot(x - a[0], y - a[1])
	rb = math.hypot(x - b[0], y - b[1])
	rc = math.hypot(x - c[0], y - c[1])
	return (x, y, max(ra, rb, rc))


_MULTIPLICATIVE_EPSILON = 1 + 1e-14

def is_in_circle(c, p):
	return c is not None and math.hypot(p[0] - c[0], p[1] - c[1]) <= c[2] * _MULTIPLICATIVE_EPSILON



def _cross_product(x0, y0, x1, y1, x2, y2):
	return (x1 - x0) * (y2 - y0) - (y1 - y0) * (x2 - x0)

n = int(input())
info = [list(map(int, input().split())) for i in range(n)]
print(make_circle(info)[2])

Solution (Python):
import math

def get_circle_center_and_radius(x1, y1, x2, y2, x3, y3):
    """
    3点を通る円の中心と半径を取得
    """
    d = 2 * ((y1 - y3) * (x1 - x2) - (y1 - y2) * (x1 - x3))
    if d == 0:
      return (0, 0), 0
    x = ((y1 - y3) * (y1 ** 2 - y2 ** 2 + x1 ** 2 - x2 ** 2) - (y1 - y2) * (y1 ** 2 - y3 ** 2 + x1 ** 2 - x3 ** 2)) / d
    y = ((x1 - x3) * (x1 ** 2 - x2 ** 2 + y1 ** 2 - y2 ** 2) - (x1 - x2) * (x1 ** 2 - x3 ** 2 + y1 ** 2 - y3 ** 2)) / -d
    r = math.sqrt((x - x1) ** 2 + (y - y1) ** 2)
    return (x, y), r
  
N = int(input())
XY = [list(map(int, input().split())) for _ in range(N)]
answer = []
if N == 2:
  xy1, xy2 = XY
  x1, y1 = xy1
  x2, y2 = xy2
  answer.append(((x2 - x1)**2 + (y2 - y1)**2)**0.5 / 2)
else:
  for i in range(N):
    xi, yi = XY[i]
    for j in range(i + 1, N):
      xj, yj = XY[j]      
      r = ((xi - xj)**2 + (yi - yj)**2)**0.5/2
      x0, y0 = ((xi + xj) / 2, (yi + yj) / 2)
      flag = 1
      for ii in range(N):
        xii, yii = XY[ii]
        if ((x0 - xii)**2 + (y0 - yii)**2)**0.5 > r:
          flag = 0
          break
      if flag:
        answer.append(r)

  for i in range(N):
    xi, yi = XY[i]
    for j in range(i + 1, N):
      xj, yj = XY[j]
      for k in range(j + 1, N):
        xk, yk = XY[k]
        xy0, r = get_circle_center_and_radius(xi, yi, xj, yj, xk, yk)
        x0, y0 = xy0
        
        if r:
          flag = 1
          for ii in range(N):
            xii, yii = XY[ii]
            if ((x0 - xii)**2 + (y0 - yii)**2)**0.5 > r:
              flag = 0
              break
          if flag:
            answer.append(r)
print(min(answer))

Solution (Python):
import numpy as np
import math

class GradientDescent:
    def __init__(self, f, df, eta=0.01, eps=1e-8):
        self.f, self.df = f, df
        self.eta = eta      
        self.eps = eps      

    def solve(self, x, count_max):
        path = []
        grad = self.df(x)
        count = 0
        while (grad**2).sum() > self.eps**2 and count < count_max:
            x = x - self.eta * grad
            grad = self.df(x)
            path.append(x)
            count += 1
        self.x = x          
        self.opt = f(x)     
        
def f(xy):
    x = xy[0]
    y = xy[1]
    ret = 0.
    for i in range(n):
        dx, dy = x - p[i][0], y - p[i][1]
        ret = max(ret, dx**2 + dy**2)
    return math.sqrt(ret)

def df(xy):
    mu = 1e-9
    x, y = xy
    return np.array([
        (f([x + mu, y]) - f([x - mu, y])) / (2. * mu),
        (f([x, y + mu]) - f([x, y - mu])) / (2. * mu)])

n = int(input())
p = []
for i in range(n):
    xy = list(map(float, input().split()))
    p.append(xy)

xy = np.array([500., 500.])
gd = GradientDescent(f, df)
for e_eta in range(2, -9, -1):
    gd.eta = 10 ** e_eta
    gd.solve(xy, 200)
    xy = gd.x

print(gd.opt)

Solution (Python):
import numpy as np
import math

class GradientDescent:
    def __init__(self, f, df, eta=0.01, eps=1e-8):
        self.f = f    
        self.df = df    
        self.eta = eta    
        self.eps = eps    
        self.path = None    

    def solve(self, x, count_max):
        path = []
        grad = self.df(x)    
        path.append(x)
        count = 0
        
        
        while (grad**2).sum() > self.eps**2 and count < count_max:
            x = x - self.eta * grad
            grad = self.df(x)
            
            path.append(x)
            count += 1
        self.path_ = np.array(path)   
        self.x_ = x    
        self.opt_ = f(x)    
        self.grad_ = grad    
        
def f(xy):
    x = xy[0]
    y = xy[1]
    ret = 0.
    for i in range(n):
        dx, dy = x - p[i][0], y - p[i][1]
        ret = max(ret, dx**2 + dy**2)
    return math.sqrt(ret)

def df(xy):
    mu = 1e-9
    x = xy[0]
    y = xy[1]
    return np.array([(f([x + mu, y]) - f([x - mu, y])) / (2. * mu), (f([x, y + mu]) - f([x, y - mu])) / (2. * mu)])

n = int(input())
p = []
for i in range(n):
    xy = list(map(float, input().split()))
    p.append(xy)

xy = np.array([500., 500.])
for e_eta in range(2, -9, -1):
    gd = GradientDescent(f, df, eta=10**e_eta)
    initial = np.array(xy)
    gd.solve(initial, 200)
    xy = gd.x_
    

print(gd.opt_)

Solution (Python):

from collections import defaultdict,deque
from heapq import heappush, heappop
from bisect import bisect_left, bisect_right
import sys, random, itertools, math
sys.setrecursionlimit(10**5)
input = sys.stdin.readline
sqrt = math.sqrt
def LI(): return list(map(int, input().split()))
def LF(): return list(map(float, input().split()))
def LI_(): return list(map(lambda x: int(x)-1, input().split()))
def II(): return int(input())
def IF(): return float(input())
def LS(): return list(map(list, input().split()))
def S(): return list(input().rstrip())
def IR(n): return [II() for _ in range(n)]
def LIR(n): return [LI() for _ in range(n)]
def FR(n): return [IF() for _ in range(n)]
def LFR(n): return [LI() for _ in range(n)]
def LIR_(n): return [LI_() for _ in range(n)]
def SR(n): return [S() for _ in range(n)]
def LSR(n): return [LS() for _ in range(n)]
mod = 1000000007
inf = 1e20


def solve():
    def fx(z):
        for i in range(n):
            if i in v:
                continue
            xi, yi = xy[i]
            x = A - xi
            y = B - yi
            if x ** 2 + y ** 2 > z ** 2:
                return False
        return True
    def ff(x1, x2, x3, y1, y2, y3):
        d = 2 * ((y1 - y3) * (x1 - x2) - (y1 - y2) * (x1 - x3))
        x = ((y1 - y3) * (y1 ** 2 - y2 ** 2 + x1 ** 2 - x2 ** 2) - (y1 - y2) * (y1 ** 2 - y3 ** 2 + x1 ** 2 - x3 ** 2)) / d
        y = ((x1 - x3) * (x1 ** 2 - x2 ** 2 + y1 ** 2 - y2 ** 2) - (x1 - x2) * (x1 ** 2 - x3 ** 2 + y1 ** 2 - y3 ** 2)) / -d
        r = math.sqrt((x - x1) ** 2 + (y - y1) ** 2)
        return x, y, r
    n = II()
    x = 0
    y = 0
    xy = LIR(n)
    ans = inf
    for i in range(n):
        for j in range(i):
            x1, y1 = xy[i]
            x2, y2 = xy[j]
            x = (x2 + x1) / 2
            y = (y2 + y1) / 2
            r = sqrt((x2 - x1) ** 2 + (y2 - y1) ** 2) / 2
            for k in range(n):
                if i != k and j != k:
                    x3, y3 = xy[k]
                    if (x3 - x) ** 2 + (y3 - y) ** 2 > r ** 2:
                        break
            else:
                ans = min(ans, r)
    
    full = itertools.combinations(range(n), 3)
    for f in full:
        a, b, c = f
        
        axy = xy[a]
        x1, y1 = axy
        bxy = xy[b]
        x2, y2 = bxy
        cxy = xy[c]
        x3, y3 = cxy
        v = [a,b,c]
        a = x1 - x2
        b = y1 - y2
        c = x2 - x3
        d = y2 - y3
        try:
            A = (d * (x1 ** 2 + y1 ** 2 - x2 ** 2 - y2 ** 2) - b * (x2 ** 2 + y2 ** 2 - x3 ** 2 - y3 ** 2)) / (2 * (a * d - b * c))
            B = (-c * (x1 ** 2 + y1 ** 2 - x2 ** 2 - y2 ** 2) + a * (x2 ** 2 + y2 ** 2 - x3 ** 2 - y3 ** 2)) / (2 * (a * d - b * c))
            r = sqrt((x1 - A) ** 2 + (y1 - B) ** 2)
            if fx(r):
                ans = min(ans, r)
        except:
            continue
    
    print(ans)
    return



if __name__ == '__main__':
    solve()

Solution (Python):
n = int(input())
p = []

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

import math
import itertools

if n == 2:
    print(math.sqrt((p[0][0]-p[1][0])**2 + (p[0][1]-p[1][1])**2)/2)
    exit()

def rad(v):
    p1 = v[0]
    p2 = v[1]
    p3 = v[2]
    x1 = p1[0]
    y1 = p1[1]
    x2 = p2[0]
    y2 = p2[1]
    x3 = p3[0]
    y3 = p3[1]
    
    px = (y1-y3)*(y1**2-y2**2+x1**2-x2**2)-(y1-y2)*(y1**2-y3**2+x1**2-x3**2)
    if (2*(y1-y3)*(x1-x2)-2*(y1-y2)*(x1-x3)) == 0:
        return None
    px = px / (2*(y1-y3)*(x1-x2)-2*(y1-y2)*(x1-x3))
    
    py = (x1-x3)*(x1**2-x2**2+y1**2-y2**2)-(x1-x2)*(x1**2-x3**2+y1**2-y3**2)
    if (2*(x1-x3)*(y1-y2)-2*(x1-x2)*(y1-y3)) == 0:
        return None
    py = py / (2*(x1-x3)*(y1-y2)-2*(x1-x2)*(y1-y3))

    return [px,py,(x1-px)**2+(y1-py)**2]

ans = 10**13

for v in itertools.combinations(p,2):
    p1 = v[0]
    p2 = v[1]
    px = (p1[0]+p2[0])/2
    py = (p1[1]+p2[1])/2
    rr = (p1[0]-p2[0])**2+(p1[1]-p2[1])**2
    rr = rr /4
    chk = True
    for i in range(n):
        if (p[i][0]-px)**2+(p[i][1]-py)**2 > rr + 10**(-6):
            chk = False
    if chk:
        ans = min(ans,rr**0.5)

for v in itertools.combinations(p,3):
    d = rad(v)
    if d == None:
        continue
    chk = True
    for i in range(n):
        if (p[i][0]-d[0])**2+(p[i][1]-d[1])**2 > d[2] + 10**(-6):
            chk = False
    if chk:
        ans = min(ans,d[2]**0.5)

print(ans)

Solution (Python):
import sys
from itertools import combinations

import numpy as np

n = int(input())
points = np.array([list(map(int, line.split())) for line in sys.stdin], dtype=np.float64)

if n == 2:
    p1, p2 = points
    print(np.sqrt(((p1 - p2) ** 2).sum()) / 2)
    exit()

ans = 0
for p1, p2, p3 in combinations(points, 3):
    A = ((p2 - p3) ** 2).sum()
    B = ((p3 - p1) ** 2).sum()
    C = ((p1 - p2) ** 2).sum()
    sa, sb, sc = sorted([A, B, C])
    if sa + sb <= sc + 1e-7:
        R2 = sc / 4
    else:
        T = A * (B + C - A)
        U = B * (C + A - B)
        W = C * (A + B - C)
        c = (T * p1 + U * p2 + W * p3) / (T + U + W)
        R2 = ((c - p1) ** 2).sum()
    ans = max(ans, R2)
print(np.sqrt(ans))

Solution (Python):
from itertools import combinations
N = int(input())
XY = [tuple(map(int, input().split())) for _ in range(N)]
l = 0
r = 2000
for _ in range(50):
    m = (l + r) / 2
    C = []
    for (x1, y1), (x2, y2) in combinations(XY, 2):
        d = ((x1-x2)**2+(y1-y2)**2)**0.5
        if d >= 2 * m:
            l = m
            break
        h = (m ** 2 - (d / 2) ** 2)**0.5
        hx = (y2-y1)/d*h
        hy = (x2-x1)/d*h
        x = (x1 + x2)/2
        y = (y1 + y2)/2
        C.append((x+hx,y-hy))
        C.append((x-hx,y+hy))
    else:
        if any(all(((x-x0)**2+(y-y0)**2)**0.5<=m+10**-8 for x0, y0 in XY) for x, y in C):
            r = m
        else:
            l = m
print(r)

Solution (Python):
import math
N, = map(int, input().split())
xs = []
ys = []
for _ in range(N):
    x, y = map(int, input().split())
    xs.append(x)
    ys.append(y)

def bb(a,b,c,d,e,f):
    if not (c-a)*(f-b)-(e-a)*(d-b):
        return None, None
    aa = a * a
    bb = b * b
    cc = c * c
    dd = d * d
    ee = e * e
    ff = f * f
    py = ((e - a) * (aa + bb - cc - dd) - (c - a) * (aa + bb - ee- ff)) / (2 * (e - a)*(b - d) - 2 * (c - a) * (b - f))
    px = (2 * (b - f) * py - aa - bb + ee + ff) / (2 * (e - a)) \
        if (c == a) else (2 * (b - d) * py - aa - bb + cc + dd) / (2 * (c - a))
    return px, py


pxs = []
pys = []
for i in range(N):
    for j in range(i+1, N):
        for k in range(j+1, N):

            px, py = bb(xs[i],ys[i],xs[j],ys[j],xs[k],ys[k])
            if not px:
                continue
            pxs.append(px)
            pys.append(py)
for i in range(N):
    for j in range(i+1, N):
        pxs.append((xs[i]+xs[j])/2)
        pys.append((ys[i]+ys[j])/2)

rr = 10**10
for px, py in zip(pxs,pys):
    r = 0
    for x, y in zip(xs, ys):
        r = max(r, (x-px)**2+(y-py)**2)
    rr = min(rr, r)
print(rr**0.5)

Solution (Python):
import sys
from math import hypot
import random










def make_circle(points):
    
    shuffled = [(float(x), float(y)) for (x, y) in points]
    random.shuffle(shuffled)

    
    c = None
    for (i, p) in enumerate(shuffled):
        if c is None or not is_in_circle(c, p):
            c = _make_circle_one_point(shuffled[: i + 1], p)
    return c



def _make_circle_one_point(points, p):
    c = (p[0], p[1], 0.0)
    for (i, q) in enumerate(points):
        if not is_in_circle(c, q):
            if c[2] == 0.0:
                c = make_diameter(p, q)
            else:
                c = _make_circle_two_points(points[: i + 1], p, q)
    return c



def _make_circle_two_points(points, p, q):
    circ = make_diameter(p, q)
    left = None
    right = None
    px, py = p
    qx, qy = q

    
    for r in points:
        if is_in_circle(circ, r):
            continue

        
        cross = _cross_product(px, py, qx, qy, r[0], r[1])
        c = make_circumcircle(p, q, r)
        if c is None:
            continue
        elif cross > 0.0 and (
                left is None or _cross_product(px, py, qx, qy, c[0], c[1]) > _cross_product(px, py, qx, qy, left[0],
                                                                                            left[1])):
            left = c
        elif cross < 0.0 and (
                right is None or _cross_product(px, py, qx, qy, c[0], c[1]) < _cross_product(px, py, qx, qy, right[0],
                                                                                             right[1])):
            right = c

    
    if left is None and right is None:
        return circ
    elif left is None:
        return right
    elif right is None:
        return left
    else:
        return left if (left[2] <= right[2]) else right


def make_diameter(a, b):
    cx = (a[0] + b[0]) / 2.0
    cy = (a[1] + b[1]) / 2.0
    r0 = hypot(cx - a[0], cy - a[1])
    r1 = hypot(cx - b[0], cy - b[1])
    return (cx, cy, max(r0, r1))


def make_circumcircle(a, b, c):
    
    ox = (min(a[0], b[0], c[0]) + max(a[0], b[0], c[0])) / 2.0
    oy = (min(a[1], b[1], c[1]) + max(a[1], b[1], c[1])) / 2.0
    ax = a[0] - ox
    ay = a[1] - oy
    bx = b[0] - ox
    by = b[1] - oy
    cx = c[0] - ox
    cy = c[1] - oy
    d = (ax * (by - cy) + bx * (cy - ay) + cx * (ay - by)) * 2.0
    if d == 0.0:
        return None
    x = ox + ((ax * ax + ay * ay) * (by - cy) + (bx * bx + by * by) * (cy - ay) + (cx * cx + cy * cy) * (ay - by)) / d
    y = oy + ((ax * ax + ay * ay) * (cx - bx) + (bx * bx + by * by) * (ax - cx) + (cx * cx + cy * cy) * (bx - ax)) / d
    ra = hypot(x - a[0], y - a[1])
    rb = hypot(x - b[0], y - b[1])
    rc = hypot(x - c[0], y - c[1])
    return (x, y, max(ra, rb, rc))


_MULTIPLICATIVE_EPSILON = 1 + 1e-14


def is_in_circle(c, p):
    return c is not None and hypot(p[0] - c[0], p[1] - c[1]) <= c[2] * _MULTIPLICATIVE_EPSILON



def _cross_product(x0, y0, x1, y1, x2, y2):
    return (x1 - x0) * (y2 - y0) - (y1 - y0) * (x2 - x0)


read = sys.stdin.read

N, *xy = map(int, read().split())
points = []
for x, y in zip(*[iter(xy)] * 2):
    points.append((x, y))

print(make_circle(points)[2])

Solution (Python):
import sys
input = sys.stdin.readline
import math

def gai(x1,y1,x2,y2,x3,y3):
    

    
    a2 = (x2 - x3) **2 + (y2 -y3) ** 2
    b2 = (x1 - x3) **2 + (y1 -y3) ** 2
    c2 = (x2 - x1) **2 + (y2 -y1) ** 2
    area = ((x2 - x1) * (y3 -y1) - (x3 -x1) * (y2- y1)) / 2

    
    if area == 0:
        xmax=max(x1,x2,x3)
        xmin=min(x1,x2,x3)
        x=xmax-xmin
        ymax=max(y1,y2,y3)
        ymin=min(y1,y2,y3)
        y=ymax-ymin
        r=(x*x+y*y)**0.5/2
        gx=(xmax+xmin)/2
        gy=(ymax+ymin)/2
        return gx,gy,r

    
    x = (a2*(b2+c2-a2)*x1 + b2*(c2+a2-b2)*x2 + c2*(a2+b2-c2)*x3) / (16*area*area)
    y = (a2*(b2+c2-a2)*y1 + b2*(c2+a2-b2)*y2 + c2*(a2+b2-c2)*y3) / (16*area*area)
    r = ((x1-x)**2 + (y1-y)**2)**0.5

    
    
    return x, y, r

def gai2(x1,y1,x2,y2):
    return (x1+x2)/2,(y1+y2)/2,((x1-x2)**2+(y1-y2)**2)**0.5/2

N=int(input())
X=[]
Y=[]
for i in range(N):
    x,y=map(int,input().split())
    X.append(x)
    Y.append(y)
if N==2:
    print(((X[1]-X[0])**2+(Y[1]-Y[0])**2)**0.5/2)
    exit()

dm=float("inf")
for i in range(N-2):
    for j in range(i+1,N-1):
        for k in range(j+1,N):
            """gx=(X[i]+X[j]+X[k])/3
            gy=(Y[i]+Y[j]+Y[k])/3
            han=((gx-X[i])**2+(gy-Y[i])**2)"""
            gx,gy,r=gai(X[i],Y[i],X[j],Y[j],X[k],Y[k])
            
            if dm>r:
                flg=True
                for l in range(N):
                    if l==i or l==j or l==k:
                        continue
                    if ((X[l]-gx)**2+(Y[l]-gy)**2)>r**2:
                        flg=False
                if flg:
                    dm=r
for i in range(N-1):
    for j in range(N):
        gx,gy,r=gai2(X[i],Y[i],X[j],Y[j])
        
        if dm>r:
            flg=True
            for l in range(N):
                if l==i or l==j:
                    continue
                if ((X[l]-gx)**2+(Y[l]-gy)**2)>r**2:
                    flg=False
            if flg:
                dm=r
print(dm)

Solution (Python):

import sys, re
from collections import deque, defaultdict, Counter
from math import ceil, sqrt, hypot, factorial, pi, sin, cos, radians, atan, degrees
from itertools import permutations, combinations, product, accumulate
from operator import itemgetter, mul
from copy import deepcopy
from string import ascii_lowercase, ascii_uppercase, digits
from fractions import gcd
from bisect import bisect
import heapq

def input(): return sys.stdin.readline().strip()
def INT(): return int(input())
def MAP(): return map(int, input().split())
def LIST(): return list(map(int, input().split()))
sys.setrecursionlimit(10 ** 9)
INF = float('inf')
mod = 10 ** 9 + 7

def circumcenter(data):
    
    (x1, y1), (x2, y2), (x3, y3) = data

    
    a2 = (x2 - x3) **2 + (y2 -y3) ** 2 
    b2 = (x1 - x3) **2 + (y1 -y3) ** 2
    c2 = (x2 - x1) **2 + (y2 -y1) ** 2
    area = ((x2 - x1) * (y3 -y1) - (x3 -x1) * (y2- y1)) / 2

    
    if area == 0:
        return -1
    
    
    x = (a2*(b2+c2-a2)*x1 + b2*(c2+a2-b2)*x2 + c2*(a2+b2-c2)*x3) / (16*area*area)
    y = (a2*(b2+c2-a2)*y1 + b2*(c2+a2-b2)*y2 + c2*(a2+b2-c2)*y3) / (16*area*area) 
    r = ((x1-x)**2 + (y1-y)**2)**0.5
    
    
    x, y, r = map(lambda x:round(x, 20) if x!=int(x) else int(x), (x, y, r))
    
    
    return x, y, r

N = INT()
xy = [LIST() for _ in range(N)]
xy.sort()

if N == 2:
    x1, y1 = xy[0]
    x2, y2 = xy[1]
    print((((x2-x1)/2)**2+((y2-y1)/2)**2)**(1/2))
    sys.exit()

ans = INF
epsilon = 10**(-6)
def dist(p1, p2):
    return ((p1[0]-p2[0])**2+(p1[1]-p2[1])**2)**(1/2)

for points in combinations(xy, 3):
    res = circumcenter(points)
    if res == -1:
        continue
    x = res[0]
    y = res[1]
    r = res[2]
    for p in xy:
        if dist(p, (x, y)) > r + epsilon:
            break
    else:
        ans = min(ans, r)

for points in combinations(xy, 2):
    x1, y1 = points[0]
    x2, y2 = points[1]
    x = (x1+x2)/2
    y = (y1+y2)/2
    r = (((x2-x1)/2)**2+((y2-y1)/2)**2)**(1/2)
    for p in xy:
        if dist(p, (x, y)) > r + epsilon:
            break
    else:
        ans = min(ans, r)
    
print(ans)

Solution (Python):
from heapq import heappush,heappop,heapify
from collections import deque,defaultdict,Counter
import itertools
from itertools import permutations,combinations,groupby
import sys
import bisect
import string
import math
import time
import random
def S_():
    return input()
def LS():
    return [i for i in input().split()]
def I():
    return int(input())
def MI():
    return map(int,input().split())
def LI():
    return [int(i) for i in input().split()]
def LI_():
    return [int(i)-1 for i in input().split()]
def NI(n):
    return [int(input()) for i in range(n)]
def NI_(n):
    return [int(input())-1 for i in range(n)]
def StoI():
    return [ord(i)-97 for i in input()]
def ItoS(nn):
    return chr(nn+97)
def LtoS(ls):
    return ''.join([chr(i+97) for i in ls])
def GI(V,E,Directed=False,index=0):
    org_inp=[]
    g=[[] for i in range(n)]
    for i in range(E):
        inp=LI()
        org_inp.append(inp)
        if index==0:
            inp[0]-=1
            inp[1]-=1
        if len(inp)==2:
            a,b=inp
            g[a].append(b)
            if not Directed:
                g[b].append(a)
        elif len(inp)==3:
            a,b,c=inp
            aa=(inp[0],inp[2])
            bb=(inp[1],inp[2])
            g[a].append(bb)
            if not Directed:
                g[b].append(aa)
    return g,org_inp
def bit_combination(k,n=2):
    rt=[]
    for tb in range(n**k):
        s=[tb//(n**bt)%n for bt in range(k)]
        rt+=[s]
    return rt
def show(*inp,end='\n'):
    if show_flg:
        print(*inp,end=end)
YN=['Yes','No']
mo=10**9+7
inf=float('inf')
l_alp=string.ascii_lowercase
u_alp=string.ascii_uppercase
ts=time.time()
sys.setrecursionlimit(10**5)
input=lambda: sys.stdin.readline().rstrip()

def ran_input():
    import random
    n=random.randint(4,16)
    rmin,rmax=1,10
    a=[random.randint(rmin,rmax) for _ in range(n)]
    return n,a

class Comb:
    def __init__(self,n,mo=10**9+7):
        self.fac=[0]*(n+1)
        self.inv=[1]*(n+1)
        self.fac[0]=1
        self.fact(n)
        for i in range(1,n+1):
            self.fac[i]=i*self.fac[i-1]%mo
            self.inv[n]*=i
            self.inv[n]%=mo
        self.inv[n]=pow(self.inv[n],mo-2,mo)
        for i in range(1,n):
            self.inv[n-i]=self.inv[n-i+1]*(n-i+1)%mo
        return
    
    def fact(self,n):
        return self.fac[n]
        
    def invf(self,n):
        return self.inv[n]

    def comb(self,x,y):
        if y<0 or y>x:
            return 0
        return self.fac[x]*self.inv[x-y]*self.inv[y]%mo

show_flg=False
show_flg=True


def check(x,y,r):
    for s,t in c:
        if (x-s)**2+(y-t)**2>r+1e-6:
            return False
    return True
    

def cir3(x1,y1,x2,y2,x3,y3):
    p = ((y1-y3)*(y1**2 -y2**2 +x1**2 -x2**2) -(y1-y2)*(y1**2 -y3**2 +x1**2 -x3**2)) / (2*(y1-y3)*(x1-x2)-2*(y1-y2)*(x1-x3))
    q = ((x1-x3)*(x1**2 -x2**2 +y1**2 -y2**2) -(x1-x2)*(x1**2 -x3**2 +y1**2 -y3**2)) / (2*(x1-x3)*(y1-y2)-2*(x1-x2)*(y1-y3))
    r=((x1-p)**2+(y1-q)**2)
    
    return p,q,r

def rot(x,y,q):
    s=x*math.cos(q)-y*math.sin(q)
    t=x*math.sin(q)+y*math.cos(q)
    return (s,t)
    
    
def cir2(x1,y1,x2,y2):
    p,q=(x1+x2)/2,(y1+y2)/2
    r=((x1-p)**2+(y1-q)**2)
    return p,q,r

ans=inf

n=I()
c=[]
CC=[]
rx=random.randint(1,10**3)
ry=random.randint(1,10**3)
for i in range(n):
    x,y=LI()
    x+=rx
    y+=ry
    CC.append((x,y))


def solve():
    ans=inf
    q=math.pi*random.randint(1,10**4)/10**4
    
    for i in range(n):
        x,y=CC[i]
        x,y=rot(x,y,q)
        c.append((x,y))
    
    
    
    for i in range(n):
        for j in range(i+1,n):
                po=c[i]+c[j]
                x1,y1,x2,y2=po
                p,q,r=cir2(x1,y1,x2,y2)
                
                if check(p,q,r):
                    ans=min(ans,r)
    
    
    for i in range(n):
        for j in range(i+1,n):
            for k in range(j+1,n):
                x1,y1,x2,y2,x3,y3=CC[i]+CC[j]+CC[k]
                a,b= (x2-x1),(y2-y1)
                C,D= (x3-x1),(y3-y1)
                if a*D==b*C:
                    continue
                x1,y1,x2,y2,x3,y3=c[i]+c[j]+c[k]
                p,q,r=cir3(x1,y1,x2,y2,x3,y3)
                
                if check(p,q,r):
                    ans=min(ans,r)

    return ans

an=inf
for rn in range(1):
    an=min(an,solve())
    

print(an**0.5)

Solution (Python):

from itertools import product
import numpy as np

def dist(p1, p2=np.zeros(2,)):
    return np.sqrt(((p1-p2)**2).sum())


def MD(P):
    c = (P[0]+P[1])/2
    r = dist(P[0], P[1])/2
    D = (c, r)
    for i in range(2,len(P)):
        if dist(D[0], P[i]) > D[1]:
            D = MDWP(P[:i], P[i])
    return D

def MDWP(P, q1):
    c = (P[0]+q1)/2
    r = dist(P[0], q1)/2
    D = (c, r)
    for i in range(1,len(P)):
        if dist(D[0], P[i]) > D[1]:
            D = MDWPs(P[:i], P[i], q1)
    return D

def MDWPs(P, q1, q2):
    c = (q1+q2)/2
    r = dist(q1, q2)/2
    D = (c, r)
    for i in range(0,len(P)):
        if dist(D[0], P[i]) > D[1]:
            D = Disk(q1, q2, P[i])
    return D

def Disk(q1, q2, q3):
    a, b = q1-q2
    c, d = q2-q3
    z = (q1**2).sum() - (q2**2).sum()
    w = (q2**2).sum() - (q3**2).sum()
    e = 2*(a*d-b*c)
    if e == 0:
        c = (q1+q2)/2
        r = dist(q1, q2)/2
        return (c, r)
    else:
        c = np.array([d*z-b*w, -c*z+a*w])/e
        r = dist(q1, c)
        return (c, r)



def solve():
    N = int(input())
    P = np.array([np.array(list(map(int, input().split()))) for _ in range(N)])
    D = MD(P)
    return str(D[1])

if __name__ == '__main__':
    print(solve())

Solution (Python):
from scipy.optimize import minimize
import numpy as np

N = int(input())
xy = [list(map(int, input().split())) for _ in range(N)]

xy = np.array(xy)


def r2(ab):
    xymab = xy - ab
    return np.square(xymab).sum(1).max()


x0 = np.zeros(2)
res = minimize(r2, x0, method="Nelder-Mead", tol=1e-6)

print(res.fun ** 0.5)

Solution (Python):
from itertools import combinations
from math import hypot

def solve():
    INF = float('inf')
    eps = 10**-6

    N = int(input())
    pts = [tuple(map(int, input().split())) for _ in range(N)]

    if N == 2:
        x1, y1 = pts[0]
        x2, y2 = pts[1]
        ans = hypot(x1-x2, y1-y2) / 2
        print(ans)
        return

    lines = {}
    for i in range(N):
        x1, y1 = pts[i]
        for j in range(i+1, N):
            x2, y2 = pts[j]
            if y1 == y2:
                a = INF
                b = (x1+x2)/2
            else:
                x0, y0 = (x1+x2)/2, (y1+y2)/2
                a = -(x2-x1)/(y2-y1)
                b = y0 - a*x0
            lines[(i, j)] = (a, b)


    ans = INF
    for i, j, k in combinations(range(N), r=3):
        a1, b1 = lines[(i, j)]
        a2, b2 = lines[(i, k)]
        if a1 == a2: continue
        if a1 == INF:
            x = b1
            y = a2*x + b2
        elif a2 == INF:
            x = b2
            y = a1*x + b1
        else:
            x = (b2-b1) / (a1-a2)
            y = (a1*b2-b1*a2) / (a1-a2)
        x1, y1 = pts[i]
        r = hypot(x1-x, y1-y)

        for xi, yi in pts:
            dist = hypot(xi-x, yi-y)
            if dist > r + eps:
                break
        else:
            ans = min(ans, r)

    for i in range(N):
        x1, y1 = pts[i]
        for j in range(i+1, N):
            x2, y2 = pts[j]
            x, y = (x1+x2)/2, (y1+y2)/2
            r = hypot(x1-x, y1-y)
            for xi, yi in pts:
                dist = hypot(xi-x, yi-y)
                if dist > r + eps:
                    break
            else:
                ans = min(ans, r)

    print(ans)


solve()

Solution (Python):
import sys

def f(p,q):
    num=0
    for ppp,qqq in xy:
        num=max(num,(p-ppp)**2+(q-qqq)**2)
    return num

inf=float("inf")
N=int(input())
xy=[list(map(int,input().split())) for i in range(N)]
ans=inf
if N==2:
    x=(xy[0][0]+xy[1][0])/2
    y=(xy[0][1]+xy[1][1])/2
    p,q=xy[0]
    num=(x-p)**2+(y-q)**2
    print(num**0.5)
    sys.exit()
for i in range(N-2):
    for j in range(i+1,N-1):
        for k in range(j+1,N):
            p,q=xy[i]
            r,s=xy[j]
            v,w=xy[k]
            a=2*v-2*r
            b=2*w-2*s
            c=v**2+w**2-r**2-s**2
            d=2*p-2*r
            e=2*q-2*s
            f=p**2+q**2-r**2-s**2
            if a*e==b*d:
                num=0
            else:
                x=(c*e-b*f)/(a*e-b*d)
                y=(c*d-a*f)/(b*d-a*e)
                num=0
                for ppp,qqq in xy:
                    num=max(num,(x-ppp)**2+(y-qqq)**2)
                ans=min(ans,num)
            x,y=(p+r)/2,(q+s)/2
            num=0
            for ppp,qqq in xy:
                num=max(num,(x-ppp)**2+(y-qqq)**2)
            ans=min(ans,num)
            x,y=(p+v)/2,(q+w)/2
            num=0
            for ppp,qqq in xy:
                num=max(num,(x-ppp)**2+(y-qqq)**2)
            ans=min(ans,num)
            x,y=(r+v)/2,(s+w)/2
            num=0
            for ppp,qqq in xy:
                num=max(num,(x-ppp)**2+(y-qqq)**2)
            ans=min(ans,num)
print(ans**0.5)