machine 一道数论题目。
题目
题目描述
一句话题意：给定 \(f _ 0(n)\)，有 \(m\) 组询问，每组询问给定 \(k,x\)，求 \((f _ 0\ast 1^k)(x)\pmod{998244353}\)。
想看详细题意可以看原题 pdf。
数据范围
对于 \(10\%\) 的数据，\(m=0\)；
对于 \(20\%\) 的数据，\(m\leq 5\)；
对于另 \(30\%\) 的数据，\(f _ 0[i]\) 全都是 \(1\)；
对于另 \(30\%\) 的数据，所有询问的 \(k\) 相等；
对于 \(100\%\) 的数据，\(n\leq 10^5,m\leq 1000,0\leq f _ 0[i]<998244353,k\leq 10^5,x\leq n\)。
时空限制
- machine.cpp/in/out
- Time limit: 4s
- Memory limit: 512MB
题解
解法 A（人口普查）
对于 \(10\%\) 的数据，\(m=0\)。
我们直接提交一个可以通过编译的空白主函数程序即可。
解法 B（小发现）
对于 \(20\%\) 的数据，\(m\leq 5\)。
我们不妨考虑直接模拟 \(\texttt{Dirichlet}\) 卷积的过程，并且加一点点的小优化（询问按 \(k\) 排序），时间复杂度为 \(\Theta(mn\ln n\max\{k\})\)，显然过不掉。
考虑 \(\texttt{Dirichlet}\) 卷积的式子。
\[f _ k(n)=\sum _ {d\mid n}f _ {k-1}(d)\]
我们不难发现，\(\texttt{Dirichlet}\) 卷积既符合交换律又符合结合律。
所以我们可以通过快速幂的方式优化卷积过程，时间复杂度为 \(\Theta(mn\ln n\log _ 2\max\{k\})\)。
解法 C（小推广）
对于另 \(30\%\) 的数据，所有询问的 \(k\) 相等。
与解法 B 相同的，我们用快速幂优化卷积过程，最后的时间复杂度为 \(\Theta(n\ln n\log _ 2k+m)\)。
解法 D（臃肿的写法）
这是 std 的解法。
对于一个数 \(n\)，我们不难分解它为
\[n=\prod _ {p _ i\in\mathbb{P} _ n}p _ i^{c _ i}\]
对于一个数 \(n\)，我们根据分解的结果定义一个数列
\[c _ 1,c _ 2,\cdots,c _ t\]
其中有
\[c _ 1\leq c _ 2\leq c _ 3\leq\cdots\leq c _ t\]
我们规定，如果两个数对应的数列相同，我们就称这两个数是乘法同构的。
不难观察到（打表即可），对于小于等于 \(10^5\)（本题 \(n\) 的范围）的数中，根据乘法同构分类，本质不同的数仅有 \(160\) 个（题解里面写的是 \(165\) 个，显然是不认真思考，抄别人代码写的题解）。
于是我们考虑用乘法同构矩阵优化 \(\texttt{Dirichlet}\) 卷积过程。
最后的时间复杂度为
\[\Theta(n\log^2 _ 2n+n\ln n+|S|^3\log _ 2k+m(|S|^2\log _ 2k+\sqrt{x}))\]
奇慢无比。
代码
std 的做法太慢了，代码还长。
#include<bits/stdc++.h>
using namespace std;
#define reg register
typedef long long ll;
#define MOD 998244353
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,100000,stdin),p1==p2)?EOF:*p1++)
static char buf[100000],*p1=buf,*p2=buf;
inline int read(void){
reg char ch=getchar();
reg int res=0;
while(ch<'0'||'9'<ch)ch=getchar();
while('0'<=ch&&ch<='9')res=10*res+ch-'0',ch=getchar();
return res;
}
inline int add(reg int a,reg int b){
reg int sum=a+b;
return sum>=MOD?sum-MOD:sum;
}
const int MAXN=100000+5;
const int MAXLOG2K=17+1;
int n,m;
int f[MAXN];
int col[MAXN];
int link[MAXN];
map<vector<int>,int> S;
bool vis[MAXN];
int tot,prime[MAXN];
int from[MAXN];
inline void Init(reg int n){
vis[1]=true;
from[1]=1;
for(reg int i=2;i<=n;++i){
if(!vis[i]){
prime[++tot]=i;
from[i]=i;
}
for(reg int j=1;j<=tot&&i*prime[j]<=n;++j){
vis[i*prime[j]]=true;
from[i*prime[j]]=prime[j];
if(i%prime[j]==0)
break;
}
}
return;
}
const int MAXCNT=160+1+5;
struct Matrix{
int n,m;
int unit[MAXCNT][MAXCNT];
inline Matrix(void){
n=m=0;
memset(unit,0,sizeof(unit));
return;
}
inline int* operator[](reg int i){
return unit[i];
}
inline Matrix operator*(const Matrix &a)const{
Matrix res;
res.n=n,res.m=a.m;
reg int r;
for(reg int i=1;i<=n;++i)
for(reg int k=1;k<=m;++k){
r=unit[i][k];
for(reg int j=1;j<=a.m;++j)
res[i][j]=add(res[i][j],1ll*r*a.unit[k][j]%MOD);
}
return res;
}
};
Matrix base[MAXLOG2K];
int main(void){
//freopen("machine.in","r",stdin);
//freopen("machine.out","w",stdout);
n=read();
for(reg int i=1;i<=n;++i)
f[i]=read();
Init(1e5);
reg int cnt=0;
col[++cnt]=1;
link[1]=1;
for(reg int i=2;i<=n;++i){
vector<int> V;
reg int tmp=i;
while(tmp!=1){
reg int p=from[tmp];
int cnt=0;
while(tmp%p==0)
tmp/=p,++cnt;
V.push_back(cnt);
}
sort(V.begin(),V.end());
if(S.count(V))
col[i]=S[V];
else{
col[i]=S[V]=++cnt;
link[i]=1;
}
}
base[0].n=base[0].m=cnt;
for(reg int i=1;i<=n;++i)
for(reg int k=1;i*k<=n;++k)
base[0][col[i]][col[i*k]]+=link[i*k];
for(reg int i=1;i<MAXLOG2K;++i)
base[i]=base[i-1]*base[i-1];
m=read();
while(m--){
static int k,x;
k=read(),x=read();
Matrix A;
A.n=1,A.m=cnt;
A[1][1]=1;
for(reg int i=0;i<MAXLOG2K;++i)
if((k>>i)&1)
A=A*base[i];
reg int ans=0;
for(reg int i=1;i*i<=x;++i)
if(x%i==0){
ans=add(ans,1ll*A[1][col[x/i]]*f[i]%MOD);
if(i*i!=x)
ans=add(ans,1ll*A[1][col[i]]*f[x/i]%MOD);
}
printf("%d\n",ans);
}
return 0;
}解法 E（人类智慧）
考虑递推式
\[
\begin{aligned}
f _ k(n)&=\sum _ {d _ 1\mid n}f _ {k-1}(d _ 1)\\
&=\sum _ {d _ 1\mid n}\sum _ {d _ 2\mid d _ 1}f _ {k-2}(d _ 2)\\
&=\sum _ {d _ 1\mid n}\sum _ {d _ 2\mid d _ 1}\sum _ {d _ 3\mid d _ 2}f _ {k-3}(d _ 3)\\
&=\sum _ {d _ 1\mid n}\sum _ {d _ 2\mid d _ 1} \cdots \sum _ {d _ {i+1}\mid d _ i}\cdots\sum _ {d _ k\mid d _ {k-1}}f _ 0(d _ k)
\end{aligned}
\]
我们不妨变换一下形式
\[
\begin{aligned}
f _ k(n)&=\sum _ {d _ 1\mid n}\sum _ {d _ 2\mid d _ 1} \cdots \sum _ {d _ {i+1}\mid d _ i}\cdots\sum _ {d _ k\mid d _ {k-1}}f _ 0(d _ k)\\
&=\sum _ {d _ k\mid d _ {k-1}}f _ 0(d _ k)[d _ k\mid d _ {k-1}\mid\cdots\mid d _ 2\mid d _ 1\mid n]\\
&=\sum _ {d _ k\mid n}f _ 0(d _ k)[d _ k\mid d _ {k-1}\mid\cdots\mid d _ 2\mid d _ 1\mid n]
\end{aligned}
\]
此时我们注意到后面的特征函数
\[
[d _ k\mid d _ {k-1}\mid\cdots\mid d _ 2\mid d _ 1\mid n]
\]
我们把上面的 \(d\) 看作分解后得到的集合，则
\[
\mathbb{P} _ {d _ k}\subseteq\mathbb{P} _ {d _ {k-1}}\subseteq\mathbb{P} _ {d _ {k-2}}\subseteq\cdots\subseteq\mathbb{P} _ {d _ 2}\subseteq\mathbb{P} _ {d _ 1}\subseteq\mathbb{P} _ n
\]
不难发现，我们本质上就是要求满足上述子集链关系的方案数，方案数与 \(f _ 0(d _ k)\) 的乘积即为对答案的贡献。
我们考虑把所有集合都消去一个 \(\mathbb{P} _ {d _ k}\)，这对答案没有影响，然后得到
\[
\begin{matrix}
\mathbb{P} _ {d _ k}&\subseteq&\mathbb{P} _ {d _ {k-1}}&\subseteq&\cdots&\subseteq&\mathbb{P} _ {d _ 1}&\subseteq&\mathbb{P} _ n\\
\downarrow&&\downarrow&&\downarrow&&\downarrow&&\downarrow\\
\varnothing&\subseteq&\mathbb{S} _ 1&\subseteq&\cdots&\subseteq&\mathbb{S} _ {k-1}&\subseteq&(\mathbb{P} _ n-\mathbb{P} _ {d _ k})
\end{matrix}
\]
\[
\varnothing\subseteq\mathbb{S} _ 1\subseteq\mathbb{S} _ 2\subseteq\cdots\subseteq\mathbb{S} _ {k-2}\subseteq\mathbb{S} _ {k-1}\subseteq(\mathbb{P} _ n-\mathbb{P} _ {d _ k})
\]
考虑子集符号象征着集合之间的差（变化），我们可以分开质因子讨论。
上述链中共有 \(k\) 个子集符号。
对于一个质因子 \(p _ i\)，它在 \((\mathbb{P} _ n-\mathbb{P} _ {d _ k})\) 中的出现次数为 \(c _ i\)，那么它单独的方案数目可以转化为：
\(c _ i\) 个相同小球放入 \(k\) 个不同盒子中，允许有的盒子空着的方案数。
用插板法可以求得为 \(\binom{k+c _ i-1}{k-1}=\binom{k+c _ i-1}{c _ i}\)。
再用乘法原理对质因子对方案数的贡献进行合并，直接乘起来即可。
最后得出答案的表达式
\[
f _ k(x)=\sum _ {d\mid x}f _ 0(d)\prod _ {c _ i\in(\mathbb{P} _ x-\mathbb{P} _ d)}\binom{k+c _ i-1}{c _ i}\pmod{998244353}
\]
这种做法快到飞起，时间复杂度为 \(\Theta(n+m\sum _ {d\mid x}\sqrt{d})\)，其中 \(\Theta(n)\) 的复杂度用于预处理阶乘及其逆元，实际上的复杂度跑不满，实际测试趋近于 \(\Theta(n+m\log _ 2m)\)。
代码
时间效率是解法 D 的 200 倍（平均一个点 0.01s），并且代码还短了一半。
#include<bits/stdc++.h>
using namespace std;
#define reg register
typedef long long ll;
#define MOD 998244353
#define getchar() (p1==p2&&(p2=(p1=buf)+fread(buf,1,100000,stdin),p1==p2)?EOF:*p1++)
static char buf[100000],*p1=buf,*p2=buf;
inline int read(void){
reg char ch=getchar();
reg int res=0;
while(ch<'0'||'9'<ch)ch=getchar();
while('0'<=ch&&ch<='9')res=10*res+ch-'0',ch=getchar();
return res;
}
const int MAXN=100000+5;
const int MAXS=500000+5;
inline int add(reg int a,reg int b){
reg int sum=a+b;
return sum>=MOD?sum-MOD:sum;
}
int n;
int f[MAXN];
int fac[MAXS],invfac[MAXS];
inline int pow(reg int x,reg int exp){
reg int res=1;
while(exp){
if(exp&1)
res=1ll*res*x%MOD;
x=1ll*x*x%MOD;
exp>>=1;
}
return res;
}
inline void Init(reg int n){
fac[1]=invfac[1]=1;
for(reg int i=2;i<=n;++i)
fac[i]=1ll*fac[i-1]*i%MOD;
invfac[n]=pow(fac[n],MOD-2);
for(reg int i=n-1;i>=2;--i)
invfac[i]=1ll*invfac[i+1]*(i+1ll)%MOD;
return;
}
inline int C(reg int n,reg int m){
return 1ll*fac[n]*invfac[m]%MOD*invfac[n-m]%MOD;
}
inline int Solve(reg int x,reg int k){
reg int res=1;
for(reg int i=2;i*i<=x;++i)
if(x%i==0){
reg int cnt=0;
while(x%i==0)
++cnt,x/=i;
res=1ll*res*C(k+cnt-1,cnt)%MOD;
}
if(x>1)
res=1ll*res*k%MOD;
return res;
}
int main(void){
Init(5e5);
n=read();
for(reg int i=1;i<=n;++i)
f[i]=read();
reg int m=read();
while(m--){
static int k,x;
k=read(),x=read();
reg int ans=0;
for(reg int i=1;i*i<=x;++i)
if(x%i==0){
ans=add(ans,1ll*f[i]*Solve(x/i,k)%MOD);
if(i*i!=x)
ans=add(ans,1ll*f[x/i]*Solve(i,k)%MOD);
}
printf("%d\n",ans);
}
return 0;
}
近期评论