庆祝!第二道数论紫题。
推式子真是太有趣了!
[prod_{i=1}^n prod_{j=1}^n frac{operatorname{lcm}(i,j)}{gcd(i,j)}
]
[= prod_{i=1}^n prod_{j=1}^n frac{i imes j}{(gcd(i,j))^2}
]
[= ( prod_{i=1}^n prod_{j=1}^n i imes j ) imes ( prod_{i=1}^n prod_{j=1}^n frac{1}{(gcd(i,j))^2})
]
先看前面的式子:
[= ( prod_{i=1}^n prod_{j=1}^n i imes j)
]
[= ( prod_{i=1}^n i )^{2 imes n}
]
[= (n!)^{2 imes n}
]
再看后面的式子:(先抛开 (-2) 次方,最后再说)
[= prod_{d=1}^n prod_{i=1}^n prod_{j=1}^n [gcd(i,j) == d]
]
[= prod_{d=1}^n d^{sum_{i=1}^{lfloor frac{n}{d}
floor} sum_{j=1}^{lfloor frac{n}{d}
floor} [gcd(i,j) == 1]}
]
先看指数:
[sum_{i=1}^{lfloor frac{n}{d}
floor} sum_{j=1}^{lfloor frac{n}{d}
floor} [gcd(i,j) == 1]
]
[= 2 imes (sum_{i=1}^{lfloor frac{n}{d}
floor} phi_i) - 1
]
(这步的依据参照 洛谷 P2568 GCD)
显然这个东西,用 欧拉筛 和 前缀和 就可以搞定!
用 (f_i) 表示 (phi_i) 的前缀和。
所以最终答案是:
[( prod_{i=1}^n prod_{j=1}^n i imes j ) imes ( prod_{i=1}^n prod_{j=1}^n frac{1}{(gcd(i,j))^2})
]
[= (n!)^{2 imes n} imes prod_{d=1}^n d^{(2 imes f_{lfloor frac{n}{d}
floor} - 1)^{-2}}
]
这样,就可以用 (O(n log n)) 解决问题。
可是,你会发现,这道题目时间限制和空间限制都有点紧,而且常数还挺大,所以不得不考虑优化。
而且, (f) 数组爆 ( exttt{long long}) 之后,似乎不能随便取模,这是个难题。
注:带个 (log) 是因为我们需要计算快速幂。
我们发现 (104857601) 是个质数。那么模质数呢,可以直接用在质数上模掉。
欧拉定理即:
[a^{phi_p} equiv 1 pmod p
]
所以,再化一步:
[= (n!)^{2 imes n} imes prod_{d=1}^n d^{((2 imes f_{lfloor frac{n}{d}
floor} - 1) \% (MOD-1))^{-2}}
]
然后呢,直接一波 欧拉筛 解决问题!
#pragma GCC optimize(2)
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll MOD=104857601;
const int N=1e6+1;
const int N1=1e5+1;
inline ll read(){char ch=getchar();int f=1;while(ch<'0' || ch>'9') {if(ch=='-') f=-f; ch=getchar();}
ll x=0;while(ch>='0' && ch<='9') x=(x<<3)+(x<<1)+ch-'0',ch=getchar();return x*f;}
int n,phi[N],prime[N1],f=0;
bool h[N]; ll k=1;
inline ll pw(ll x,ll y) {
ll ans=1; while(y) {
if(y&1) ans=(ans*x)%MOD;
x=(x*x)%MOD; y>>=1;
} return ans;
} //快速幂模板
inline void Euler() {
phi[1]=1; h[1]=1;
for(int i=2;i<=n;i++) {
k=(ll(k*i)%MOD); //顺便处理个阶乘
if(!h[i]) prime[++f]=i,phi[i]=i-1;
for(int j=1;j<=f && i*prime[j]<=n;j++) {
h[i*prime[j]]=1;
if(i%prime[j]==0) {
phi[i*prime[j]]=phi[i]*prime[j];
break;
}
phi[i*prime[j]]=phi[i]*(prime[j]-1);
}
}
} //欧拉筛模板
int main(){
n=read(); Euler();
// printf("%lld
",k);
k=pw(k,n<<1); ll ans=1; //第一部分的值
// printf("%lld
",k);
// for(int i=1;i<=n;i++) printf("%d ",phi[i]);
// putchar('
');
for(int i=1;i<=n;i++) phi[i]=(phi[i]<<1)+phi[i-1]%(MOD-1); //计算2倍的时候直接模掉
for(int i=2;i<=n;i++) ans=(ans*pw(i,phi[n/i]-1))%MOD; //第二部分
printf("%lld
",(k*pw(ll(ans*ans)%MOD,MOD-2))%MOD); //然后 -2 次方,就是2次方的倒数,而倒数已经转换为逆元
return 0;
}