题目描述
一位科学家正在尝试制造一种非常大的晶体,具体来说是一种大的碳晶体。他认为,既然钻石是碳的晶体并且非常珍贵,那么从长远来看,他的新碳晶体也会像钻石一样珍贵。他晶体中的原子无法自然结合在一起,因此他希望在晶体中心施加一个强大的力,吸引所有碳原子并使它们保持在一起。
钻石晶体中的碳原子可以看作是放置在一个立方体中。科学家也希望将他的晶体碳原子放置在一个N×N×NN \times N \times NN×N×N的立方体中,其中NNN为偶数。如果该立方体的中心是(0,0,0)(0, 0, 0)(0,0,0)且所有边均平行于xyxyxy、yzyzyz或xzxzxz平面,则所有原子都将放置在三维整数坐标上。因此,如果(x,y,z)(x, y, z)(x,y,z)是N×N×NN \times N \times NN×N×N立方体中某个原子的坐标,则xxx、yyy和zzz为整数且满足(−N/2≤x,y,z≤N/2)(-N/2 \le x, y, z \le N/2)(−N/2≤x,y,z≤N/2)。由于中心的强大引力会吸引所有原子,因此原子的放置方式需确保没有任何原子位于中心和另一个原子之间。例如,如果坐标(2,2,2)(2, 2, 2)(2,2,2)处有一个原子,则不应在坐标(1,1,1)(1, 1, 1)(1,1,1)处放置原子,因为(1,1,1)(1, 1, 1)(1,1,1)处的原子会阻挡(2,2,2)(2, 2, 2)(2,2,2)与中心之间的引力。同样地,如果(1,1,1)(1, 1, 1)(1,1,1)处有原子,则(2,2,2)(2, 2, 2)(2,2,2)处不应有原子。
给定要制造晶体的立方体尺寸(边长),你的任务是找出在上述约束下可以放置的最大原子数。
输入格式
输入文件最多包含303030行。每行包含一个偶数整数NNN(0<N≤2000000 < N \le 2000000<N≤200000),表示科学家计划放置原子的立方体的边长。输入以一行NNN的值为000结束。
输出格式
对于除最后一行外的每一行输入,输出一行。该行应包含输出序号(格式为Crystal i:)和一个整数,表示可以放置的最大原子数。
样例输入
4 2 0样例输出
Crystal 1: 98 Crystal 2: 26数学基础:莫比乌斯函数与莫比乌斯反演
一、莫比乌斯函数
莫比乌斯函数μ(n)\mu(n)μ(n)是定义在正整数上的函数:
μ(n)={1如果 n=1(−1)k如果 n 是 k 个不同质数的乘积0如果 n 被一个质数的平方整除 \mu(n) = \begin{cases} 1 & \text{如果 } n = 1 \\ (-1)^k & \text{如果 } n \text{ 是 } k \text{ 个不同质数的乘积} \\ 0 & \text{如果 } n \text{ 被一个质数的平方整除} \end{cases}μ(n)=⎩⎨⎧1(−1)k0如果n=1如果n是k个不同质数的乘积如果n被一个质数的平方整除
重要性质:
- 积性函数:如果gcd(a,b)=1\gcd(a, b) = 1gcd(a,b)=1,则μ(ab)=μ(a)μ(b)\mu(ab) = \mu(a)\mu(b)μ(ab)=μ(a)μ(b)
- 求和性质:
∑d∣nμ(d)={1如果 n=10如果 n>1 \sum_{d \mid n} \mu(d) = \begin{cases} 1 & \text{如果 } n = 1 \\ 0 & \text{如果 } n > 1 \end{cases}d∣n∑μ(d)={10如果n=1如果n>1
二、莫比乌斯反演
设f(n)f(n)f(n)和g(n)g(n)g(n)是定义在正整数上的两个函数。
第一形式(约数和形式):
如果f(n)=∑d∣ng(d)f(n) = \sum_{d \mid n} g(d)f(n)=∑d∣ng(d),那么g(n)=∑d∣nμ(d)f(nd)g(n) = \sum_{d \mid n} \mu(d) f\left(\frac{n}{d}\right)g(n)=∑d∣nμ(d)f(dn)
第二形式(倍数和形式):
如果f(n)=∑n∣dg(d)f(n) = \sum_{n \mid d} g(d)f(n)=∑n∣dg(d),那么g(n)=∑n∣dμ(dn)f(d)g(n) = \sum_{n \mid d} \mu\left(\frac{d}{n}\right) f(d)g(n)=∑n∣dμ(nd)f(d)
莫比乌斯反演可以看作是"容斥原理"的数论形式。当我们知道一个"包含所有因子"的函数f(n)f(n)f(n)时,可以用莫比乌斯函数"筛出"我们真正关心的函数g(n)g(n)g(n)。
题目分析与解题思路
1. 问题转化
题目要求在一个边长为偶数NNN的立方体网格中放置尽可能多的点(原子),使得从原点(0,0,0)(0,0,0)(0,0,0)(即立方体中心)到任意一个被放置的点的线段上,没有其他被放置的点(整数点)存在。
换句话说,所有被放置的点必须是从原点可见的,即从原点到该点的线段上不存在其他整数点(除了端点)。
在三维整数网格中,一个点(x,y,z)(x, y, z)(x,y,z)从原点可见的充要条件是:
gcd(∣x∣,∣y∣,∣z∣)=1 \gcd(|x|, |y|, |z|) = 1gcd(∣x∣,∣y∣,∣z∣)=1
理由:如果gcd(∣x∣,∣y∣,∣z∣)=d>1\gcd(|x|,|y|,|z|) = d > 1gcd(∣x∣,∣y∣,∣z∣)=d>1,那么点(xd,yd,zd)(\frac{x}{d}, \frac{y}{d}, \frac{z}{d})(dx,dy,dz)也在该线段上,并且更靠近原点,从而原点与该点之间存在其他整数点,违反了规则。
因此,我们的目标是:在坐标范围[−M,M][-M, M][−M,M]内(其中M=N/2M = N/2M=N/2),统计所有满足gcd(∣x∣,∣y∣,∣z∣)=1\gcd(|x|,|y|,|z|) = 1gcd(∣x∣,∣y∣,∣z∣)=1的整数点(x,y,z)(x, y, z)(x,y,z)的数量。注意,原点(0,0,0)(0,0,0)(0,0,0)的gcd\gcdgcd定义为000,我们需要排除原点。
2. 数学建模与推导
设M=N/2M = N/2M=N/2,坐标范围[−M,M][-M, M][−M,M]。我们定义:
- g(d)g(d)g(d)= 满足gcd(∣x∣,∣y∣,∣z∣)=d\gcd(|x|,|y|,|z|) = dgcd(∣x∣,∣y∣,∣z∣)=d的点数(d≥1d \ge 1d≥1)
- f(d)f(d)f(d)= 满足d∣gcd(∣x∣,∣y∣,∣z∣)d \mid \gcd(|x|,|y|,|z|)d∣gcd(∣x∣,∣y∣,∣z∣)的点数(即三个坐标都是ddd的倍数)
显然有:
f(d)=∑k≥1g(k⋅d) f(d) = \sum_{k \ge 1} g(k \cdot d)f(d)=k≥1∑g(k⋅d)
这是因为如果gcd\gcdgcd是ddd的倍数,那么它可以是d,2d,3d,…d, 2d, 3d, \dotsd,2d,3d,…。
使用第二形式的莫比乌斯反演,令n=1n = 1n=1:
g(1)=∑d≥1μ(d)f(d) g(1) = \sum_{d \ge 1} \mu(d) f(d)g(1)=d≥1∑μ(d)f(d)
这里g(1)g(1)g(1)就是我们需要的可见点数(gcd=1\gcd = 1gcd=1)。
3. 计算f(d)f(d)f(d)
对于给定的ddd,f(d)f(d)f(d)是三个坐标都是ddd的倍数的点数。
在范围[−M,M][-M, M][−M,M]中,xxx是ddd的倍数的值有:
- 000
- ±d,±2d,…,±kd\pm d, \pm 2d, \dots, \pm kd±d,±2d,…,±kd,其中k=⌊M/d⌋k = \lfloor M/d \rfloork=⌊M/d⌋
所以每个坐标有2k+12k + 12k+1个可能值。三个坐标独立,总共有(2k+1)3(2k + 1)^3(2k+1)3种组合。
但原点(0,0,0)(0,0,0)(0,0,0)的gcd\gcdgcd是000,不属于gcd=d≥1\gcd = d \ge 1gcd=d≥1的情况,所以排除原点:
f(d)=(2⋅⌊M/d⌋+1)3−1 f(d) = (2 \cdot \lfloor M/d \rfloor + 1)^3 - 1f(d)=(2⋅⌊M/d⌋+1)3−1
4. 最终公式
将f(d)f(d)f(d)代入反演公式,得到:
可见点数=∑d=1Mμ(d)⋅[(2⋅⌊M/d⌋+1)3−1] \text{可见点数} = \sum_{d=1}^{M} \mu(d) \cdot \left[ (2 \cdot \lfloor M/d \rfloor + 1)^3 - 1 \right]可见点数=d=1∑Mμ(d)⋅[(2⋅⌊M/d⌋+1)3−1]
其中M=N/2M = N/2M=N/2。
5. 验证样例
对于N=4N = 4N=4(M=2M = 2M=2):
- d=1d = 1d=1:μ(1)=1\mu(1) = 1μ(1)=1,⌊2/1⌋=2\lfloor 2/1 \rfloor = 2⌊2/1⌋=2,(2⋅2+1)3−1=53−1=124(2 \cdot 2 + 1)^3 - 1 = 5^3 - 1 = 124(2⋅2+1)3−1=53−1=124
- d=2d = 2d=2:μ(2)=−1\mu(2) = -1μ(2)=−1,⌊2/2⌋=1\lfloor 2/2 \rfloor = 1⌊2/2⌋=1,(2⋅1+1)3−1=33−1=26(2 \cdot 1 + 1)^3 - 1 = 3^3 - 1 = 26(2⋅1+1)3−1=33−1=26
- d=3d = 3d=3:μ(3)=−1\mu(3) = -1μ(3)=−1,⌊2/3⌋=0\lfloor 2/3 \rfloor = 0⌊2/3⌋=0,(2⋅0+1)3−1=0(2 \cdot 0 + 1)^3 - 1 = 0(2⋅0+1)3−1=0
- d=4d = 4d=4:μ(4)=0\mu(4) = 0μ(4)=0,⌊2/4⌋=0\lfloor 2/4 \rfloor = 0⌊2/4⌋=0,(2⋅0+1)3−1=0(2 \cdot 0 + 1)^3 - 1 = 0(2⋅0+1)3−1=0
总和:124−26=98124 - 26 = 98124−26=98,与样例输出一致。
对于N=2N = 2N=2(M=1M = 1M=1):
- d=1d = 1d=1:μ(1)=1\mu(1) = 1μ(1)=1,⌊1/1⌋=1\lfloor 1/1 \rfloor = 1⌊1/1⌋=1,(2⋅1+1)3−1=33−1=26(2 \cdot 1 + 1)^3 - 1 = 3^3 - 1 = 26(2⋅1+1)3−1=33−1=26
- d=2d = 2d=2:μ(2)=−1\mu(2) = -1μ(2)=−1,⌊1/2⌋=0\lfloor 1/2 \rfloor = 0⌊1/2⌋=0,(2⋅0+1)3−1=0(2 \cdot 0 + 1)^3 - 1 = 0(2⋅0+1)3−1=0
总和:262626,与样例输出一致。
算法设计与实现
1. 算法步骤
- 预处理莫比乌斯函数:使用线性筛法计算μ(1)\mu(1)μ(1)到μ(Mmax)\mu(M_{\max})μ(Mmax),其中Mmax=100000M_{\max} = 100000Mmax=100000(因为N≤200000N \le 200000N≤200000,所以M≤100000M \le 100000M≤100000)。
- 处理每个查询:
- 读入NNN(偶数),计算M=N/2M = N/2M=N/2。
- 初始化答案ans=0ans = 0ans=0。
- 对ddd从111到MMM循环,累加μ(d)⋅[(2⋅⌊M/d⌋+1)3−1]\mu(d) \cdot \left[ (2 \cdot \lfloor M/d \rfloor + 1)^3 - 1 \right]μ(d)⋅[(2⋅⌊M/d⌋+1)3−1]。
- 输出
Crystal i: ans。
2. 时间复杂度分析
- 预处理莫比乌斯函数:O(Mmax)O(M_{\max})O(Mmax),其中Mmax=100000M_{\max} = 100000Mmax=100000。
- 每个查询:O(M)O(M)O(M),其中M≤100000M \le 100000M≤100000。
- 总操作数:最多303030个查询,约3×1063 \times 10^63×106次运算,在合理范围内。
3. 空间复杂度分析
需要存储莫比乌斯函数数组,大小为O(Mmax)O(M_{\max})O(Mmax),即100001100001100001个整数,空间充足。
代码实现
// Make a Crystal// UVa ID: 11014// Verdict: Accepted// Submission Date: 2025-12-17// UVa Run Time: 0.000s//// 版权所有(C)2025,邱秋。metaphysis # yeah dot net#include<bits/stdc++.h>usingnamespacestd;typedeflonglongLL;constintMAXM=100000;intmu[MAXM+5];boolisPrime[MAXM+5];vector<int>primes;// 预处理莫比乌斯函数:使用线性筛法voidsieve(){fill(isPrime,isPrime+MAXM+1,true);mu[1]=1;for(inti=2;i<=MAXM;++i){if(isPrime[i]){primes.push_back(i);mu[i]=-1;// 质数的莫比乌斯函数值为 -1}for(intp:primes){if(i*p>MAXM)break;isPrime[i*p]=false;if(i%p==0){mu[i*p]=0;// 有平方因子break;}else{mu[i*p]=-mu[i];// 积性函数性质}}}}intmain(){sieve();intcaseNo=1;intN;while(scanf("%d",&N)==1&&N!=0){LL M=N/2;LL ans=0;for(LL d=1;d<=M;++d){LL t=M/d;// floor(M/d)LL term=(2*t+1);term=term*term*term-1;// (2t+1)^3 - 1ans+=mu[d]*term;}printf("Crystal %d: %lld\n",caseNo++,ans);}return0;}代码说明
线性筛法计算莫比乌斯函数:
- 初始化所有数为质数,μ(1)=1\mu(1) = 1μ(1)=1。
- 遍历iii从222到MAXMMAXMMAXM:
- 如果iii是质数,加入质数表,μ(i)=−1\mu(i) = -1μ(i)=−1。
- 用当前质数表筛去合数i×pi \times pi×p:
- 如果iii能被ppp整除,则i×pi \times pi×p有平方因子,μ(i×p)=0\mu(i \times p) = 0μ(i×p)=0。
- 否则,μ(i×p)=−μ(i)\mu(i \times p) = -\mu(i)μ(i×p)=−μ(i)(积性函数性质)。
主循环:
- 读取每个NNN,直到N=0N = 0N=0结束。
- 计算M=N/2M = N/2M=N/2。
- 对ddd从111到MMM累加贡献。
- 输出结果,注意使用
%lld格式输出long long类型。
注意事项:
- 使用
long long类型存储中间结果和答案,避免溢出。 - ddd循环上界为MMM,因为当d>Md > Md>M时,⌊M/d⌋=0\lfloor M/d \rfloor = 0⌊M/d⌋=0,贡献为000。
- 使用
算法优化思考
虽然当前算法已经可以通过题目测试,但还可以进一步优化:
整除分块优化:计算∑d=1Mμ(d)⋅[(2⋅⌊M/d⌋+1)3−1]\sum_{d=1}^{M} \mu(d) \cdot \left[ (2 \cdot \lfloor M/d \rfloor + 1)^3 - 1 \right]∑d=1Mμ(d)⋅[(2⋅⌊M/d⌋+1)3−1]时,⌊M/d⌋\lfloor M/d \rfloor⌊M/d⌋的值在连续区间内相同,可以使用整除分块将复杂度从O(M)O(M)O(M)降为O(M)O(\sqrt{M})O(M)。
预处理前缀和:可以预处理莫比乌斯函数的前缀和,结合整除分块进一步优化。
记忆化:对于重复的MMM值,可以缓存计算结果。
但对于本题M≤100000M \le 100000M≤100000且最多303030个查询的情况,当前O(M)O(M)O(M)算法已经足够高效。
总结
本题的核心在于将几何约束转化为数论条件:从原点可见的点等价于gcd(∣x∣,∣y∣,∣z∣)=1\gcd(|x|,|y|,|z|) = 1gcd(∣x∣,∣y∣,∣z∣)=1。通过引入莫比乌斯函数和莫比乌斯反演,我们避免了复杂的容斥计数,得到了简洁高效的数学公式。算法实现主要分为两部分:
- 预处理莫比乌斯函数(线性筛法)
- 对每个查询计算和式
掌握莫比乌斯反演这一工具,对于解决类似的数论计数问题非常有帮助,它能够将复杂的容斥过程转化为简洁的数学表达式,是算法竞赛中处理gcd\gcdgcd相关计数问题的利器。