洛谷 P5221 Product 题解

原题链接

庆祝!第二道数论紫题。

推式子真是太有趣了!

[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;
}
原文地址:https://www.cnblogs.com/bifanwen/p/12509625.html