• 周二. 3月 28th, 2023

5G编程聚合网

5G时代下一个聚合的编程学习网

热门标签

洛谷 P4831 Scarlet loves WenHuaKe

admin

11月 28, 2021

一题多解

这题首先可以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)));
}
“`]

发表回复

您的电子邮箱地址不会被公开。 必填项已用*标注