P9060 [Ynoi2002] Goedel Machine
前言
- 循环数瞎搞见祖宗

- 准备推荐这题了
\(solve\)
-
我们把哥德尔机放一边 -
考虑一个数 \(m\)
-
如果在一段区间内有 \(k\) 个这个数的倍数,那共有 \(2^k-1\) 个子集的 \(\gcd\) 是 \(m\)(去掉空集),总共会产生 \(m^{2^k-1}\) 的贡献
-
那么我们考虑统计贡献,由于我们要把这些贡献乘起来,对于一个数 \(m\),我们进行拆分,可以证明可以拆成若干个质数相乘来表示,反正最后都是乘起来没啥影响,所以我们就将一个任意数的贡献简化成对于质数的贡献,这就很好做了
-
然后你被卡亖了 -
我们貌似没法很好的做区间的东西,所以我们直接上莫队,由于我们知道一个数可以拆成质因数,我们就直接对这个数分解质因数,然后加入删除可以求逆
-
不过你发现一个问题,就是你多多少少要在莫队复杂度上乘上一个不可忽略的复杂度,但是时限就一秒,不太好过去

- 考虑一下怎么办
- 我们发现枚举质因数影响了时间复杂度,如果我们直接令跑莫队的数都是质数,那就很好办了
- 如何使跑莫队时数是质数呢
- 我们知道一个数在试除法求质因数时枚举上界是 \(\sqrt n\),我们循其源,即若这个数不为质数,那最大的约数一定不超过 \(\sqrt n\)
- 这启发我们根号分治,由于上面的影响是值域的,对值域根号分治
- 我们设序列中最大的数为 \(mx\),那阈值为 \(\sqrt {mx}\),我们直接用线性筛算出在这个阈值内的所有质数,对每个数除去这些质数,剩下的就是大于阈值的质数或 \(1\),大于的数直接跑莫队,其他都是可以被阈值内的质数表示的,这块我们还要处理(顺便提一嘴,可被阈值内质数表示的数的贡献可以和上面找大于阈值的质数一起算)
- 我们开始算阈值内的那一部分
- 上面说一个数的贡献是什么样的了,我们枚举每一个质数 \(g\) 时,先预处理出其的 \(n\) 次方,方便后面算贡献
- 然后枚举 \(g\) 的幂 \(i\),我们算一个前缀,\(s_j\) 表示从 \(1\) 到 \(j\) 有几个满足 \(a_j \mid i\),通过差分统计贡献,最后把 \(a_i\) 中所有 \(g\) 全部除去,方便莫队和后面的计算
点击查看代码
void solve1()
{
lim=sqrt(mx);
Prime(lim);
for (int i=1;i<=mx;i++)
inv[i]=Inv(i);
for (int i=1;i<=pr[0];i++)
{
num[0]=pr[i];
for (int j=1;j<=n;j++)
num[j]=dou(num[j-1]);
for (int j=pr[i];j<=mx;j*=pr[i])
{
for (int k=1;k<=n;k++)
s[k]=s[k-1]+(!(a[k]%j));
for (int k=1;k<=m;k++)
(ans[k]*=num[s[q[k].r]-s[q[k].l-1]]%p*inv[pr[i]]%p)%=p;
}
for (int j=1;j<=n;j++)
while (!(a[j]%pr[i]))
a[j]/=pr[i];
}
}
- 然后是莫队,需要统计出现次数,开桶记就行,同样预处理次方
- 其实你可以将撤销减去的贡献和加入增加的贡献分开求,最后再统一算贡献和逆,好写
- 代码是分块写的方便当时调试,也清晰好看
点击查看代码
#include <bits/stdc++.h>
#define int long long
using namespace std;
constexpr int maxn=1e5+10,p=998244353;
int read()
{
int x=0,f=1;
char ch=getchar();
while (ch<'0' || ch>'9')
{
if (ch == '-') f=-1;
ch=getchar();
}
while (ch>='0' && ch<='9')
{
x=(x<<1)+(x<<3)+ch-'0';
ch=getchar();
}
return x*f;
}
int n,m;
int mx,lim;
int a[maxn];
struct _ {int l,r,id;}q[maxn];
int ans[maxn];
namespace preprocessing
{
int pr[maxn];
int is_no_p[maxn]={1,1};
void Prime(int n)
{
for (int i=2;i<=n;i++)
{
if (!is_no_p[i]) pr[++pr[0]]=i;
for (int j=1;j<=pr[0] && i*pr[j]<=n;j++)
{
is_no_p[i*pr[j]]=1;
if (i%pr[j] == 0) break;
}
}
}
int dou(int x) {return x*x%p;}
int power(int x,int y)
{
int res=1;
while (y)
{
if (y&1) (res*=x)%=p;
(x*=x)%=p;
y>>=1;
}
return res;
}
int Inv(int x) {return power(x,p-2);}
int inv[maxn];
int num[maxn],s[maxn];
void solve1()
{
lim=sqrt(mx);
Prime(lim);
for (int i=1;i<=mx;i++)
inv[i]=Inv(i);
for (int i=1;i<=pr[0];i++)
{
num[0]=pr[i];
for (int j=1;j<=n;j++)
num[j]=dou(num[j-1]);
for (int j=pr[i];j<=mx;j*=pr[i])
{
for (int k=1;k<=n;k++)
s[k]=s[k-1]+(!(a[k]%j));
for (int k=1;k<=m;k++)
(ans[k]*=num[s[q[k].r]-s[q[k].l-1]]%p*inv[pr[i]]%p)%=p;
}
for (int j=1;j<=n;j++)
while (!(a[j]%pr[i]))
a[j]/=pr[i];
}
}
}using namespace preprocessing;
namespace MO
{
int L[maxn],R[maxn],pos[maxn],cnt;
void init()
{
cnt=sqrt(n);
for (int i=1;i<=cnt;i++)
L[i]=R[i-1]+1,R[i]=cnt*i;
if (R[cnt]<n)
cnt++,L[cnt]=R[cnt-1]+1,R[cnt]=n;
for (int i=1;i<=cnt;i++)
for (int j=L[i];j<=R[i];j++)
pos[j]=i;
}
int t[maxn];
int g[maxn<<2],bel[maxn];
int col[maxn];
int as=1,den=1;
void add(int x)
{
if (x == 1) return;
(as*=g[bel[x]+col[x]])%=p;
col[x]++;
}
void del(int x)
{
if (x == 1) return;
col[x]--;
(den*=g[bel[x]+col[x]])%=p;
}
void solve2()
{
for (int i=1;i<=n;i++)
t[a[i]]++;
for (int i=2;i<=mx;i++)
{
if (!t[i]) continue;
bel[i]=++g[0];
g[g[0]]=i;
t[i]--;
while (t[i])
g[0]++,g[g[0]]=dou(g[g[0]-1]),t[i]--;
}
init();
sort(q+1,q+1+m,[](_ &a,_ &b)
{
return pos[a.l] == pos[b.l] ?
(pos[a.l]&1) ? a.r<b.r : a.r>b.r :
pos[a.l]<pos[b.l];
});
int l=1,r=0;
for (int i=1;i<=m;i++)
{
while (l>q[i].l) add(a[--l]);
while (r<q[i].r) add(a[++r]);
while (l<q[i].l) del(a[l++]);
while (r>q[i].r) del(a[r--]);
(ans[q[i].id]*=as*Inv(den)%p)%=p;
}
}
}using namespace MO;
signed main()
{
n=read(),m=read();
for (int i=1;i<=n;i++)
a[i]=read(),mx=max(mx,a[i]);
for (int i=1;i<=m;i++)
{
int l=read(),r=read();
q[i]={l,r,i};
ans[i]=1;
}
solve1();
solve2();
for (int i=1;i<=m;i++)
printf("%lld\n",ans[i]);
return 0;
}
后话
$\mathscr{msjing}$

浙公网安备 33010602011771号