蜂后所在单元格的中心为原点,可以看出各单位格的中心坐标是$(3a,\sqrt 3b)$或$(3a+\frac{3}{2},\sqrt 3b+\frac{\sqrt 3}{2})$。
前者到原点的距离是$\sqrt{9a^2+3b^2}$,后者到原点的距离是$\sqrt{9a^2+3b^2+9a+3b+3}$。
即$3a^2+b^2=L^2/3$或 $3a^2+3a+b^2+b+1=L^2/3$,将后者两侧同乘以4可得
$$3(2a+1)^2+(2b+1)^2=\frac{4L^2}{3}.$$
之后可以运用第143题等题中提及的有理点的方法来分析,不过此处直接使用此文的结论有
$$B(L)=6d(\frac{L^2}{3},1,3)$$
其中$d(L^2/3,1,3)$表示$L^2/3$的模3余1的约数的数量。
而$450=6\times 5^2\times3$,因此只需考虑
$$\frac{L^2}{3}=p_1^{a_1}p_2^{a_2}p_3^{a_3}m,$$
其中$m$是没有模3余1型素因子的数,$p_1,p_2,p_3$都是模3余1型的素数 且
$$(a_1+1)(a_2+1)(a_3+1)=75,$$
(注意$a_1,a_2,a_3$可以为0)。
因此先用筛法枚举素数并找出模3余1型的素数列表,在此表的基础上再用筛法和前缀和对每个n找出不超过n且没有模3余1型素因子的数的数量,最后用深度优先搜索$p_1,p_2,p_3$和$a_1,a_2,a_3$、以及m的组合。
最终结果是$58065134$。
注:各单元格的中心点也可以在复平面上表示为$a+b\omega$,其中$\omega=e^{2\pi i/3}$是三阶单位根,$a,b\in\mathbb Z$。如此则本问题转化Eisenstein整数的分解,即求$\frac{4L^2}{3}=z \bar z$,其中$z\in\mathbb Z(\omega)$是Eisenstein整数。
注2:因用到math.prod函数,以下代码为Python 3。
import math
L=5*(10**11)
target=int(L/(7*7*13*13*math.sqrt(3)))
d=[0]*target
n=2
while n < target:
for i in range(n*n,target,n):
d[i]=d[i]+1
n=n+1
while n < target and d[n]>0:
n=n+1
primes1=[k for k in range(2,target) if d[k]==0 and k%3==1]
bound=int(L/(7*7*13*13*19*math.sqrt(3)))
d=[0]+[1]*bound
for p in primes1:
for q in range(p,bound+1,p):
d[q]=0
for i in range(1,bound+1):
d[i]+=d[i-1]
def g(a,k=0,m=1):
if a:
r=0
i=k
while i <= len(primes1)-len(a) and m*math.prod(primes1[i+x]**y for x,y in enumerate(a)) <= L:
r+=g(a[1:],i+1,m*(primes1[i]**a[0]))
i+=1
return r
else:
return d[L//(m*3)]+d[int(L/(m*math.sqrt(3)))]
f=lambda n,a:g(a) if n==1 else sum(f(n//i,a+[i//2]) for i in [3,5,15,25,75] if n%i==0)
print(f(75,[]))
|