1. 程式人生 > 其它 >「ZJOI2015」地震後的幻想鄉

「ZJOI2015」地震後的幻想鄉

題目

點這裡看題目。

分析

首先,設原圖的最小生成樹的邊集為 \(T\),則容易得到:

\[\begin{aligned} E(\max_{x\in T}e_x) &=\int_{0}^1P(t<\max_{x\in T}e_x)\mathrm dt \end{aligned} \]

而可以發現 \(P(t<\max_{x\in T}e_x)\) 恰好為加入了 \(\le t\) 的邊之後,圖仍然不連通的概率。因此我們可以直接列舉一個邊集 \(E'\subseteq E\),使得加入了 \(E'\) 後圖不連通,期望則可以被表述為:

\[\begin{aligned} \int_{0}^1P(t<\max_{x\in T}e_x)\mathrm dt &=\sum_{E'\subseteq E}\int_{0}^1t^{|E'|}(1-t)^{m-|E'|}\mathrm d t\\ &=\sum_{E'\subseteq E}\int_{0}^1\sum_{j=0}^{m-|E'|}(-1)^j\binom{m-|E'|}{j}t^{j+|E'|}\mathrm d t\\ &=\sum_{j=0}^m\sum_{k=0}^{m-j}(-1)^j\binom{m-k}{j}\left(\int_0^1t^{j+k}\mathrm d t\right)\left(\sum_{E'\subseteq E}[|E'|=k]\right)\\ &=\sum_{j=0}^m\sum_{k=0}^{m-j}(-1)^j\binom{m-k}{j}\frac{1}{j+k+1}\left(\sum_{E'\subseteq E}[|E'|=k]\right) \end{aligned} \]

所以,我們嘗試計算出 \(\sum_{E'\subseteq E}[|E'|=k]\)

,整個問題也便迎刃而解了。注意 \(E'\) 還需要保證不連通,這其實頗有一點麻煩——我們反之保證連通。設 \(f_{S,k}\) 表示在 \(S\subseteq V\)匯出子圖中,能夠保證 \(S\) 內點連通且大小恰好為 \(k\) 的邊集的數量。設 \(S\) 的匯出子圖的邊集為 \(E(S)\)。正難則反,我們可以想到容斥不連通的情況:

\[f_{S,k}=\binom{|E(S)|}{k}-\sum_{T\subseteq S,T\neq\varnothing}[\min T=\min S]\sum_{j=0}^kf_{T,j}\binom{|E(S\setminus T)|}{k-j} \]

如果直接暴力 DP 是 \(O(3^nm^2)\)

的。不過注意到,第二維可以看作是生成函式的指數,我們相當於在進行生成函式運算。而我們又知道,最終 \(f_V\) 的生成函式的指數一定是 \(\le m\) 的,因此我們可以預先插入 \(m+1\) 個點值,最後再做拉格朗日插值得到每一項係數

最終我們得到了一個 \(O(3^nm+m^2)\) 的演算法。


實現還需要一點技巧:

為了方便,由於 \(\binom{45}{22}\approx4.1\times 10^{12}\) 並不是很大,long long 也還裝得下,我們可以選一個略大於 \(\binom{45}{22}\) 的質數,在它的模域下做運算,這樣就簡單了許多。

最後,由於我們需要輸出浮點數結果,運算精度也是我們需要考慮的問題。注意到,答案最終一定是在 \((0,1)\)

範圍內,因此我們可以輸出它 \(\bmod 1\) 之後的結果。這樣的話,像 \(\binom{m-k}{j}\frac{1}{j+k+1}\bmod 1\) 這樣的運算,就可以轉化為 \(\frac{1}{j+k+1}\left(\binom{m-k}{j}\bmod (j+k+1)\right)\)。我們可以直接將分子對分母取模,最後算的時候再計算 \(ans \bmod 1\) 即可。

如果不是不調整精度過不了完全圖,我會卡嗎我?

小結:

  1. 仍然需要注意常見的運算技巧。這裡主要的處理期望的技巧,還是把貢獻滾到字首上,從而轉化為概率

  2. 注意一下圖上的子集 DP 的方法。一般來說,我們會傾向於通過容斥計算,在無向圖的連通性和有向圖的強連通中都可以這麼做;另一方面,強連通圖中也有用耳分解的演算法。

  3. 注意實現中的小技巧。尤其是處理精度這一塊的,平時接觸不多更要注意。

    一般來說,我們會把浮點數壓縮到一個較小的範圍內,從而保證精度。這裡進行 \(\bmod 1\) 及後續運算就是為了達成這一點。

    你也可以說,我是發現所有整數部分都被抵消了才想到了這個方法

程式碼

#include <cmath>
#include <cstdio>

#define rep( i, a, b ) for( int i = (a) ; i <= (b) ; i ++ )
#define per( i, a, b ) for( int i = (a) ; i >= (b) ; i -- )

typedef long long LL;

const LL mod = 4116715363889;
const int MAXN = 105, MAXS = ( 1 << 10 ) + 5;

template<typename _T>
void read( _T &x ) {
    x = 0; char s = getchar(); bool f = false;
    while( s < '0' || '9' < s ) { f = s == '-', s = getchar(); }
    while( '0' <= s && s <= '9' ) { x = ( x << 3 ) + ( x << 1 ) + ( s - '0' ), s = getchar(); }
    if( f ) x = -x;
}

template<typename _T>
void write( _T x ) {
    if( x < 0 ) putchar( '-' ), x = -x;
    if( 9 < x ) write( x / 10 );
    putchar( x % 10 + '0' );
}

int fin[MAXN];

LL f[MAXS], g[MAXS];
int cnt[MAXS];

LL C[MAXN][MAXN];

LL ways[MAXN], coe[MAXN];

int N, M;

inline LL Mul( const LL a, const LL b ) {
    return ( a * b - ( LL ) ( ( long double ) a / mod * b ) * mod + mod ) % mod;
}

inline LL Qkpow( LL, LL );
inline LL Inv( const LL a ) { return Qkpow( a, mod - 2 ); }
inline LL Sub( LL x, const LL v ) { return ( x -= v ) < 0 ? x + mod : x; }
inline LL Add( LL x, const LL v ) { return ( x += v ) >= mod ? x - mod : x; }

inline LL Qkpow( LL base, LL indx ) {
    LL ret = 1;
    while( indx ) {
        if( indx & 1 ) ret = Mul( ret, base );
        base = Mul( base, base ), indx >>= 1;
    }
    return ret;
}

LL Query( const int x ) {
    for( int S = 0 ; S < ( 1 << N ) ; S ++ ) {
        f[S] = 0;
        per( i, M, 0 ) f[S] = Add( Mul( f[S], x ), C[cnt[S]][i] );
        g[S] = f[S];
        for( int T = ( S - 1 ) & S ; T ; T = ( T - 1 ) & S )
            if( ( T & ( - T ) ) == ( S & ( - S ) ) )
                g[S] = Sub( g[S], Mul( g[T], f[S ^ T] ) );
    }
    return g[( 1 << N ) - 1];
}

inline long double Mod1( const long double x ) { 
    return x - floorl( x );
}

int main() {
    read( N ), read( M );
    rep( i, 1, M ) {
        int u, v; read( u ), read( v );
        cnt[( 1 << ( u - 1 ) ) | ( 1 << ( v - 1 ) )] ++;
    }
    rep( i, 0, M ) {
        C[i][0] = C[i][i] = 1;
        rep( j, 1, i - 1 ) C[i][j] = Add( C[i - 1][j], C[i - 1][j - 1] );
    }
    for( int i = 0 ; i < N ; i ++ )
        for( int S = 0 ; S < ( 1 << N ) ; S ++ )
            if( ! ( S >> i & 1 ) ) cnt[S | ( 1 << i )] += cnt[S];
    coe[0] = 1;
    rep( i, 1, M + 1 )
        per( j, M + 1, 0 ) {
            coe[j] = Mul( coe[j], mod - i );
            if( j ) coe[j] = Add( coe[j], coe[j - 1] );
        }
    rep( i, 1, M + 1 ) {
        LL inv = Inv( i ), tmp = 1;
        rep( j, 0, M + 1 ) {
            if( j ) coe[j] = Sub( coe[j], coe[j - 1] );
            coe[j] = Mul( coe[j], mod - inv );
        }
        rep( j, 1, M + 1 ) if( i ^ j )
            tmp = Mul( tmp, Sub( i, j ) );
        tmp = Mul( Query( i ), Inv( tmp ) );
        rep( j, 0, M ) ways[j] = Add( ways[j], Mul( tmp, coe[j] ) );
        per( j, M + 1, 0 ) {
            coe[j] = Mul( coe[j], mod - i );
            if( j ) coe[j] = Add( coe[j], coe[j - 1] );
        }
    }
    long double ans = 0;
    rep( i, 0, M ) 
        rep( j, 0, M - i ) {
            int cur = i + j + 1;
            int res = 1ll * ( C[M - i][j] % cur ) * ( Sub( C[M][i], ways[i] ) % cur ) % cur;
            fin[cur] = j & 1 ? ( ( fin[cur] - res ) % cur + cur ) % cur : ( fin[cur] + res ) % cur;
        }
    rep( i, 1, M + 1 ) ans = Mod1( ans + ( long double ) fin[i] / i );
    printf( "%.6Lf\n", ans );
    return 0;
}