无平方因子数

不能被除了1以外的任何完全平方数整除的数称为无平方因子数，例如1, 2, 3, 5, 6, 7, 10, 11都是无平方因子数，而4, 8, 9, 12都不是。

在小于

2^{50}

的数中，有多少个是无平方因子数？

本题难度:

解答

显然素数都是无平方因子数，因此直接统计无平方因子数是不可取的，因为

2^{50}

以内的素数总数就不容易有效地计算，从而本题的主要策略是利用容斥原理计算有平方因子数的总数。

设p是不超过

2^{25}

的一个素数，那么在

2^{50}

以内共有

a_{p} = ⌊ 2^{50} / p^{2} ⌋

个数是

p^{2}

的倍数，因此有平方因子数的总数不超过

\sum_{p < 2^{25}} a_{p}

。

在上述和式中既是p的倍数又是q的倍数（

p, q

是不同的素数）的数计算了两次，因此考虑不超过

2^{25}

，且可以写成两个互异素数乘积的数

m = p q

，需要从上一段的和式中排除

\sum_{m < 2^{25}} ⌊ 2^{50} / m^{2} ⌋

。

接下去再用同样的方式考虑三个互异素数的乘积、四个互异素数的乘积。。。以此类推。

为方便书写，引入以下莫比乌斯函数：

μ (n) = {\begin{cases} 1 & n = 1, \\ 0 & n 不是无平方因子数, \\ (- 1)^{k} & n 是 k 个互异素数的乘积, \end{cases}

则按上述分析，用筛法生成莫比乌斯函数的值并计算即得结果为

\sum_{m = 1}^{2^{25}} ⌊ \frac{2^{50}}{m^{2}} ⌋ \cdot μ (m) = 684465067343069.

注：实现时要注意先作整除再乘负数。

target=(1 << 25)+1
d=[i for i in range(target)]
mob=[0 for i in range(target)]
mob[1]=1
n=2
while n < target:
    for i in range(n,target,n):
        if d[i]>1:
            if i%(n*n)==0:
                d[i]=1
                mob[i]=0
            else:
                d[i]/=n
                mob[i]+=1
        if d[i]==1 and mob[i]>0:
            mob[i]=-1 if mob[i]%2==1 else 1
    while n < target and d[n] < n:
        n+=1

print sum(mob[n]*((1 << 50)/(n*n)) for n in range(1,(1 << 25)+1))