2018秦皇岛ccpc-camp Steins;Gate (原根+FFT)

因为给定的模数P保证是素数,所以P一定有原根.
根据原根的性质,若(g)(P)的原根,则(g^k)能够生成([1,P-1])中所有的数,这样的k一共有P-2个.
(a_i*a_j(mod P)=a_k) 就可以转化为(g^i*g^j(mod P) = g^{i+j}(mod P)=g^k).
问题转化为了求有多少对有序的<i,j>满足 ((i+j)(mod (P-1)) = k).
求出原根后,对([1,P-1])中的每个数编号, 统计每个编号出现的次数,然后FFT求卷积
要特判0,因为原根不会生成0.所以用总的有序对数-其他不含0的有序对数得到含0的有序对,这是0的答案;超过P-1的数肯定没有符合的有序对.

#include <bits/stdc++.h>
using namespace std;
typedef long long LL;
const int maxn = 2e5+5;
bool notpri[maxn];
int pri[maxn],zyz[maxn];
typedef long long LL;

void pre(int N){
    notpri[1]=1;
    for(int i=2;i<N;++i){
        if(!notpri[i]){
            pri[++pri[0]]=i;
        }
        for(int j=1;j<=pri[0] && (LL)i*(LL)pri[j]<N;++j){
            notpri[i*pri[j]]=1;
            if(i%pri[j]==0){
                break;
            }
        }
    }
}
LL Quick_Pow(LL x,LL p,LL mod){
    if(!p){
        return 1ll;
    }
    LL res=Quick_Pow(x,p>>1,mod);
    res=res*res%mod;
    if((p&1ll)==1ll){
        res=(x%mod*res)%mod;
    }
    return res;
}
int FindRoot(int x){/*求素奇数的最小原根,倘若x不是奇数,但是也有原根的话,将质
因子分解改成对phi(x)即可。倘若要求多个原根,直接接着暴力验证即可*/
    int tmp=x-1;
    for(int i=1;tmp && i<=pri[0];++i){
        if(tmp%pri[i]==0){
            zyz[++zyz[0]]=pri[i];
            while(tmp%pri[i]==0){
                tmp/=pri[i];
            }
        }
    }
    for(int g=2;g<=x-1;++g){
        bool flag=1;
        for(int i=1;i<=zyz[0];++i){
            if(Quick_Pow((LL)g,(LL)((x-1)/zyz[i]),(LL)x)==1){
                flag=0;
                break;
            }
        }
        if(flag){
            return g;
        }
    }
    return 0;
}

const int MAXN = 4e5 + 10;
const double PI = acos(-1.0);
struct Complex{
    double x, y;
    inline Complex operator+(const Complex b) const {
        return (Complex){x +b.x,y + b.y};
    }
    inline Complex operator-(const Complex b) const {
        return (Complex){x -b.x,y - b.y};
    }
    inline Complex operator*(const Complex b) const {
        return (Complex){x *b.x -y * b.y,x * b.y + y * b.x};
    }
} va[MAXN * 2 + MAXN / 2], vb[MAXN * 2 + MAXN / 2];
int lenth = 1, rev[MAXN * 2 + MAXN / 2];
int N, M;   // f 和 g 的数量
    //f g和 的系数
    // 卷积结果
    // 大数乘积
int f[MAXN],g[MAXN];
vector<LL> conv;
vector<LL> multi;
//f g
void init()
{
    int tim = 0;
    lenth = 1;
    conv.clear(), multi.clear();
    memset(va, 0, sizeof va);
    memset(vb, 0, sizeof vb);
    while (lenth <= N + M - 2)
        lenth <<= 1, tim++;
    for (int i = 0; i < lenth; i++)
        rev[i] = (rev[i >> 1] >> 1) + ((i & 1) << (tim - 1));
}
void FFT(Complex *A, const int fla)
{
    for (int i = 0; i < lenth; i++){
        if (i < rev[i]){
            swap(A[i], A[rev[i]]);
        }
    }
    for (int i = 1; i < lenth; i <<= 1){
        const Complex w = (Complex){cos(PI / i), fla * sin(PI / i)};
        for (int j = 0; j < lenth; j += (i << 1)){
            Complex K = (Complex){1, 0};
            for (int k = 0; k < i; k++, K = K * w){
                const Complex x = A[j + k], y = K * A[j + k + i];
                A[j + k] = x + y;
                A[j + k + i] = x - y;
            }
        }
    }
}
void getConv(){             //求多项式
    init();
    for (int i = 0; i < N; i++)
        va[i].x = f[i];
    for (int i = 0; i < M; i++)
        vb[i].x = g[i];
    FFT(va, 1), FFT(vb, 1);
    for (int i = 0; i < lenth; i++)
        va[i] = va[i] * vb[i];
    FFT(va, -1);
    for (int i = 0; i <= N + M - 2; i++)
        conv.push_back((LL)(va[i].x / lenth + 0.5));
}

LL vz[MAXN];
int id[MAXN];
int cnt[MAXN];
int rnk[MAXN];
LL ans[MAXN];

int main()
{
    #ifndef ONLINE_JUDGE
        freopen("in.txt","r",stdin);
        freopen("out.txt","w",stdout);
    #endif
    pre(maxn-1);
    int n,P ;
    scanf("%d %d",&n, &P);
    int rt = FindRoot(P);
    memset(id,0,sizeof(id));


    LL idx = 1;
    for(int i=0;i<P-1;++i){         //一共只有[0,P-2],P-1个id
        id[idx] = i;
        rnk[i]= idx;
        idx = idx*rt%P;
    }

    for(int i=1;i<=n;++i){
        LL tmp;
        scanf("%lld",&tmp);
        vz[i] = tmp;
        tmp%=P;
        if(tmp==0) continue;        //卷积中不考虑0的贡献
        cnt[id[tmp]]++;
    }

    N = M = P;
    for(int i=0;i<P-1;++i){
        f[i] = g[i] = cnt[i];
    }

    getConv();
    int sz = conv.size();

    for(int i=0;i<sz;++i){
        LL tmp = conv[i];
        ans[rnk[i%(P-1)]] += tmp;
    }

    LL tot = (LL)n*n;           //全部的枚举可能-不选0之外的组合 = 包含0的组合
    for(int i=1;i<P;++i){
        tot -= ans[i];
    }
    ans[0] = tot;
    for(int i=1;i<=n;++i){
        if(vz[i]>=P) printf("0
");
        else{
            printf("%lld
",ans[vz[i]]);
        }
    }
    return 0;
}

原文地址:https://www.cnblogs.com/xiuwenli/p/9735804.html