luogu 3790 文艺数学题 - 矩阵树定理 - 容斥原理

题目传送门

　　戳我来传送

题目大意

　　给定一个图，问它的所有生成树的边权的最大公约数之和。

　　可以考虑计算边权的最大公约数为$i$的生成树的个数$f(i)$，最后累加答案。

　　然后考虑这样的生成树的个数怎么求，根据某个经典套路，我们可以容斥。

　　因为可以求出边权的最大公约数为$i$的倍数的生成树的个数$F(i)$，所以减去它的倍数的$f$就是$f(i)$了。

　　但是这么做会T掉。

　　可以用$O(W\log W)$的时间内预处理出为边权$i$的倍数的边数有多少条。然后高消前判断一下边数是否大于等于$n - 1$。

　　复杂度的话考虑超过 $\sqrt{W}$ 的边数至多有 $O(m)$ 条，因此复杂度为 $O((\sqrt{W}+\frac{m}{n - 1})n^3)$

Code

  1 /**
  2  * luogu
  3  * Problem#3790
  4  * Accepted
  5  * Time: 6016ms
  6  * Memory: 9402k
  7  */
  8 #include <iostream>
  9 #include <cstring>
 10 #include <cstdio>
 11 using namespace std;
 12 typedef bool boolean;
 13 #define ll long long
 14 
 15 const int N = 65, M = 1e9 + 7;
 16 
 17 void exgcd(int a, int b, int& d, int& x, int& y) {
 18     if(!b)
 19         d = a, x = 1, y = 0;
 20     else {
 21         exgcd(b, a % b, d, y, x);
 22         y -= (a / b) * x;        
 23     }
 24 }
 25 
 26 int inv(int a, int n) {
 27     int d, x, y;
 28     exgcd(a, n, d, x, y);
 29     return (x < 0) ? (x + n) : (x);
 30 }
 31 
 32 typedef class Matrix {
 33     public:
 34         int a[N][N];
 35 
 36         void reset() {
 37             memset(a, 0, sizeof(a));
 38         }
 39 
 40         int det(int n) {
 41             int pro = 1;
 42             for (int i = 1, e; i <= n; i++) {
 43                 e = 0;
 44                 for (int j = i; j <= n && !e; j++)
 45                     if (a[j][i])
 46                         e = j;
 47                 if (e == 0)    return 0;
 48                 if (e != i)    
 49                     for (int j = 1; j <= n; j++)
 50                         swap(a[e][j], a[i][j]);
 51                 for (int j = 1, x, y; j <= n; j++) {
 52                     if (j == i)    continue;
 53                     x = a[i][i], y = a[j][i], pro = pro * 1ll * x % M;
 54                     for (int k = 1; k <= n; k++) {
 55                         a[j][k] = (a[j][k] * 1ll * x - a[i][k] * 1ll * y) % M;
 56                         if (a[j][k] < 0)    a[j][k] += M;
 57                     }
 58 
 59                 }
 60             }
 61             int rt = inv(pro, M);
 62             for (int i = 1; i <= n; i++)
 63                 rt = rt * 1ll * a[i][i] % M;
 64             return rt;
 65         }
 66 }Matrix;
 67 
 68 typedef class Edge {
 69     public:
 70         int u, v, w;
 71 }Edge;
 72 
 73 const int Val = 1e6 + 1;
 74 
 75 int n, m;
 76 Edge* es;
 77 Matrix a;
 78 int ce[Val], f[Val];
 79 
 80 inline void init() {
 81     scanf("%d%d", &n, &m);
 82     es = new Edge[(m + 1)];
 83     for (int i = 1; i <= m; i++)
 84         scanf("%d%d%d", &es[i].u, &es[i].v, &es[i].w), ce[es[i].w]++;
 85 }
 86 
 87 int calc(int g) {
 88     if (ce[g] < n - 1)    return 0;
 89     a.reset();
 90     for (int i = 1, u, v; i <= m; i++)
 91         if (!(es[i].w % g)) {
 92             u = es[i].u - 1, v = es[i].v - 1;
 93             a.a[u][u]++, a.a[v][v]++;
 94             a.a[u][v]--, a.a[v][u]--;
 95             if (a.a[u][v] < 0)    a.a[u][v] += M;
 96             if (a.a[v][u] < 0)    a.a[v][u] += M;
 97         }
 98     int rt = a.det(n - 1);
 99 //    cerr << rt << endl;
100     return rt;
101 }
102 
103 int res = 0;
104 inline void solve() {
105     for (int i = 1; i < Val; i++)
106         for (int j = (i << 1); j < Val; j += i)
107             ce[i] += ce[j];
108     for (int i = Val - 1; i; i--) {
109         f[i] = calc(i);
110         for (int j = (i << 1); j < Val; j += i) {
111             f[i] -= f[j];
112             if (f[i] < 0)
113                 f[i] += M;
114         }
115         res = (res + f[i] * 1ll * i) % M;
116     }
117     printf("%d", res);
118 }
119 
120 int main() {
121     init();
122     solve();
123     return 0;
124 }

 

posted @ 2018-03-02 11:55  阿波罗2003  阅读(217)  评论(0)  编辑  收藏  举报