【GDOI2018模拟7.9】相逢是问候

Description

Solution

首先我们要知道一个叫做扩展欧拉定理的东西：

cxmodp≡cxmodφ(p)+φ(p)∗[x>=φ(p)]modp

如果我们要求

ccx

那么就是

ccxmodp≡ccxmodφ(p)+φ(p)∗[x>=φ(p)]modp

然后再把上面的

cxmodφ(p)展开

ccxmodp≡ccxmodφ(φ(p))+[x>=φ(φ(p))]∗(φ(φ(p)))modφ(p)+φ(p)∗[x>=φ(p)]modp

我们可以发现，迭代越多次，

φ就重复的越多。

推论：预计迭代log(p)次就可以变成1

证明：偶数的

φ肯定减一半，奇数就姑且让它减一，最后一定在log此变成1

可以发现，最多不超过26次。当

φ迭代出1的时候，很显然，就变成了一个固定的值，这个很好证明。

那么，我们可以预处理一个f[1~26][i]表示迭代这么多次的时候的值是多少，这个可以nlog^3来做（包括快速幂的复杂度）

然后直接用线段树做就好了，维护一个最小值，修改的时候，如果最小值是<=26的话，那么就暴力修改，每个点最多只会被修改26次。

预处理超时！！！

我们可以优化一下快速幂的复杂度。
因为快速幂的底数是相同的，所以我们可以把幂拆成一半。
前15个和后15个，然后预处理一下，求快速幂的时候直接拼起来就好了。

注意一个小细节

预处理的时候，处理cx的次方的时候，要先判断x的原值是否大于φ(p)，然后x再模原本要模的数。所以在快速幂的预处理的时候，还要开一个数组去储存它是否>φ(p)。

Code

#include<iostream>
#include<stdio.h>
#include<string.h>
#include<algorithm>
#include<math.h>
#define fo(i,a,b) for(i=a;i<=b;i++)
using namespace std;
typedef long long ll;
const int maxn=5e4+7;
ll i,j,k,l,n,m,ans,c,mo,x,y,z,lim,sh;
ll t[maxn*4],a[maxn],f[maxn][30],pp[30],d[maxn*4],mi[maxn*4];
ll b[30][maxn*2],e[30][maxn*2],er[31],u[30][maxn*2],v[30][maxn*2];
ll w[30];
ll phi(ll x){
    ll i,j=x;
    for(i=2;i<=x;i++){
        if(x%i==0)j=j*(i-1)/i;
        while(x%i==0)x/=i;
    }
    if(x!=1)j=j*(x-1)/x;
    return j;
}
ll qsm(ll x,ll y,ll mo){
    ll z=1;
    for(;y;y/=2,x=x*x%mo)if(y&1)z=z*x%mo;
    return z;
}
ll qsm1(ll y,ll o){
    ll p=y/er[15],q=y%er[15],z;
    z=b[o][q]*e[o][p]%pp[o];
    sh=1;if(!p&&!u[o][q])sh=0;
    return z;
}
ll qsm2(ll x,ll y){
    ll z=1;
    for(;y;y/=2,x=x*x)if(y&1)z=z*x;
    return z;
}
void gao(int x,int p){
    if(x==27)return;
    pp[x]=p;
    gao(x+1,phi(p));
}
ll dfs(int x,int y){
    if(y==lim+1){
        sh=1;if(x<pp[y-1])sh=0;
        return x;
    }
    ll o=dfs(x,y+1),q;
    q=qsm1(o%pp[y]+pp[y]*sh,y-1);
    return q;
}
ll dfs1(int x,int y){
    if(y==lim+1)return x;
    ll o=dfs1(x,y+1),q;
    q=qsm2(c,o);
    return q;
}
void merge(int x){
    t[x]=(t[x*2]+t[x*2+1])%mo;
    mi[x]=min(mi[x*2],mi[x*2+1]);
}
void build(int x,int l,int r){
    if(l==r){
        t[x]=f[l][0];
        return;
    }
    int mid=(l+r)/2;
    build(x*2,l,mid);build(x*2+1,mid+1,r);
    merge(x);
}
void gao(int x,int l,int r){
    if(mi[x]>=26)return;
    if(l==r){d[l]++;mi[x]++;t[x]=f[l][d[l]];return;}
    int mid=(l+r)/2;
    gao(x*2,l,mid),gao(x*2+1,mid+1,r);
    merge(x);
}
void change(int x,int l,int r,int y,int z){
    if(l==y&&r==z){gao(x,l,r);mi[x]++;return;}
    int mid=(l+r)/2;
    if(z<=mid)change(x*2,l,mid,y,z);
    else if(y>mid)change(x*2+1,mid+1,r,y,z);
    else change(x*2,l,mid,y,mid),change(x*2+1,mid+1,r,mid+1,z);
    merge(x);
}
ll find(int x,int l,int r,int y,int z){
    if(l==y&&r==z)return t[x];
    int mid=(l+r)/2;
    if(z<=mid)return find(x*2,l,mid,y,z);
    else if(y>mid)return find(x*2+1,mid+1,r,y,z);
    else return (find(x*2,l,mid,y,mid)+find(x*2+1,mid+1,r,mid+1,z))%mo;
    merge(x);
}
int main(){
//  freopen("fan.in","r",stdin);
//  freopen("fan.out","w",stdout);
    scanf("%d%d%d%d",&n,&m,&mo,&c);
    gao(0,mo);
    er[0]=1;fo(i,1,30)er[i]=er[i-1]*2;
    fo(j,0,26){
        b[j][0]=e[j][0]=1;if(b[j][0]>=pp[j])u[j][0]=1;
        if(!pp[j])break;
        fo(i,1,er[15]-1){
            b[j][i]=b[j][i-1]*c%pp[j],e[j][i]=qsm(b[j][i],er[15],pp[j]);
            if(b[j][i-1]*c>=pp[j]||u[j][i-1])u[j][i]=1;
        }
    }
    fo(lim,0,26)w[lim]=dfs(c,1);
    fo(i,1,n){
        scanf("%d",&a[i]);
        fo(lim,0,26)
        f[i][lim]=dfs(a[i],1);
    }
    build(1,1,n);
    while(m--){
        scanf("%d%d%d",&z,&x,&y);
        if(!z)change(1,1,n,x,y);
        else printf("%lld\n",find(1,1,n,x,y));
    }
}