C. 在表格裡造序列

静谧幽蓝發表於2024-08-19

題意

對於每一對滿足 \(1\le i,j\le n\) 的 \((i,j)\)，計算有多少個長度為 \(m\) 的序列，權值在 \([1,n]\) 之間且 \(\gcd(a_1,a_2,...,a_m)=\gcd(i,j)\)。答案對 \(998244353\) 取模。

思路

方法：莫比烏斯反演+杜教篩

不會莫比烏斯反演？ 出門右轉：OI-wiki。
不會杜教篩？ 出門右轉：OI-wiki。

現在進入正文部分。

首先，我們考慮他的子問題：計算有多少個長度為 \(m\) 的序列，權值在 \([1,n]\) 之間且 \(\gcd(a_1,a_2,...,a_m)=k\)。

顯然的，\(a_i\) 一定是 \(k\) 的倍數。則我們給所有數同時除 \(k\)，現在就權值就變成了 \(\left[1,\left\lfloor\frac nk\right\rfloor\right]\)，限制也變成了 \(\gcd(a_1,a_2,...,a_m)=1\)。

我們設 \(f_i\) 表示 \(\gcd(a_1,a_2,...,a_m)=i\) 的方案數，\(g_i\) 表示僅僅滿足 \(i|a_j\) 的方案數。這也意味著 \(g_i\) 計算的是 \(i|\gcd(a_1,a_2,...,a_m)\) 的方案數。

那麼顯然的，\(g_k=\left\lfloor\frac nk\right\rfloor^m\)。而且，\(f\) 和 \(g\) 還滿足以下關係式：

\[g_k=\sum_{k|d}f_d \]

那麼根據莫比烏斯反演的式子，則有：

\[f_k=\sum_{k|d}\mu\left(\frac dk\right)g(d) \]

我們令 \(d=ik\)，則：

\[f_k=\sum_{i=1}^{\left\lfloor n/k\right\rfloor}\mu(i)g(ik)=\sum_{i=1}^{\left\lfloor n/k\right\rfloor}\mu(i)\left\lfloor\frac n{ik}\right\rfloor^m \]

現在我們知道了子問題的答案，那麼總問題的呢？

答案就很顯然了！就是這個式子：

\[\sum_{i=1}^n\sum_{j=1}^nf_{\gcd(i,j)} \]

當然，這個式子肯定是過不了的。於是我們考慮推式子。

\[\begin{align*} &\sum_{i=1}^n\sum_{j=1}^nf_{\gcd(i,j)}\\ =&\sum_{k=1}^n\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=k]f_k\\ =&\sum_{k=1}^nf_k\sum_{i=1}^{\left\lfloor n/k\right\rfloor}\sum_{j=1}^{\left\lfloor n/k\right\rfloor}[\gcd(i,j)=1]\\ \end{align*} \]

現在又要用到狄利克雷卷積的基本式：\(\epsilon=\mu*1\)。也就是說：

\[\sum_{d|n}\mu(d)=[n=1] \]

所以我們的式子又變成了這個樣子：

\[\begin{align*} &\sum_{k=1}^nf_k\sum_{i=1}^{\left\lfloor n/k\right\rfloor}\sum_{j=1}^{\left\lfloor n/k\right\rfloor}[\gcd(i,j)=1]\\ =&\sum_{k=1}^nf_k\sum_{i=1}^{\left\lfloor n/k\right\rfloor}\sum_{j=1}^{\left\lfloor n/k\right\rfloor}\sum_{d|\gcd(i,j)}\mu(d)\\ =&\sum_{k=1}^nf_k\sum_{i=1}^{\left\lfloor n/k\right\rfloor}\sum_{j=1}^{\left\lfloor n/k\right\rfloor}\sum_{d|i,d|j}\mu(d)\\ =&\sum_{k=1}^nf_k\sum_{d=1}^{\left\lfloor n/k\right\rfloor}\sum_{i=1}^{\left\lfloor n/kd\right\rfloor}\sum_{j=1}^{\left\lfloor n/kd\right\rfloor}\mu(d)\\ =&\sum_{k=1}^nf_k\sum_{d=1}^{\left\lfloor n/k\right\rfloor}\mu(d)\left\lfloor\frac{n}{kd}\right\rfloor^2\\ =&\sum_{k=1}^n\left(\sum_{i=1}^{\left\lfloor n/k\right\rfloor}\mu(i)\left\lfloor\frac n{ik}\right\rfloor^m\right)\left(\sum_{d=1}^{\left\lfloor n/k\right\rfloor}\mu(d)\left\lfloor\frac{n}{kd}\right\rfloor^2\right)\\ \end{align*} \]

發現了嗎？ 兩個式子竟長的如此相似！

於是我們考慮裡面的式子可以用整除分塊做。而觀察到每一個除法都帶有 \(\left\lfloor\frac nk\right\rfloor\) 的因式，於是也可以對 \(k\) 做整除分塊。

因為 \(n\) 有 \(10^9\)，於是要用杜教篩算 \(\mu\) 函式。

複雜度 \(O(n^{0.75})\)。

程式碼

#include <bits/stdc++.h>
#pragma GCC optimize(2)

using namespace std;

#define map unordered_map
#define int long long

#define fileio(x,y) freopen(x,"r",stdin),freopen(y,"w",stdout)
#define tup tuple<int,int,int>
#define pii pair<int,int>
#define bit(x) bitset<x>
#define pb emplace_back
#define mt make_tuple 
#define mp make_pair

const int   mod=998244353;
const int   maxn=2e6+10;

map<int,int>    mu,pp;
bit(maxn)       vis;
int             summu[maxn];
int             pri[maxn];
int             n,m,cnt;

int power(int a,int b,int p) { int tar=1; for(; b; b>>=1,a=1ll*a*a%p) if(b&1) tar=1ll*tar*a%p; return tar; }

int calcMu(int x)
{
    if(x<maxn) return summu[x];
    if(mu[x]) return mu[x];
    int val=1;
    for(int l=2,r; l<=x; l=r+1)
    {
        r=(x/(x/l));
        (val-=calcMu(x/l)*(r-l+1)%mod)%=mod;
    }
    return (mu[x]=val);
}

void init()
{
    summu[1]=1;
    for(int i=2; i<maxn; i++)
    {
        if(!vis[i]) { pri[++cnt]=i; summu[i]=-1; }
        for(int j=1; j<=cnt&&i*pri[j]<maxn; j++)
        {
            vis[i*pri[j]]=true;
            if(i%pri[j]==0) { summu[i*pri[j]]=0; break; }
            summu[i*pri[j]]=-summu[i];
        }
    }
    for(int i=1; i<maxn; i++) (summu[i]+=summu[i-1])%=mod;
    return;
}

void work()
{
    /* Code */
    cin>>n>>m;
    // n=1e9,m=1e18;
    int tar=0; m%=(mod-1);
    for(int l=1,r; l<=n; l=r+1) { r=(n/(n/l)); pp[n/l]=power(n/l,m,mod); }
    for(int l=1,r; l<=n; l=r+1)
    {
        r=(n/(n/l));
        // initialization
        int temp=n/l,tar1=0,tar2=0;
        for(int p=1,q; p<=temp; p=q+1)
        {
            q=(temp/(temp/p));
            int qmu=calcMu(q)-calcMu(p-1);
            (tar1+=pp[temp/p]*qmu%mod)%=mod;
            (tar2+=qmu*(temp/p)%mod*(temp/p)%mod)%=mod;
        }
        tar+=(r-l+1)*tar1%mod*tar2%mod;
    }
    cout<<(tar%mod+mod)%mod<<'\n';
    return;
}

signed main()
{
    // fileio(".in",".out");
    ios::sync_with_stdio(false);
    cin.tie(0);
    int t,stTime=clock();
    init();
    t=1;
    while(t--) work();
    // cerr<<"Time : "<<(double)(clock()-stTime)/CLOCKS_PER_SEC<<'\n';
    return 0;
} // Texas yyds !!!

相關文章