标签：style blog http color os io for art
这是前几天多校的题目，高了好久突然听旁边的大神推出来说是可以用高斯消元，一直喊着赶快敲模板，对于从来没有接触过高斯消元的我来说根本就是一头雾水，无赖之下这几天做DP，正好又做到了这个题，没办法得从头开始看，后来在网上找了别人的高斯消元的模板后发现其实也还是很好理解，就是先构造一个增广矩阵，然后化行阶梯形，最后迭代求解
首先有一个介绍高斯消元感觉过于详细的博客http://blog.csdn.net/tsaid/article/details/7329301
首先看一下这个题怎么构造这个增广矩阵，我们把所有可以达到的分数组合作为一个点，再考虑它与其他点所连的边，例如：
(300, 200)<--(1-p)----(300, 300)-----(p)----->(350, 300)
a b c
我们就可以理解为有一条b--->a的权值为（1-p）的边，有一条b---->c的权值为p的边，那这样首先构造状态转移方程：
DP[b] = 1 + (1-p) * DP[a] + p * DP[c]
变形后得到:
DP[b] - p * DP[c] - (1-p) * DP[a] = 1;
这就是我们的增广矩阵的系数了，对于方程的一般形式Ax = B，可以理解为第b个方程中的变元的系数为A[b][b] = 1, A[b][c] = p, A[b][a] = (1-p),B[b] = 1
这样就构造出了一个（A，B）的一个增广矩阵，保存在a中
然后就是高斯消元化行阶梯行，看看代码很好理解，就只有两个操作r1<-->r2，交换两行，r2 = r2 - r1 * a (其中a为一个常系数)
第一次写Gauss 代码比较戳，可以根据网上详细的介绍结合代码看，很容易懂的
void gauss() { int col = 0; for(int k=0;k<cnt && col < cnt;k++, col ++) { double Max = fabs(a[k][col]); int Maxr = k; for(int r = k + 1; r < cnt; r ++) if( fabs(a[r][col]) - Max > eps ) Max = fabs( a[Maxr = r][col] ); if(fabs(Max) < eps) { k --; continue; } for(int c = col;c<=cnt; c ++ ) SWAP(a[Maxr][c], a[k][c] ); for(int r = k + 1; r < cnt; r ++ ) if( fabs(a[r][col]) > eps ) { double tmp = a[r][col] / a[k][col]; for(int c = col; c <= cnt; c ++ ) a[r][c] -= tmp * a[k][c]; } } }
化行阶梯行后就是迭代求解了，由于这里一定有解，所以不需要判断无解的情况，直接算就是了
1 for(int r=cnt-1;r>=0;r--) 2 { 3 for(int c = cnt-1;c>r;c--) 4 a[r][cnt] -= a[r][c] * ans[c]; 5 ans[r] = a[r][cnt] / a[r][r]; 6 }
最后的总复杂度就是O(n^3)n接近200
另外还有一个复杂度O(n)的解法的思路（N=20）http://www.cnblogs.com/gj-Acit/p/3888390.html
1 #include <map> 2 #include <set> 3 #include <stack> 4 #include <queue> 5 #include <cmath> 6 #include <ctime> 7 #include <vector> 8 #include <cstdio> 9 #include <cctype> 10 #include <cstring> 11 #include <cstdlib> 12 #include <iostream> 13 #include <algorithm> 14 using namespace std; 15 #define INF ((LL)100000000000000000) 16 #define inf (-((LL)1<<40)) 17 #define lson k<<1, L, mid 18 #define rson k<<1|1, mid+1, R 19 #define mem0(a) memset(a,0,sizeof(a)) 20 #define mem1(a) memset(a,-1,sizeof(a)) 21 #define mem(a, b) memset(a, b, sizeof(a)) 22 #define FOPENIN(IN) freopen(IN, "r", stdin) 23 #define FOPENOUT(OUT) freopen(OUT, "w", stdout) 24 template<class T> T CMP_MIN ( T a, T b ) { return a < b; } 25 template<class T> T CMP_MAX ( T a, T b ) { return a > b; } 26 template<class T> T MAX ( T a, T b ) { return a > b ? a : b; } 27 template<class T> T MIN ( T a, T b ) { return a < b ? a : b; } 28 template<class T> T GCD ( T a, T b ) { return b ? GCD ( b, a % b ) : a; } 29 template<class T> T LCM ( T a, T b ) { return a / GCD ( a, b ) * b; } 30 template<class T> void SWAP( T& a, T& b ) { T t = a; a = b; b = t; } 31 32 //typedef __int64 LL; 33 typedef long long LL; 34 const int MAXN = 255; 35 const int MAXM = 110000; 36 const double eps = 1e-12; 37 int dx[4] = { -1, 0, 0, 1}; 38 int dy[4] = {0, -1, 1, 0}; 39 40 double p; 41 int id[30][30], cnt, N = 20; 42 double G[400][400], a[300][300], ans[300]; 43 struct NODE 44 { 45 int u, d; 46 NODE(){} 47 NODE(int _u, int _d):u(_u), d(_d){} 48 }; 49 50 void preInit() 51 { 52 cnt = 0; 53 for(int i=0;i<=N-1;i++) 54 { 55 for(int j = 0; j <= i; j ++ ) 56 { 57 id[i][j] = cnt++; 58 } 59 } 60 id[N][N-1] = cnt ++; 61 } 62 63 void init() 64 { 65 preInit(); 66 mem0(G);mem0(a); 67 for(int i=0;i<=N-1;i++) 68 { 69 for(int j=0;j<=i;j++) 70 { 71 int nx = MAX(i, j + 1), ny = MIN(i, j + 1); 72 G[id[i][j]][id[nx][ny]] = p; 73 nx = i; ny = (j- 2) >= 0 ? j-2 : 0; 74 G[id[i][j]][id[nx][ny]] = 1-p; 75 } 76 } 77 for ( int i = 0; i < cnt; i++ ) 78 { 79 a[i][i] = a[i][cnt] = 1.0; 80 if ( i == cnt-1 ) { a[i][cnt] = 0; } 81 for ( int j = 0; j < cnt; j++ ) if ( fabs ( G[i][j] ) > eps ) 82 a[i][j] -= G[i][j]; 83 } 84 } 85 86 87 void gauss() 88 { 89 int col = 0; 90 for(int k=0;k<cnt && col < cnt;k++, col ++) 91 { 92 double Max = fabs(a[k][col]); 93 int Maxr = k; 94 for(int r = k + 1; r < cnt; r ++) 95 if( fabs(a[r][col]) - Max > eps ) 96 Max = fabs( a[Maxr = r][col] ); 97 if(fabs(Max) < eps) { k --; continue; } 98 for(int c = col;c<=cnt; c ++ ) 99 SWAP(a[Maxr][c], a[k][c] ); 100 for(int r = k + 1; r < cnt; r ++ ) if( fabs(a[r][col]) > eps ) 101 { 102 double tmp = a[r][col] / a[k][col]; 103 for(int c = col; c <= cnt; c ++ ) 104 a[r][c] -= tmp * a[k][c]; 105 } 106 } 107 for(int r=cnt-1;r>=0;r--) 108 { 109 for(int c = cnt-1;c>r;c--) 110 a[r][cnt] -= a[r][c] * ans[c]; 111 ans[r] = a[r][cnt] / a[r][r]; 112 } 113 } 114 115 int main() 116 { 117 // FOPENIN ( "in.txt" ); 118 // FOPENOUT("out.txt"); 119 while ( ~scanf ( "%lf", &p ) ) 120 { 121 init(); 122 gauss(); 123 printf("%.9lf\n", ans[0]); 124 } 125 return 0; 126 }
HDU 4870 Rating（高斯消元 ）,布布扣,bubuko.com
标签：style blog http color os io for art
原文地址：http://www.cnblogs.com/gj-Acit/p/3888382.html
