POJ.3532 Resistance

1.Problem
http://poj.org/problem?id=3532

2.Idea
See this blog.

3.Source

 typedef vector<double> vec;  
 typedef vector<vec> mat;  
 double resistor[MAX_N][MAX_N];    
 int N, M;  
 vec gauss(const mat&A, const vec&b) {  
   int n = A.size();//!!!!  
   mat B(n, vec(n + 1));  
   for (int i = 0; i < n; i++)  
     for (int j = 0; j < n; j++)B[i][j] = A[i][j];  
   for (int i = 0; i < n; i++)B[i][n] = b[i];  
   for (int i = 0; i < n; i++) {  
     int pivot = i;  
     for (int j = i; j < n; j++) {  
       if (abs(B[j][i]) > abs(B[pivot][i])) {//!!!  
         pivot = j;  
       }  
     }  
     swap(B[i], B[pivot]);  
     if (abs(B[i][i]) < EPS)return vec(); 
     for (int j = i + 1; j <= n; j++) {  
       B[i][j] /= B[i][i];  
     }  
     for (int j = 0; j < n; j++) {  
       if (i != j) {  
         for (int k = i + 1; k <= n; k++) {  
           B[j][k] -= B[j][i] * B[i][k];  
         }  
       }  
     }  
   }  
   vec x(n);  
   for (int i = 0; i < n; i++) {  
     x[i] = B[i][n];  
   }  
   return x;  
 }  
 int main() {  
   while (scanf("%d%d", &N, &M) != EOF) {  
     memset(resistor, 0, sizeof(resistor));  
     for (int i = 0; i < M; i++) {  
       int from, to;  
       double R;  
       scanf("%d%d%lf", &from, &to, &R);  
       if (R == 0)continue;  
       from--, to--;  
       resistor[from][to] += 1 / R;  
       resistor[to][from] += 1 / R;  
     }  
     for (int i = 0; i < N; i++) {  
       for (int j = 0; j < N; j++) {  
         resistor[i][j] = 1.0 / resistor[i][j];  
       }  
     }  
     mat A(N, vec(N, 0));  
     vec b(N, 0);  
     b[0] = 1.0;  
     b[N - 1] = 0.0;  
     A[0][0] = 1, A[N - 1][N - 1] = 1;  
     for (int i = 1; i < N - 1; i++) {  
       for (int j = 0; j < N; j++) {  
         if (resistor[i][j] > 0) {  
           double I = 1.0 / resistor[i][j];  
           A[i][i] -= I;  
           A[i][j] += I;  
         }  
       }  
     }  
     vec voltage = gauss(A, b);  
     double current = 0;  
     for (int i = 0; i < N; i++) {  
       if (resistor[0][i] > 0) {  
         current += (voltage[0] - voltage[i]) / resistor[0][i];  
       }  
     }  
     printf("%.2f\n", 1.0 / current);  
   }  
   return 0;  
 }