题目描述：给出一个长度为 (n) 的序列 (a_i)，对于每一个 (pin [0,n)cap )，求有多少个长为 (n) 的整数序列 (i_j) 满足：(i_jin [1,n])，且在十进制不进位加法下 (sum_{j=1}^na_{i_j}=p)。答案对 (2^{58}) 取模。
数据范围：(nle 10^5,0le x_i<10^5)
实际上就相当于对于一个幂级数 (F)，在 (5) 维长度为 (10) 的循环卷积意义下，求 (F^n)。
于是考虑如何做 dft 和 idft。这个跟平常的 ntt 差别实际上在于：
- 模 (2^{58}) 意义下除以 (10^5)？
- (omega_{10}) 怎么办？
首先看第一个问题，我们把模数改为 (2^{64})，也就是 ULL 自然溢出。(5^5) 就可以算逆元，(2^5) 显然就可以直接除掉。
关于第二个问题，众所周知 (omega_{10}) 是无理数，这个模数还这么丑，于是就只能暴力扩域。设 (x=omega_5)，(omega_{10}=-x^3)，则可以将用到的数都表示为 (x) 的多项式。
大家都知道余数定理：设 (F,G) 为两个多项式，(G(a)=0)，则 ((Fmod G)(a)=F(a))。我们取 (G) 为五次分圆多项式 (1+x+x^2+x^3+x^4)。于是在系数模 (2^{64}) ，整体模 (G) 的意义下，得到的三次多项式将 (omega_5) 代入就得到答案。
剩下的问题在于，如何求出最后的答案。我们知道答案都为整数，然而我们得到的是 (a_0+a_1omega_5+a_2omega_5^2+a_3omega_5^3)。因为它是整数，所以可以设它为 (a_0-y)，则 (omega_5) 为 (y+a_1x+a_2x^2+a_3x^3=0) 的根。根据分圆多项式在有理数域上的不可约性，(G) 是被 ((x-omega_5)) 整除的，次数最小的，首1的有理多项式，所以 (G|(y+a_1x+a_2x^2+a_3x^3))，所以 (a_1=a_2=a_3=y=0)，答案即为 (a_0)。
#include<bits/stdc++.h>
#define Rint register int
using namespace std;
typedef unsigned long long LL;
const int N = 100000, po10[5] = {1, 10, 100, 1000, 10000};
const LL inv = 6723469279985657373ull, mod = (1ull << 58) - 1;
template<typename T>
inline void read(T &x){
int ch = getchar(); x = 0; bool f = false;
for(;ch < '0' || ch > '9';ch = getchar()) f |= ch == '-';
for(;ch >= '0' && ch <= '9';ch = getchar()) x = x * 10 + ch - '0';
if(f) x = -x;
}
struct comp {
LL a[4];
comp(){memset(a, 0, sizeof a);}
comp operator = (const comp &o){memcpy(a, o.a, sizeof a); return *this;}
comp operator * (LL val) const {
comp res;
for(Rint i = 0;i < 4;++ i) res.a[i] = a[i] * val;
return res;
}
comp operator + (const comp &o) const {
comp res;
for(Rint i = 0;i < 4;++ i) res.a[i] = a[i] + o.a[i];
return res;
}
comp operator - (const comp &o) const {
comp res;
for(Rint i = 0;i < 4;++ i) res.a[i] = a[i] - o.a[i];
return res;
}
comp operator * (const comp &o) const {
comp res; static LL x[7]; memset(x, 0, sizeof x);
for(Rint i = 0;i < 4;++ i)
for(Rint j = 0;j < 4;++ j) x[i + j] += a[i] * o.a[j];
for(Rint i = 2;~i;-- i){
for(Rint j = 0;j < 4;++ j) x[i + j] -= x[i + 4];
x[i + 4] = 0;
}
res.a[0] = x[0]; res.a[1] = x[1]; res.a[2] = x[2]; res.a[3] = x[3];
return res;
}
} w10[11], A[N];
int n;
inline comp ksm(comp a, int b){
comp res; res.a[0] = 1;
for(;b;b >>= 1, a = a * a) if(b & 1) res = res * a;
return res;
}
void dft(comp *A){
static comp B[10];
for(Rint d = 0;d < 5;++ d)
for(Rint i = 0;i < N;++ i) if(!((i / po10[d]) % 10)){
memset(B, 0, sizeof B);
for(Rint j = 0;j < 10;++ j)
for(Rint k = 0;k < 10;++ k) B[j] = B[j] + A[i + k * po10[d]] * w10[j * k % 10];
for(Rint j = 0;j < 10;++ j) A[i + j * po10[d]] = B[j];
}
}
void idft(comp *A){
static comp B[10];
for(Rint d = 0;d < 5;++ d)
for(Rint i = 0;i < N;++ i) if(!((i / po10[d]) % 10)){
memset(B, 0, sizeof B);
for(Rint j = 0;j < 10;++ j)
for(Rint k = 0;k < 10;++ k) B[j] = B[j] + A[i + k * po10[d]] * w10[10 - j * k % 10];
for(Rint j = 0;j < 10;++ j) A[i + j * po10[d]] = B[j];
}
for(Rint i = 0;i < N;++ i) A[i] = A[i] * inv;
}
int main(){
read(n); w10[0].a[0] = 1; w10[1].a[3] = -1;
for(Rint i = 2;i <= 10;++ i) w10[i] = w10[i - 1] * w10[1];
for(Rint i = 1, x;i <= n;++ i){read(x); ++ A[x].a[0];}
dft(A); for(Rint i = 0;i < N;++ i) A[i] = ksm(A[i], n); idft(A);
for(Rint i = 0;i < n;++ i) cout << ((A[i].a[0] >> 5) & mod) << endl;
}