下标模意义下的多项式乘法及其应用

已知两个数组\(a,b\)，求一个数组\(c\)，满足\(c_i=\sum_{j+k\equiv i(MOD\ K)}a_jb_k(MOD\ M)\)。这里我们把这个东西称为下标模意义下的多项式乘法。那么这个东西怎么做呢？

先说结论：如果MOD M意义下存在K次单位根，那么把平时用NTT做多项式乘法时的原根G统一换成这个K次单位根就可以了(不能跑\(O(nlogn)\)的NTT，因为把数组长度扩充到2的倍数会导致出错)。

复习一下单位根的定义

这是为什么呢？

举个例子，当K=17，M=998244353的时候。此时MOD M意义下是存在K次单位根的，它的值是503044，令其为\(w\)。假设现在有一个多项式\(a\)，且项数不止17个，令它对\(x^{17}\)取模之后的多项式为\(b\)，\(b\)只有\(17\)项，且\(b_i=a_i+a_{i+17}+a_{i+34}\cdots\)。以503044作为单位根的值对\(a,b\)分别做一次DFT(系数转点值)，得到\(c,d\)。由于\(w^{17}=1\)，所以\(c,d\)都有长度为17的循环节，并且容易发现\(c=d\)。这时候截取\(c\)的前17项做IDFT，就可以得到\(b\)，也就是下标取模之后的多项式。这说明只要对一个多项式先DFT，取前17项再IDFT，就可以完成一次下标的取模。

那么两个多项式相乘也是一样的，最后把点值相乘的结果IDFT回去的时候只取前17项就行了。

但是注意到\(O(nlogn)\)的NTT不能这么搞，因为NTT会把数组长度先扩充到2的倍数，但是这个问题里数组长度不是17就不对了。

放一张毛啸2016年候选队论文的截图，仅供参考：

例题

n 个点 m 条带权无向边，权值是 0 到 16。问多少种方案选择一些边（不选重边），使得图联通，且边权和模 17 正好为 x。对每个 0≤x≤16 都求答案。\(2 \le n \le 17,1\le m\le 10^5\)。

解法

前置知识：

1. 二维FFT(高维FFT)。有两个矩阵\(a,b\)，现在要求一个新的矩阵\(c\)满足\(c_{i,j}=\sum_{u+v=i,x+y=j}a_{u,x}b_{v,y}\)。做法


2. 子集exp和子集ln。一个子集幂级数指的是一个数组\(f_s\)，下标\(s\)是一个"子集"，可以看成是\([0,2^n-1]\)中的一个整数，用二进制的方式表示n个元素中的哪些在子集中。子集exp和ln则类似EGF的exp和ln。子集exp的组合意义：如果\(g=e^f\)，且\(f_s\)表示将子集\(s\)中所有的元素进行操作A的方案是，则\(g_s\)表示把\(s\)进一步分成若干子集，每个子集进行A操作的总方案数。子集ln则是子集exp的逆运算，用来把g变成f的。关于子集exp和ln的实现方法和证明具体见这里。

首先令\(f_{i,s}\)表示子集\(s\)内部的边全部乱选，且不能选重边，和MOD 17为i的方案数。\(s\)内部的边指的是两段都在\(s\)中的边。同样定义\(g_{i,s}\)表示只选\(s\)中的边，图连通，没有重边且和MOD 17为i的方案数。观察得到转移式：令\(l\)表示\(s\)中的最低的1，则\(g_{i,s}=f_{i,s}-\sum_{t\subsetneq s,l\in t,j+k=i(MOD\ 17)}g_{j,t}\cdot f_{k,s \ xor \ t}\)，也就是枚举l所在的连通块。直接写这玩意儿能做到\(O(3^n\cdot17^2)\)。

如果所有边的权值都为0，那就可以把\(f,g\)的第一维去掉，只留下子集那一维。根据子集exp的组合意义，很容易发现\(f=e^g\)，\(g=ln(f)\)。根据上面链接里的子集ln求法，直接\(O(2^n\cdot n^2)\)求\(f\)的ln即可。\(f\)本身可以通过简单的dp求出，复杂度是\(O(2^n\cdot n\cdot 17^2)\)。在所有边权都为0的时候则可以直接\(O(2^n\cdot n)\)求出。

边权值不全为0的情况则需要加上第一维。注意到第一维是一个上面说的"下标模意义下的多项式乘法"。可以根据二维FFT的方法依葫芦画瓢，先把第一维以503044位单位根进行DFT，然后把第二维子集ln，最后再把第一维IDFT。这样套娃的正确性我不会证明，感性理解。

\(n\le 17\)，所以把n用17表示了。总复杂度\(O(17^3\cdot 2^n)\)。

#include
#define rep(i,n) for(int i=0;i#define repn(i,n) for(int i=1;i<=n;++i)
#define LL long long
#define pii pair
#define fi first
#define se second
#define mpr make_pair
#define pb push_back
void fileio()
{
#ifdef LGS
freopen("in.txt","r",stdin);
freopen("out.txt","w",stdout);
#endif
}
void termin()
{
#ifdef LGS
std::cout<<"\n\nPROGRAM TERMINATED";
#endif
exit(0);
}
using namespace std;
const LL MOD=998244353;
const LL w=503044;
LL qpow(LL x,LL a)
{
LL res=x,ret=1;
while(a>0)
{
if(a&1) (ret*=res)%=MOD;
a>>=1;
(res*=res)%=MOD;
}
return ret;
}
LL n,m,es[20][20][20],f[20][140000],g[20][140000],dp[20][20],modv[20][20],inv17,inv[20];
//f：乱选(但不选重边)；g：连通
vector DFT(vector v)
{
LL cur=1;
vector ret;
rep(i,v.size())
{
LL curr=1,res=0;
rep(j,v.size())
{
(res+=v[j]*curr)%=MOD;
(curr*=cur)%=MOD;
}
(cur*=w)%=MOD;
ret.pb(res);
}
return ret;
}
vector IDFT(vector v)
{
LL cur=1,ww=qpow(w,MOD-2);
vector ret;
rep(i,v.size())
{
LL curr=1,res=0;
rep(j,v.size())
{
(res+=v[j]*curr)%=MOD;
(curr*=cur)%=MOD;
}
(cur*=ww)%=MOD;
ret.pb(res);
}
rep(i,ret.size()) (ret[i]*=inv17)%=MOD;
return ret;
}
void add(LL &x,LL y){x+=y;if(x>=MOD) x-=MOD;}
LL cur[140000],h[20][140000],fr[20],to[20];
void subsetLn()
{
rep(i,18) rep(j,1< rep(i,1< //第三维FMT
rep(i,n+1)
{
for(int j=n-1;j>=0;--j) rep(k,1< add(h[i][k|(1< }
//第二维ln
rep(i,1< {
rep(j,n+1) fr[j]=h[j][i],to[j]=0;
repn(j,n)
{
to[j]=fr[j];
LL sum=0;
rep(k,j) (sum+=to[k]*fr[j-k]%MOD*k)%=MOD;
(sum*=inv[j])%=MOD;
add(to[j],MOD-sum);
}
rep(j,n+1) h[j][i]=to[j];
}
//第三维IFMT
rep(i,n+1)
{
rep(j,n) rep(k,1< add(h[i][k|(1< }
rep(i,1<}
int main()
{
fileio();
//freopen("roads.in","r",stdin);
//freopen("roads.out","w",stdout);
inv17=qpow(17,MOD-2);
rep(i,19) inv[i]=qpow(i,MOD-2);
rep(i,18) rep(j,18) modv[i][j]=(i+j)%17;
cin>>n>>m;
LL x,y,z;
rep(i,m)
{
scanf("%lld%lld%lld",&x,&y,&z);--x;--y;
if(x>y) swap(x,y);
++es[x][y][z];
}
rep(i,n) for(int j=i+1;j rep(msk,1<0)
{
int lowbit=msk&(-msk),lft=msk^lowbit,id=__builtin_ctz(msk);
if(lft==0) f[0][msk]=1;
else
{
vector nds;rep(j,n) if(lft&(1< rep(i,18) rep(j,18) dp[i][j]=0;
rep(i,17) dp[0][i]=f[i][lft];
rep(i,nds.size())
{
int ii=nds[i];
rep(j,17) rep(k,17) (dp[i+1][modv[j][k]]+=dp[i][j]*es[id][ii][k])%=MOD;
}
rep(i,17) f[i][msk]=dp[nds.size()][i];
}
}
//第一维DFT
rep(i,1< {
vector vv;
rep(j,17) vv.pb(f[j][i]);
vv=DFT(vv);
rep(j,17) g[j][i]=vv[j];
}
//子集ln
rep(i,17)
{
rep(j,1< subsetLn();
rep(j,1< }
//第二维IDFT
rep(i,1< {
vector vv;
rep(j,17) vv.pb(g[j][i]);
vv=IDFT(vv);
rep(j,17) g[j][i]=vv[j];
}
rep(i,17) printf("%lld ",g[i][(1< puts("");
termin();
}