我的数学 SE 问题确定一对旋转六边形格子中的重合点是否最接近原点?解释了我遇到的问题。我不会在这里详细再现整个事情,但这张图片以图形方式说明了问题。
- 右侧:定义闭合点 (5, 6) 和开圆点 (5, 3) 的重合点(见红色箭头)是该重合点阵中距原点最近的六个点之一。
- 左侧:由 (6, 5) 和 (5, 3) (红色箭头)类似定义的重合点远不是离原点最近的重合点。
唯一的区别是我们通过将闭合点阵旋转大约 6 度将 (6, 5) 替换为 (5, 6)。
我的问题的答案建议将众所周知的欧几里得算法的修改应用于该问题。
它现在应用于两个复数Eisenstein 整数,而不是两个实整数,这是复平面中点的六边形分布。有关更多信息,请参阅:
为了尝试从那个答案到 Python 中实现这个非常非常好的和清晰描述的分步方法,我定义了一个Point
类并为该类提供了以下方法:.subtract(P)
, .multiply(P),
.divround(P) ,
.prod(P) and
.sub(P)` 尽我所能复制该答案中描述的内容。
请注意,一个Point
对象P
由两个整数 P.p
定义,P.q
并且所有方法都对这些整数进行操作并返回Point
. 大多数操作是加法和乘法,唯一暂时进入浮点数的是.divround(P)
.
下面是一个包含 300 多个 (p, q), (k, l) 的测试程序,由函数测试euclid_2021(p, q, k, l)
。
到目前为止,不太严格的测试表明该算法给出了我认为正确的答案,至少它似乎与我之前编写的一个非常长、丑陋、看起来很疯狂的脚本非常吻合,但我没有已知的正确方法来检查两者。
但是在这里,我只询问一个明显的无限循环,该方法仅针对以下 300 多个案例中的一个陷入困境。
Group_1
并Group_3
快速运行并返回表示它是否最接近的布尔值。
Group_2
只要我省略第二种情况(9、2、2、6),一切都很好。但是对于这种情况,循环会无限期地运行。
关于为什么这个案例失败的任何数学想法?
当仅测试“坏”情况(9、2、2、6)并插入以下打印行时:
def euclid_2021_with_print(p, q, k, l):
A = Point(p, q)
B = Point(k, l)
if A.abs < B.abs:
A, B = B, A
while B.nonzero:
D = A.divround(B)
B, A = A.subtract(D.multiply(B)), B
print((D.p, D.q), (A.p, A.q), (B.p, B.q)) ### Print inside while loop
return A.isone
该算法只是在两种状态之间来回循环:
(0, 0) (4, -3) (-2, 3)
(0, 0) (-2, 3) (4, -3)
(0, 0) (4, -3) (-2, 3)
(0, 0) (-2, 3) (4, -3)
(0, 0) (4, -3) (-2, 3)
(0, 0) (-2, 3) (4, -3)
第 1、2、3 组的真/假输出,但跳过“坏”情况:
import numpy as np
import matplotlib.pyplot as plt
class Point():
ones = ((1, 0), (-1, 0), (0, 1), (0, -1), (-1, 1), (1, -1))
u = 0.5 * (1 + np.sqrt(3) * 1j)
def __init__(self, p, q):
self.p = int(p)
self.q = int(q)
self.pq = (self.p, self.q)
self.nonzero = not ((self.p == 0) and (self.q == 0))
self.isone = self.pq in self.ones
self.abs = np.sqrt(self.p**2 + self.q**2)
self.xy = self.p + self.u * self.q
self.x, self.y = np.real(self.xy), np.imag(self.xy)
def subtract(self, P):
return Point(self.p - P.p, self.q - P.q)
def multiply(self, P):
return Point(self.p * P.p - self.q * P.q,
self.p * P.q + self.q * P.p + self.q * P.q)
def divround(self, P):
(x, y), (z, t) = self.pq, P.pq
ptop, qtop = x * (z + t) + y*t, (y*z - x*t)
bot = z**2 + z*t + t**2
pfrac, qfrac = ptop/bot, qtop/bot
pint, qint = [int(thing) for thing in (pfrac, qfrac)]
four = ((pint, qint), (pint, qint+1), (pint+1, qint), (pint+1, qint+1))
ndsqs = []
for (p, q) in four:
dp, dq = pfrac - p, qfrac - q
ndsqs.append( -1 * (dp**2 + dp * dq + dq**2))
return Point(*four[np.argmax(ndsqs)])
def prod(self, B):
a, b, c, d = self.x, self.y, B.x, B.y
return Point(a*c - b*d, a*d + b*c + b*d)
def sub(self, B):
a, b, c, d = self.x, self.y, B.x, B.y
return Point(a-c, b-d)
def euclid_2021(p, q, k, l):
A = Point(p, q)
B = Point(k, l)
if A.abs < B.abs:
A, B = B, A
while B.nonzero:
D = A.divround(B)
B, A = A.subtract(D.multiply(B)), B
return A.isone
group_1 = [[4, 0, 3, 0], [4, 0, 0, 3], [0, 4, 3, 0], [0, 4, 0, 3],
[5, 0, 3, 1], [5, 0, 1, 3], [0, 5, 3, 1], [0, 5, 1, 3],
[6, 0, 4, 1], [6, 0, 1, 4], [0, 6, 4, 1], [0, 6, 1, 4],
[6, 0, 3, 2], [6, 0, 2, 3], [0, 6, 3, 2], [0, 6, 2, 3],
[7, 0, 5, 0], [7, 0, 0, 5],
[5, 3, 5, 0], [5, 3, 0, 5], [3, 5, 5, 0], [3, 5, 0, 5],
[0, 7, 5, 0], [0, 7, 0, 5], [7, 0, 4, 2], [7, 0, 2, 4],
[5, 3, 4, 2], [5, 3, 2, 4], [3, 5, 4, 2], [3, 5, 2, 4],
[0, 7, 4, 2], [0, 7, 2, 4],
[7, 0, 3, 3],
[5, 3, 3, 3], [3, 5, 3, 3],
[0, 7, 3, 3],
[8, 0, 6, 0], [8, 0, 0, 6], [0, 8, 6, 0], [0, 8, 0, 6],
[8, 0, 4, 3], [8, 0, 3, 4], [0, 8, 4, 3], [0, 8, 3, 4],
[9, 0, 6, 1], [9, 0, 1, 6], [0, 9, 6, 1], [0, 9, 1, 6],
[10, 0, 7, 1], [10, 0, 1, 7], [0, 10, 7, 1], [0, 10, 1, 7],
[10, 0, 6, 2], [10, 0, 2, 6], [0, 10, 6, 2], [0, 10, 2, 6],
[11, 0, 8, 0], [11, 0, 0, 8], [0, 11, 8, 0], [0, 11, 0, 8],
[11, 0, 7, 2], [11, 0, 2, 7], [0, 11, 7, 2], [0, 11, 2, 7],
[11, 0, 6, 3], [11, 0, 3, 6], [0, 11, 6, 3], [0, 11, 3, 6],
[12, 0, 9, 0], [12, 0, 0, 9], [0, 12, 9, 0], [0, 12, 0, 9],
[12, 0, 8, 1], [12, 0, 1, 8], [0, 12, 8, 1], [0, 12, 1, 8],
[12, 0, 8, 2], [12, 0, 2, 8], [0, 12, 8, 2], [0, 12, 2, 8],
[12, 0, 7, 3], [12, 0, 3, 7], [0, 12, 7, 3], [0, 12, 3, 7],
[12, 0, 6, 4], [12, 0, 4, 6], [0, 12, 6, 4], [0, 12, 4, 6],
[12, 0, 5, 5], [0, 12, 5, 5],
[2, 1, 2, 0], [2, 1, 0, 2], [1, 2, 2, 0], [1, 2, 0, 2],
[3, 1, 2, 1], [3, 1, 1, 2], [1, 3, 2, 1], [1, 3, 1, 2],
[4, 1, 2, 2], [1, 4, 2, 2],
[5, 1, 4, 0], [5, 1, 0, 4], [1, 5, 4, 0], [1, 5, 0, 4],
[6, 1, 5, 0], [6, 1, 0, 5], [1, 6, 5, 0], [1, 6, 0, 5],
[7, 1, 5, 1], [7, 1, 1, 5], [1, 7, 5, 1], [1, 7, 1, 5],
[8, 1, 6, 1], [8, 1, 1, 6], [1, 8, 6, 1], [1, 8, 1, 6],
[8, 1, 5, 2], [8, 1, 2, 5], [1, 8, 5, 2], [1, 8, 2, 5],
[8, 1, 4, 3], [8, 1, 3, 4], [1, 8, 4, 3], [1, 8, 3, 4],
[9, 1, 7, 0],
[9, 1, 5, 3], [9, 1, 3, 5],
[9, 1, 0, 7],
[6, 5, 7, 0],
[6, 5, 5, 3], [6, 5, 3, 5],
[6, 5, 0, 7], [5, 6, 7, 0],
[5, 6, 5, 3], [5, 6, 3, 5],
[5, 6, 0, 7],
[1, 9, 7, 0],
[1, 9, 5, 3], [1, 9, 3, 5],
[1, 9, 0, 7],
[9, 1, 6, 2], [9, 1, 2, 6],
[6, 5, 6, 2], [6, 5, 2, 6],
[5, 6, 6, 2], [5, 6, 2, 6],
[1, 9, 6, 2], [1, 9, 2, 6],
[9, 1, 4, 4],
[6, 5, 4, 4], [5, 6, 4, 4],
[1, 9, 4, 4],
[10, 1, 8, 0], [10, 1, 0, 8], [1, 10, 8, 0], [1, 10, 0, 8],
[10, 1, 7, 1], [10, 1, 1, 7], [1, 10, 7, 1], [1, 10, 1, 7],
[10, 1, 6, 3], [10, 1, 3, 6], [1, 10, 6, 3], [1, 10, 3, 6],
[10, 1, 5, 4], [10, 1, 4, 5], [1, 10, 5, 4], [1, 10, 4, 5],
[11, 1, 8, 1], [11, 1, 1, 8], [1, 11, 8, 1], [1, 11, 1, 8],
[11, 1, 6, 4], [11, 1, 4, 6], [1, 11, 6, 4], [1, 11, 4, 6],
[11, 1, 5, 5], [1, 11, 5, 5],
[2, 2, 2, 1], [2, 2, 1, 2],
[4, 2, 4, 0], [4, 2, 0, 4], [2, 4, 4, 0], [2, 4, 0, 4],
[5, 2, 4, 1], [5, 2, 1, 4], [2, 5, 4, 1], [2, 5, 1, 4],
[6, 2, 4, 2], [6, 2, 2, 4], [2, 6, 4, 2], [2, 6, 2, 4],
[6, 2, 3, 3], [2, 6, 3, 3],
[7, 2, 6, 0], [7, 2, 0, 6], [2, 7, 6, 0], [2, 7, 0, 6],
[7, 2, 5, 2], [7, 2, 2, 5], [2, 7, 5, 2], [2, 7, 2, 5],
[7, 2, 4, 3], [7, 2, 3, 4], [2, 7, 4, 3], [2, 7, 3, 4],
[8, 2, 7, 0], [8, 2, 5, 3], [8, 2, 3, 5], [8, 2, 0, 7],
[2, 8, 7, 0],
[2, 8, 5, 3], [2, 8, 3, 5],
[2, 8, 0, 7],
[8, 2, 6, 1], [8, 2, 1, 6], [2, 8, 6, 1], [2, 8, 1, 6],
[8, 2, 4, 4], [2, 8, 4, 4],
[9, 2, 7, 1], [9, 2, 1, 7], [2, 9, 7, 1], [2, 9, 1, 7]]
group_2 = [[9, 2, 6, 2], [9, 2, 2, 6], [2, 9, 6, 2], [2, 9, 2, 6]]
group_3 = [[9, 2, 5, 4], [9, 2, 4, 5], [2, 9, 5, 4], [2, 9, 4, 5],
[10, 2, 8, 0], [10, 2, 0, 8], [2, 10, 8, 0], [2, 10, 0, 8],
[10, 2, 8, 1], [10, 2, 1, 8], [2, 10, 8, 1], [2, 10, 1, 8],
[10, 2, 7, 2], [10, 2, 2, 7], [2, 10, 7, 2], [2, 10, 2, 7],
[10, 2, 6, 3], [10, 2, 3, 6], [2, 10, 6, 3], [2, 10, 3, 6],
[4, 3, 4, 1], [4, 3, 1, 4], [3, 4, 4, 1], [3, 4, 1, 4],
[4, 3, 3, 2], [4, 3, 2, 3], [3, 4, 3, 2], [3, 4, 2, 3],
[6, 3, 6, 0], [6, 3, 0, 6], [3, 6, 6, 0], [3, 6, 0, 6],
[6, 3, 4, 3], [6, 3, 3, 4], [3, 6, 4, 3], [3, 6, 3, 4],
[7, 3, 6, 1], [7, 3, 1, 6], [3, 7, 6, 1], [3, 7, 1, 6],
[8, 3, 7, 0], [8, 3, 5, 3], [8, 3, 3, 5], [8, 3, 0, 7],
[3, 8, 7, 0], [3, 8, 5, 3], [3, 8, 3, 5], [3, 8, 0, 7],
[8, 3, 7, 1], [8, 3, 1, 7], [3, 8, 7, 1], [3, 8, 1, 7],
[8, 3, 6, 2], [8, 3, 2, 6], [3, 8, 6, 2], [3, 8, 2, 6],
[9, 3, 8, 0], [9, 3, 0, 8], [3, 9, 8, 0], [3, 9, 0, 8],
[9, 3, 7, 2], [9, 3, 2, 7], [3, 9, 7, 2], [3, 9, 2, 7],
[9, 3, 6, 3], [9, 3, 3, 6], [3, 9, 6, 3], [3, 9, 3, 6],
[9, 3, 5, 4], [9, 3, 4, 5], [3, 9, 5, 4], [3, 9, 4, 5],
[4, 4, 5, 0], [4, 4, 0, 5],
[4, 4, 4, 2], [4, 4, 2, 4],
[4, 4, 3, 3],
[5, 4, 6, 0], [5, 4, 0, 6], [4, 5, 6, 0], [4, 5, 0, 6],
[5, 4, 5, 1], [5, 4, 1, 5], [4, 5, 5, 1], [4, 5, 1, 5],
[6, 4, 6, 1], [6, 4, 1, 6], [4, 6, 6, 1], [4, 6, 1, 6],
[6, 4, 5, 2], [6, 4, 2, 5], [4, 6, 5, 2], [4, 6, 2, 5],
[7, 4, 7, 0],
[7, 4, 5, 3], [7, 4, 3, 5],
[7, 4, 0, 7],
[4, 7, 7, 0],
[4, 7, 5, 3], [4, 7, 3, 5],
[4, 7, 0, 7],
[7, 4, 6, 2], [7, 4, 2, 6], [4, 7, 6, 2], [4, 7, 2, 6],
[7, 4, 4, 4], [4, 7, 4, 4],
[8, 4, 8, 0], [8, 4, 0, 8], [4, 8, 8, 0], [4, 8, 0, 8],
[8, 4, 7, 1], [8, 4, 1, 7], [4, 8, 7, 1], [4, 8, 1, 7],
[8, 4, 6, 3], [8, 4, 3, 6], [4, 8, 6, 3], [4, 8, 3, 6],
[8, 4, 5, 4], [8, 4, 4, 5], [4, 8, 5, 4], [4, 8, 4, 5],
[5, 5, 6, 1], [5, 5, 1, 6],
[5, 5, 5, 2], [5, 5, 2, 5],
[7, 5, 8, 0], [7, 5, 0, 8], [5, 7, 8, 0], [5, 7, 0, 8],
[7, 5, 7, 1], [7, 5, 1, 7], [5, 7, 7, 1], [5, 7, 1, 7],
[7, 5, 6, 3], [7, 5, 3, 6], [5, 7, 6, 3], [5, 7, 3, 6],
[7, 5, 5, 4], [7, 5, 4, 5], [5, 7, 5, 4], [5, 7, 4, 5],
[6, 6, 7, 1], [6, 6, 1, 7],
[6, 6, 6, 3], [6, 6, 3, 6],
[6, 6, 5, 4], [6, 6, 4, 5]]
G1 = [euclid_2021(*pqkl) for pqkl in group_1]
G2 = [euclid_2021(*pqkl) for pqkl in [group_2[i] for i in (0, 2, 3)]]
G3 = [euclid_2021(*pqkl) for pqkl in group_3]
plt.plot(G1)
plt.plot(G2, '-r', linewidth=6)
plt.plot(G3)
plt.ylim(-0.1, 1.1)
plt.show()