BZOJ 4174 tty的求助 莫比乌斯反演

题目大意：求$\sum_{n=1}^N\sum_{m=1}^M\sum_{k=0}^{m-1}\lfloor\frac{nk+x}m\rfloor\ mod\ 998244353$

假设$n$和$m$都已经确定了，现在要求这坨玩应：
$\sum_{k=0}^{m-1}\lfloor\frac{nk+x}m\rfloor$
$=\sum_{k=0}^{m-1}(\lfloor\frac{nk\%m+x}m\rfloor+\frac{nk-nk\%m}m)$
$=\sum_{k=0}^{m-1}(\lfloor\frac{nk\%m+x}m\rfloor+\frac{nk}m-\frac{nk\%m}m)$

我们一项一项考虑

令$d=\gcd(n,m)$，那么有

$\sum_{k=0}^{m-1}\lfloor\frac{nk\%m+x}m\rfloor$
$=d*\sum_{k=0}^{\frac md-1}\lfloor\frac{kd+x}m\rfloor$
$=d*(\frac md*\frac{x-x\%m}m+\sum_{k=0}^{\frac md-1}\lfloor\frac{kd+x\%m}m\rfloor)$
$=d*(\frac md*\frac{x-x\%m}m+\sum_{k=0}^{\frac md-1}[kd+x\%m\geq m])$
$=d*(\frac{x-x\%m}d+\lfloor\frac{x\%m}d\rfloor)$
$=d*\lfloor\frac xd\rfloor$

$\sum_{k=0}^{m-1}\frac{nk}m=\frac nm*\frac{m*(m-1)}2=\frac{n*m-n}2$

$\sum_{k=0}^{m-1}\frac{nk\%m}m=d*\sum_{k=0}^{\frac md-1}\frac{kd}m=\frac{d^2}m*\frac{(\frac md-1)*\frac md}2=\frac{m-d}2$

故答案为
$\sum_{n=1}^N\sum_{m=1}^M(d*\lfloor\frac xd\rfloor+\frac{n*m-n}2-\frac{m-d}2)$
$=\frac12*\sum_{n=1}^N\sum_{m=1}^M(2*d*\lfloor\frac xd\rfloor+d+n*m-n-m)$
$=\frac12*(S(N)*S(M)-S(N)*m-S(M)*n+\sum_{d=1}^{min(N,M)}(d+2*d*\lfloor\frac xd\rfloor)\sum_{k=1}^{min(\lfloor\frac Nd\rfloor,\lfloor\frac Md\rfloor)}\mu(k)*\lfloor\frac N{d*k}\rfloor*\lfloor\frac M{d*k}\rfloor)$

其中$S(n)=\frac{n*(n+1)}2$

然后$O(nlogn)$枚举$d$和$k$即可

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define M 500500
#define MOD 998244353
using namespace std;
int n,m,x;
long long ans;
int mu[M];
int prime[M],tot;
bool not_prime[M];
void Linear_Shaker()
{
    int i,j;
    mu[1]=1;
    for(i=2;i<=500000;i++)
    {
        if(!not_prime[i])
        {
            prime[++tot]=i;
            mu[i]=MOD-1;
        }
        for(j=1;prime[j]*i<=500000;j++)
        {
            not_prime[prime[j]*i]=true;
            if(i%prime[j]==0)
            {
                mu[prime[j]*i]=0;
                break;
            }
            mu[prime[j]*i]=(MOD-mu[i])%MOD;
        }
    }
}
long long Sum(long long n)
{
    return (n*(n+1)>>1)%MOD;
}
int main()
{
    int i,j;
    cin>>n>>m>>x;
    Linear_Shaker();
    ans=((Sum(n)*Sum(m)-Sum(n)*m-Sum(m)*n)%MOD+MOD)%MOD;
    if(n>m) swap(n,m);
    for(i=1;i<=n;i++)
    {
        long long temp=i+x/i*i*2;
        for(j=1;j*i<=n;j++)
            (ans+=temp*mu[j]%MOD*(n/i/j)%MOD*(m/i/j)%MOD)%=MOD;
    }
    cout<<(ans*(MOD+1>>1)%MOD)<<endl;
    return 0;
}