Problem:


Problem D: Circle and Points

You are given N points in the xy-plane. You have a
circle of radius one and move it on the xy-plane, so as to
enclose as many of the points as possible. Find how many points can be
simultaneously enclosed at the maximum. A point is considered
enclosed by a circle when it is inside or on the circle.

Fig 1. Circle and Points

Input

The input consists of a series of data sets, followed by a single line
only containing a single character '0', which indicates the end of the
input.  Each data set begins with a line containing an integer N,
which indicates the number of points in the data set. It is
followed by N lines describing the coordinates of the
points. Each of the N lines has two decimal fractions X
and Y, describing the x- and y-coordinates of a
point, respectively. They are given with five digits after the decimal
point.

You may assume 1 <= N <= 300,  0.0 <= X <= 10.0, and 0.0
<= Y <= 10.0. No two points are closer than 0.0001.  No two
points in a data set are approximately at a distance of 2.0. More
precisely, for any two points in a data set, the distance d
between the two never satisfies 1.9999 <= d <= 2.0001. Finally,
no three points in a data set are simultaneously very close to a
single circle of radius one. More precisely,
let P1, P2, and P3 be any three points in a
data set, and d1, d2, and d3
the distances from an arbitrarily selected point in the xy-plane to each of them
respectively. Then it never simultaneously holds that 0.9999 <=
di <= 1.0001 (i = 1, 2, 3).

Output

For each data set, print a single line containing the maximum number
of points in the data set that can be simultaneously enclosed by a
circle of radius one. No other characters including leading and
trailing spaces should be printed.

Sample Input

3
6.47634 7.69628
5.16828 4.79915
6.69533 6.20378
6
7.15296 4.08328
6.50827 2.69466
5.91219 3.86661
5.29853 4.16097
6.10838 3.46039
6.34060 2.41599
8
7.90650 4.01746
4.10998 4.18354
4.67289 4.01887
6.33885 4.28388
4.98106 3.82728
5.12379 5.16473
7.84664 4.67693
4.02776 3.87990
20
6.65128 5.47490
6.42743 6.26189
6.35864 4.61611
6.59020 4.54228
4.43967 5.70059
4.38226 5.70536
5.50755 6.18163
7.41971 6.13668
6.71936 3.04496
5.61832 4.23857
5.99424 4.29328
5.60961 4.32998
6.82242 5.79683
5.44693 3.82724
6.70906 3.65736
7.89087 5.68000
6.23300 4.59530
5.92401 4.92329
6.24168 3.81389
6.22671 3.62210
0

Output for the Sample Input

2
5
5
11



Sample input:
3
6.47634 7.69628
5.16828 4.79915
6.69533 6.20378
6
7.15296 4.08328
6.50827 2.69466
5.91219 3.86661
5.29853 4.16097
6.10838 3.46039
6.34060 2.41599
8
7.90650 4.01746
4.10998 4.18354
4.67289 4.01887
6.33885 4.28388
4.98106 3.82728
5.12379 5.16473
7.84664 4.67693
4.02776 3.87990
20
6.65128 5.47490
6.42743 6.26189
6.35864 4.61611
6.59020 4.54228
4.43967 5.70059
4.38226 5.70536
5.50755 6.18163
7.41971 6.13668
6.71936 3.04496
5.61832 4.23857
5.99424 4.29328
5.60961 4.32998
6.82242 5.79683
5.44693 3.82724
6.70906 3.65736
7.89087 5.68000
6.23300 4.59530
5.92401 4.92329
6.24168 3.81389
6.22671 3.62210
0


Sample output:
2
5
5
11

Submission status: Wrong Answer

Solution code:
while True:
    n = int(raw_input())
    if n == 0:
        break
    pts = [complex(*[float(_) for _ in raw_input().split()]) for i in xrange(n)]
    ans = 1
    for i, p in enumerate(pts):
        for j, q in enumerate(pts[i + 1:], i + 1):
            pq = abs(q - p)
            if pq > 2.0:
                continue
            m = (p + q) / 2.0
            mo = (4.0 - pq ** 2) / 4.0 * (q - p) / pq * 1j
            o1 = m + mo
            o2 = m - mo
            count1 = count2 = 2
            for k, r in enumerate(pts):
                if k == i or k == j:
                    continue
                if abs(r - o1) <= 1.0:
                    count1 += 1
                if abs(r - o2) <= 1.0:
                    count2 += 1
            ans = max(ans, count1, count2)
    print ans

Submission status: Wrong Answer

Solution code:
#!/usr/bin/env python

from collections import deque
import itertools as it
import sys
import math

sys.setrecursionlimit(1000000)

def func():
    while True:
        N = input()
        if N == 0:
            break
        pos = []
        for loop in range(N):
            S = map(float, raw_input().split())
            pos.append(S)
        ans = 0
        for i in range(N):
            for j in range(N):
                v = [pos[j][0] - pos[i][0], pos[j][1] - pos[i][1]]
                norm = math.sqrt(v[0] * v[0] + v[1] * v[1])
                if norm > 2 or norm == 0:
                    continue
                vn = [v[0] / norm, v[1] / norm]
                norm /= 2
                norm2 = math.sqrt(1 - norm * norm)
                v90 = [-vn[1], vn[0]]
                cent = [pos[i][0] + v[0] / 2 + v90[0] * norm2, 
                pos[i][1] + v[1] / 2 + v90[1] * norm2]
                S = 2
                for k in range(N):
                    if i == k or j == k:
                        continue
                    x = cent[0] - pos[k][0]
                    y = cent[1] - pos[k][1]
                    if x * x + y * y <= 1:
                        S += 1
                ans = max(S, ans)
        if ans == 80:
            return 0
        print ans

func()

Submission status: Wrong Answer

Solution code:
#!/usr/bin/env python

from collections import deque
import itertools as it
import sys
import math

sys.setrecursionlimit(1000000)

def func():
    while True:
        N = input()
        if N == 0:
            break
        R = range(N)
        pos = []
        for loop in R:
            S = map(float, raw_input().split())
            pos.append(S)
        ans = 0
        for i in R:
            for j in R:
                v = [pos[j][0] - pos[i][0], pos[j][1] - pos[i][1]]
                norm = math.sqrt(v[0] * v[0] + v[1] * v[1])
                if norm > 2 or norm == 0:
                    continue
                vn = [v[0] / norm, v[1] / norm]
                norm /= 2
                norm2 = math.sqrt(1 - norm * norm)
                v90 = [-vn[1], vn[0]]
                cent = [pos[i][0] + v[0] / 2 + v90[0] * norm2, 
                pos[i][1] + v[1] / 2 + v90[1] * norm2]
                S = 2
                for k in R:
                    if i == k or j == k:
                        continue
                    x = cent[0] - pos[k][0]
                    y = cent[1] - pos[k][1]
                    if x * x + y * y <= 1:
                        S += 1
                ans = max(S, ans)
        if ans == 80:
            return 0
        print max(1, ans)

func()

Submission status: Wrong Answer

Solution code:
#!/usr/bin/env python

from collections import deque
import itertools as it
import sys
import math

sys.setrecursionlimit(1000000)

def func():
    lst = []
    while True:
        S = input()
        lst.append(S)
        if S == '0':
            break
    return 0
    lst = list(reversed(lst))
    while True:
        N = int(lst.pop())
        if N == 0:
            break
        R = range(N)
        pos = []
        for loop in R:
            S = list(map(float, lst.pop().split()))
            pos.append(S)
        ans = 0
        for i in R:
            for j in R:
                v = [pos[j][0] - pos[i][0], pos[j][1] - pos[i][1]]
                norm = math.sqrt(v[0] * v[0] + v[1] * v[1])
                if norm > 2 or norm == 0:
                    continue
                vn = [v[0] / norm, v[1] / norm]
                norm /= 2
                norm2 = math.sqrt(1 - norm * norm)
                v90 = [-vn[1], vn[0]]
                cent = [pos[i][0] + v[0] / 2 + v90[0] * norm2, 
                pos[i][1] + v[1] / 2 + v90[1] * norm2]
                S = 2
                for k in R:
                    if i == k or j == k:
                        continue
                    x = cent[0] - pos[k][0]
                    y = cent[1] - pos[k][1]
                    if x * x + y * y <= 1:
                        S += 1
                ans = max(S, ans)
        print(max(1, ans))

if __name__ == '__main__':
    func()

Submission status: Wrong Answer

Solution code:
while 1:
	n = input()
	if n == 0: break
	pxy = [map(float,raw_input().split()) for i in range(n)]
	ans = 0
	for y in range(100):
		for x in range(100):
			count = 0
			for px,py in pxy:
				if (x/10.0-px)**2+(y/10.0-py)**2 <= 1.0:
					count += 1
			ans = max(ans,count)
	print ans

Submission status: Wrong Answer

Solution code:
import random
ans = 0
while 1:
	n = input()
	if n == 0: break
	xy = [map(float,raw_input().split()) for i in range(n)]
	for i in range(10000):
		cx,cy = 10*random.random(),10*random.random()
		cur = 0
		for x,y in xy:
			if (cx-x)**2 + (cy-y)**2 <= 1:
				cur += 1
		ans = max(ans,cur)
	print ans

Submission status: Wrong Answer

Solution code:
import random
ans = 0
while 1:
	n = input()
	if n == 0: break
	xy = [map(float,raw_input().split()) for i in range(n)]
	for r in range(10000):
		cx,cy = 10*random.random(),10*random.random()
		cur = 0
		for i in range(n):
			x,y = xy[i]
			if (cx-x)**2 + (cy-y)**2 <= 1:
				cur += 1
			if n - (i+1) + cur < ans: break
		ans = max(ans,cur)
	print ans

Submission status: Accepted

Solution code:
from math import sqrt
from bisect import bisect_left
def circle_center(x1, y1, x2, y2):
    xd = x2 - x1; yd = y2 - y1
    d = xd**2 + yd**2
    k = sqrt((4.0 - d) / d) / 2.0
    xc = (x1 + x2) / 2.0
    yc = (y1 + y2) / 2.0
    return [[xc - k*yd, yc + k*xd], [xc + k*yd, yc - k*xd]]
while 1:
    n = int(raw_input())
    if n==0: break
    p = [map(float, raw_input().split()) for i in xrange(n)]
    p.sort()
    prev = 0
    ans = 1
    for i in xrange(n):
        bx, by = p[i]
        while bx - p[prev][0] >= 2.0: prev += 1
        for j in xrange(i+1, n):
            cx, cy = p[j]
            if cx - bx >= 2.0: break
            if (bx - cx)**2 + (by - cy)**2 <= 4.0:
                for ex, ey in circle_center(bx, by, cx, cy):
                    count = 2
                    for k in xrange(prev, n):
                        if k==i or k==j: continue
                        dx, dy = p[k]
                        if dx - bx >= 2.0: break
                        if (ex - dx)**2 + (ey - dy)**2 <= 1.0:
                            count += 1
                    ans = max(ans, count)
    print ans

Submission status: Accepted

Solution code:
from math import sqrt
from bisect import bisect_left
def circle_center(x1, y1, x2, y2):
    xd = x2 - x1; yd = y2 - y1
    d = xd**2 + yd**2
    k = sqrt((4.0 - d) / d) / 2.0
    xc = (x1 + x2) / 2.0; yc = (y1 + y2) / 2.0
    xd *= k; yd *= k
    return [[xc - yd, yc + xd], [xc + yd, yc - xd]]
while 1:
    n = int(raw_input())
    if n==0: break
    p = sorted(map(float, raw_input().split()) for i in xrange(n))
    prev = 0
    ans = 1
    for i in xrange(n):
        bx, by = p[i]
        while bx - p[prev][0] >= 2.0: prev += 1
        for j in xrange(i+1, n):
            cx, cy = p[j]
            if cx - bx >= 2.0: break
            if (bx - cx)**2 + (by - cy)**2 <= 4.0:
                for ex, ey in circle_center(bx, by, cx, cy):
                    count = 2
                    for k in xrange(prev, n):
                        if k==i or k==j: continue
                        dx, dy = p[k]
                        if dx - ex >= 1.0: break
                        if (ex - dx)**2 + (ey - dy)**2 <= 1.0:
                            count += 1
                    ans = max(ans, count)
    print ans

Submission status: Accepted

Solution code:
from math import sqrt
from bisect import bisect_left
def circle_center(x1, y1, x2, y2):
    xd = x2 - x1; yd = y2 - y1
    d = xd**2 + yd**2
    k = sqrt((4.0 - d) / d) / 2.0
    xc = (x1 + x2) / 2.0; yc = (y1 + y2) / 2.0
    xd *= k; yd *= k
    return [[xc - yd, yc + xd], [xc + yd, yc - xd]]
while 1:
    n = int(raw_input())
    if n==0: break
    p = sorted(map(float, raw_input().split()) for i in xrange(n))
    px = [e[0] for e in p]; py = [e[1] for e in p]
    prev = 0
    ans = 1
    for i in xrange(n):
        bx = px[i]; by = py[i]
        while bx - px[prev] >= 2.0: prev += 1
        for j in xrange(i+1, n):
            cx = px[j]; cy = py[j]
            if cx - bx >= 2.0: break
            if (bx - cx)**2 + (by - cy)**2 <= 4.0:
                for ex, ey in circle_center(bx, by, cx, cy):
                    count = 2
                    for k in xrange(prev, n):
                        if k==i or k==j: continue
                        dx = px[k]; dy = py[k]
                        if dx - ex >= 1.0: break
                        if (ex - dx)**2 + (ey - dy)**2 <= 1.0:
                            count += 1
                    ans = max(ans, count)
    print ans

Submission status: Accepted

Solution code:
from math import sqrt
def circle_center(x1, y1, x2, y2):
    xd = x2 - x1; yd = y2 - y1
    d = xd**2 + yd**2
    k = sqrt((4.0 - d) / d) / 2.0
    xc = (x1 + x2) / 2.0; yc = (y1 + y2) / 2.0
    xd *= k; yd *= k
    return [[xc - yd, yc + xd], [xc + yd, yc - xd]]
while 1:
    n = int(raw_input())
    if n==0: break
    p = sorted(map(float, raw_input().split()) for i in xrange(n))
    prev = 0
    ans = 1
    for i in xrange(n):
        bx, by = p[i]
        while bx - p[prev][0] >= 2.0: prev += 1
        for j in xrange(i+1, n):
            cx, cy = p[j]
            if cx - bx >= 2.0: break
            if abs(by - cy) <= 2.0 and (bx - cx)**2 + (by - cy)**2 <= 4.0:
                for ex, ey in circle_center(bx, by, cx, cy):
                    count = 2
                    for k in xrange(prev, n):
                        if k==i or k==j: continue
                        dx, dy = p[k]
                        if ex - dx >= 1.0: continue
                        if dx - ex >= 1.0: break
                        if abs(ey - dy) <= 1.0 and (ex - dx)**2 + (ey - dy)**2 <= 1.0:
                            count += 1
                    ans = max(ans, count)
    print ans

Submission status: Accepted

Solution code:
from math import sqrt
def circle_center(x1, y1, x2, y2):
    xd = x2 - x1; yd = y2 - y1
    k = sqrt(1.0/(xd**2 + yd**2) - 0.25)
    xc = (x1 + x2) / 2.0; yc = (y1 + y2) / 2.0
    xd *= k; yd *= k
    return [[xc - yd, yc + xd], [xc + yd, yc - xd]]
while 1:
    n = int(raw_input())
    if n==0: break
    p = sorted(map(float, raw_input().split()) for i in xrange(n))
    prev = 0; pvx = p[0][0]
    ans = 1
    for i in xrange(n):
        bx, by = p[i]
        bxm2 = bx - 2.0; bxp2 = bx + 2.0
        bym2 = by - 2.0; byp2 = by + 2.0
        while bxm2 >= pvx:
            prev += 1; pvx = p[prev][0]
        for j in xrange(i+1, n):
            cx, cy = p[j]
            if cx >= bxp2: break
            if bym2 <= cy <= byp2 and (bx - cx)**2 + (by - cy)**2 <= 4.0:
                for ex, ey in circle_center(bx, by, cx, cy):
                    count = 2
                    exm1 = ex - 1.0; exp1 = ex + 1.0
                    eym1 = ey - 1.0; eyp1 = ey + 1.0
                    q = prev
                    while exm1 >= p[q][0]: q += 1
                    for k in xrange(q, n):
                        if k==i or k==j: continue
                        dx, dy = p[k]
                        if dx >= exp1: break
                        if eym1 <= dy <= eyp1 and (ex - dx)**2 + (ey - dy)**2 <= 1.0:
                            count += 1
                    ans = max(ans, count)
    print ans

Submission status: Accepted

Solution code:
import math


def find_point(x1, y1, x2, y2):
    # 二頂点(x1,y1),(x2,y2)を円周上にもつ時の半径1の円の中心
    if (x1-x2)**2+(y1-y2)**2 > 4:
        return False
    mx, my = (x1+x2)/2, (y1+y2)/2
    L = math.sqrt((y2-y1)**2+(x2-x1)**2)
    k = math.sqrt(1-((y2-y1)**2+(x2-x1)**2)/4)
    # (mx,my)+-(k/L)(-y2+y1,x2-x1)
    X1, Y1 = mx+(k/L)*(-y2+y1), my+(k/L)*(x2-x1)
    X2, Y2 = mx-(k/L)*(-y2+y1), my-(k/L)*(x2-x1)
    return ((X1, Y1), (X2, Y2))


def solve():
    N = int(input())
    if N == 0:
        return
    que = [tuple(map(float, input().split())) for _ in range(N)]
    ans = 1
    for i in range(N):
        for j in range(i+1, N):
            A = find_point(que[i][0], que[i][1], que[j][0], que[j][1])
            if A != False:
                a, b = A
                suba = 0
                subb = 0
                for k in range(N):
                    if math.hypot(que[k][0]-a[0], que[k][1]-a[1]) < 1+10**(-7):
                        suba += 1
                    if math.hypot(que[k][0]-b[0], que[k][1]-b[1]) < 1+10**(-7):
                        subb += 1
                ans = max(ans, suba, subb)
    print(ans)
    return solve()


def main():
    solve()


if __name__ == "__main__":
    main()

Submission status: Accepted

Solution code:

import cmath
def counter(z,a):
    count= 0
    for i in a:
        if abs(i-z)<=1.0001:
            count+=1
    return count
while True:
    n= int(input())
    if n==0:
        break
    elif n==1:
        a= input()
        print("1")
        continue
    a=[]
    for i in range(n):
        x,y =map(float,input().split())
        a.append(complex(x,y))
    ans=1

    for i in range(0,n):
        for j in range(i+1,n):
            p= (a[j]-a[i])/2
            q= (a[j]+a[i])/2
            if abs(p)> 1:
                continue
            r= ((1-abs(p)**2)**0.5)/abs(p) * p*complex(0,1)
                
            z1= q+r
            z2= q-r

            ans= max(ans,counter(z1,a),counter(z2,a))
    print(ans)

Submission status: Accepted

Solution code:
import cmath
EPS = 1e-5

while True:
	n = int(input())
	if n == 0:
		break
	d = []
	for _ in range(n):
		x, y = list(map(float, input().split()))
		d.append(complex(x, y))

	res = 1
	for i in range(n):
		o1 = d[i]
		cnt0 = 1
		for j in range(n):
			if j == i:
				continue
			o2 = d[j]
			if abs(o1 - o2) + EPS <= 1:
				cnt0 += 1

		res = max(res, cnt0)
		for j in range(i+1, n):
			o2 = d[j]
			if abs(o1 - o2) - EPS >= 2:
				continue
			else:
				co = (abs(o1 - o2)**2) / (2 * abs(o1 - o2))
				si = cmath.sqrt(1 - co**2)
				p1 = o1 + (o2 - o1) * (co + 1j * si) / abs(o1 - o2)
				p2 = o1 + (o2 - o1) * (co - 1j * si) / abs(o1 - o2)

			cnt1, cnt2 = 2, 2
			for l in range(n):
				if l == i or l == j:
					continue
				else:
					o3 = d[l]
					if abs(p1 - o3) + EPS <= 1:
						cnt1 += 1
					if abs(p2 - o3) + EPS <= 1:
						cnt2 += 1
			res = max(res, cnt1, cnt2)
	print(res)

Submission status: Accepted

Solution code:
import math
def count(a, b, P):
    if abs(a-b) > 2:
        return 0
    if a == b:
        return 0
    m = (a+b)/2
    c = m + (m-a)*(math.sqrt(1-abs(m-a)**2))*1j/abs(m-a)
    cnt = 0
    for x in P:
        if abs(x-c) <=1.0001:
            cnt += 1
    return cnt

while True:
    N = int(input())
    if N == 0:
        break
    P = []
    for i in range(N):
        x, y = map(float, input().split())
        P.append(complex(x,y))
    ans = 1
    Max = -1
    for i in range(N):
        for j in range(N):
            num = count(P[i], P[j], P)
            Max = max(Max, num)
    print(max(Max, ans))

Submission status: Accepted

Solution code:
#複素数でやってみる
from itertools import combinations
import math

while True:
    N=int(input())
    if N==0:
        break
    mp=[]
    for _ in range(N):
        a,b=map(float,input().split())
        mp.append(complex(a,b))
    center_lst=[]
    for p,q in combinations(mp,2):
        mid=(p+q)/2
        if abs(p-q)>2:
            continue
        d_vec=complex(-(p-q).imag,(p-q).real)/abs(p-q)
        t=1-abs(p-mid)**2
        center_vec1=mid+d_vec*math.sqrt(t)
        center_vec2=mid-d_vec*math.sqrt(t)
        center_lst.append(center_vec1)
        center_lst.append(center_vec2)
    
    ans=1
    for center in center_lst:
        tmp=0
        for x in mp:
            x_center=x-center
            if abs(x_center)<=1+1e-7:
                tmp+=1
        if tmp>ans:
            ans=tmp
    print(ans)

Submission status: Accepted

Solution code:
#import matplotlib.pyplot as plt
#import matplotlib.patches as patches
from itertools import combinations
printed_later = []
while True:
    n = int(input())
    if(n == 0):
        break
    dots = []
    for _ in range(n):
        x, y = map(float, input().split())
        dots.append(complex(x, y))
    inside_max = 1
    for p, q in combinations(dots, 2):
        abs_p_q = abs(p - q)
        if abs_p_q > 2:
            continue
        for pl_mn_j in [1j, -1j]:
            center = (p + q) / 2 + (p - q) * pl_mn_j / abs_p_q * (1 - (abs_p_q / 2) ** 2) ** 0.5
            dots_inside_num = 0
            for dot in dots:
                #plt.gca().add_patch(patches.Circle(xy=(center.real, center.imag), radius=1, fill=False, ec='r'))  # debug
                if abs(dot - center) < 1 or dot == p or dot == q:
                    dots_inside_num += 1
            if dots_inside_num > inside_max:
                inside_max = dots_inside_num
    printed_later.append(inside_max)
    #plt.scatter([z.real for z in dots], [z.imag for z in dots])
    #plt.show()
for line in printed_later:
    print(line)

Submission status: Accepted

Solution code:
from math import sqrt

def center(a,b):
    m = (a+b)/2
    d = abs(a-m)
    vec= complex((b-a).imag,(a-b).real)
    length = sqrt(1-d**2)
    c1,c2 = m+vec*(length/(2*d)),m-vec*(length/(2*d))
    return c1,c2



def check(c):
    cnt = 0
    for loc in location:
        if abs(c-loc) <= 1.0+EPS:
            cnt += 1
    return cnt 

while True:
    n = int(input())
    if n==0:
        break
    location = []
    ans = 1
    EPS = 1e-7
    for _ in range(n):
        a,b = map(float ,input().split())
        location.append(complex(a,b))

    for i in range(n):
        for j in range(i+1,n):
            if abs(location[i]-location[j]) >=1.0*2:
                continue
            c1,c2 = center(location[i],location[j])
            ans = max(ans,check(c1),check(c2))
            
    print(ans)