python画图:小圆覆盖大圆问题

python画图:小圆覆盖大圆问题,第1张

小圆覆盖大圆

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()

结果:

欢迎分享,转载请注明来源:内存溢出

原文地址: http://outofmemory.cn/langs/718625.html

(0)
打赏 微信扫一扫 微信扫一扫 支付宝扫一扫 支付宝扫一扫
上一篇 2022-04-25
下一篇 2022-04-25

发表评论

登录后才能评论

评论列表(0条)

保存