MATLAB
实现:https://blog.csdn.net/qq_41845823/article/details/120324970
用至少多少个半径为r的小圆覆盖半径为R的大圆
是的任意两个小圆圆心距离为根号3r,且任意三个圆之间交于一点
python实现:
import numpy as np
import matplotlib.pyplot as plt
def circle1(p,r,theta):
"""
p为列表,原点
r为半径
"""
X = p[0] + r * np.cos(theta)
Y = p[1] + r * np.sin(theta)
plt.plot(X, Y)
# ==========================================
# 圆的基本信息
# 1.圆半径
R = 1000
r = 60.5
test_r1 = r * np.sqrt(3)
theta = np.arange(0, 2*np.pi, 0.01)
# 2.圆心坐标
a, b = (0., 0.)
# ==========================================
X = a + R * np.cos(theta)
Y = b + R * np.sin(theta)
plt.plot(X, Y)
X = a + r * np.cos(theta)
Y = b + r * np.sin(theta)
plt.plot(X, Y)
p = []
clo = 0
clo_old = 0
clo_a = 0
# 奇数列
for j in range(1, 12, 2):
clo += 2
clo_a = 0
for i in range(1, 20, 2):
p.append([j*r*3/2, i*test_r1/2])
p.append([j*r*3/2, -i*test_r1/2])
p.append([-j*r * 3 / 2, i * test_r1 / 2])
p.append([-j*r * 3 / 2, -i * test_r1 / 2])
clo_a +=2
# 偶数列
for i in range(0, 12, 2):
clo += 2
clo_old = 0
for j in range(11):
p.append([i*r*3/2, j*test_r1])
p.append([i*r*3/2, -j*test_r1])
p.append([-i * r * 3 / 2, j * test_r1])
p.append([-i * r * 3 / 2, -j * test_r1])
clo_old += 2
print("总列数:%d" % (clo-1))
print("偶数列:%d" % (clo_old-1))
print("奇数列:%d" % clo_a)
for p1 in p:
circle1(p1, r, theta)
# ==========================================
plt.show()
最终结果为:
修改一下,去掉多余小圆
import numpy as np
import matplotlib.pyplot as plt
def circle1(p,r,theta):
X = p[0] + r * np.cos(theta)
Y = p[1] + r * np.sin(theta)
plt.plot(X, Y)
# ==========================================
# 圆的基本信息
# 1.圆半径
R = 1000
r = 60.5
test_r1 = r * np.sqrt(3)
theta = np.arange(0, 2*np.pi, 0.01)
# 2.圆心坐标
a, b = (0., 0.)
# ==========================================
X = a + R * np.cos(theta)
Y = b + R * np.sin(theta)
plt.plot(X, Y)
X = a + r * np.cos(theta)
Y = b + r * np.sin(theta)
plt.plot(X, Y)
p = []
clo = 0
clo_old = 0
clo_a = 0
# 奇数列
for j in range(1, 12, 2):
clo += 2
clo_a = 0
for i in range(1, 20, 2):
p.append([j*r*3/2, i*test_r1/2])
p.append([j*r*3/2, -i*test_r1/2])
p.append([-j*r * 3 / 2, i * test_r1 / 2])
p.append([-j*r * 3 / 2, -i * test_r1 / 2])
clo_a +=2
# 偶数列
for i in range(0, 12, 2):
clo += 2
clo_old = 0
for j in range(11):
p.append([i*r*3/2, j*test_r1])
p.append([i*r*3/2, -j*test_r1])
p.append([-i * r * 3 / 2, j * test_r1])
p.append([-i * r * 3 / 2, -j * test_r1])
clo_old += 2
count = 0
for p1 in p:
if p1[0]*p1[0]+p1[1]*p1[1] <= 1000000:
count += 1
circle1(p1, r, theta)
elif p1[0]*p1[0]+p1[1]*p1[1] <= (1000+60.5)*(1000+60.5):
count += 1
circle1(p1,r,theta)
# ==========================================
plt.show()
结果如图:
综合各个元素最终修改版
import numpy as np
import matplotlib.pyplot as plt
import math
# 求两圆交点
def insec(p1, r1, p2, r2):
x = p1[0]
y = p1[1]
R = r1
a = p2[0]
b = p2[1]
S = r2
c = (a - x) ** 2 + (b - y) ** 2
d = np.sqrt(c)
if (d > (R + S)) or d < (abs(R - S)):
# print("Two circles have no intersection")
return None, None
elif d.all == 0:
# print("Two circles have same center!")
return None, None
else:
A = (R ** 2 - S ** 2 + d ** 2) / (2 * d)
h = np.sqrt(R ** 2 - A ** 2)
x2 = x + A * (a - x) / d
y2 = y + A * (b - y) / d
x3 = round(x2 - h * (b - y) / d, 2)
y3 = round(y2 + h * (a - x) / d, 2)
x4 = round(x2 + h * (b - y) / d, 2)
y4 = round(y2 - h * (a - x) / d, 2)
c1 = [x3, y3]
c2 = [x4, y4]
return c1, c2
# 判断一点是否在(0,0)大圆内
def in_circle(p, r):
n = p[0] ** 2 + p[1] ** 2
if n < r ** 2:
return True
else:
return False
# 画圆
def circle1(p, r, theta):
X = p[0] + r * np.cos(theta)
Y = p[1] + r * np.sin(theta)
plt.plot(X, Y)
# 判断重复画第几象限的圆
def re_circle(p, r, theta):
if p[0] == 0:
# (0,0) 位置画一个
if p[1] == 0:
circle1([p[0], p[1]], r, theta)
# print("(%.3f,%.3f)" % (p[0], p[1]))
return 1
# 0列之重复画上下两个
else:
circle1([p[0], p[1]], r, theta)
# print("(%.3f,%.3f)" % (p[0], p[1]))
circle1([p[0], -p[1]], r, theta)
# print("(%.3f,%.3f)" % (p[0], -p[1]))
return 2
# y=0的小圆只需要画x周上的两个
elif y1[i] == 0:
circle1([p[0], p[1]], r, theta)
# print("(%.3f,%.3f)" % (p[0], p[1]))
circle1([-p[0], p[1]], r, theta)
# print("(%.3f,%.3f)" % (-p[0], p[1]))
return 2
# 第一象限内的需要画四个
else:
circle1([p[0], p[1]], r, theta)
# print("(%.3f,%.3f)" % (p[0], p[1]))
circle1([p[0], -p[1]], r, theta)
# print("(%.3f,%.3f)" % (p[0], -p[1]))
circle1([-p[0], p[1]], r, theta)
# print("(%.3f,%.3f)" % (-p[0], p[1]))
circle1([-p[0], -p[1]], r, theta)
# print("(%.3f,%.3f)" % (-p[0], -p[1]))
return 4
# ==========================================
# 圆的基本信息
# 1.圆半径
R = 1000
r = 30.5
test_r1 = r * np.sqrt(3)
theta = np.arange(0, 2 * np.pi, 0.01)
n = math.ceil(R/(3*r/2))+1
# 2.圆心坐标
# ==========================================
# 画大圆
fig = plt.figure()
axes = fig.add_subplot()
axes.axis('equal')
X = R * np.cos(theta)
Y = R * np.sin(theta)
plt.plot(X, Y)
# 每列x,y坐标
x1 = [0.] * n
y1 = [0.] * n
N = [0] *n
# 存取每一列初始圆的圆心
for i in range(n):
x1[i] = i * r * 3 / 2
if i % 2:
y1[i] = test_r1 / 2
else:
y1[i] = 0
# count存放小圆数量
count = 0
for i in range(n):
c = insec([x1[i], y1[i]], r, [0, 0], R)
# 当x坐标表示的小圆和大圆没有交点但圆心在大圆内,或者小圆与大圆有焦点的时候则考虑该列
# 否则从该列开始就不需要进行循环
if (c == (None, None) and x1[i]**2 + y1[i]**2 < 1000**2) or c != (None, None):
for j in range(n):
c = insec([x1[i], y1[i]], r, [0, 0], R)
if c == (None, None):
if x1[i]**2 +y1[i]**2 < 1000**2:
N[i] += re_circle([x1[i], y1[i]], r, theta)
else:
count += N[i]
if i != 0:
N[i] /= 2
y1[i] = temp
break
else:
# # 求(x1[i],y1[i])圆与下面的一个圆的交点
# c = insec([x1[i], y1[i]], r, [x1[i], ], R)
N[i] += re_circle([x1[i], y1[i]], r, theta)
temp = y1[i]
y1[i] += test_r1
else:
break
# ==========================================
plt.axhline(0, -1111, 1111)
plt.axvline(0, -1111, 1111)
print(count)
# for i in range(len(x1)):
# print("(%.3f,%.3f)" % (x1[i], y1[i]))
# for i in range(len(N)):
# print("%d" % int(N[i]))
plt.show()
结果:
欢迎分享,转载请注明来源:内存溢出
评论列表(0条)