P5071 此时此刻的光辉
题面
题解
小清新题。
首先,这个题要求你会 Pollard-Rho,如果不会可参考 OI Wiki 中相关文章。
一眼想到莫队,然后发现直接维护的话复杂度带 \(\log\),无法通过。
然后你发现可以阈值分治,我们将 \(1000\) 以内的质数(共 \(168\) 个)用前缀和维护出现次数,然后发现每个 \(a_i\) 中只有至多 \(2\) 个大于 \(1000\) 的质因数。
这样,复杂度变为 \(O(n\sqrt n+n\sqrt[4]{a_i})\),其中 \(O(n\sqrt[4]{a_i})\) 带较大常数,可以压线通过。
Code
#include<bits/stdc++.h>
using namespace std;
const int p=19260817;
namespace pm{
int d[]={2,3,5,7,11,13};
int pw(int a,int b,int p){
int res=1;
while(b){
if(b&1) res=(long long)(res)*a%p;
a=(long long)(a)*a%p;
b>>=1;
}
return res;
}
bool is_prime(int x){
for(int i=0;i<6;i++){
int u=d[i];
if(x==u) return true;
if(x<u) return false;
if(pw(u,x-1,x)!=1) return false;
int pos=x-1;
while(!(pos&1)) pos>>=1;
int now=pw(u,pos,x);
while(pos!=x-1){
if(now==x-1) break;
if(now==1) break;
if((long long)(now)*now%x==1) return false;
now=(long long)(now)*now%x;
pos<<=1;
}
}
return true;
}
vector<int>ans;
mt19937 rnd(6527);
void fact(int u){
if(u==1) return;
if(is_prime(u)){
ans.push_back(u);
return;
}
if(u==4){
ans.push_back(2); ans.push_back(2); return;
}
if(u==6){
ans.push_back(2); ans.push_back(3); return;
}
if(u==8){
ans.push_back(2); ans.push_back(2); ans.push_back(2); return;
}
if(u==9){
ans.push_back(3); ans.push_back(3); return;
}
while(1){
int c=rnd()%u;
int now1=1,now2=1;
int cnt=0;
int mul=1;
while(1){
int new1=((long long)(now1)*now1+c)%u;
int new2=((long long)(now2)*now2+c)%u;
new2=((long long)(new2)*new2+c)%u;
cnt++; mul=1ll*mul*abs(new1-new2)%u;
if(!(cnt&15)){
int h;
if((h=__gcd(mul,u))!=1){
if(h==u) break;
fact(h);
fact(u/h);
return;
}
mul=1;
}
now1=new1; now2=new2;
}
}
}
}
int n,q;
int a[100005];
int b1[100005],b2[100005];
int ps[200005];
int nh[200005];
int prm[170];
int ny[200005];
vector<int>fact(int u){
vector<int>ans;
for(int j=1;j<=168;j++){
int i=prm[j];
while(u%i==0){
u/=i; ans.push_back(i);
}
}
if(u!=1){
pm::ans.clear(); pm::fact(u); for(auto d:pm::ans){
ans.push_back(d);
}
sort(ans.begin(),ans.end());
}
return ans;
}
struct qry{
int l,r,ans,id;
};
qry s[100005];
bool cmp1(qry x,qry y){
if(x.l/400!=y.l/400) return x.l<y.l;
if(x.l/400%2) return x.r>y.r;
return x.r<y.r;
}
bool cmp2(qry x,qry y){
return x.id<y.id;
}
int sum[100005][170];
int tmp[1005];
void sol1(){
int pos=0;
for(int i=2;i<=1000;i++){
if(pm::is_prime(i)) {
prm[++pos]=i; tmp[i]=pos;
}
}
pos=0;
for(int i=1;i<=n;i++){
vector<int>g=fact(a[i]);
if(!g.empty()&&g.back()>1000){
int u=g.back(); g.pop_back();
b1[i]=u; ps[++pos]=u;
}
if(!g.empty()&&g.back()>1000){
int u=g.back(); g.pop_back();
b2[i]=u; ps[++pos]=u;
}
for(auto u:g){
sum[i][tmp[u]]++;
}
}
sort(ps+1,ps+1+pos); pos=unique(ps+1,ps+1+pos)-ps-1;
for(int i=1;i<=n;i++){
if(b1[i]) b1[i]=lower_bound(ps+1,ps+1+pos,b1[i])-ps;
if(b2[i]) b2[i]=lower_bound(ps+1,ps+1+pos,b2[i])-ps;
}
for(int i=1;i<=n;i++){
for(int j=1;j<=168;j++){
sum[i][j]+=sum[i-1][j];
}
}
for(int i=1;i<=q;i++){
for(int j=1;j<=168;j++){
s[i].ans=1ll*s[i].ans*(sum[s[i].r][j]-sum[s[i].l-1][j]+1)%p;
}
}
}
int ans;
void add(int id){
if(b1[id]){
ans=1ll*ans*ny[nh[b1[id]]]%p;
nh[b1[id]]++;
ans=1ll*ans*nh[b1[id]]%p;
}
if(b2[id]){
ans=1ll*ans*ny[nh[b2[id]]]%p;
nh[b2[id]]++;
ans=1ll*ans*nh[b2[id]]%p;
}
}
void del(int id){
if(b1[id]){
ans=1ll*ans*ny[nh[b1[id]]]%p;
nh[b1[id]]--;
ans=1ll*ans*nh[b1[id]]%p;
}
if(b2[id]){
ans=1ll*ans*ny[nh[b2[id]]]%p;
nh[b2[id]]--;
ans=1ll*ans*nh[b2[id]]%p;
}
}
void sol2(){
for(int i=1;i<=n*2;i++) ny[i]=pm::pw(i,p-2,p);
ans=1;
int nl=1,nr=0;
sort(s+1,s+1+q,cmp1);
for(int i=1;i<=q;i++){
while(nl>s[i].l) add(--nl);
while(nr<s[i].r) add(++nr);
while(nl<s[i].l) del(nl++);
while(nr>s[i].r) del(nr--);
s[i].ans=1ll*s[i].ans*ans%p;
}
}
int main(){
ios::sync_with_stdio(false); cin.tie(0); cout.tie(0);
cin>>n>>q;
for(int i=1;i<=n;i++) cin>>a[i];
for(int i=1;i<=n*2;i++) nh[i]++;
for(int i=1;i<=q;i++){
s[i].id=i; s[i].ans=1;
cin>>s[i].l>>s[i].r;
}
sol1();
sol2();
sort(s+1,s+1+q,cmp2);
for(int i=1;i<=q;i++) cout<<s[i].ans<<'\n';
return 0;
}