【原根+FFT】牛客国庆Day1 I Steins;Gate

最新推荐文章于 2021-03-25 13:16:01 发布

Toooooocold

最新推荐文章于 2021-03-25 13:16:01 发布

阅读量371

点赞数

分类专栏：数论-FFT

本文链接：https://blog.csdn.net/qq_32506797/article/details/82951628

版权

数论-FFT 专栏收录该内容

2 篇文章 0 订阅

订阅专栏

该博客介绍了如何利用原根和快速傅里叶变换(FFT)来解决一个数值计算问题：在N个数a1, a2, ...,aN中，找出满足a[i]*a[j] ≡ a[k] (mod P)的有序对(i, j)，其中P为质数。通过将乘法转换为加法，利用FFT加速计算过程。" 115641126,10429134,VLAN配置：基于IP地址的网络通信分析,"['网络配置', '交换机', 'VLAN划分', '网络通信']

摘要由CSDN通过智能技术生成

$S o u c e :$ 牛客国庆集训派对Day1
$P r o b l e m :$
$N = 1 e 5 ， a 1, a 2, . . ., a N$ ，对于每一个 $a k$ ，求有多少个有序二元组 $(i, j)$ 满足 $a [i] * a [j] = a [k] (m o d P)$ ，其中 $P$ 为一给定质数。
$I d e a :$
1到P-1都可以用P的某一原根的幂次进行表示，那么乘法转化为加法，用FFT加速。
$C o d e :$

#include<bits/stdc++.h>
using namespace std;

#define I inline
#define pb push_back
typedef double DB;
typedef long long LL;
typedef vector<int> VI;

const int N = 2e5+10;

I int qpow(int a, int b, int P) {
    int ret = 1;
    while(b) {
        if(b&1) ret = 1LL*ret*a%P;
        a = 1LL*a*a%P; b >>= 1;
    }
    return ret;
}

bool isn[N]; int p[N];
I void init() {
    int pnt = 0; isn[1] = 1;
    for(int i = 2; i < N; i++) {
        if(!isn[i]) p[pnt++] = i;
        for(int j = 0; j<pnt && i*p[j]<N; j++) {
            isn[i*p[j]] = 1;
            if(i%p[j] == 0) break;
        }
    }
}

I int findroot(int P) {
    static VI a; a.clear();
    int t = P-1;
    for(int i = 0; p[i]*p[i] <= t; i++) if(t%p[i] == 0) {
        a.pb((P-1)/p[i]);
        while(t%p[i] == 0) t /= p[i];
    }
    if(t > 1) a.pb((P-1)/t);
    for(int g = 2; g < P; g++) {
        bool f = 1;
        for(int &t:a) if(qpow(g, t, P) == 1) { f = 0; break; }
        if(f) return g;
    }
    return -1;
}

I namespace FFT {
    const DB PI = acos(-1.0);
    struct CP {
        DB x,y;
        CP(DB x=0, DB y=0):x(x),y(y){}
        CP operator-(const CP &b) const { return CP(x-b.x, y-b.y); }
        CP operator+(const CP &b) const { return CP(x+b.x, y+b.y); }
        CP operator*(const CP &b) const { return CP(x*b.x-y*b.y, x*b.y+y*b.x); }
    };

    int rev[N<<2];
    I void FFT(CP *A, int n, int op) {
        for(int i = 0; i < n; i++) if(i < rev[i]) swap(A[i], A[rev[i]]);
        for(int i = 1; i < n; i<<=1) {
            CP wn(cos(PI/i), op*sin(PI/i));
            for(int j = 0, R = i<<1; j < n; j += R) {
                CP w(1,0);
                for(int k = 0; k < i; k++, w = w*wn) {
                    CP x = A[j+k], y = w*A[j+i+k];
                    A[j+k] = x+y; A[j+i+k] = x-y;
                }
            }
        }
        if(op == -1) for(int i = 0; i < n; i++) A[i].x /= n;
    }

    CP A[N<<2], B[N<<2];
    I void mul(LL *a, int la, LL *b, int lb) {
        int len = 1; while(len < la+lb) len <<= 1;
        for(int i = 0; i < len; i++) {
            rev[i] = ((rev[i>>1]>>1)|(len>>(i&1)))&(len-1);
            A[i] = a[i], B[i] = b[i];
        }
        FFT(A, len, 1); FFT(B, len, 1);
        for(int i = 0; i < len; i++) A[i] = A[i]*B[i];
        FFT(A, len, -1);
        for(int i = 0; i < len; i++) a[i] = (LL)(A[i].x+0.5);
    }
}

int id[N];
LL a[N<<2], arr[N];

I void work() {
    int n, P; scanf("%d%d", &n, &P);
    int g = findroot(P), t = g;
    for(int i = 1; i <= P-1; i++) { id[t] = i; t = 1LL*t*g%P; }
    LL cnt = 0;
    for(int i = 1; i <= n; i++) {
        int x; scanf("%d", &x);
        arr[i] = x; x %= P;
        if(x == 0) cnt++; else a[id[x]]++;
    }
    mul(a, P, a, P);
    for(int i = P; i <= 2*P-2; i++) a[i-P+1] += a[i];
    cnt = 2*cnt*(n-cnt)+cnt*cnt;
    for(int i = 1; i <= n; i++) {
        if(arr[i] >= P) { puts("0"); continue; }
        arr[i] %= P;
        if(arr[i] == 0) printf("%lld\n", cnt);
        else printf("%lld\n", a[id[arr[i]]]);
    }
}

int main() {
    init(); work();
    return 0;
}