【大意】
令 (A_{i, j}) 表示长度为 (L) 的字符串，只用前 (k(2leq kleq 26)) 个小写字母，且 (a) 字母出现次数模 (n) 为 (i) ，(b) 字母出现次数模 (n) 为 (j) 的方案数
输出矩阵 (A_{0cdots (n-1), 0cdots (n-1)})
【分析】
单位根反演题，现场过了一堆人我十分吃惊
一直认为是冷门算法
(displaystyle A_{i, j}=sum_{p+q+r=L}dbinom{L}{p, q, r}[pmod n=i]cdot [qmod n=j]cdot (k-2)^r)
考虑单位根反演的公式：
(displaystyle [kmod n=0]={1over n}sum_{i=0}^{n-1}omega_n^{ki})
代入得到：
(displaystyle A_{i, j}=sum_{p+q+r=L}dbinom{L}{p, q, r}{1over n}sum_{a=0}^{n-1}omega_n^{(p-i)a}{1over n}sum_{b=0}^{n-1}omega_n^{(q-j)b}cdot (k-2)^r)
调换求和顺序得到
(displaystyle A_{i, j}={1over n^2}sum_{a=0}^{n-1}omega_n^{ai}sum_{b=0}^{n-1}omega_n^{bj}sum_{p+q+r=L}dbinom{L}{p, q, r}(omega_n^a)^pcdot (omega_n^b)^qcdot (k-2)^r)
后面是个多项式定理的形式，合并得到：
(displaystyle A_{i, j}={1over n^2}sum_{a=0}^{n-1}omega_n^{ai}sum_{b=0}^{n-1}omega_n^{bj}(omega_n^a+omega_n^b+k-2)^L)
考虑到 (nmid (p-1))
故取 (p) 的原根 (g)，即可令 (displaystyle omega_nequiv g^{p-1over n}pmod p)
如果只求解单个，则复杂度没问题，但求和多个时，复杂度是不够的
考虑预处理 (displaystyle W_{a, b}=(omega_n^a+omega_n^b+k-2)^L)
计算所有 (W_{a, b}) 的复杂度为 (O(n^2log p))
之后，我们预处理 (displaystyle S_a=sum_{b=0}^{n-1}omega_n^{bj}(omega_n^a+omega_n^b+k-2)^L=sum_{b=0}^{n-1}omega_n^{bj}cdot W_{a, b})
计算 (S_{0cdots (n-1)}) 的复杂度为 (O(n^2))
最后，我们求解所有 (A_{i, j}) 时，由 (displaystyle A_{i, j}={1over n^2}sum_{a=0}^{n-1}omega_n^{ai}sum_{b=0}^{n-1}omega_n^{bj}(omega_n^a+omega_n^b+k-2)^L={1over n^2}sum_{a=0}^{n-1}omega_n^{ai}cdot S_a)
求解对所有的复杂度优化为 (O(n^3))
故总复杂度为 (O(n^3))
【代码】
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef pair<ll, ll> pii;
typedef double db;
#define fi first
#define se second
const int MOD=1e9+9, g=13;
ll k, l, n, omega[512], powit[512][512], t[512][512];
inline ll fpow(ll a,ll x) { ll ans=1; for(;x;x>>=1,a=a*a%MOD) if(x&1) ans=ans*a%MOD; return ans; }
inline int w(int val) { return (val%=n)>=0?omega[val]:omega[n+val]; }
inline int add(int a, int b) { return (a+=b)>=MOD?a-MOD:a; }
inline void pre(){
cin>>k>>l>>n;
omega[0]=1;
omega[1]=fpow(g, (MOD-1)/n);
for(int i=2;i<n;++i) omega[i]=omega[i-1]*omega[1]%MOD;
for(int p=0;p<n;++p)
for(int q=0;q<n;++q)
powit[p][q]=fpow(k-2+w(-p)+w(-q), l);
for(int a=0;a<n;++a)
for(int j=0;j<n;++j){
int res=0;
for(int b=0;b<n;++b)
res=add(res, w(b*j)*powit[a][b]%MOD);
t[a][j]=res;
}
}
int main(){
ios::sync_with_stdio(0);
cin.tie(0); cout.tie(0);
pre();
ll C=fpow(n, MOD+MOD-4);
for(int i=0;i<n;++i){
for(int j=0;j<n;++j){
int res=0;
for(int a=0;a<n;++a)
res=add(res, w(a*i)*t[a][j]%MOD);
cout<<res*C%MOD<<" ";
}
cout<<"
";
}
cout.flush();
return 0;
}