洛谷5577

2019-11-08 洛谷 FFT, 扩域 0 评论

洛谷5577

Description

个不超过位的进制数，每次随意取出一些数 （注意可以一个数都不选），然后把这些数进行进制下不进位加法。

求最终的每个数字被表示出来的方案数对取模

Input

第一行三个整数意义为题目中所述。

第二行个整数，之间用空格隔开，表示纸条上的数字。

（注意 , 纸条上的数均为进制数）

Output

共行，你需要在第行输出最终得到数字的方案数模的值。

Sample Input

1 2	3 5 1 1 2 3

Sample Output

Constraints

时时

Solution

答案是形如的所有项的系数并且此处卷积的定义是：两个下标乘起来得到的下标是这两个下标的进制下不进位加法得到的值。

于是就有了个的优质做法：每次加入一个数就枚举所有下标并更新它乘上到达的位置。可以获得 20pts 的好成绩

如果这是普通的多项式乘法那么直接对于每个数 DFT 一遍算出许多点值再对应乘起来再把这些点值 IDFT 回去就是答案。

然后你大概感性理解意会一下题目中的「进制下不进位加法」的卷积其实也可以类似的做：把每个下标看成一个维的向量且其中每一位（就是进制下的每一位），那么一个下标和的卷积就是向量中对应的每一位在模意义下的循环卷积。

循环卷积的例子：模意义下的就变成了

我们尝试对于一个如何进行维下的 DFT

回想一下一般的一维 DFT 含义：对于第项系数它对 DFT 的结果中位置的贡献的系数是

二维的话位置对位置的贡献就是直接先对每一行 DFT 再对每一列 DFT 就行

那么似乎可以直接套到多维上去：对于一个向量，其对结果中位置的贡献是

回到题目。

因为模意义下不存在和的整数解，但我们又需要表示出和，以下内容均在扩域后的意义下进行

即把每个「数」看成是

易证这个域里的「数」满足加法和乘法的封闭性

结合 20pts 做法，注意到每次乘上的只有两项，其中它对所有位置都会造成的贡献；而对于每个位置会造成的贡献

那我们对于个数可以每次 DFT 得到个形如的的点值，一共求遍对于每一位把它对应的个乘起来最后再 IDFT 回去就是答案

这样做复杂度可能是的无疑是过不了的。

但是发现的范围很大但却小得多似乎指引我们按照值域做。

仔细一想对于某一位来说，似乎并不需要求出所有的个分别是多少，只需要对于每个知道有多少个就行了（只需要个数，至于它们具体属于个中的哪几个是无关紧要的）

那问题就转化为有个数，每一项对于位置的贡献都是上面那个很长的式子。求对于每个位置的每个，在这个数中有多少个数对于他的贡献是

这东西看上去异常的难以下手，出题人的标程里写的好像是什么奇奇怪怪的折半引理消去引理的奇淫技巧~~看不懂~~

讨论区里又有一位 djq 大爷说是什么单位根反演，但（道理我都懂~~板子我都抄过~~，哪里有单位根反演啊/晕）

然后我不知道怎么搞直接把每个值对应下标加一然后做一遍 DFT ，令人惊奇的发现它！就！对！了！

（我写到这里一看上面我写的，就想问自己，这**为什么不对啊）

但这里给读者讲一下~~（也给将来的我自己）~~（看懂的就直接跳过好了）：

大概是对应下标加一后位置的系数就是这个出现的次数，对于每一个它对于每个下标的影响都是个定值。那么有个，贡献就是。那么我们直接把这个一起做 DFT 肯定是没有问题的（我们每次是加一（而没有改变后面的的项的系数），所以最后乘出来对于每个位置的贡献确实就是有多少个算多少个）~~感觉跟没讲一样，有啥问题评论区问吧~~

现在我们知道了对于一个位置和有多少个数最后对其的贡献是（记为。那么这些数点积时都会是的形式，那么最后点积出来位置的值就是

最后再做一遍 IDFT 就能得到最终的答案了，总复杂度大概是左右，如果快速幂没有预处理（见实现细节）可能还要多个

好你就这样高高兴兴的写完然后输出每个答案 ( ) 中的就发现样例都过不了？？？

你一看题解发现要输出才行？？？而且写了就直接能过？？？

如果你做过 cf1103E 的话可能也会对洛谷上唯一的题解里最终输出的是 a[i].a[0]-a[i].a[4] 而不是 a[i].a[0] 感到疑惑（然而题解并没有提到这个 magic number 是哪来的（以及他把写成为啥是对的都没有解释，像我等蒟蒻看的真是难受））

这里就方便一下我等蒟蒻略作证明（不需要的可以跳过看代码了）：

先说结论：可以证得：

（可以像那篇题解一样「带进去检验」，然而对于我等来说真不容易）

我的证法是：你可以发现最终这个东西一定是个整数，并且无论我们把中的哪一个带入答案都会是整数（这个我确实不会证，应该和 FFT 最后会是不带虚部的实数是同一个原因）

那么我们在复平面中画出所有的的幂

我们最后要得到是一个不带虚部的整数

有

而

所以得到

又因为，由上面的第二个式子得到，又由第一个式子得到于是就有

实际上不严谨的证明大概想一下上下对应的系数一定相同才能使得最后被抵消完落在轴上成为不带虚部的整数，左右两边如果相同的话就是是个有理数 ~~感觉是对的~~

知道这个后求值就很简单了：当然也可以直接用刚才三角函数的式子：

也差不多，不知道其他人怎么证的我还是画了个图：

和刚才类似的，我们可以得到

，我们得到以及

这里并不能保证上下对应位置的系数相等了（即不能保证），但由可以得到得到

答案就是（题解里好像是所以应该可以证明？感觉我还不会证，欢迎交流。~~反正打表出来是对的那我们就当它是对的~~ ）

一些实现上的细节

求在的时候是注意不要对应系数直接赋值成了
的时候快速幂可能被卡常，可以预处理出的幂来加速
也可以用表示（）不用分类讨论
DFT 时的单位根只有一项有值，可以把的乘法优化成 ~~我太懒了就没写~~

Code

#include<bits/stdc++.h>
#define FIO "P5577"
#define ll long long

char buf[1<<20]; int bufl,bufr;
#define getch ((bufl^bufr||(bufl=0,bufr=fread(buf,1,1<<20,stdin)))?buf[bufl++]:EOF)
template <class T>inline void read (T &x,const int &bas=10) {
	T f=1;
	x=0;
	char ch=getch;
	for (; !isdigit (ch)&&ch!='-'; ch=getch);
	if (ch=='-')f=-1,ch=getch;
	for (; isdigit (ch); ch=getch)x=x*bas+ch-'0';
	x*=f;
}

const int MOD=998244353,V=6,N=pow(5,7)+5;

inline int add(int a,const int &b){a+=b;return a>=MOD?a-MOD:a;}
inline int sub(int a,const int &b){a-=b;return a<		0?a+MOD:a;}
inline int mul(const int &a,const int &b){return 1ll*a*b%MOD;}
inline int& inc(int &a,const int &b){return a=add(a,b);}
inline int& dec(int &a,const int &b){return a=sub(a,b);}
inline int& pro(int &a,const int &b){return a=mul(a,b);}
inline int qpow(int a,int b){int c=1;for(;b;b>>=1,pro(a,a))if(b&1)pro(c,a);return c;}

int n,k,m;

struct Com{
  int a[V];
  inline void clear(){a[0]=a[1]=a[2]=a[3]=a[4]=a[5]=0;}
  inline Com(){clear();}
  inline int& operator [](int x){return a[x];}
  inline const int& operator [](const int &x)const{return a[x];}
  inline Com operator +(const Com &b)const{Com c;for(int i=0;i<k;i++)c[i]=add(a[i],b[i]);return c;}
  inline Com operator -(const Com &b)const{Com c;for(int i=0;i<k;i++)c[i]=sub(a[i],b[i]);return c;}
  inline Com operator *(const int &b)const{Com c;for(int i=0;i<k;i++)c[i]=mul(a[i],b);return c;}
  inline Com& operator +=(const Com &b){return *this=*this+b;}
  inline Com& operator *=(const int &b){return *this=*this*b;}
  inline Com operator *(const Com &b)const{
    static ll r[V<<1];
    const ll Moc=8ll*MOD*MOD;
    for(int i=0;i<k;i++)if(a[i])for(int j=0;j<k;j++)if(b[j]){
      r[i+j]+=1ll*a[i]*b[j];
      if(r[i+j]>=Moc)r[i+j]-=Moc;
    }
    Com c;
    for(int i=0;i<k;i++)c[i]=(r[i]+r[i+k])%MOD,r[i]=r[i+k]=0;
    return c;
  }
  inline Com& operator *=(const Com &b){return *this=*this*b;}
  inline void out(){for(int i=0;i<k;i++)printf("%d%c",a[i],i^(k-1)?' ':'\n');}
  inline int val(){
    return sub(add(a[0],a[1]),add(a[2],a[3]));
  }
}a[N],w[2][V],b[N];

inline void pre(){
  w[0][0][0]=w[1][0][0]=1;
  w[0][1][1]=w[1][1][k-1]=1;
  for(int i=2;i<k;i++)
    w[0][i]=w[0][i-1]*w[0][1],w[1][i]=w[1][i-1]*w[1][1];
}

//inline Com qpow(Com a,int b){Com c;c[0]=1;for(;b;b>>=1,a*=a)if(b&1)c*=a;return c;}
const int B=1000;
Com pw0[B+5],pw1[B+5];
inline Com qpow(Com a,int b){return pw1[b/B]*pw0[b%B];}

inline void getpow(){
  for(int i=0;i<n;i++)b[i][0]=1;
  for(int i=0;i<k;i++){
    Com omega;
    //!!! i=1 gg 
    omega[0]=1;omega[i]++;
    pw0[0][0]=pw1[0][0]=1;
    for(int i=1;i<=B;i++)pw0[i]=pw0[i-1]*omega;
    pw1[1]=pw0[B];
    for(int i=1;i<=B;i++)pw1[i]=pw1[i-1]*pw1[1];
    for(int id=0;id<n;id++)b[id]*=qpow(omega,a[id][i]);
  }
  for(int i=0;i<n;i++)a[i]=b[i];
}

int pw[V];
inline void fft(int opt){
  pw[0]=1;for(int i=1;i<m;i++)pw[i]=pw[i-1]*k;
  for(int i=0;i<m;i++){
    for(int j=0;j<n;j++)
    if(!((j/pw[i])%k)){
      static Com f[V];
      for(int x=0;x<k;x++)
        for(int y=0;y<k;y++)
          f[y]+=a[j+pw[i]*x]*w[opt][x*y%k];
      for(int x=0;x<k;x++)
        a[j+pw[i]*x]=f[x],f[x].clear();
    }
  }
  if(opt)for(int i=0,inv=qpow(n,MOD-2);i<n;i++)a[i]*=inv;
}

int main(){
	freopen(FIO".in","r",stdin);
	freopen(FIO".out","w",stdout);
 read(n);read(k);read(m);
  for(int i=1,x;i<=n;i++)read(x,k),a[x][0]++;
  n=qpow(k,m);
  pre();

  fft(0);

  getpow();

  fft(1);

  for(int i=0;i<n;i++)printf("%d\n",a[i].val());

	return 0;
}

本文链接： https://20021123.xyz/2019/11/08/洛谷5577/

版权声明： 本博客所有文章除特别声明外，均采用 CC BY 4.0 CN协议许可协议。转载请注明出处！

不知不觉Acmer & Student

Never Settle