莫比乌斯反演（洛谷P5221）

题意：
给一个n（1<=n<=1000000）;
求
$$\prod _{i=1}^n\prod _{j=1}^n\frac{lcm(i,j)}{gcd(i,j)}(mod104857601)$$
题解：
公式推演：
${\prod }_{i=1}^n{\prod }_{j=1}^n\frac{lcm(i,j)}{gcd(i,j)}$

$={\prod }_{i=1}^n{\prod }_{j=1}^n\frac{i\ast j}{gcd(i,j{)}^2}$

$=({\prod }_{i=1}^n{\prod }_{j=1}^{n}i\ast j)({\prod }_{i=1}^n{\prod }_{j=1}^{n}gcd(i,j{)}^{-2})$

$=({\prod }_{i=1}^n{i}^n\ast n!)({\prod }_{i=1}^n{\prod }_{j=1}^{n}gcd(i,j{)}^{-2})$

$=(n!{)}^{2n}\ast ({\prod }_{i=1}^n{\prod }_{j=1}^{n}gcd(i,j{)}^{-2})$

$=(n!{)}^{2n}\ast ({\prod }_{i=1}^n{\prod }_{j=1}^{n}gcd(i,j{)}^{-2})$

把${\prod }_{i=1}^n{\prod }_{j=1}^{n}gcd(i,j)$拿出来单独看

${\prod }_{i=1}^n{\prod }_{j=1}^{n}gcd(i,j)$

$={\prod }_{d=1}^n{\prod }_{i=1}^n{\prod }_{j=1}^n[gcd(i,j)==d]$

$={\prod }_{d=1}^n{d}^{{\sum }_{i=1}^n{\sum }_{j=1}^n[gcd(i,j)==d]}$

$={\prod }_{d=1}^n{d}^{{\sum }_{i=1}^{\lfloor \frac{n}{d}\rfloor }{\sum }_{j=1}^{\lfloor \frac{n}{d}\rfloor }[gcd(i,j)==1]}$

令$sum[x]={\sum }_{i=1}^x\varphi (i)$

由莫比乌斯反演可以得到：
$$(n!{)}^{2n}\ast (\prod _{d=1}^n{d}^{2sum\lfloor \frac{n}{d}\rfloor -1}{)}^{-2}$$

边上那个-2用欧拉降幂整掉，这题只给7MB空间，很容易被卡，用原地做法，分块也可以，分块确实可以更快点，但是得多开点空间，我这里就不用分块了
code：

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const int N=1e6+5;
const int mod=104857601;
int prime[80000],phi[N],tot=0;
ll fac=1;
inline int read() {
    int x=0,w=0;char ch=getchar();
    while(!isdigit(ch))w|=ch=='-',ch=getchar();
    while(isdigit(ch))x=(x<<1)+(x<<3)+(ch^48),ch=getchar();
    return w?-x:x;
}
void init(int maxn){
    phi[1]=1;
    for(int i=2;i<maxn;i++){
        fac=fac*i%mod;
        if(!phi[i]){
            prime[tot++]=i;
            phi[i]=i-1;
        }
        for(int j=0;j<tot&&prime[j]*i<maxn;j++){
            if(i%prime[j]==0){
                phi[i*prime[j]]=phi[i]*prime[j];
                break;
            }
            phi[i*prime[j]]=phi[i]*(prime[j]-1);
        }
    }
    for(int i=1;i<maxn;i++)phi[i]=(phi[i-1]+phi[i])%(mod-1);
}
ll fast(ll x,ll y=mod-1){
    ll ans=1;
    while(y){
        if(y&1)ans=ans*x%mod;
        x=x*x%mod;
        y>>=1;
    }
    return ans;
}
int main()
{
    int n=read();
    init(n+1);
    ll ans1,ans2=1;
    ans1=fast(fac,2*n);
    for(int i=1;i<=n;i++){
        ll power=((2-4*phi[n/i])%(mod-1)+mod-1)%(mod-1);
        ans2=ans2*fast(i,power)%mod;
    }
    printf("%lld",ans1*ans2%mod);
    return 0;
}