洛谷 P4831 Scarlet loves WenHuaKe

一题多解

这题首先可以DP
直接三个类型暴力转移。
(f[i][j])表示无限制填(i)行(j)列的方案数。
(f1[i][j])表示有一列已经有棋子的方案数。
(f2[i][j])表示有两列有棋子的方案数。
转移一下。复杂度(O(n*(m-n)))

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int M=998244353;
const int N=100010;
int fac[N],invfac[N];
int f[N];
int n,m;
int C(int x,int y){
    return (ll)fac[x]*invfac[y]%M*invfac[x-y]%M;
}
int fp(int x,int y){
    int ret=1;
    for (; y; y>>=1,x=(ll)x*x%M)
    if (y&1) ret=(ll)ret*x%M;
    return ret;
}
int sqr(int x){
    return (ll)x*x%M;
}
int MUL(int x,int y){
    return (ll)x*y%M;
}
int F(int,int);
//int debug;
int f__k[N][20];
int F1(int n,int m){
    //++debug;
    //cerr<<"F_!"<<n<<" "<<m<<" "<<n-m<<endl;
    if (n>m) return 0;
    if (n==0) return 1;
    if (m<2) return 0;
    int &ret=f__k[n][m-n];
    if (ret) return ret-1;
    ret=(MUL(MUL(n,m-1),F1(n-1,m-1))+F(n,m-1))%M+1;
    return ret-1;
}
int mm[N][20];
int F2(int n,int m){
    //getchar();
    //cerr<<m-n<<endl;
    if (n>m) return 0;
    if (n==0) return 1;
    if (m<2) return 0;
    int &ret=mm[n][m-n];
    if (ret) return ret-1;
    //cerr<<"FF"<<n<<" "<<m<<endl;
    ret=F(n,m-2);
    //cerr<<"ORG"<<n<<" "<<m<<" "<<ret<<endl;
    if (n>=2){
      (ret+=MUL(2,MUL(MUL(C(n,2),F2(n-2,m-2)),MUL(m-2,m-3))))%=M;
      (ret+=MUL(2,MUL(MUL(C(n,2),F(n-2,m-3)),m-2)))%=M;
    }
    //cerr<<"ORGG"<<n<<" "<<m<<" "<<ret<<endl;
      if (n>=1)
      (ret+=MUL(2,MUL(F1(n-1,m-2),MUL(C(n,1),m-2))))%=M;
      //cerr<<"OGGGG"<<n<<" "<<m<<" "<<ret<<endl;
      if (n>=1)
      (ret+=MUL(C(n,1),F(n-1,m-2)))%=M;
      //cerr<<"FFRET"<<n<<" "<<m<<" "<<ret<<endl;
      ++ret;
      //++debug;
      return ret-1;
}
int mp[N][20];
int F(int n,int m){
    //cerr<<"F"<<n<<" "<<m<<endl;
    if (n>m) return 0;
    if (n==0) return 1;
    if (m<2) return 0;
    //cerr<<m-n<<endl;
    int &ret=mp[n][m-n];
    if (ret) return ret-1;
    //cerr<<"F"<<n<<" "<<m<<" "<<mp.size()<<endl;
    ret=MUL(C(m,2),F2(n-1,m))+1;
    //if (debug%300==0) cerr<<"FFFF"<<n<<" "<<m<<" "<<ret<<" "<<debug<<endl;
    return ret-1;
}
int FF[2010][2010];
int main(){
    scanf("%d%d",&n,&m);
    if (n<=2000){
    FF[0][0]=1;
    for (int i=0; i<=n; ++i){
        for (int j=0; j<=m; ++j)
        if ((i*2-j)%2==0){
            int k=(i*2-j)/2;
            //cerr<<i<<" "<<j<<" "<<k<<" "<<f[i][j][k]<<endl; 
            int z0=m-k-j,c=FF[i][j];
            //cerr<<i<<" "<<j<<" "<<z0<<endl;
            if (j>=1&&z0>=1) (FF[i+1][j]+=(ll)j*z0%M*c%M)%=M;
            if (z0>=2) (FF[i+1][j+2]+=((ll)z0*(z0-1)/2)%M*c%M)%=M;
            if (j>=2) (FF[i+1][j-2]+=((ll)j*(j-1)/2)%M*c%M)%=M;
        }
    }
    int ans=0;
    for (int i=0; i<=m; ++i){
        //cerr<<i<<" "<<f[n][i]<<endl;
        (ans+=FF[n][i])%=M;
    }
    cout<<ans;
    return 0;
    }
    fac[0]=1;
    for (int i=1; i<=m; ++i) fac[i]=(ll)fac[i-1]*i%M;
    invfac[m]=fp(fac[m],M-2);
    for (int i=m-1; i>=0; --i) invfac[i]=(ll)invfac[i+1]*(i+1)%M;
    //cerr<<invfac[2]<<" "<<invfac[3]<<endl;
    cout<<F(n,m)<<endl;
    //cerr<<debug<<endl;
}

还可以用题解的反演，没有实现过。
接下来就是重头戏了。
还可以用生成函数和二分图计数，
将原题抽象成(n+m)个点的图。
假设在第(i)行，第(j)列放了一个棋子，
可视为将连接了(i)号点和(n+j)号点的无向边。

由题目条件，可知(le n)的点，度为(2)。
$ > n$ 的点，度(le 2)。
那么思考一下，图的形态就是若干个环，和(m-n)条起点(>n)，终点(>n)的链。
由题目条件，可知(le n)的点，度为(2)。
$ > n$ 的点，度$ le 2$。

那么思考一下，图的形态就是若干个环，和(m-n)条起点(>n)，终点(>n)的链。
写出链的函数表示该链有(x)个大于(n)的点，(x-1)个小于等于(n)的点

[f(0)=0 ]

[f(1)=1$$（单点合法） $$f(x)=frac{x!*(x-1)!}{2} (x>1)]

写出环的函数表示该环有(x)个大于(n)的点，(x)个小于等于(n)的点

[g(0)=0 ]

[g(1)=0$$（长为$2$的环不存在） $$g(x)=frac{x!*(x-1)!}{2} (x>1)]

虽然有了这个东西，但是直接写到普通型生成函数里还是不行的。

将链写进

[F(x)= sum_{i=0}^{infty} frac{f(x)}{x!*(x-1)!} x^i ]

将环写进

[G(x)= sum_{i=0}^{infty} frac{g(x)}{(x!)^2} x^i ]

这样

[F^2(x)=sum_{i=0}^{infty} sum_{j=0}^i frac{f(j)*f(i-j)}{j!*(j-1)!*(i-j)!*(i-j-1)!} x^i ]

[F^2(x)=sum_{i=0}^{infty} sum_{j=0}^i frac{f(j)*f(i-j)* binom{i}{j} * binom{i-2}{j-1}}{i!*(i-2)!} x^i ]

[G^2(x)=sum_{i=0}^{infty} sum_{j=0}^i frac{g(j)*g(i-j)}{(j!)^2*((i-j)!)^2} x^i ]

[G^2(x)=sum_{i=0}^{infty} sum_{j=0}^i frac{g(j)*g(i-j)* binom{i}{j}^2}{(i!)^2} x^i ]

相当的优美，链的卷积考虑了二分图上两边选择的组合数，环的卷积也是。
对于链

[ binom{i}{j}$$为在$>n $的$i$个点中选$j$个 $$ binom{i-2}{j-1}$$为在$ le n $的$i-2$个点中选$j-1$个 环类似。 那么既然卷积有组合意义，那么答案为 $F(x)^{m-n}$表示把链搞在一起。（这里对链的顺序会算重） $G(x)^i$表示把$i$个环搞在一起。（这里对环也会算重） 会算重，是因为把计算的是排列，而要求的是组合。 考虑$sum_{i=0}^{infty} frac{G^i(x)}{i!}$ 那么意义就是环的组合的上面定义的那一种生成函数。 即$e^{G(x)}$ $F^{m-n}(x)/(m-n)!$再取系数后就是链的组合。 那么答案就只需要划分集合后乘上对应的链和环即可。 可以不写多项式exp，利用$G(x)$的特殊性质。 复杂度为$O((m-n)*m log (m))$ 如果多项式快速幂，$O(m log(m))$ 暴力exp版 ```cpp // luogu-judger-enable-o2 #include <bits/stdc++.h> using namespace std; #define FAST typedef long long ll; const int M=998244353,G=3; struct poly{ vector<int> g; __attribute((always_inline)) size_t size() const { return g.size(); } __attribute((always_inline)) poly(){ } __attribute((always_inline)) poly(size_t x):g(x){ } __attribute((always_inline)) void resize(size_t x){ g.resize(x); } __attribute((always_inline)) poly cat(size_t x){ g.resize(x); return *this; } __attribute((always_inline)) int& operator [](int x){ return g[x]; } __attribute((always_inline)) int operator [](int x) const { return g[x]; } void operator =(int x){ g.resize(1); g[0]=x; } friend istream& operator >>(istream &is,poly &A){ int n; is>>n; A.resize(n); for (int i=0; i<n; ++i) is>>A[i]; return is; } friend ostream& operator <<(ostream &os,const poly &A){ for (size_t i=0; i<A.size(); ++i) os<<A[i]<<" "; return os; } }; __attribute((always_inline)) int fp(int x,int y){ int ret=1; for (; y; y>>=1,x=(ll)x*x%M) if (y&1) ret=(ll)ret*x%M; return ret; } __attribute((always_inline)) int D(int x){ return x>=M?x-M:x; } __attribute((always_inline)) int U(int x){ return x<0?x+M:x; } const int P=270000; int a[P],w[P],b[P],rev[P]; void fft(int *a,int n){ for (int i=0; i<n; ++i) if (i>rev[i]) swap(a[i],a[rev[i]]); for (int i=1; i<n; i<<=1){ w[0]=1; w[1]=fp(G,(M-1)/(i<<1)); for (int j=2; j<i; ++j) w[j]=(ll)w[j-1]*w[1]%M; for (int j=0; j<n; j+=(i<<1)) for (int k=j; k<j+i; ++k){ int x=a[k],y=(ll)a[k+i]*w[k-j]%M; a[k]=x+y>=M?x+y-M:x+y; a[k+i]=x-y<0?x-y+M:x-y; // a[k]=D(x+y); // a[k+i]=U(x-y); } } } void print(const poly &A){ for (size_t i=0; i<A.size(); ++i) cout<<A[i]<<" "; cout<<endl; } poly operator *(const poly &A,const int &x){ poly B(A.size()); for (size_t i=0; i<A.size(); ++i) B[i]=(ll)A[i]*x%M; return B; } poly operator +(const poly &A,const int &x){//x>=0 x<M poly B=A; B[0]=D(B[0]+x); return B; } poly operator -(const poly &A,const int &x){//x>=0 x<M poly B=A; B[0]=U(B[0]-x); return B; } poly operator *(const int &x,const poly &A){ poly B(A.size()); for (size_t i=0; i<A.size(); ++i) B[i]=(ll)A[i]*x%M; return B; } poly operator +(const int &x,const poly &A){//x>=0 x<M poly B=A; B[0]=D(B[0]+x); return B; } poly operator -(const int &x,const poly &A){//x>=0 x<M poly B=A; B[0]=U(B[0]-x); return B; } poly operator *(const poly &A,const poly &B){ size_t u,bit=0; for (u=1; u<A.size()+B.size()-1; u<<=1,++bit); --bit; for (size_t i=0; i<A.size(); ++i) a[i]=A[i]; memset(a+A.size(),0,sizeof(*a)*(u-A.size())); for (size_t i=0; i<B.size(); ++i) b[i]=B[i]; memset(b+B.size(),0,sizeof(*b)*(u-B.size())); for (size_t i=0; i<u; ++i) rev[i]=rev[i>>1]>>1|((i&1)<<bit); fft(a,u); fft(b,u); for (size_t i=0; i<u; ++i) a[i]=(ll)a[i]*b[i]%M; reverse(a+1,a+u); fft(a,u); int ni=fp(u,M-2); for (size_t i=0; i<u; ++i) a[i]=(ll)a[i]*ni%M; poly ret(A.size()+B.size()-1); for (size_t i=0; i<ret.size(); ++i) ret[i]=a[i]; return ret; } poly cat(const poly &_,size_t x){//_=cat(_,x) is equal to _.resize(x) poly ret(x); if (_.size()>=x) for (size_t i=0; i<x; ++i) ret[i]=_[i]; else for (size_t i=0; i<_.size(); ++i) ret[i]=_[i]; return ret; } poly pow(const poly &A,const int &x){//A^x has limit size_t u,bit=0; for (u=1; u<(A.size()*x<<1)-1; u<<=1,++bit); --bit; for (size_t i=0; i<A.size(); ++i) a[i]=A[i]; memset(a+A.size(),0,sizeof(*a)*(u-A.size())); for (size_t i=0; i<u; ++i) rev[i]=rev[i>>1]>>1|((i&1)<<bit); fft(a,u); for (size_t i=0; i<u; ++i) a[i]=fp(a[i],x); reverse(a+1,a+u); fft(a,u); int ni=fp(u,M-2); for (size_t i=0; i<u; ++i) a[i]=(ll)a[i]*ni%M; poly ret((A.size()*x<<1)-1); for (size_t i=0; i<ret.size(); ++i) ret[i]=a[i]; return ret; } poly pow(const poly &A,const int &x,const int &C){ poly ret(1); ret[0]=1; poly B=cat(A,C); #ifdef FAST for (int y=x; y; y>>=1,B=(B*B).cat(C)) if (y&1) ret=(ret*B).cat(C); #else for (int y=x; y; y>>=1,B=cat(B*B,C)) if (y&1) ret=cat(ret*B,C); #endif return ret; } poly operator +(const poly &A,const poly &B){ if (A.size()>=B.size()){ poly ret=A; for (size_t i=0; i<B.size(); ++i) ret[i]=D(ret[i]+B[i]); return ret; } poly ret=B; for (size_t i=0; i<A.size(); ++i) ret[i]=D(ret[i]+A[i]); return ret; } poly operator -(const poly &A,const poly &B){ poly ret=A; if (ret.size()<B.size()) ret.resize(B.size()); for (size_t i=0; i<B.size(); ++i) ret[i]=U(ret[i]-B[i]); return ret; } poly getrev(const poly &A){ size_t u; for (u=1; u<A.size()+A.size()-1; u<<=1); poly B(1); B[0]=fp(A[0],M-2); #ifdef FAST for (size_t i=2; i<=u; i<<=1) B=B+B-(cat(A,i>>1)*(B*B).cat(i>>1)).cat(i); #else for (size_t i=2; i<=u; i<<=1) B=B+B-cat(cat(A,i>>1)*cat(B*B,i>>1),i); #endif return B; } poly diff(const poly &A){ poly B(A.size()-1); for (size_t i=0; i<B.size(); ++i) B[i]=(ll)A[i+1]*(i+1)%M; return B; } #ifdef FAST int inv[P]; #endif poly inte(const poly &A){//slow poly B(A.size()+1); #ifdef FAST //prework inv for (size_t i=1; i<B.size(); ++i) B[i]=(ll)A[i-1]*inv[i]%M; #else for (size_t i=1; i<B.size(); ++i) B[i]=(ll)A[i-1]*fp(i,M-2)%M; #endif return B; } poly ln(const poly &A){//A[0]=1 #ifdef FAST return inte(diff(A)*getrev(A).cat(A.size())); #else return inte(diff(A)*cat(getrev(A),A.size())); #endif } poly exp(const poly &A){//A[0]=0 size_t u; for (u=1; u<A.size()+A.size()-1; u<<=1); poly B(1); B[0]=1; #ifdef FAST for (size_t i=2; i<=u; i<<=1) B=B+(B*(cat(A,i>>1)-ln(B).cat(i>>1))).cat(i); #else for (size_t i=2; i<=u; i<<=1) B=B+cat(B*(cat(A,i>>1)-cat(ln(B),i>>1)),i); #endif return B; } poly sqrt(const poly &A){//??? size_t u; for (u=1; u<A.size()+A.size()-1; u<<=1); poly B(1); B[0]=sqrt(A[0]); #ifdef FAST for (size_t i=2; i<=u; i<<=1) B=((cat(A,i>>1)+(B*B).cat(i>>1))*(getrev(B+B)).cat(i>>1)).cat(i); #else for (size_t i=2; i<=u; i<<=1) B=cat((cat(A,i>>1)+cat(B*B,i>>1))*cat(getrev(B+B),i>>1),i); #endif return B; } namespace mymath{ const int N=100010; int rec[N]; int HALF(int x){ return x&1?(x+M)>>1:x>>1; } int ADD(int x,int y){ return (x+=y)>=M?x-M:x; } int SUB(int x,int y){ return (x-=y)<0?x+M:x; } int MUL(int x,int y){ return (ll)x*y%M; } int qp(int x,int y=M-2){ int ret=1; for (; y; y>>=1,x=MUL(x,x)) if (y&1) ret=MUL(ret,x); return ret; } int sqr(int x){ return (ll)x*x%M; } int frac(int x){ if (x==0) return 1; if (rec[x]) return rec[x]; return rec[x]=MUL(frac(x-1),x); } int choose(int x,int y){ if (x<y) return 0; return MUL(frac(x),qp(MUL(frac(y),frac(x-y)))); } } int main(){ #ifdef FAST ios::sync_with_stdio(0); inv[0]=inv[1]=1; for (int i=2; i<P; ++i) inv[i]=(ll)inv[M%i]*(M-M/i)%M; #endif int n,m; cin>>n>>m; //cerr<<n<<" "<<m<<endl; poly b,a; b.resize(1); b[0]=1; a.resize(m+1); a[0]=0; a[1]=1; for (int i=2; i<=m; ++i) a[i]=mymath::HALF(1); for (int i=1; i<=m-n; ++i){ b=b*a; b.cat(m+1); } b.cat(m+1); poly f; f.resize(m+1); f[0]=0; f[1]=0; for (int i=2; i<=m; ++i) f[i]=mymath::qp(2*i); f=exp(f); f.resize(m+1); for (int i=m-n; i<=m; ++i) b[i]=mymath::MUL(mymath::choose(m,i),mymath::MUL(b[i],mymath::MUL(mymath::frac(i),mymath::frac(i-m+n)))); for (int i=0; i<=m; ++i) f[i]=mymath::MUL(mymath::choose(n,i),mymath::MUL(f[i],mymath::sqr(mymath::frac(i)))); b=(b*f).cat(m+1); cout<<mymath::MUL(b[m],mymath::qp(mymath::frac(m-n))); } ``` 非暴力exp版 ```cpp // luogu-judger-enable-o2 // luogu-judger-enable-o2 #include <bits/stdc++.h> using namespace std; #define FAST typedef long long ll; const int M=998244353,G=3; struct poly{ vector<int> g; __attribute((always_inline)) size_t size() const { return g.size(); } __attribute((always_inline)) poly(){ } template<class T> __attribute((always_inline)) poly(const T &x):g(x){ } __attribute((always_inline)) void resize(size_t x){ g.resize(x); } __attribute((always_inline)) poly cat(size_t x){ g.resize(x); return *this; } __attribute((always_inline)) int& operator [](int x){ return g[x]; } __attribute((always_inline)) int operator [](int x) const { return g[x]; } void operator =(int x){ g.resize(1); g[0]=x; } friend istream& operator >>(istream &is,poly &A){ int n; is>>n; A.resize(n); for (int i=0; i<n; ++i) is>>A[i]; return is; } friend ostream& operator <<(ostream &os,const poly &A){ for (size_t i=0; i<A.size(); ++i) os<<A[i]<<" "; return os; } }; __attribute((always_inline)) int fp(int x,int y){ int ret=1; for (; y; y>>=1,x=(ll)x*x%M) if (y&1) ret=(ll)ret*x%M; return ret; } __attribute((always_inline)) int D(int x){ return x>=M?x-M:x; } __attribute((always_inline)) int U(int x){ return x<0?x+M:x; } __attribute((always_inline)) int ADD(int x,int y){ return (x+=y)>=M?x-M:x; } __attribute((always_inline)) int SUB(int x,int y){ return (x-=y)<0?x+M:x; } const int P=270000; int a[P],w[P],b[P],rev[P]; void fft(int *a,int n){ for (int i=0; i<n; ++i) if (i>rev[i]) swap(a[i],a[rev[i]]); for (int i=1; i<n; i<<=1){ w[0]=1; w[1]=fp(G,(M-1)/(i<<1)); for (int j=2; j<i; ++j) w[j]=(ll)w[j-1]*w[1]%M; for (int j=0; j<n; j+=(i<<1)) for (int k=j; k<j+i; ++k){ int x=a[k],y=(ll)a[k+i]*w[k-j]%M; a[k]=ADD(x,y); a[k+i]=SUB(x,y); // a[k]=D(x+y); // a[k+i]=U(x-y); } } } void print(const poly &A){ for (size_t i=0; i<A.size(); ++i) cout<<A[i]<<" "; cout<<endl; } poly operator *(const poly &A,const int &x){ poly B(A.size()); for (size_t i=0; i<A.size(); ++i) B[i]=(ll)A[i]*x%M; return B; } poly operator +(const poly &A,const int &x){//x>=0 x<M poly B=A; B[0]=D(B[0]+x); return B; } poly operator -(const poly &A,const int &x){//x>=0 x<M poly B=A; B[0]=U(B[0]-x); return B; } poly operator *(const int &x,const poly &A){ poly B(A.size()); for (size_t i=0; i<A.size(); ++i) B[i]=(ll)A[i]*x%M; return B; } poly operator +(const int &x,const poly &A){//x>=0 x<M poly B=A; B[0]=D(B[0]+x); return B; } poly operator -(const int &x,const poly &A){//x>=0 x<M poly B=A; B[0]=U(B[0]-x); return B; } poly operator *(const poly &A,const poly &B){ size_t u,bit=0; for (u=1; u<A.size()+B.size()-1; u<<=1,++bit); --bit; for (size_t i=0; i<A.size(); ++i) a[i]=A[i]; memset(a+A.size(),0,sizeof(*a)*(u-A.size())); for (size_t i=0; i<B.size(); ++i) b[i]=B[i]; memset(b+B.size(),0,sizeof(*b)*(u-B.size())); for (size_t i=0; i<u; ++i) rev[i]=rev[i>>1]>>1|((i&1)<<bit); fft(a,u); fft(b,u); for (size_t i=0; i<u; ++i) a[i]=(ll)a[i]*b[i]%M; reverse(a+1,a+u); fft(a,u); int ni=fp(u,M-2); for (size_t i=0; i<u; ++i) a[i]=(ll)a[i]*ni%M; poly ret(A.size()+B.size()-1); for (size_t i=0; i<ret.size(); ++i) ret[i]=a[i]; return ret; } poly cat(const poly &_,size_t x){//_=cat(_,x) is equal to _.resize(x) poly ret(x); if (_.size()>=x) for (size_t i=0; i<x; ++i) ret[i]=_[i]; else for (size_t i=0; i<_.size(); ++i) ret[i]=_[i]; return ret; } poly pow(const poly &A,const int &x){//A^x has limit size_t u,bit=0; for (u=1; u<(A.size()*x<<1)-1; u<<=1,++bit); --bit; for (size_t i=0; i<A.size(); ++i) a[i]=A[i]; memset(a+A.size(),0,sizeof(*a)*(u-A.size())); for (size_t i=0; i<u; ++i) rev[i]=rev[i>>1]>>1|((i&1)<<bit); fft(a,u); for (size_t i=0; i<u; ++i) a[i]=fp(a[i],x); reverse(a+1,a+u); fft(a,u); int ni=fp(u,M-2); for (size_t i=0; i<u; ++i) a[i]=(ll)a[i]*ni%M; poly ret((A.size()*x<<1)-1); for (size_t i=0; i<ret.size(); ++i) ret[i]=a[i]; return ret; } poly pow(const poly &A,const int &x,const int &C){ poly ret(1); ret[0]=1; poly B=cat(A,C); #ifdef FAST for (int y=x; y; y>>=1,B=(B*B).cat(C)) if (y&1) ret=(ret*B).cat(C); #else for (int y=x; y; y>>=1,B=cat(B*B,C)) if (y&1) ret=cat(ret*B,C); #endif return ret; } poly operator +(const poly &A,const poly &B){ if (A.size()>=B.size()){ poly ret=A; for (size_t i=0; i<B.size(); ++i) ret[i]=D(ret[i]+B[i]); return ret; } poly ret=B; for (size_t i=0; i<A.size(); ++i) ret[i]=D(ret[i]+A[i]); return ret; } poly operator -(const poly &A,const poly &B){ poly ret=A; if (ret.size()<B.size()) ret.resize(B.size()); for (size_t i=0; i<B.size(); ++i) ret[i]=U(ret[i]-B[i]); return ret; } poly getrev(const poly &A){ size_t u; for (u=1; u<A.size()+A.size()-1; u<<=1); poly B(1); B[0]=fp(A[0],M-2); #ifdef FAST for (size_t i=2; i<=u; i<<=1) B=B+B-(cat(A,i>>1)*(B*B).cat(i>>1)).cat(i); #else for (size_t i=2; i<=u; i<<=1) B=B+B-cat(cat(A,i>>1)*cat(B*B,i>>1),i); #endif return B; } poly diff(const poly &A){ poly B(A.size()-1); for (size_t i=0; i<B.size(); ++i) B[i]=(ll)A[i+1]*(i+1)%M; return B; } #ifdef FAST int inv[P]; #endif poly inte(const poly &A){//slow poly B(A.size()+1); #ifdef FAST //prework inv for (size_t i=1; i<B.size(); ++i) B[i]=(ll)A[i-1]*inv[i]%M; #else for (size_t i=1; i<B.size(); ++i) B[i]=(ll)A[i-1]*fp(i,M-2)%M; #endif return B; } poly ln(const poly &A){//A[0]=1 #ifdef FAST return inte(diff(A)*getrev(A).cat(A.size())); #else return inte(diff(A)*cat(getrev(A),A.size())); #endif } poly exp(const poly &A){//A[0]=0 size_t u; for (u=1; u<A.size()+A.size()-1; u<<=1); poly B(1); B[0]=1; #ifdef FAST for (size_t i=2; i<=u; i<<=1) B=B+(B*(cat(A,i>>1)-ln(B).cat(i>>1))).cat(i); #else for (size_t i=2; i<=u; i<<=1) B=B+cat(B*(cat(A,i>>1)-cat(ln(B),i>>1)),i); #endif return B; } poly sqrt(const poly &A){//??? size_t u; for (u=1; u<A.size()+A.size()-1; u<<=1); poly B(1); B[0]=sqrt(A[0]); #ifdef FAST for (size_t i=2; i<=u; i<<=1) B=((cat(A,i>>1)+(B*B).cat(i>>1))*(getrev(B+B)).cat(i>>1)).cat(i); #else for (size_t i=2; i<=u; i<<=1) B=cat((cat(A,i>>1)+cat(B*B,i>>1))*cat(getrev(B+B),i>>1),i); #endif return B; } namespace mymath{ const int N=100010; int rec[N]; int HALF(int x){ return x&1?(x+M)>>1:x>>1; } int ADD(int x,int y){ return (x+=y)>=M?x-M:x; } int SUB(int x,int y){ return (x-=y)<0?x+M:x; } int MUL(int x,int y){ return (ll)x*y%M; } int qp(int x,int y=M-2){ int ret=1; for (; y; y>>=1,x=MUL(x,x)) if (y&1) ret=MUL(ret,x); return ret; } int sqr(int x){ return (ll)x*x%M; } int frac(int x){ if (x==0) return 1; if (rec[x]) return rec[x]; return rec[x]=MUL(frac(x-1),x); } int choose(int x,int y){ if (x<y) return 0; return MUL(frac(x),qp(MUL(frac(y),frac(x-y)))); } } int n,m; void doit(poly &f){ int now=1; f.resize(m+1); for (int i=0; i<=m; ++i){ f[i]=now; now=mymath::MUL(now,mymath::HALF(1)); now=mymath::MUL(now,inv[i+1]); now=mymath::SUB(0,now); } poly g(m+1); int fz=1,fm=1; for (int i=0; i<=m; ++i){ g[i]=mymath::MUL(fz,fm); fz=mymath::MUL(fz,2*i+1); fm=mymath::MUL(fm,mymath::HALF(1)); fm=mymath::MUL(fm,inv[i+1]); } f=(f*g).cat(m+1); } int main(){ #ifdef FAST ios::sync_with_stdio(0); inv[0]=inv[1]=1; for (int i=2; i<P; ++i) inv[i]=(ll)inv[M%i]*(M-M/i)%M; #endif cin>>n>>m; //cerr<<n<<" "<<m<<endl; poly b,a; b.resize(1); b[0]=1; a.resize(m+1); a[0]=0; a[1]=1; for (int i=2; i<=m; ++i) a[i]=mymath::HALF(1); for (int i=1; i<=m-n; ++i){ b=b*a; b.cat(m+1); } b.cat(m+1); poly f; doit(f); for (int i=m-n; i<=m; ++i) b[i]=mymath::MUL(mymath::choose(m,i),mymath::MUL(b[i],mymath::MUL(mymath::frac(i),mymath::frac(i-m+n)))); for (int i=0; i<=m; ++i) f[i]=mymath::MUL(mymath::choose(n,i),mymath::MUL(f[i],mymath::sqr(mymath::frac(i)))); b=(b*f).cat(m+1); cout<<mymath::MUL(b[m],mymath::qp(mymath::frac(m-n))); } ```]