P4240 毒瘤之神的考验【莫比乌斯反演+各种数学技巧】

毒瘤题

题意:

T T T 次询问,每次给定 n , m n,m n,m ,求: ( i = 1 n j = 1 m φ ( i j ) )    m o d    998244353 (\sum_{i=1}^{n}{\sum_{j=1}^{m}{\varphi(ij)}})\;mod\;998244353 (i=1nj=1mφ(ij))mod998244353 的值。

数据范围: 1 T 1 0 4 , 1 n , m 1 0 5 1\leq T \leq 10^4,1\leq n,m\leq10^5 1T104,1n,m105

解题过程:

以下过程中 n = m i n ( n , m ) n=min(n,m) n=min(n,m)
首先要知道如下性质:

φ ( a b ) = φ ( a ) φ ( b ) g c d ( a , b ) φ ( g c d ( a , b ) ) \varphi(ab)=\frac{\varphi(a)\varphi(b)*gcd(a,b)}{\varphi(gcd(a,b))} φ(ab)=φ(gcd(a,b))φ(a)φ(b)gcd(a,b)

证明

φ ( a ) φ ( b ) = a <munder> p a </munder> p 1 p    b <munder> p b </munder> p 1 p \varphi(a)\varphi(b)=a\prod_{p|a}{\frac{p-1}{p}\;b\prod_{p|b}{\frac{p-1}{p}}} φ(a)φ(b)=apapp1bpbpp1
= a b <munder> p a b </munder> p 1 p <munder> p g c d ( a , b ) </munder> p 1 p =ab\prod_{p|ab}{\frac{p-1}{p}\prod_{p|gcd(a,b)}{\frac{p-1}{p}}} =abpabpp1pgcd(a,b)pp1
两边同时乘以 g c d ( a , b ) gcd(a,b) gcd(a,b),有:
φ ( a )    φ ( b )    g c d ( a , b ) = φ ( a b )    φ ( g c d ( a b ) ) \varphi(a)\;\varphi(b)\;gcd(a,b)=\varphi(ab)\;\varphi(gcd(ab)) φ(a)φ(b)gcd(a,b)=φ(ab)φ(gcd(ab))
即证。

接下来就可以进行神奇的化简过程:

<munderover> i = 1 n </munderover> <munderover> j = 1 m </munderover> φ ( i j ) \sum_{i=1}^{n}{\sum_{j=1}^{m}{\varphi(ij)}} i=1nj=1mφ(ij)
= <munderover> i = 1 n </munderover> <munderover> j = 1 m </munderover> φ ( i ) φ ( j ) g c d ( i , j ) φ ( g c d ( i , j ) ) =\sum_{i=1}^{n}{\sum_{j=1}^{m}{\varphi(i)\varphi(j)\frac{gcd(i,j)}{\varphi(gcd(i,j))}}} =i=1nj=1mφ(i)φ(j)φ(gcd(i,j))gcd(i,j)
= <munderover> d = 1 n </munderover> <munderover> i = 1 n </munderover> <munderover> j = 1 m </munderover> φ ( i ) φ ( j ) d φ ( d ) [ g c d ( i , j ) = d ]    ( ) =\sum_{d=1}^{n}{\sum_{i=1}^{n}{\sum_{j=1}^{m}{\varphi(i)\varphi(j)\frac{d}{\varphi(d)}[gcd(i,j)=d]}}}\;(改变枚举顺序) =d=1ni=1nj=1mφ(i)φ(j)φ(d)d[gcd(i,j)=d]()
= <munderover> d = 1 n </munderover> <munderover> i = 1 n d </munderover> <munderover> j = 1 m d </munderover> φ ( i d ) φ ( j d ) d φ ( d ) [ g c d ( i , j ) = 1 ] =\sum_{d=1}^{n}{\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}{\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}{\varphi(id)\varphi(jd)\frac{d}{\varphi(d)}[gcd(i,j)=1]}}} =d=1ni=1dnj=1dmφ(id)φ(jd)φ(d)d[gcd(i,j)=1]
= <munderover> d = 1 n </munderover> <munderover> i = 1 n d </munderover> <munderover> j = 1 m d </munderover> φ ( i d ) φ ( j d ) d φ ( d ) <munder> k g c d ( i , j ) </munder> μ ( k ) =\sum_{d=1}^{n}{\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}{\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}{\varphi(id)\varphi(jd)\frac{d}{\varphi(d)}\sum_{k|gcd(i,j)}{\mu(k)}}}} =d=1ni=1dnj=1dmφ(id)φ(jd)φ(d)dkgcd(i,j)μ(k)
= <munderover> d = 1 n </munderover> <munderover> i = 1 n d </munderover> <munderover> j = 1 m d </munderover> φ ( i d ) φ ( j d ) d φ ( d ) <munder> k i , k j </munder> μ ( k ) =\sum_{d=1}^{n}{\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}{\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}{\varphi(id)\varphi(jd)\frac{d}{\varphi(d)}\sum_{k|i,k|j}{\mu(k)}}}} =d=1ni=1dnj=1dmφ(id)φ(jd)φ(d)dki,kjμ(k)
= <munderover> d = 1 n </munderover> d φ ( d ) <munderover> i = 1 n d </munderover> <munderover> j = 1 m d </munderover> φ ( i d ) φ ( j d ) <munder> k i , k j </munder> μ ( k ) =\sum_{d=1}^{n}{\frac{d}{\varphi(d)}\sum_{i=1}^{\lfloor \frac{n}{d}\rfloor}{\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}{\varphi(id)\varphi(jd)\sum_{k|i,k|j}{\mu(k)}}}} =d=1nφ(d)di=1dnj=1dmφ(id)φ(jd)ki,kjμ(k)
= <munderover> d = 1 n </munderover> d φ ( d ) <munderover> k = 1 n d </munderover> μ ( k ) <munderover> i = 1 n d k </munderover> <munderover> j = 1 m d k </munderover> φ ( i d k ) φ ( j d k ) =\sum_{d=1}^{n}{\frac{d}{\varphi(d)}\sum_{k=1}^{\lfloor\frac{n}{d}\rfloor}{\mu(k)\sum_{i=1}^{\lfloor\frac{n}{dk}\rfloor}{\sum_{j=1}^{\lfloor\frac{m}{dk}\rfloor}{\varphi(idk)\varphi(jdk)}}}} =d=1nφ(d)dk=1dnμ(k)i=1dknj=1dkmφ(idk)φ(jdk)
T = d k T=dk T=dk,
= <munderover> T = 1 n </munderover> <munder> k T </munder> k φ ( k ) μ ( T k ) <munderover> i = 1 n T </munderover> φ ( i T ) <munderover> j = 1 m T </munderover> φ ( j T ) =\sum_{T=1}^{n}{\sum_{k|T}{\frac{k}{\varphi(k)}\mu(\frac{T}{k})}\sum_{i=1}^{\lfloor\frac{n}{T}\rfloor}{\varphi(iT)\sum_{j=1}^{\lfloor\frac{m}{T}\rfloor}{\varphi(jT)}}} =T=1nkTφ(k)kμ(kT)i=1Tnφ(iT)j=1Tmφ(jT)

F ( x ) = k x k φ ( k ) μ ( x k )    F(x)=\sum_{k|x}{\frac{k}{\varphi(k)}\mu(\frac{x}{k})}\; F(x)=kxφ(k)kμ(kx), 显然可以通过先枚举因子再枚举倍数的方法(即埃氏筛的原理)在 O ( n l o g n ) O(nlogn) O(nlogn) 的时间内预处理出来。

G [ x ] [ y ] = i = 1 y φ ( i x ) G[x][y]=\sum_{i=1}^{y}{\varphi(ix)} G[x][y]=i=1yφ(ix) ,显然有 G [ x ] [ y ] = G [ x ] [ y 1 ] + φ ( y x ) G[x][y]=G[x][y-1]+\varphi(yx) G[x][y]=G[x][y1]+φ(yx) ,同时对于一个 x x x ,有 y m x y\leq\lfloor\frac{m}{x}\rfloor yxm ,因此也可以利用递推式在 O ( n l o g n ) O(nlogn) O(nlogn) 下预处理出。

S [ y ] [ z ] [ x ] = x = 1 n k x k φ ( k ) μ ( x k ) i = 1 y φ ( i x ) j = 1 z φ ( j x ) S[y][z][x]=\sum_{x=1}^{n}{\sum_{k|x}{\frac{k}{\varphi(k)}\mu(\frac{x}{k})}\sum_{i=1}^{y}{\varphi(ix)\sum_{j=1}^{z}{\varphi(jx)}}} S[y][z][x]=x=1nkxφ(k)kμ(kx)i=1yφ(ix)j=1zφ(jx)
有递推式: S [ y ] [ z ] [ x ] = S [ y ] [ z ] [ x 1 ] + F [ x ] G [ x ] [ y ] G [ x ] [ z ] S[y][z][x]=S[y][z][x-1]+F[x]*G[x][y]*G[x][z] S[y][z][x]=S[y][z][x1]+F[x]G[x][y]G[x][z]
这个有三维,如果直接算的话,时间上肯定会爆炸,而且空间上也不能实现。
根据数论分块的理论,对 x x x 进行分块讨论,当 n x = m x \lfloor\frac{n}{x}\rfloor=\lfloor\frac{m}{x}\rfloor xn=xm 时,对应 x x x 的区间为 [ l , r ] [l,r] [l,r] ,那么该段区间上的答案为 S [ n / l ] [ m / l ] [ r ] S [ n / l ] [ m / l ] [ l 1 ] S[n/l][m/l][r]-S[n/l][m/l][l-1] S[n/l][m/l][r]S[n/l][m/l][l1]。但是无法预处理出全部的 S S S,我们只能在 n l \lfloor\frac{n}{l}\rfloor ln 比较小,即 x x x 比较大的时候利用预处理的结果。当 x x x 比较小的时候,可以利用递推式直接计算。但是大和小之间的界限是什么?这时,我们可以人为对 y y y z z z 设置范围。当 0 y , z m a x n 0\leq y,z \leq maxn 0y,zmaxn 时,用分块;否则,用递推暴力。
预处理的那部分的时间复杂度小于 O ( n m a x n 2 ) O(n*maxn^2) O(nmaxn2)
这样,总的时间复杂度为: O ( n l o g n + n m a x n 2 + T ( n + n / m a x n ) ) O(nlogn+n*maxn^2+T(\sqrt n+n/maxn)) O(nlogn+nmaxn2+T(n +n/maxn)),实际上要更小。

代码实现:由于 G G G S S S 都比较大,且每一维的长度都不相等,所以用 v e c t o r vector vector会比较省空间。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=998244353;
const int N=1e5+5;
const int maxn=31;
vector<int>prime;
int mu[N],phi[N];
ll inv[N],F[N];
vector<ll>G[N],S[maxn][maxn];
bool vis[N];
void read(int &x)
{
    x=0;
    int f=1;
    char c=getchar();
    while(!isdigit(c))
    {
        if(c=='-')
            f=-1;
        c=getchar();
    }
    while(isdigit(c))
    {
        x=(x<<3)+(x<<1)+c-'0';
        c=getchar();
    }
    x*=f;
}
void init()
{
    mu[1]=1;
    phi[1]=1;
    inv[1]=1;
    for(int i=2;i<=N-5;i++)//预处理欧拉函数和莫比乌斯函数
    {
        if(!vis[i])
        {
            prime.push_back(i);
            mu[i]=-1;
            phi[i]=i-1;
        }
        for(int j=0;j<prime.size()&&i*prime[j]<=N-5;j++)
        {
            vis[i*prime[j]]=1;
            if(i%prime[j]==0)
            {
                mu[i*prime[j]]=0;
                phi[i*prime[j]]=phi[i]*prime[j];
                break;
            }
            else
            {
                mu[i*prime[j]]=-mu[i];
                phi[i*prime[j]]=phi[i]*phi[prime[j]];
            }
        }
    }
    for(int i=2;i<=N-5;i++)//递推逆元
        inv[i]=1LL*(mod-mod/i)*inv[mod%i]%mod;
    for(int i=1;i<=N-5;i++)//预处理函数:F(x)
        for(int j=1;j*i<=N-5;j++)
            F[i*j]=(F[i*j]+1LL*mu[j]*i%mod*inv[phi[i]]%mod)%mod;
    for(int i=1;i<=N-5;i++)//预处理函数G[x][y]
    {
        G[i].push_back(0);//G[i][0]=0;
        for(int j=1;j<=(N-5)/i;j++)
            G[i].push_back(G[i][j-1]+1LL*phi[i*j]%mod);//G[i][j]=G[i][j-1]+phi[i*j];
    }
    for(int j=1;j<maxn;j++)
    {
        for(int k=1;k<maxn;k++)
        {
            S[j][k].push_back(0);
            for(int i=1;i<=(N-5)/max(j,k);i++)//除最大
                S[j][k].push_back((S[j][k][i-1]+F[i]*G[i][j]%mod*G[i][k]%mod)%mod);
        }
    }
}
ll solve(int n,int m)
{
    ll ans=0;
    if(n>m)
        swap(n,m);
    for(int i=1;i<=m/(maxn-1);i++)//用m除以界限
        ans=(ans+F[i]*G[i][n/i]%mod*G[i][m/i]%mod)%mod;
    for(int l=m/(maxn-1)+1,r=0;l<=n;l=r+1)
    {
        r=min(n,min(n/(n/l),m/(m/l)));
        ans=(ans+(S[n/l][m/l][r]-S[n/l][m/l][l-1]+mod)%mod)%mod;
    }
    return ans;
}
int main()
{
    int t,n,m;
    init();
    read(t);//cout<<"t= "<<t<<endl;
    while(t--)
    {
        read(n);
        read(m);
        printf("%lld\n",solve(n,m));
    }
    return 0;
}

参考博客1
参考博客2

全部评论

相关推荐

点赞 收藏 评论
分享
牛客网
牛客企业服务