[Bluestein套路+多项式ln]2022“杭电杯”中国大学生算法设计超级联赛（7）1010【Connectivity of Erdős-Rényi Graph】题解

ZigZagK

2022年8月12日 15:17

HDU

查看标签

0 条评论

题目概述

HDU7218

解题报告

这道题看起来和城市规划很像，但是由于要考虑概率，所以卷积式子和求方案数的式子不一致。

定义 $f_i$ 表示 $i$ 个点成为一个连通块的概率，考虑两个连通块的概率：

$$ f_if_j(1-p)^{ij}{i+j\choose i}\to g_{i+j} $$

考虑指数型生成函数我们就可以把 ${i+j\choose j}$ 去掉，但是 $f_if_j(1-p)^{ij}\to g_{i+j}$ 我们还是没有办法进行卷积。

这时候利用Bluestein算法里面的一个套路：$ij={i+j\choose 2}-{i\choose 2}-{j\choose 2}$ ，就可以把 $ij$ 转化成卷积的形式：

$$ {f_{i}\over (1-p)^{i\choose 2}}{f_j\over (1-p)^{j\choose 2}}\to{g_{i+j}\over (1-p)^{i+j\choose 2}} $$

令 $f'_i$ 表示处理过之后的概率，$F(x)$ 是 $f'_i$ 的指数型生成函数，那么类似城市规划那道题，$F^i(x)\over i!$ 就是 $i$ 个连通块的概率。

令 $g_i$ 表示 $i$ 个点无向图的概率，显然 $g_i=1$ 。令 $g'_i$ 表示处理过之后的 $i$ 个点无向图的概率，$G(x)$ 表示 $g'_i$ 的指数型生成函数，则 $G(x)=e^{F(x)}$ 。

所以我们就知道了 $F(x)$ ，考虑枚举块数，统计答案：

$$ (1-p)^{n\choose 2}[{x^n\over n!}]\sum_{i=1}^{n}i\cdot {F^i(x)\over i!}\\ =(1-p)^{n\choose 2}[{x^n\over n!}]\sum_{i=1}^{n}\cdot {F^i(x)\over (i-1)!}\\ =(1-p)^{n\choose 2}[{x^n\over n!}]F(x)\sum_{i=0}^{n}\cdot {F^i(x)\over i!}\\ =(1-p)^{n\choose 2}[{x^n\over n!}]F(x)e^{F(x)} $$

多项式exp太慢了，为了避免就用 $G(x)$ 代替，那么答案就是 $(1-p)^{n\choose 2}[{x^n\over n!}]G(x)\ln G(x)$ 。

另一种做法是考虑容斥，这也是城市规划的另一种做法。

定义 $f_i$ 表示 $i$ 个点成为一个连通块的概率，$g_i$ 表示 $i$ 个点的概率（ $g_i=1$ ），那么：

$$ f_0=0,f_i=1-\sum_{j=0}^{i-1}{i-1\choose j-1}(1-p)^{j(i-j)}f_jg_{i-j} $$

即枚举第一个点所在连通块的大小，和上面同理，我们利用Bluestein套路，令 $f'_i={f_i\over (1-p)^{i\choose 2}(i-1)!}$ ，则：

$$ f'_i={1\over (1-p)^{i\choose 2}(i-1)!}-\sum_{j=0}^{i-1}f'_jg'_{i-j} $$

令 $F_1(x)$ ，$G(x)$ 为 $f',g'$ 对应的生成函数，$H(x)$ 为前面那项的生成函数，则：

$$ F_1(x)=H(x)-F_1(x)G(x)\\ F_1(x)={H(x)\over 1+G(x)} $$

考虑答案，由于期望可以拆分，我们考虑统计连通块大小为 $i$ 时的贡献，最后求和即可：

$$ ans_n=\sum_{i=0}^{n}{n\choose i}(1-p)^{i(n-i)}f_ig_{n-i} $$

这个式子也可以用套路处理，卷积一遍就可以求出所有的 $ans_n$ 。

示例程序

多项式ln

#include<cstdio>
#include<cctype>
#include<algorithm>
using namespace std;
typedef long long LL;
const int maxn=500000,maxt=1<<20,MOD=998244353;

int te,Q,A,B,P,a[maxn+5];
int pw[maxn+5],INV[maxn+5],fac[maxn+5];
int wn[maxt+5],tem[maxt+5],G[maxt+5],H[maxt+5];

#define EOLN(x) ((x)==10 || (x)==13 || (x)==EOF)
inline char readc(){
    static char buf[100000],*l=buf,*r=buf;
    return l==r && (r=(l=buf)+fread(buf,1,100000,stdin),l==r)?EOF:*l++;
}
template<typename T> int readi(T &x){
    T tot=0;char ch=readc(),lst='+';
    while (!isdigit(ch)) {if (ch==EOF) return EOF;lst=ch;ch=readc();}
    while (isdigit(ch)) tot=(tot<<3)+(tot<<1)+(ch^48),ch=readc();
    lst=='-'?x=-tot:x=tot;return EOLN(ch);
}
struct fastO{
    int si;char buf[100000];
    fastO() {si=0;}
    void putc(char ch){
        if (si==100000) fwrite(buf,1,si,stdout),si=0;
        buf[si++]=ch;
    }
    ~fastO() {fwrite(buf,1,si,stdout);}
}fo;
template<typename T> void writei(T x,char ch='\n'){
    int len=0,buf[100];
    if (x<0) fo.putc('-'),x=-x;
    do buf[len++]=x%10,x/=10; while (x);
    while (len) fo.putc(buf[--len]+48);
    fo.putc(ch);
}
inline int ADD(int x,int y) {return x+y>=MOD?x+y-MOD:x+y;}
inline int MUL(int x,int y) {return (LL)x*y%MOD;}
int Pow(int w,LL b) {int s;for (s=1;b;b>>=1,w=MUL(w,w)) if (b&1) s=MUL(s,w);return s;}
void NTTPre(){
    int x=Pow(3,(MOD-1)/maxt);
    wn[maxt/2]=1;
    for (int i=maxt/2+1;i<maxt;i++) wn[i]=MUL(wn[i-1],x);
    for (int i=maxt/2-1;i>=0;i--) wn[i]=wn[i<<1];
}
void Make(int n){
    INV[0]=INV[1]=1;for (int i=2;i<=n;i++) INV[i]=MUL(MOD-MOD/i,INV[MOD%i]);
    fac[0]=1;for (int i=1;i<=n;i++) fac[i]=MUL(fac[i-1],i),INV[i]=MUL(INV[i-1],INV[i]);
}
void NTT(int *a,int n,int f){
    if (f>0){
        for (int k=n>>1;k;k>>=1)
            for (int i=0;i<n;i+=k<<1)
                for (int j=0;j<k;j++){
                    int x=a[i+j],y=a[i+j+k];
                    a[i+j+k]=MUL(x+MOD-y,wn[k+j]);
                    a[i+j]=ADD(x,y);
                }
    } else {
        for (int k=1;k<n;k<<=1)
            for (int i=0;i<n;i+=k<<1)
                for (int j=0;j<k;j++){
                    int x=a[i+j],y=MUL(a[i+j+k],wn[k+j]);
                    a[i+j+k]=ADD(x,MOD-y);
                    a[i+j]=ADD(a[i+j],y);
                }
        for (int i=0,INV=MOD-(MOD-1)/n;i<n;i++) a[i]=MUL(a[i],INV);
        reverse(a+1,a+n);
    }
}
void Inte(int *F,int *a,int n,int K){ // F=integral(a)
    for (int i=n+K;i>=K;i--) F[i]=MUL(a[i-K],MUL(INV[i],fac[i-K]));
    for (int i=0;i<K;i++) F[i]=0;
}
void Deri(int *F,int *a,int n,int K){ // F=a'
    for (int i=0;i<=n-K;i++) F[i]=MUL(a[i+K],MUL(fac[i+K],INV[i]));
    for (int i=n-K+1;i<=n;i++) F[i]=0;
}
void Inv(int *F,int *a,int n){ // F=1/a
    if (n==1) {F[0]=Pow(a[0],MOD-2);return;}
    Inv(F,a,n>>1);
    for (int i=0;i<n;i++) tem[i]=a[i],tem[i+n]=F[i+n]=0;
    NTT(tem,n<<1,1);NTT(F,n<<1,1);
    for (int i=0;i<(n<<1);i++) tem[i]=MUL(F[i],2+MOD-MUL(tem[i],F[i]));
    NTT(tem,n<<1,-1);for (int i=0;i<n;i++) F[i]=tem[i],F[i+n]=0;
}
void Ln(int *F,int *a,int n){ // F=ln(a)
    static int A[maxt+5],B[maxt+5];
    Deri(A,a,n-1,1);for (int i=0;i<n;i++) A[i+n]=0;
    Inv(B,a,n);NTT(A,n<<1,1);NTT(B,n<<1,1);
    for (int i=0;i<(n<<1);i++) B[i]=MUL(A[i],B[i]);
    NTT(B,n<<1,-1);Inte(F,B,n-2,1);
}
int main(){
    Make(maxn);NTTPre();
    for (readi(te);te;te--){
        readi(Q);readi(A);readi(B);P=MUL(A,Pow(B,MOD-2));P=ADD(1,MOD-P);
        int inv=Pow(P,MOD-2),N=0;
        for (int i=1;i<=Q;i++) readi(a[i]),N=max(N,a[i]);
        int t;for (t=1;t<=N;t<<=1);
        for (int i=0;i<=N;i++) G[i]=MUL(INV[i],Pow(inv,(LL)i*(i-1)/2));
        for (int i=N+1;i<t;i++) G[i]=0;
        Ln(H,G,t);
        for (int i=t;i<(t<<1);i++) H[i]=G[i]=0;
        NTT(G,t<<1,1);NTT(H,t<<1,1);
        for (int i=0;i<(t<<1);i++) G[i]=MUL(G[i],H[i]);
        NTT(G,t<<1,-1);
        for (int i=1;i<=Q;i++){
            int val=MUL(fac[a[i]],Pow(P,(LL)a[i]*(a[i]-1)/2));
            writei(MUL(G[a[i]],val),i<Q?' ':'\n');
        }
    }
    return 0;
}

多项式求逆

#include<cstdio>
#include<cctype>
#include<algorithm>
using namespace std;
typedef long long LL;
const int maxn=500000,maxt=1<<20,MOD=998244353;

int te,Q,A,B,P,a[maxn+5];
int pw[maxn+5],INV[maxn+5],fac[maxn+5];
int wn[maxt+5],tem[maxt+5],G[maxt+5],H[maxt+5],F[maxt+5];

#define EOLN(x) ((x)==10 || (x)==13 || (x)==EOF)
inline char readc(){
    static char buf[100000],*l=buf,*r=buf;
    return l==r && (r=(l=buf)+fread(buf,1,100000,stdin),l==r)?EOF:*l++;
}
template<typename T> int readi(T &x){
    T tot=0;char ch=readc(),lst='+';
    while (!isdigit(ch)) {if (ch==EOF) return EOF;lst=ch;ch=readc();}
    while (isdigit(ch)) tot=(tot<<3)+(tot<<1)+(ch^48),ch=readc();
    lst=='-'?x=-tot:x=tot;return EOLN(ch);
}
struct fastO{
    int si;char buf[100000];
    fastO() {si=0;}
    void putc(char ch){
        if (si==100000) fwrite(buf,1,si,stdout),si=0;
        buf[si++]=ch;
    }
    ~fastO() {fwrite(buf,1,si,stdout);}
}fo;
template<typename T> void writei(T x,char ch='\n'){
    int len=0,buf[100];
    if (x<0) fo.putc('-'),x=-x;
    do buf[len++]=x%10,x/=10; while (x);
    while (len) fo.putc(buf[--len]+48);
    fo.putc(ch);
}
inline int ADD(int x,int y) {return x+y>=MOD?x+y-MOD:x+y;}
inline int MUL(int x,int y) {return (LL)x*y%MOD;}
int Pow(int w,LL b) {int s;for (s=1;b;b>>=1,w=MUL(w,w)) if (b&1) s=MUL(s,w);return s;}
void NTTPre(){
    int x=Pow(3,(MOD-1)/maxt);
    wn[maxt/2]=1;
    for (int i=maxt/2+1;i<maxt;i++) wn[i]=MUL(wn[i-1],x);
    for (int i=maxt/2-1;i>=0;i--) wn[i]=wn[i<<1];
}
void Make(int n){
    INV[0]=INV[1]=1;for (int i=2;i<=n;i++) INV[i]=MUL(MOD-MOD/i,INV[MOD%i]);
    fac[0]=1;for (int i=1;i<=n;i++) fac[i]=MUL(fac[i-1],i),INV[i]=MUL(INV[i-1],INV[i]);
}
void NTT(int *a,int n,int f){
    if (f>0){
        for (int k=n>>1;k;k>>=1)
            for (int i=0;i<n;i+=k<<1)
                for (int j=0;j<k;j++){
                    int x=a[i+j],y=a[i+j+k];
                    a[i+j+k]=MUL(x+MOD-y,wn[k+j]);
                    a[i+j]=ADD(x,y);
                }
    } else {
        for (int k=1;k<n;k<<=1)
            for (int i=0;i<n;i+=k<<1)
                for (int j=0;j<k;j++){
                    int x=a[i+j],y=MUL(a[i+j+k],wn[k+j]);
                    a[i+j+k]=ADD(x,MOD-y);
                    a[i+j]=ADD(a[i+j],y);
                }
        for (int i=0,INV=MOD-(MOD-1)/n;i<n;i++) a[i]=MUL(a[i],INV);
        reverse(a+1,a+n);
    }
}
void Inte(int *F,int *a,int n,int K){ // F=integral(a)
    for (int i=n+K;i>=K;i--) F[i]=MUL(a[i-K],MUL(INV[i],fac[i-K]));
    for (int i=0;i<K;i++) F[i]=0;
}
void Deri(int *F,int *a,int n,int K){ // F=a'
    for (int i=0;i<=n-K;i++) F[i]=MUL(a[i+K],MUL(fac[i+K],INV[i]));
    for (int i=n-K+1;i<=n;i++) F[i]=0;
}
void Inv(int *F,int *a,int n){ // F=1/a
    if (n==1) {F[0]=Pow(a[0],MOD-2);return;}
    Inv(F,a,n>>1);
    for (int i=0;i<n;i++) tem[i]=a[i],tem[i+n]=F[i+n]=0;
    NTT(tem,n<<1,1);NTT(F,n<<1,1);
    for (int i=0;i<(n<<1);i++) tem[i]=MUL(F[i],2+MOD-MUL(tem[i],F[i]));
    NTT(tem,n<<1,-1);for (int i=0;i<n;i++) F[i]=tem[i],F[i+n]=0;
}
int main(){
    Make(maxn);NTTPre();
    for (readi(te);te;te--){
        readi(Q);readi(A);readi(B);P=MUL(A,Pow(B,MOD-2));P=ADD(1,MOD-P);
        int inv=Pow(P,MOD-2),N=0;
        for (int i=1;i<=Q;i++) readi(a[i]),N=max(N,a[i]);
        int t;for (t=1;t<=N;t<<=1);
        for (int i=1;i<=N;i++) G[i]=MUL(INV[i],Pow(inv,(LL)i*(i-1)/2));
        for (int i=1;i<=N;i++) H[i]=MUL(INV[i-1],Pow(inv,(LL)i*(i-1)/2));
        G[0]=1;H[0]=0;for (int i=N+1;i<(t<<1);i++) G[i]=H[i]=0;
        Inv(F,G,t);
        NTT(F,t<<1,1);NTT(H,t<<1,1);
        for (int i=0;i<(t<<1);i++) F[i]=MUL(F[i],H[i]);
        NTT(F,t<<1,-1);
        G[0]=1;F[0]=0;
        for (int i=1;i<=N;i++) F[i]=MUL(F[i],MUL(fac[i-1],INV[i]));
        for (int i=N+1;i<(t<<1);i++) F[i]=0;
        NTT(F,t<<1,1);NTT(G,t<<1,1);
        for (int i=0;i<(t<<1);i++) F[i]=MUL(F[i],G[i]);
        NTT(F,t<<1,-1);
        for (int i=1;i<=Q;i++){
            int val=MUL(fac[a[i]],Pow(P,(LL)a[i]*(a[i]-1)/2));
            writei(MUL(F[a[i]],val),i<Q?' ':'\n');
        }
    }
    return 0;
}

本文链接：https://zigzagk.top/2022/08/12/HDU7218

请不要发毫无意义或内容不文明的评论。与本文无关评论请发留言板！

邮箱格式错误

网址请用http://或https://开头

支持Markdown和LaTeX数学公式

Emoji QwQ 贴吧 QQ Orz

😀 😁 😂 🤣 😃 😄 😅 😆 😉 😊 😋 😎 😍 😘 🥰 😗 😙 😚 🙂 🤗 🤩 🤔 🤨 😐 😑 😶 🙄 😏 😣 😥 😮 🤐 😯 😪 😫 😴 😌 😛 😜 😝 🤤 😒 😓 😔 😕 🙃 🤑 😲 ☹️ 🙁 😖 😞 😟 😤 😢 😭 😦 😧 😨 😩 😬 😰 😱 🥵 🥶 😳 🤪 😵 😡 😠 🤬 😷 🤒 🤕 🤢 🤮 🤧 😇 🤠 🥳 🥴 🥺 🤥 🤫 🤭 🧐 🤓 👻 🤲 👐 🙌 👏 🤝 👍 👎 👊 ✊ 🤛 🤜 🤞 ✌️ 🤟 🤘 👌 👈 👉 👆 👇 ✋ 🤚 🖐 🖖 👋 🤙 💪 🦵 🦶 🖕 ✍️ 🙏 👀 👴 🐴 🍋