ZJT's Blog

大変に気分がいい~

FFT相关的一些优化技巧

2018-4-17 17:132018-4-17 17:13
FFT卡常数

大概就是myy论文里的那些东西吧,闲着没事整理一下。

记FFT变换的矩阵为$F$(即单位根的范德蒙德矩阵),将FFT视作$F$乘上待变换的行向量:$FFT(a)=Fa$。我们记$\overline F$为$F$的复共轭矩阵(矩阵的每个元素都共轭),$a'$为$a$除了第$0$项外都翻转后的结果(本文的所有矩阵、向量的下标均从$0$开始),即第$a_i$与$a_{n-i}$交换($n$为$a$的长度)。显然$a'$可以写成$a$乘上一个翻转矩阵$T$,即$a'=aT$。

  • 性质 1:$\overline Fa=(Fa)'$

这个根据单位根的性质$\omega^{-k}=\omega^{n-k}$可以轻易得到。

  • 性质 2:$F^{-1}=\overline F/n$

这其实就是IFFT最后一步干的事情。。

技巧 1:快速求两个实多项式的点值(以及插值)

现在有两个实向量$a,b$,我们想计算$Fa,Fb$。按照一般的做法,直接跑两遍FFT就行了。其实这个过程可以做到只用一次FFT。

我们构造复向量$a+bi$,考虑把这个向量FFT: $$F(a+bi)=Fa+F(bi)=Fa+iFb$$

我们把这整个东西取复共轭,看看会发生什么: $$ \begin{aligned} \overline{F(a+bi)}&=\overline {Fa}+\overline {iFb}\\ &=\overline Fa-i\overline Fb\\ &=(Fa)'-i(Fb)'\\ \end{aligned} $$

两边同时乘上$T$,再联立之前FFT的结果: $$ Fa-iFb=(\overline{F(a+bi)})'\\ Fa+iFb=F(a+bi) $$

知道这两个东西之后,可以直接解出$Fa,Fb$。也就是说,我们直接算出$a+bi$的FFT,把它复共轭之后翻转系数(从第$1$项开始翻转),之后就可以解出$Fa$和$Fb$了。

同理,我们也可以同时做两个实多项式点值的插值。把这整个过程反过来做一遍即可。

技巧 2:4次FFT的任意模数多项式乘法

这其实是上面那个东西的一个应用。现在有两个整系数多项式$A(x),B(x)$,我们要算它们循环卷积的结果(如果要做普通卷积,直接把次数扩大一倍即可),所有系数都对一个$10^9$级别的整数$P$取模(这个整数可以什么性质都没有)。我们考虑选取整数$M\approx\sqrt P$,我们把$A,B$的系数对$M$作带余除法,得到$A(x)=A_0(x)+MA_1(x),B(x)=B_0(x)+MB_1(x)$。这里$A_0,A_1,B_0,B_1$的系数都是$M$级别的,他们之间乘起来之后系数是$O(M^2n)=O(Pn)$级别的,$n$不太大的话(大概$5*10^5$以内)可以认为是在double的精度范围内,可以直接用double做乘法后取模。

我们用上面提到的那个技巧,用$2$次FFT算出$A_0,A_1,B_0,B_1$的FFT点值$a_0,a_1,b_0,b_1$,之后根据$$\begin{aligned}A(x)B(x)&=A_0(x)B_0(x)\\&+M(A_0(x)B_1(x)+A_1(x)B_0(x))\\&+M^2A_1(x)B_1(x)\end{aligned}$$我们可以再次使用上面那个技巧,用$2$次FFT算出$a_0\otimes b_0,a_0\otimes b_1,a_1\otimes b_0,a_1\otimes b_1$的IFFT插值,就可以直接进行计算了。

总共$4$次FFT。

技巧 3:2次FFT的实多项式乘法

这个用技巧$1$也可以做到,但较技巧$1$更容易实现。
有两个实多项式$A(x),B(x)$,我们要计算它们的循环卷积。

我们可以用两次FFT计算$(A(x)+iB(x))^2=A^2(x)-B^2(x)+2iA(x)B(x)$。

直接取虚部并除以2即可。

查看详细内容

扩展拉格朗日反演

2018-3-4 9:432018-4-6 16:8
拉格朗日反演FFT生成函数

突然发现拉格朗日反演其实还有扩展的形式。。。

定义$F^{-1}$为$F$的复合逆:

$$ F^{-1}(F(x))=x $$

$$[x^n]G(F^{-1}(x))=\frac1n[x^{-1}] ( \frac{dG(x)}{dx}\frac1{F^n(x)})$$

(这个好像叫Lagrange–Bürmann formula

注意一般$G$的表达式都同时含有$x$和$F^{-1}(x)$,我们把$F(x)$代入$x$后就能得到$G(x)$的表达式,进而求出$\frac {dG(x)}{dx}$。

感觉还是蛮有用的(大雾)

查看详细内容

[自选题 #120] Mike学OI

2018-1-4 20:442018-1-4 20:44
集训队作业自选题生成函数FFT

传送门

多项式多点插值+多点求值的模板题。

不过之前一直以为插值是$O(n\log^3n)$的,直到这题solution给了一个两个log的做法。插值的瓶颈在于求出$g_i=\prod_{i\ne j}(x_i-x_j)$,我们可以把他表示成$f(x)=\frac{\prod{x-x_j}}{x-x_i}$在$x=x_i$处的极限。根据洛必达法则,我们对分式上下求导,得到$\lim_{x\to x_i}f(x)=(\prod(x-x_j))'|_{x=x_i}$。也就是说,我们先分治算出$\prod(x-x_j)$,然后再把这个多项式求导,并在所有$x_i$处多点求值一下,就能得到所有的$g_i$。得到$g_i$之后就好办了,代到拉格朗日插值的式子里面,随便分治一下就出来了。

好久没写这种东西了,代码写的很乱。。。感觉要调整一下写法了。

#include <bits/stdc++.h>
#define LL long long
#define MAXN 270000
#define MAXW 262144
#define B 600
using namespace std;

const LL P=998244353;
LL w[MAXN];

LL getPow(LL x,LL y){
    LL res=1;
    while(y){
        if(y&1) res=res*x%P;
        x=x*x%P;
        y>>=1;
    }
    return res;
}

struct Poly{
    LL a[MAXN];
    int len;
    inline LL &operator[](int i){ return a[i]; }
    void reset(int _len){ len=_len; memset(a,0,(len+1)*sizeof(LL)); }
};

void FFT(Poly &a,int len,int flag){
    static int rev[MAXN],revlen;
    if(revlen!=len){
        revlen=len;
        for(int i=1;i<len;i++) rev[i]=rev[i>>1]>>1|((i&1)?(len>>1):0);
    }
    for(int i=0;i<len;i++) if(i<rev[i]) swap(a.a[i],a.a[rev[i]]);
    for(int l=2;l<=len;l<<=1){
        int l2=l>>1;
        for(int i=0;i<len;i+=l){
            for(int j=0;j<l2;j++){
                LL t1=a.a[i+j],t2=a.a[i+j+l2]*w[MAXW/l*j];
                a.a[i+j]=(t1+t2)%P;
                a.a[i+j+l2]=(t1-t2)%P;
            }
        }
    }
    if(flag==-1){
        LL invn=getPow(len,P-2);
        for(int i=0;i<len;i++) a.a[i]=a.a[i]*invn%P;
        for(int i=1;i<len;i++) if(i<len-i) swap(a.a[i],a.a[len-i]);
    }
}

void getMul(Poly &a,Poly &b,Poly &c){
    int sizew=1;
    while(sizew<=a.len+b.len) sizew<<=1;
    sizew<<=1;
    static Poly t1,t2;
    t1.reset(sizew); t2.reset(sizew);
    for(int i=0;i<=a.len;i++) t1[i]=a[i];
    for(int i=0;i<=b.len;i++) t2[i]=b[i];
    FFT(t1,sizew,1); FFT(t2,sizew,1);
    for(int i=0;i<sizew;i++) t1[i]=t1[i]*t2[i]%P;
    FFT(t1,sizew,-1);
    c.len=a.len+b.len;
    for(int i=0;i<=c.len;i++) c[i]=t1[i];
}

void getInv(Poly &b,Poly &a,int len=-1){
    if(len==-1) for(len=1;len<=a.len;len<<=1);
    if(len==1){
        b.len=len;
        b[0]=getPow(a[0],P-2);
        return;
    }
    static Poly t1,t2;
    getInv(t1,a,len>>1);
    for(int i=len>>1;i<len<<1;i++) t1[i]=0;
    t2.reset((len<<1)-1);
    for(int i=0;i<len;i++) t2[i]=a[i];
    FFT(t1,len<<1,1); FFT(t2,len<<1,1);
    for(int i=0;i<len<<1;i++) t1[i]=(2*t1[i]-t2[i]*t1[i]%P*t1[i])%P;
    FFT(t1,len<<1,-1);
    b.len=len-1;
    for(int i=0;i<len;i++) b[i]=t1[i];
}

void reverse(Poly &a){
    for(int i=0;i<=a.len;i++) if(i<a.len-i) swap(a[i],a[a.len-i]);
}

void getMod(Poly &a,Poly &b,Poly &c){
    if(a.len<b.len){
        c.len=a.len;
        for(int i=0;i<=a.len;i++) c[i]=a[i];
        return;
    }
    static Poly t1,t2,t3,t4;
    t1.len=a.len; for(int i=0;i<=a.len;i++) t1[i]=a[i];
    t2.len=b.len; for(int i=0;i<=b.len;i++) t2[i]=b[i];
    reverse(t1); reverse(t2);
    while(t2.len<a.len-b.len) t2[++t2.len]=0;
    t2.len=a.len-b.len;
    getInv(t4,t2);
    while(t1.len<a.len-b.len) t1[++t1.len]=0;
    while(t4.len<a.len-b.len) t4[++t4.len]=0;
    t1.len=t4.len=a.len-b.len;
    getMul(t1,t4,t3); t3.len=a.len-b.len;
    t2.len=b.len; getMul(t3,t2,t2);
    int lenc=b.len-1;
    c.len=a.len;
    for(int i=0;i<=a.len;i++) c[i]=t1[i]-t2[i];
    reverse(c);
    c.len=lenc;
}

LL *xp[MAXN],*fp[MAXN];

void getXP(int id,LL *x,int len){
    if(len==1){
        xp[id]=new LL[2];
        xp[id][0]=-x[0];
        xp[id][1]=1;
        fp[id]=new LL[2];
        return;
    }
    int mid=len>>1;
    getXP(id<<1,x,mid);
    getXP(id<<1|1,x+mid,len-mid);
    static Poly t1,t2;
    t1.len=mid; t2.len=len-mid;
    for(int i=0;i<=mid;i++) t1[i]=xp[id<<1][i];
    for(int i=0;i<=len-mid;i++) t2[i]=xp[id<<1|1][i];
    getMul(t1,t2,t1);
    xp[id]=new LL[len+1];
    fp[id]=new LL[len+1];
    for(int i=0;i<=len;i++) xp[id][i]=t1[i];
}

void getV(int id,LL *x,LL *y,int len){
    if(len<=B){
        for(int i=0;i<len;i++){
            LL res=0,temp=1;
            for(int j=0;j<=len;j++)
                res=(res+temp*fp[id][j])%P,temp=temp*x[i]%P;
            y[i]=res;
        }
        return;
    }
    int mid=len>>1;
    static Poly tp,t1,t2;
    t1.len=len;
    for(int i=0;i<=len;i++) t1[i]=fp[id][i];
    tp.len=mid;
    for(int i=0;i<=tp.len;i++) tp[i]=xp[id<<1][i];
    getMod(t1,tp,t2);
    for(int i=0;i<=mid;i++) fp[id<<1][i]=t2[i];
    tp.len=len-mid;
    for(int i=0;i<=tp.len;i++) tp[i]=xp[id<<1|1][i];
    getMod(t1,tp,t2);
    for(int i=0;i<=mid;i++) fp[id<<1|1][i]=t2[i];
    getV(id<<1,x,y,mid);
    getV(id<<1|1,x+mid,y+mid,len-mid);
}

void getV(Poly &a,LL *x,LL *y,int len){
    getXP(1,x,len);
    static Poly t1;
    t1.len=len;
    for(int i=0;i<=len;i++) t1[i]=xp[1][i];
    getMod(a,t1,t1);
    for(int i=0;i<=t1.len;i++) fp[1][i]=t1[i];
    for(int i=t1.len+1;i<=len;i++) fp[1][i]=0;
    getV(1,x,y,len);
}

void getPoly(int id,LL *x,LL *y,int len){
    if(len<=B){
        for(int i=0;i<=len;i++) fp[id][i]=0;
        for(int i=0;i<len;i++){
            LL last=0;
            for(int j=len;j;j--){
                LL temp=(xp[id][j]+last)%P;
                last=temp*x[i]%P;
                fp[id][j-1]=(fp[id][j-1]+temp*y[i])%P;
            }
        }
        return;
    }
    int mid=len>>1;
    getPoly(id<<1,x,y,mid);
    getPoly(id<<1|1,x+mid,y+mid,len-mid);
    static Poly t1,t2;
    t1.len=mid; for(int i=0;i<=t1.len;i++) t1[i]=fp[id<<1][i];
    t2.len=len-mid; for(int i=0;i<=t2.len;i++) t2[i]=xp[id<<1|1][i];
    getMul(t1,t2,t1);
    for(int i=0;i<=len;i++) fp[id][i]=t1[i];
    t1.len=len-mid; for(int i=0;i<=t1.len;i++) t1[i]=fp[id<<1|1][i];
    t2.len=mid; for(int i=0;i<=t2.len;i++) t2[i]=xp[id<<1][i];
    getMul(t1,t2,t1);
    for(int i=0;i<=len;i++) fp[id][i]=(fp[id][i]+t1[i])%P;
}

void getPoly(LL *x,LL *y,Poly &f,int len){
    getXP(1,x,len);
    static Poly py;
    static LL yv[MAXN];
    py.len=len;
    for(int i=0;i<=len;i++) py[i]=xp[1][i];
    for(int i=0;i<len;i++) fp[1][i]=py[i+1]*(i+1)%P;
    fp[1][len]=0;
    getV(1,x,yv,len);
    for(int i=0;i<len;i++) yv[i]=y[i]*getPow(yv[i],P-2)%P;
    getPoly(1,x,yv,len);
    for(int i=0;i<=len;i++) f[i]=fp[1][i];
    f.len=len;
}

void init(){
    w[1]=getPow(3,(P-1)/MAXW);
    w[0]=1;
    for(int i=2;i<MAXW;i++) w[i]=w[i-1]*w[1]%P;
}

void gao(){
    int n,m;
    static LL xv[MAXN],yv[MAXN];
    static Poly f;
    scanf("%d",&n);
    for(int i=0;i<n;i++) scanf("%lld%lld",xv+i,yv+i);
    getPoly(xv,yv,f,n);
    scanf("%d",&m);
    for(int i=0;i<m;i++) scanf("%lld",xv+i);
    getV(f,xv,yv,m);
    for(int i=0;i<m;i++) printf("%lld ",(yv[i]%P+P)%P);
}

int main(){
#ifdef DEBUG
    freopen("120.in","r",stdin);
#endif
    init();
    gao();
}

查看详细内容

[自选题 #124] Path

2018-1-4 19:272018-1-4 19:27
集训队作业自选题FFTHook length定理

传送门

这个题要用到一个Hook length定理,这个定理是专门用来求这一类填数问题的方案数的。有了这个定理之后就好做了。我们考虑第$i$行,这一行每个元素的hook length就是一个递减数列中间删掉一些数,我们可以通过一个卷积来算出少了哪些数。这样就能求出所有hook length的乘积,然后直接套用Hook length定理就行了。

#include <bits/stdc++.h>
#define MAXN 4200000
#define LL long long
using namespace std;

const int MAXW=2097152;
const LL P=1004535809,G=3;

LL getPow(LL x,LL y){
    LL res=1;
    while(y){
        if(y&1) res=res*x%P;
        x=x*x%P;
        y>>=1;
    }
    return res;
}

int n,m,len,sizew;
int p[MAXN];
LL fac[MAXN],invfac[MAXN],w[MAXN];
LL t1[MAXN],t2[MAXN];
LL a[MAXN],b[MAXN];

void FFT(LL *a,int len,int flag){
    static int rev[MAXN],revlen;
    if(revlen!=len){
        revlen=len;
        for(int i=1;i<len;i++) rev[i]=rev[i>>1]>>1|((i&1)?(len>>1):0);
    }
    for(int i=0;i<len;i++) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int l=2;l<=len;l<<=1){
        int l2=l>>1;
        for(int i=0;i<len;i+=l){
            for(int j=0;j<l2;j++){
                LL t1=a[i+j],t2=a[i+j+l2]*w[MAXW/l*j];
                a[i+j]=(t1+t2)%P;
                a[i+j+l2]=(t1-t2)%P;
            }
        }
    }
    if(flag==-1){
        LL invn=getPow(len,P-2);
        for(int i=0;i<len;i++) a[i]=a[i]*invn%P;
        for(int i=1;i<len;i++) if(i<len-i) swap(a[i],a[len-i]);
    }
}

int main(){
#ifdef DEBUG
    freopen("124.in","r",stdin);
#endif
    scanf("%d",&n);
    for(int i=1;i<=n;i++){
        scanf("%d",p+i);
        if(p[i]>m) m=p[i];
    }
    len=(n+m)*2;
    for(sizew=1;sizew<=len;sizew<<=1);
    fac[0]=1;
    for(int i=1;i<=len;i++) fac[i]=fac[i-1]*i%P;
    invfac[len]=getPow(fac[len],P-2);
    for(int i=len-1;i>=0;i--) invfac[i]=invfac[i+1]*(i+1)%P;
    w[1]=getPow(G,(P-1)/MAXW);
    w[0]=1;
    for(int i=2;i<MAXW;i++) w[i]=w[i-1]*w[1]%P;
    for(int i=1;i<=n;i++) t1[(p[i]-i+sizew)%sizew]++;
    for(int i=1;i<n;i++)
        t2[(i+1-(p[i+1]+1)+sizew)%sizew]++;
    FFT(t1,sizew,1); FFT(t2,sizew,1);
    for(int i=0;i<sizew;i++) t1[i]=t1[i]*t2[i]%P;
    FFT(t1,sizew,-1);
    LL ans=1;
    int last=n;
    for(int i=n;i>=1;i--){
        ans=ans*fac[p[i]+n-i]%P;//*invfac[last-i]%P;
        if(p[i]!=p[i-1]){
            int len=last-i+1;
            last=i-1;
        }
    }
    for(int i=0;i<=n+m;i++)
        ans=ans*getPow(invfac[i+1]*fac[i]%P,(t1[i]+P)%P)%P;
    ans=getPow(ans,P-2);
    for(int i=1;i<=n;i++) ans=ans*fac[p[i]]%P;
    printf("%lld\n",ans);
    return 0;
}

查看详细内容

[自选题 #127] Ball

2018-1-2 12:212018-1-2 12:22
集训队作业自选题生成函数FFT线性递推

传送门

考虑DP:$f_{i,j}$表示$i$个球分$j$组的方案数。显然$f_{i,j}=f_{i-1,j}+f_{i-1,j-1}+f_{i-2,j-1}$。

如果我们把$f_i$看成一个$x$多项式,则$f_i(x)=(x+1)f_{i-1}(x)+xf_{i-2}(x)$。这是个常系数线性递推,特征方程为$\lambda^2-(x+1)\lambda-x=0$,解出特征根为$\lambda_1=\frac{x+1+\sqrt{x^2+6x+1}}{2},\lambda_2=\frac{x+1-\sqrt{x^2+6x+1}}{2}$。

我们钦定$f_{-1}(x)=0$,此时转移方程在所有正整数处都成立。我们设$f_i(x)=c_1(x)\lambda_1^{i+1}(x)+c_2(x)\lambda_2^{i+1}(x)$,把$f_{-1},f_0$代进去,得到$c_1=\frac 1{\lambda_1-\lambda_2},c_2=\frac 1{\lambda_2-\lambda_1}$。所以$f_n$就能表示成:$$f_n(x)=\frac {\lambda_1^{n+1}(x)-\lambda_2^{n+1}(x)}{\lambda_1(x)-\lambda_2(x)}$$

直接开根把$\lambda$算出来之后$\ln+\exp$求幂再求逆就行了。

由于$m\leq n$(当$m>n$时后面那些显然就是$0$),而$\lambda_2$的常数项是$0$,所以$\lambda_2^{n+1}(x)\equiv 0\pmod {x^{m+1}}$,这意味着分子后面那项根本没有贡献,可以直接忽略,从而减小常数。

#include <bits/stdc++.h>
#define MAXN 1100000
#define MAXW 1048576
#define LL long long 
using namespace std;

const LL P=998244353;

int n,m,sizew;
LL w[MAXW],inv_v[MAXN];

LL getPow(LL x,LL y){
    LL res=1;
    while(y){
        if(y&1) res=res*x%P;
        x=x*x%P;
        y>>=1;
    }
    return res;
}

void init(){
    w[0]=1;
    w[1]=getPow(3,(P-1)/MAXW);
    for(int i=2;i<MAXW;i++) w[i]=w[i-1]*w[1]%P;
    inv_v[1]=1;
    for(int i=2;i<MAXN;i++)
        inv_v[i]=-inv_v[P%i]*(P/i)%P;
}

void FFT(LL *a,int n,int flag){
    static int rev[MAXN],revn;
    if(revn!=n){
        revn=n;
        for(int i=1;i<n;i++) rev[i]=rev[i>>1]>>1|((i&1)?(n>>1):0);
    }
    for(int i=0;i<n;i++) if(i<rev[i]) swap(a[i],a[rev[i]]);
    for(int l=2;l<=n;l<<=1){
        int l2=l>>1;
        for(int i=0;i<n;i+=l){
            for(int j=0;j<l2;j++){
                LL t1=a[i+j],t2=a[i+j+l2]*w[MAXW/l*j];
                a[i+j]=(t1+t2)%P;
                a[i+j+l2]=(t1-t2)%P;
            }
        }
    }
    if(flag==-1){
        LL invn=getPow(n,P-2);
        for(int i=0;i<n;i++) a[i]=a[i]*invn%P;
        for(int i=1;i<n;i++) if(i<n-i) swap(a[i],a[n-i]);
    }
}

void getInv(LL *b,LL *a,int len){
    if(len==1){
        b[0]=getPow(a[0],P-2);
        return;
    }
    static LL t1[MAXN],t2[MAXN];
    getInv(t1,a,len>>1);
    for(int i=len>>1;i<len<<1;i++) t1[i]=t2[i]=0;
    for(int i=0;i<len;i++) t2[i]=a[i];
    FFT(t1,len<<1,1); FFT(t2,len<<1,1);
    for(int i=0;i<len<<1;i++) t1[i]=(2*t1[i]-t2[i]*t1[i]%P*t1[i])%P;
    FFT(t1,len<<1,-1);
    for(int i=0;i<len;i++) b[i]=t1[i];
}

void getLn(LL *b,LL *a,int len){
    static LL t1[MAXN],t2[MAXN];
    for(int i=0;i<len;i++) t1[i]=a[i+1]*(i+1)%P,t2[i]=a[i];
    t1[len-1]=0;
    getInv(t2,t2,len);
    for(int i=0;i<len;i++) t1[i+len]=t2[i+len]=0;
    FFT(t1,len<<1,1); FFT(t2,len<<1,1);
    for(int i=0;i<len<<1;i++) t1[i]=t1[i]*t2[i]%P;
    FFT(t1,len<<1,-1);
    for(int i=1;i<len;i++) b[i]=t1[i-1]*inv_v[i]%P;
    b[0]=0;
}

void getExp(LL *b,LL *a,int len){
    if(len==1){
        b[0]=1;
        return;
    }
    static LL t1[MAXN],t2[MAXN],t3[MAXN];
    getExp(t1,a,len>>1);
    for(int i=0;i<len>>1;i++) t3[i]=t1[i],t1[i+(len>>1)]=t3[i+(len>>1)]=0;
    getLn(t2,t1,len);
    for(int i=0;i<len>>1;i++) t2[i]=t2[i+(len>>1)]-a[i+(len>>1)],t2[i+(len>>1)]=0;
    FFT(t3,len,1); FFT(t2,len,1);
    for(int i=0;i<len;i++) t3[i]=t3[i]*t2[i]%P;
    FFT(t3,len,-1);
    for(int i=0;i<len>>1;i++) b[i]=t1[i],b[i+(len>>1)]=-t3[i];
}

void getPow(LL *b,LL *a,int len,LL k){
    static LL t1[MAXN];
    getLn(t1,a,len);
    for(int i=0;i<len;i++) t1[i]=t1[i]*k%P;
    getExp(b,t1,len);
}

void solve(){
    for(sizew=1;sizew<=n;sizew<<=1);
    static LL t1[MAXN],t2[MAXN];
    t1[0]=1; t1[1]=6; t1[2]=1;
    getPow(t1,t1,sizew,inv_v[2]);
    for(int i=0;i<sizew;i++) t2[i]=t1[i];
    t2[0]=(t2[0]+1)%P; t2[1]=(t2[1]+1)%P;
    for(int i=0;i<sizew;i++) t2[i]=t2[i]*inv_v[2]%P;
    getPow(t2,t2,sizew,m+1);
    getInv(t1,t1,sizew);
    for(int i=0;i<sizew;i++) t1[i+sizew]=t2[i+sizew]=0;
    FFT(t1,sizew<<1,1); FFT(t2,sizew<<1,1);
    for(int i=0;i<sizew<<1;i++) t1[i]=t1[i]*t2[i]%P;
    FFT(t1,sizew<<1,-1);
    for(int i=1;i<=n;i++) printf("%lld ",(t1[i]+P)%P);
}

int main(){
#ifdef DEBUG
    freopen("127.in","r",stdin);
#endif
    scanf("%d%d",&m,&n);
    int _n=n;
    if(n>m) n=m;
    init();
    solve();
    for(int i=1;i<=_n-n;i++) printf("0 ");
    return 0;
}

查看详细内容