CodeForces 568 B.Symmetric and Transitive(dp+生成函数+FFT)

202 篇文章 1 订阅
25 篇文章 0 订阅

Description

找出 n n 个点的二元关系中满足对称性和传递性但不满足自反性的个数

Input

一个整数n(1n4000)

Output

输出满足条件的二元关系数,结果模 109+7 10 9 + 7

Sample Input

1

Sample Output

1

Solution1

把二元关系看作边,那么一个连通块里的点满足自反对称传递,只有存在单点集的图才对应一个合法的二元关系,枚举从 1 1 n枚举单点的个数,那么剩下的点可以任意划分,即贝尔数,下面考虑求出 i i 个点分出j个连通块的方案数 dp[i][j] d p [ i ] [ j ] ,考虑第 i i 个元素,要么新开一个连通块,要么并到之前的j个连通块中,故有转移方程 dp[i][j]+=dp[i1][j1]+jdp[i1][j] d p [ i ] [ j ] + = d p [ i − 1 ] [ j − 1 ] + j ⋅ d p [ i − 1 ] [ j ] ,令 f[i]=j=1idp[i][j] f [ i ] = ∑ j = 1 i d p [ i ] [ j ] ,则 f[i] f [ i ] 表示将 i i 个点分成若干连通块的方案数,答案即为i=1nCnif[i],时间复杂度 O(n2) O ( n 2 )

Code1

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#include<ctime>
using namespace std;
typedef long long ll;
typedef pair<int,int>P;
const int INF=0x3f3f3f3f,maxn=4005;
#define mod 1000000007
int n,fact[maxn],inv[maxn],dp[maxn][maxn],f[maxn];
void init(int n=4000)
{
    fact[0]=1;
    for(int i=1;i<=n;i++)fact[i]=(ll)i*fact[i-1]%mod;
    inv[1]=1;
    for(int i=2;i<=n;i++)inv[i]=mod-(ll)(mod/i)*inv[mod%i]%mod;
    inv[0]=1;
    for(int i=1;i<=n;i++)inv[i]=(ll)inv[i-1]*inv[i]%mod;
    dp[0][0]=f[0]=1;
    for(int i=1;i<=n;i++)
        for(int j=1;j<=i;j++)
        {
            dp[i][j]=(dp[i-1][j-1]+(ll)j*dp[i-1][j]%mod)%mod;
            f[i]=(f[i]+dp[i][j])%mod;
        }
}
int C(int n,int m)
{
    return (ll)fact[n]*inv[m]%mod*inv[n-m]%mod;
}
int main()
{
    init();
    while(~scanf("%d",&n))
    {
        int ans=0;
        for(int i=0;i<n;i++)ans=(ans+(ll)C(n,i)*f[i]%mod)%mod;
        printf("%d\n",ans);
    }
    return 0;
}

Solution2

n n 个点及它们之间的二元关系看作一张无向图,问题转化为求没有自环且满足传递性的无向图个数,考虑每个连通块:规模为1的连通块有两种:单点集和自环,规模超过 1 1 的连通块只有一种,即完全图

先不考虑点的编号,则一个连通块的生成函数为A(x)=2x1!+x22!+x33!+...=ex+x1,若干连通块组成的图的指数型生成函数即为 eA(x)=eex+x1 e A ( x ) = e e x + x − 1 ​ ;若干不是自环的连通块组成的图的生成函数即为 A(x)x A ( x ) − x ​ ,故其指数型生成函数为 eex1 e e x − 1 ​ 。进而答案的生成函数为 B(x)=eex+x1eex1=(ex1)eex1=f(x)ef(x) B ( x ) = e e x + x − 1 − e e x − 1 = ( e x − 1 ) ⋅ e e x − 1 = f ( x ) ⋅ e f ( x ) ​ ,故只要对 f(x) f ( x ) ​ 求一下多项式 exp e x p ​ 然后和自身乘一下即可,此时考虑节点编号,答案即为 n![xn]B(x) n ! ⋅ [ x n ] B ( x ) ​ ,时间复杂度 O(nlogn) O ( n l o g n )

Code2

#include<cstdio>
#include<iostream>
#include<cstring>
#include<algorithm>
#include<cmath>
#include<vector>
#include<queue>
#include<map>
#include<set>
#include<ctime>
using namespace std;
typedef long long ll;
#define maxfft 8192+5
#define maxn 4005
#define mod 1000000007
const double pi=acos(-1.0);
struct cp
{
    double a,b;
    cp operator +(const cp &o)const {return (cp){a+o.a,b+o.b};}
    cp operator -(const cp &o)const {return (cp){a-o.a,b-o.b};}
    cp operator *(const cp &o)const {return (cp){a*o.a-b*o.b,b*o.a+a*o.b};}
    cp operator *(const double &o)const {return (cp){a*o,b*o};}
    cp operator !() const{return (cp){a,-b};}
}w[maxfft];
int pos[maxfft];
void fft_init(int len)
{
    int j=0;
    while((1<<j)<len)j++;
    j--;
    for(int i=0;i<len;i++)
        pos[i]=pos[i>>1]>>1|((i&1)<<j);
}
void fft(cp *x,int len,int sta)
{
    for(int i=0;i<len;i++)
        if(i<pos[i])swap(x[i],x[pos[i]]);
    w[0]=(cp){1,0};
    for(unsigned i=2;i<=len;i<<=1)
    {
        cp g=(cp){cos(2*pi/i),sin(2*pi/i)*sta};
        for(int j=i>>1;j>=0;j-=2)w[j]=w[j>>1];
        for(int j=1;j<i>>1;j+=2)w[j]=w[j-1]*g;
        for(int j=0;j<len;j+=i)
        {
            cp *a=x+j,*b=a+(i>>1);
            for(int l=0;l<i>>1;l++)
            {
                cp o=b[l]*w[l];
                b[l]=a[l]-o;
                a[l]=a[l]+o;
            }
        }
    }
    if(sta==-1)for(int i=0;i<len;i++)x[i].a/=len,x[i].b/=len;
}
cp x[maxfft],y[maxfft],z[maxfft];
int temp[maxfft];
void FFT(int *a,int *b,int n,int m,int *c)
{
    if(n<=100&&m<=100||min(n,m)<=5)
    {
        for(int i=0;i<n+m-1;i++)temp[i]=0;
        for(int i=0;i<n;i++)
            for(int j=0;j<m;j++)
            {
                temp[i+j]+=(ll)a[i]*b[j]%mod;
                if(temp[i+j]>=mod)temp[i+j]-=mod;
            }
        for(int i=0;i<n+m-1;i++)c[i]=temp[i];
        return ;
    }
    int len=1;
    while(len<n+m)len<<=1;
    fft_init(len);
    for(int i=0;i<len;i++)
    {
        int aa=i<n?a[i]:0,bb=i<m?b[i]:0;
        x[i]=(cp){(aa>>15),(aa&32767)},y[i]=(cp){(bb>>15),(bb&32767)};
    }
    fft(x,len,1),fft(y,len,1);
    for(int i=0;i<len;i++)
    {
        int j=len-1&len-i;
        z[i]=((x[i]+!x[j])*(y[i]-!y[j])+(x[i]-!x[j])*(y[i]+!y[j]))*(cp){0,-0.25};
    }
    fft(z,len,-1);
    for(int i=0;i<n+m-1;i++)
    {
        ll ta=(ll)(z[i].a+0.5)%mod;
        ta=(ta<<15)%mod;
        c[i]=ta;
    }
    for(int i=0;i<len;i++)
    {
        int j=len-1&len-i;
        z[i]=(x[i]-!x[j])*(y[i]-!y[j])*(cp){-0.25,0}+(x[i]+!x[j])*(y[i]+!y[j])*(cp){0,0.25};
    }
    fft(z,len,-1);
    for(int i=0;i<n+m-1;i++)
    {
        ll ta=(ll)(z[i].a+0.5)%mod,tb=(ll)(z[i].b+0.5)%mod;
        ta=(ta+(tb<<30))%mod;
        c[i]=(c[i]+ta)%mod;
    }
}
int Pow(int a,int b)
{
    int ans=1;
    while(b)
    {
        if(b&1)ans=(ll)ans*a%mod;
        a=(ll)a*a%mod;
        b>>=1;
    }
    return ans;
} 
int inv[maxfft],sum[maxfft];
void init(int n=8192)
{
    inv[1]=1;
    for(int i=2;i<=n;i++)inv[i]=mod-(ll)(mod/i)*inv[mod%i]%mod;
    sum[0]=1;
    for(int i=1;i<=n;i++)sum[i]=(ll)sum[i-1]*inv[i]%mod;
}
int temp1[maxfft],temp2[maxfft],temp3[maxfft];
void Poly_Inv(int *poly,int n,int *ans)
{
    ans[0]=Pow(poly[0],mod-2);//inv[poly[0]];
    for(int i=2;i<=n;i<<=1)
    {
        FFT(poly,ans,i,i/2,temp1);
        FFT(ans,temp1+i/2,i/2,i/2,temp1);
        for(int j=0;j<i/2;j++)ans[j+i/2]=temp1[j]==0?0:mod-temp1[j];
    }
}
void Poly_Log(int *poly,int n,int *ans)
{
    Poly_Inv(poly,n,temp2);
    for(int i=0;i<n-1;i++)ans[i]=(ll)poly[i+1]*(i+1)%mod;
    FFT(ans,temp2,n-1,n,ans);
    for(int i=n-1;i>0;i--)ans[i]=(ll)ans[i-1]*inv[i]%mod;
    ans[0]=0;
}
void Poly_Exp(int *poly,int n,int *ans)
{
    if(n==1)
    {
        ans[0]=1;
        return ;
    }
    Poly_Exp(poly,n/2,ans);
    Poly_Log(ans,n,temp3);
    for(int i=0;i<n;i++)
    {
        temp3[i]=poly[i]-temp3[i];
        if(temp3[i]<0)temp3[i]+=mod;    
    }
    temp3[0]++;
    if(temp3[0]==mod)temp3[0]=0;
    FFT(ans,temp3,n,n,ans);
    for(int i=n;i<2*n;i++)ans[i]=0;
}

int n,f[maxfft],g[maxfft];
int main()
{
    init();
    int n;
    while(~scanf("%d",&n))
    {
        memset(f,0,sizeof(f));
        int m=1;
        while(m<=n)m<<=1;
        for(int i=1;i<m;i++)f[i]=sum[i];
        Poly_Exp(f,m,g);
        FFT(f,g,m,m,f);
        int ans=f[n];
        for(int i=2;i<=n;i++)ans=(ll)ans*i%mod;
        printf("%d\n",ans);
    }
    return 0;
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值