我意识到圆不必以点为中心，因此计算通过两点任意组合的所有圆，包括以每个点为中心的圆。
然后我找到每个圆覆盖的点并使用贪心算法找到一个最小的圆集来覆盖所有点，但同样，它可能不是 最小的圆集，但很容易计算.
from collections import namedtuple
from itertools import product
from math import sqrt
from pprint import pprint as pp
Pt = namedtuple('Pt', 'x, y')
Cir = namedtuple('Cir', 'x, y, r')
def circles_from_p1p2r(p1, p2, r):
'Following explanation at http://mathforum.org/library/drmath/view/53027.html'
(x1, y1), (x2, y2) = p1, p2
if p1 == p2:
#raise ValueError('coincident points gives infinite number of Circles')
return None, None
# delta x, delta y between points
dx, dy = x2 - x1, y2 - y1
# dist between points
q = sqrt(dx**2 + dy**2)
if q > 2.0*r:
#raise ValueError('separation of points > diameter')
return None, None
# halfway point
x3, y3 = (x1+x2)/2, (y1+y2)/2
# distance along the mirror line
d = sqrt(r**2-(q/2)**2)
# One answer
c1 = Cir(x = x3 - d*dy/q,
y = y3 + d*dx/q,
r = abs(r))
# The other answer
c2 = Cir(x = x3 + d*dy/q,
y = y3 - d*dx/q,
r = abs(r))
return c1, c2
def covers(c, pt):
return (c.x - pt.x)**2 + (c.y - pt.y)**2 <= c.r**2
if __name__ == '__main__':
for r, points in [(3, [Pt(*i) for i in [(1, 3), (0, 2), (4, 5), (2, 4), (0, 3)]]),
(2, [Pt(*i) for i in [(1, 3), (0, 2), (4, 5), (2, 4), (0, 3)]]),
(3, [Pt(*i) for i in [(-5, 5), (-4, 4), (3, 2), (1, -1), (-3, 2), (4, -2), (6, -6)]])]:
n, p = len(points), points
# All circles between two points (which can both be the same point)
circles = set(sum([[c1, c2]
for c1, c2 in [circles_from_p1p2r(p1, p2, r) for p1, p2 in product(p, p)]
if c1 is not None], []))
# points covered by each circle
coverage = {c: {pt for pt in points if covers(c, pt)}
for c in circles}
# Ignore all but one of circles covering points covered in whole by other circles
#print('\nwas considering %i circles' % len(coverage))
items = sorted(coverage.items(), key=lambda keyval:len(keyval[1]))
for i, (ci, coveri) in enumerate(items):
for j in range(i+1, len(items)):
cj, coverj = items[j]
if not coverj - coveri:
coverage[cj] = {}
coverage = {key: val for key, val in coverage.items() if val}
#print('Reduced to %i circles for consideration' % len(coverage))
# Greedy coverage choice
chosen, covered = [], set()
while len(covered) < n:
_, nxt_circle, nxt_cov = max((len(pts - covered), c, pts)
for c, pts in coverage.items())
delta = nxt_cov - covered
covered |= nxt_cov
chosen.append([nxt_circle, delta])
# Output
print('\n%i points' % n)
pp(points)
print('A minimum of circles of radius %g to cover the points (And the extra points they covered)' % r)
pp(chosen)
显示三个运行的输出是：
5 points
[Pt(x=1, y=3), Pt(x=0, y=2), Pt(x=4, y=5), Pt(x=2, y=4), Pt(x=0, y=3)]
A minimum of circles of radius 3 to cover the points (And the extra points they covered)
[[Cir(x=2.958039891549808, y=2.5, r=3),
{Pt(x=4, y=5), Pt(x=0, y=3), Pt(x=1, y=3), Pt(x=0, y=2), Pt(x=2, y=4)}]]
5 points
[Pt(x=1, y=3), Pt(x=0, y=2), Pt(x=4, y=5), Pt(x=2, y=4), Pt(x=0, y=3)]
A minimum of circles of radius 2 to cover the points (And the extra points they covered)
[[Cir(x=1.9364916731037085, y=2.5, r=2),
{Pt(x=0, y=3), Pt(x=1, y=3), Pt(x=0, y=2), Pt(x=2, y=4)}],
[Cir(x=4, y=5, r=2), {Pt(x=4, y=5)}]]
7 points
[Pt(x=-5, y=5),
Pt(x=-4, y=4),
Pt(x=3, y=2),
Pt(x=1, y=-1),
Pt(x=-3, y=2),
Pt(x=4, y=-2),
Pt(x=6, y=-6)]
A minimum of circles of radius 3 to cover the points (And the extra points they covered)
[[Cir(x=3.9951865152835286, y=-0.8301243435223524, r=3),
{Pt(x=3, y=2), Pt(x=1, y=-1), Pt(x=4, y=-2)}],
[Cir(x=-2.0048134847164714, y=4.830124343522352, r=3),
{Pt(x=-4, y=4), Pt(x=-3, y=2), Pt(x=-5, y=5)}],
[Cir(x=6.7888543819998315, y=-3.1055728090000843, r=3), {Pt(x=6, y=-6)}]]