【CF1103E】Radix sum

foreverlasting

2019-02-21 08:22:40

Solution

[同步发表在我的博客里哦](https://foreverlasting1202.github.io/) 多项式。 神题。 首先拿到题目,看懂题目后,你就会知道这题要怎么做。显然的卷积。但问题在于它给出了一个神仙模数$2^{58}$。然后我就不会了,因此当时就咕了这道题。不过后来[这位神仙](https://www.cnblogs.com/Mr-Spade/p/10390667.html)给出了一个中文题解,研究了一下,突然找到了做法。不过还是好神仙啊! 我们考虑$DFT$和$IDFT$我们要做什么。 $1$、找到单位根,这里是$\omega_{10}$。 $2$、找到$10^5$的逆元,这里是在$2^{58}$下的逆元。 首先看向第二个问题,这里看起来比较好解决。首先$5$在$2^{58}$下是有逆元的,所以只用考虑$2^5$的逆元。这里就有一步神仙转换,即发现$x\ast 2^5$在$2^{64}$下的结果可以直接计算,假设结果为$y$。 那我们就有$\frac {y}{2}=x\ast 2^4$在模$2^{63}$下。$\frac {y}{2^5}=x$在模$2^{59}$下。而模数是$2^{58}$,这没关系,再取个模就可以了。 现在再去处理第一个问题。我们可以发现的是$\omega_{10}^5$在模$2^{64}$下是有意义的,且答案为$-1$。并且有$\omega_{10}^2=\omega_{5}$。于是得到了$\omega_{5}^5-1=0$。 这有什么用呢?于是我们把模$x^5-1$的多项式环作为一种数据类型,以此来做乘法等操作。(这里我和$Mr-Spade$的思路有点不同) 这里先讲一下怎么做乘法,再证明一下正确性。 做乘法就是把两个最高次为$x^4$的多项式相乘,然后模$x^5-1$即可。于是$x^8$就变成$+x^3$,$x^7$就变成$+x^2$以此类推,用个数组记录一下,然后像$FFT$原理那样处理一下单位根就好了。 最后考虑正确性。其实在模$x^5-1$下直接输出常数项是错误的。我们按照$DFT$的正确性,$a_1+a_2\ast x+a_3\ast x^2+a_4\ast x^3+a_5\ast x^4$应该是个整数,但事实上原方程$x^5-1=0$的根中$\omega_5\omega_5^2\omega_5^3\omega_5^4$和$1$都可能令这个多项式为整数。那怎么处理呢?试根显然太麻烦,我们考虑证明$a_2=a_3=a_4=a_5$。首先说明为什么要证明这个。因为如果这个成立,我们只要再模$x-1$就可以消除其他系数,只留下常数项。 这个怎么证明呢。首先我们要知道$1$不是解。因为如果$1$是解,$\omega_{10}=\pm 1$。这个显然错误的。 我们令$a_1+a_2\ast x+a_3\ast x^2+a_4\ast x^3+a_5\ast x^4=y$,则$y$为整数。我们假设有四个根分别为$x_1,x_2,x_3,x_4$,移项因式分解得到 $k\ast \prod_{i=1}^4 (x-x_i)=0$ 拆掉,就会得到 $a_5=k$ $a_4=-k\ast \sum_{i=1}^4 x_i$ $a_3=k\ast \sum_{i=1}^4 \sum_{j=i+1}^4 x_i\ast x_j$ $a_2=-k\ast \sum_{i=1}^4 \sum_{j=i+1}^4 \sum_{l=j+1}^4 x_i\ast x_j\ast x_l$ 当你去将上述四个根带入检验的时候,根据$x^5-1=0$这条方程,你会得到$a_5=a_4=a_3=a_2=k$。因此就得知了正确性。 这道神仙题最终还是做完了。 code: ```cpp //2019.2.20 by ljz #include<bits/stdc++.h> using namespace std; #define res register int #define LL long long #define inf 0x3f3f3f3f #define eps 1e-10 inline int read() { res s=0; bool w=0; char ch=getchar(); while(ch<'0'||ch>'9') { if(ch=='-')w=1; ch=getchar(); } while(ch>='0'&&ch<='9')s=s*10+ch-'0',ch=getchar(); return w?-s:s; } inline void _swap(res &x,res &y) { x^=y^=x^=y; } inline int _abs(const res &x) { return x>0?x:-x; } inline int _max(const res &x,const res &y) { return x>y?x:y; } inline int _min(const res &x,const res &y) { return x<y?x:y; } #define ull unsigned LL const int pw[]={1,10,100,1000,10000,100000,1000000}; const ull INV5=14757395258967641293ull; const ull kcz=1ull<<58; namespace MAIN { int n; struct cp{ ull a[5]; cp() {a[0]=a[1]=a[2]=a[3]=a[4]=0;} }w10[10],a[100000+10]; inline cp operator + (const cp &x,const cp &y){ cp z; for(res i=0;i<=4;i++)z.a[i]=x.a[i]+y.a[i]; return z; } inline cp operator - (const cp &x,const cp &y){ cp z; for(res i=0;i<=4;i++)z.a[i]=x.a[i]-y.a[i]; return z; } inline cp operator * (const cp &x,const ull &y){ cp z; for(res i=0;i<=4;i++)z.a[i]=x.a[i]*y; return z; } inline cp operator * (const cp &x,const cp &y){ ull b[9]; for(res i=0;i<=8;i++)b[i]=0; for(res i=0;i<=4;i++) for(res j=0;j<=4;j++)b[i+j]+=x.a[i]*y.a[j]; for(res i=3;~i;i--)b[i]+=b[i+5],b[i+5]=0; cp z; for(res i=0;i<=4;i++)z.a[i]=b[i]; return z; } inline void pre(){ w10[0].a[0]=1,w10[1].a[3]=-1; for(res i=2;i<10;i++)w10[i]=w10[i-1]*w10[1]; } inline void FFT(cp *A,const res &type,const res &lim){ for(res mid=0;mid<5;mid++){ for(res i=0;i<lim;i++) if((i/pw[mid])%10==0){ cp b[10]; for(res j=0;j<10;j++) for(res k=0;k<10;k++) b[j]=b[j]+A[i+k*pw[mid]]*w10[(type*j*k%10+10)%10]; for(res j=0;j<10;j++) A[i+j*pw[mid]]=b[j]; } } if(type==-1){ ull INV=INV5*INV5*INV5*INV5*INV5; for(res i=0;i<lim;i++)A[i]=A[i]*INV; } } inline cp qpow(cp x,res y){ cp ret=w10[0]; while(y){ if(y&1)ret=ret*x; x=x*x,y>>=1; } return ret; } inline void MAIN(){ pre(); n=read(); for(res i=1;i<=n;i++)a[read()].a[0]++; FFT(a,1,100000); for(res i=0;i<100000;i++)a[i]=qpow(a[i],n); FFT(a,-1,100000); for(res i=0;i<n;i++)printf("%llu\n",((a[i].a[0]-a[i].a[4])/32)%kcz); } } int main() { MAIN::MAIN(); return 0; } ```