wmq的A×B Problem

wmq的A×B Problem

题目链接：http://oj.xjtuacm.com/problem/13/

题目大意：$T$组数据，每组给出$n$个数$a_i$及一个素数$m$，求这$n$个数两两相乘模$m$余$k$有多少个（$0leqslant k < m$）.

数论+FFT

原根的概念

设$n geqslant 1$，$(a,n)=1$，使得$a^d equiv 1(mod n)$成立的最小的正整数$d$，被称为$a$对模$n$的阶，记做$delta_n(a)$.

当$delta_n(a)=varphi(n)$时，称$a$为模$n$的一个原根.

一些性质

1.质数必存在原根（有缘再证）；

2.若$(a,n)=1$，$a^d equiv 1 (mod n)$，则$delta_n(a)|d$.

    证明：设$d=kdelta_n(a)+r$，其中$0 leqslant r < delta_n(a)$，则有

            $a^d equiv a^{kdelta_n(a)+r} equiv a^r equiv 1(mod n)$.

            由$delta_n(a)$的最小性得，$r=0$.

            $ herefore delta_n(a)|d$.

2.若$(a,n)=1$，$a^k equiv a^l (mod n)$，则$k equiv l(mod delta_n(a))$.

    证明：不妨设$k geqslant l$，则有$a^{k-l} equiv 1 (mod n)$.

           $ecause delta_n(a)|k-l$.

           $ herefore k equiv l(mod delta_n(a))$.

3.若$(a,n)=1$，则$a^0,a^1,...,a^{delta_n(a)-1}$模$n$两两不同余。特别地，当$a$为$n$的原根时，这$varphi(n)$个数是模$n$的一个缩系.

    证明：假设存在$0 leqslant k,l < delta_n(a)$，使得$a^k equiv a^l (mod n)$.

            由2.知，$k equiv l (mod delta_n(a))$.

            $ herefore a^0,a^1,...a^{delta_n(a)-1}$两两不同余。

　　　　 当$a$是模$n$的原根时，$delta_n(a) = varphi(n)$，故这$varphi(n)$个数是模$n$的一个缩系.

---

因为题目中$m$为质数，故必存在原根，设其中一个原根为$r$。

将模$m$为零的数单独处理：$ans[0]=mod[0] imes (n-mod[0])+mod[0] imes (mod[0]-1)/2$;

将所有的模m不为零的数$a_i$写成$r^x$的形式，从而构造出一个关于$r$的多项式：$S=sum_{i=0}^{varphi(m)-1}A_ir^i$.

考虑任意两个数的乘积$r^i imes r^j = r^{i+j}$，将多项式$S$作平方，删去多余项（直接相乘包含了$a_i imes a_i$的结果），即为答案.

复杂度为$O(varphi(m) imes lg(varphi(m)) )$.

代码如下：

#include <cstdio>
#include <cmath>
#include <iostream>
#include <cstring>
#define N 60010
using namespace std;
typedef long long ll;
ll T,n,m,t,r,len,phi;
ll ans[N],mod[N],rt[N],temp[4*N];
const double pi=acos(-1.0);
struct Complex{
    double r,i;
    Complex(double r=0,double i=0):r(r),i(i){};
    Complex operator+(const Complex &rhs){return Complex(r+rhs.r,i+rhs.i);}
    Complex operator-(const Complex &rhs){return Complex(r-rhs.r,i-rhs.i);}
    Complex operator*(const Complex &rhs){return Complex(r*rhs.r-i*rhs.i,i*rhs.r+r*rhs.i);}
}a[4*N],b[4*N],c[4*N];
void sincos(double theta,double &p0,double &p1){
    p0=sin(theta);p1=cos(theta);
}
void fft_main(Complex P[], ll n, ll oper){
    for(ll i=1,j=0;i<n-1;i++){
        for(ll s=n;j^=s>>=1,~j&s;);
        if(i<j)swap(P[i],P[j]);
    }
    Complex unit_p0;
    for(ll d=0;(1<<d)<n;d++){
        ll m=1<<d,m2=m*2;
        double p0=pi/m*oper;
        sincos(p0,unit_p0.i,unit_p0.r);
        for(ll i=0;i<n;i+=m2){
            Complex unit=1;
            for(ll j=0;j<m;j++){
                Complex &P1=P[i+j+m],&P2=P[i+j];
                Complex t=unit*P1;
                P1=P2-t; P2=P2+t;
                unit=unit*unit_p0;
            }
        }
    }if(oper==-1)for(ll i=0;i<len;i++)P[i].r/=len;
}
void fft(Complex a[],Complex b[],ll len){
    fft_main(a,len,1);fft_main(b,len,1);
    for(ll i=0;i<len;++i)c[i]=a[i]*b[i];
    fft_main(c,len,-1);
}
void fft_init(){
    len=1;
    while(len<2*phi)len<<=1;
    ll i=0;
    for(ll x=1;i<phi;++i){
        rt[i]=mod[x];
        a[i].r=b[i].r=rt[i];
        a[i].i=b[i].i=0;
        x=(x*r)%m;
    }
    for(;i<len;++i){
        a[i].r=b[i].r=0;
        a[i].i=b[i].i=0;
    }
}
ll powmod(ll a,ll n,ll m){
    ll r=1,t=a;
    while(n){
        if(n&1)r=(r*t)%m;
        t=(t*t)%m;
        n>>=1;
    }
    return r;
}
ll Root(ll m){
    for(ll i=2,j;i<m;++i){
        for(j=2;j<phi;++j)
            if(powmod(i,j,m)==1)
                break;
        if(j==phi)return i;
    }
    return 0;
}
ll in(){
    ll res=0,flag=0,ch;
    if((ch=getchar())=='-')flag=1;
    else if('0'<=ch&&ch<='9')res=ch-'0';
    while('0'<=(ch=getchar())&&ch<='9')res=res*10+ch-'0';
    return flag?-res:res;
}
void out(ll x){
    if(x>9)out(x/10);
    putchar(x%10+'0');
}
void solve(){
    for(ll i=0;i<len;++i){
        temp[i]=c[i].r+0.5;
        temp[i]/=2;
    }
    for(ll i=0;i<phi;++i){
        temp[2*i]-=rt[i]*rt[i]/2;
        temp[2*i]+=rt[i]*(rt[i]-1)/2;
    }
    for(ll i=0,x=1;i<len;++i){
        ans[x]+=temp[i];
        x=(x*r)%m;
    }
    ans[0]=mod[0]*(n-mod[0])+mod[0]*(mod[0]-1)/2;
}
void init(){
    memset(temp,0,sizeof(temp));
    n=in();m=in();
    for(ll i=0;i<m;++i)mod[i]=0;
    for(ll i=0;i<m;++i)ans[i]=0;
    for(ll i=0;i<n;++i){
        t=in();
        mod[t%m]++;
    }
}
int main(void){
    T=in();
    while(T--){
        init();
        phi=m-1;
        r=Root(m);
        fft_init();
        fft(a,b,len);
        solve();
        for(ll i=0;i<m;++i){
            out(ans[i]);
            puts("");
        }
    }
}