1. 程式人生 > >【bzoj 3534】 [SDOI2014] 重建

【bzoj 3534】 [SDOI2014] 重建

題意:
  給一個圖,每條邊有出現概率,求這個圖恰好為一棵樹的概率。
 
解法:
  考慮Kirchhoff矩陣的意義:K[G]=D[G]A[G]=B[G]B[G]T,之所以能夠進行生成樹計數是對於其伴隨矩陣在計數n1條邊的集合時,當n1條邊中存在環就會產生線性組合而導致行列式為零,否則恰好對角線上均為伴隨矩陣中所賦的值,使得det(Bi,j)2就為1。
  因此我們可以令鄰接矩陣中存下邊權,使得可以用Kirchhoff矩陣計數生成樹中邊權的乘積。
  考慮到答案為Tu,v(p(u,v)[(u,v)T])((1p(u,v))[u,vT]),我們可以令邊權為p

(u,v)/(1p(u,v)),求出行列式然後再乘回分母。
  但是注意到p(u,v)==1時分母為0,此時我們理論上可以用一個足夠小的ϵ來表示0,考慮到double的表示法,取ϵ=220就足夠了,這樣可以避免過大/過小數運算所導致的精度誤差,而原先的概率的精度(兩位小數)也是可以接受的。
  時間複雜度O(n3),不知道為什麼出題人出這麼小,可能是怕資料都是0然後被水過去吧。

/*
    I will chase the meteor for you, a thousand times over.
    Please wait for me, until I fade forever.
    Just 'coz GEOTCBRL.
*/
#include <bits/stdc++.h> using namespace std; #define fore(i,u) for (int i = head[u] ; i ; i = nxt[i]) #define rep(i,a,b) for (int i = a , _ = b ; i <= _ ; i ++) #define per(i,a,b) for (int i = a , _ = b ; i >= _ ; i --) #define For(i,a,b) for (int i = a , _ = b ; i < _ ; i ++) #define Dwn(i,a,b) for (int i = ((int) a) - 1 , _ = (b) ; i >= _ ; i --)
#define fors(s0,s) for (int s0 = (s) , _S = s ; s0 ; s0 = (s0 - 1) & _S) #define foreach(it,s) for (__typeof(s.begin()) it = s.begin(); it != s.end(); it ++) #define mp make_pair #define pb push_back #define pii pair<int , int> #define fir first #define sec second #define MS(x,a) memset(x , (a) , sizeof (x)) #define gprintf(...) fprintf(stderr , __VA_ARGS__) #define gout(x) std::cerr << #x << "=" << x << std::endl #define gout1(a,i) std::cerr << #a << '[' << i << "]=" << a[i] << std::endl #define gout2(a,i,j) std::cerr << #a << '[' << i << "][" << j << "]=" << a[i][j] << std::endl #define garr(a,l,r,tp) rep (__it , l , r) gprintf(tp"%c" , a[__it] , " \n"[__it == _]) template <class T> inline void upmax(T&a , T b) { if (a < b) a = b ; } template <class T> inline void upmin(T&a , T b) { if (a > b) a = b ; } typedef long long ll; const int maxn = 57; const int maxm = 200007; const int mod = 1000000007; const int inf = 0x7fffffff; const double eps = pow(2 , -19); typedef int arr[maxn]; typedef int adj[maxm]; inline int fcmp(double a , double b) { if (fabs(a - b) <= eps) return 0; if (a < b - eps) return -1; return 1; } inline int add(int a , int b) { return ((ll) a + b) % mod ; } inline int mul(int a , int b) { return ((ll) a * b) % mod ; } inline int dec(int a , int b) { return add(a , mod - b % mod) ; } inline int Pow(int a , int b) { int t = 1; while (b) { if (b & 1) t = mul(t , a); a = mul(a , a) , b >>= 1; } return t; } #define rd RD<int> #define rdll RD<long long> template <typename Type> inline Type RD() { Type x = 0; int flag = 0; char c = getchar(); while (!isdigit(c) && c != '-') c = getchar(); (c == '-') ? (flag = 1) : (x = c - '0'); while (isdigit(c = getchar())) x = x * 10 + c - '0'; return flag ? -x : x; } inline char rdch() { char c = getchar(); while (!isalpha(c)) c = getchar(); return c; } // beginning int n; double to_mul , a[maxn][maxn]; void input() { n = rd(); to_mul = 1; rep (i , 1 , n) rep (j , 1 , n) { scanf("%lf" , &a[i][j]); double b = 1 - a[i][j]; upmax(b , eps); if (i < j) to_mul *= b; if (i != j) a[i][j] = - a[i][j]; a[i][j] /= b; } rep (i , 1 , n) rep (j , i + 1 , n) a[i][i] -= a[i][j] , a[j][j] -= a[i][j]; } inline double det() { double ret = 1; For (i , 1 , n) { int j = i; while (j < n && fabs(a[j][i]) < eps) j ++; if (j == n) return 0; if (i != j) { ret *= -1; For (k , 1 , n) swap(a[i][k] , a[j][k]); } For (j , i + 1 , n) if (fabs(a[j][i]) >= eps) { double t = a[j][i] / a[i][i]; For (k , i , n) a[j][k] -= t * a[i][k]; } } For (i , 1 , n) ret *= a[i][i]; return fabs(ret); } void solve() { double ans = det(); ans *= to_mul; printf("%.10lf\n" , ans); } int main() { #ifndef ONLINE_JUDGE freopen("data.txt" , "r" , stdin); #endif input(); solve(); return 0; }