[Medium-Hard] P6624 [省选联考 2020 A 卷] 作业题

xiezheyuan2024/12/22约 1286 字大约 6 分钟

whk 归来的第一篇题解！

简要题意

定义一个 $n$ 个点的树的权值，为其各边的边权和乘上各边的边权的 $\gcd$ 。

给定一个 $n$ 个点 $m$ 条边的无向图，边有边权。你需要计算其所有生成树的权值总和。答案对 $p=998,244,353$ 取模。

$1\leq n\leq 30,1\leq m\leq \frac{1}{2}n(n-1)$ ，值域 $[1,152501]$ 。

思路

首先看到 $\gcd$ ，考虑利用欧拉反演：

\begin{aligned} &\left(\sum_{i=1}^{n-1}w_i\right)\cdot \gcd(w_1,\cdots,w_{n-1})\\ &=\left(\sum_{i=1}^{n-1}w_i\right)\sum_{d\mid w_1,\cdots,d\mid w_{n-1}} \varphi(d)\\ &=\sum_{d}\varphi(d) [d\mid w_1\land\cdots\land d\mid w_{n-1}]\sum_{i=1}^{n-1}w_i \end{aligned}

预处理 $\varphi$ ，枚举 $d$ ，将所有边权为 $d$ 的倍数的边组成一个新图 $G_d$ ，现在的问题就变成了求 $G_d$ 的所有生成树的边权和的和。可以想到 Matrix-Tree 定理，不过美中不足的是，Matrix-Tree 定理计算的是所有生成树的边权积的和，而不是所有生成树的边权和的和。这就很难办。

这个时候你打开了题解区，看到了一种非常神秘的方法：将边权从 $\mathbb{Z_{p}}$ 换成一元一次多项式 $\{ax+b\mid a,b\in \mathbb{Z_{p}}\}$ 。

这样做有什么好处呢？根据多项式的基本性质，我们有：

$(ax+b)\pm (cx+d)=(a\pm c)x+(b\pm d)$ 。
$(ax+b)(cx+d)\equiv (ad+bc)x+bd\pmod{x^2}$ 。

注意所谓“模 $x^2$ 意义下”的实际意义是截断 $x^2$ 及更高次项。

你发现如果这样定义，乘积的结果的一次项系数会带有加法，不难想到直接令 $b,d=1$ ，这样就可以实现类似 $a+c$ 的效果了。

现在还有一个遗留问题，如果用朴素的初等行变换化为上三角矩阵的方法计算行列式，必须使用到除法，那么 $(ax+b)\div (cx+d)$ 是什么呢？考虑这个问题，不妨考虑 $(ax+b)^{-1}$ 是什么。

令 $(ax+b)^{-1}\equiv (cx+d)\pmod{x^2}$ ，则有 $(ax+b)(cx+d)\equiv (ad+bc)x+bd\equiv 1\pmod{x^2}$ ，得到 $ad+bc=0,bd=1$ 。

故 $d=b^{-1},c=-ab^{-2}$ 。得到逆元的结果后，计算除法只需要简单地乘上逆元即可。

按照 Matrix-Tree 定理的方法建立 Laplacian 矩阵（度数矩阵和邻接矩阵都是一元多项式，且常数项为 $1$ ），然后正常初等行变换即可，最后得到的行列式的一次项系数就是答案。

时间复杂度 $O(n^3V\log p)$ ，精细的实现可以做到 $O(n^3V)$ ，其中 $V$ 为值域。实际上这只是一个非常宽松的上界（我们可以将 $G_d$ 没有边的情况进行特判），因而可以通过本题。

代码

#include <bits/stdc++.h>
using namespace std;

const int N = 35, M = ((30 * 29) >> 1) + 5, V = 152501 + 5;

constexpr int mod = 998244353;
int Add(int x, int y){ return (x + y) >= mod ? (x + y - mod) : (x + y); }
int Sub(int x, int y){ return (x - y) < 0 ? (x - y + mod) : (x - y); }
int Mul(int x, int y){ return 1ll * x * y % mod; }

int fastpow(int x, int y){
    int ret = 1;
    for(;y;y>>=1,x=Mul(x, x)){ if(y & 1) ret = Mul(ret, x); }
    return ret;
}

int n, m;

struct edge{
    int u, v, w;
} e[M];

int pri[V], phi[V], tot;
bool vis[V];

void sieve(int n){
    phi[1] = 1;
    for(int i=2;i<=n;i++){
        if(!vis[i]) pri[++tot] = i, phi[i] = i - 1;
        for(int j=1;j<=tot&&1ll*i*pri[j]<=n;j++){
            vis[i * pri[j]] = 1, phi[i * pri[j]] = Mul(phi[i], phi[pri[j]]);
            if(!(i % pri[j])){
                phi[i * pri[j]] = Mul(phi[i], pri[j]);
                break;
            }
        }
    }
}

struct unary{
    int a, b; // ax + b
    unary(int _a = 0, int _b = 0) : a(_a), b(_b) {}
    unary operator+(const unary& rhs) const { return {Add(a, rhs.a), Add(b, rhs.b)}; }
    unary operator-(const unary& rhs) const { return {Sub(a, rhs.a), Sub(b, rhs.b)}; }
    unary operator*(const unary& rhs) const { return {Add(Mul(a, rhs.b), Mul(b, rhs.a)), Mul(b, rhs.b)}; }
    unary inv() const {
        int invb = fastpow(b, mod - 2);
        return {Mul(mod - a, Mul(invb, invb)), invb};
    }
    unary operator/(const unary& rhs) const { return (*this) * rhs.inv(); }
    unary operator+=(const unary& rhs){ return *this = *this + rhs; }
};

unary deg[N][N], a[N][N];

unary det(int n){
    for(int i=1;i<=n;i++){
        for(int j=i+1;j<=n;j++){
            unary v = a[j][i] / a[i][i];
            for(int k=i;k<=n;k++) a[j][k] = a[j][k] - v * a[i][k];
        }
    }
    unary ret = a[1][1];
    for(int i=2;i<=n;i++) ret = ret * a[i][i];
    return ret;
}

signed main(){
    ios::sync_with_stdio(false);
    cin.tie(0); cout.tie(0);
    cin >> n >> m;
    int mx = 152501; sieve(mx);
    for(int i=1;i<=m;i++) cin >> e[i].u >> e[i].v >> e[i].w;
    int ans = 0;
    for(int d=1;d<=mx;d++){
        for(int i=1;i<=n;i++) deg[i][i] = {0, 0};
        for(int i=1;i<=n;i++){
            for(int j=1;j<=n;j++) a[i][j] = {0, 0};
        }
        int tot = 0;
        for(int i=1;i<=m;i++){
            if(e[i].w % d) continue;
            tot++;
            deg[e[i].u][e[i].u] += {e[i].w, 1}, deg[e[i].v][e[i].v] += {e[i].w, 1};
            a[e[i].u][e[i].v] += {e[i].w, 1}, a[e[i].v][e[i].u] += {e[i].w, 1};
        }
        if(!tot) continue;
        for(int i=1;i<=n;i++){
            for(int j=1;j<=n;j++) a[i][j] = deg[i][j] - a[i][j];
        }
        ans = Add(ans, Mul(det(n - 1).a, phi[d]));
    }
    cout << ans << '\n';
    return 0;
}

// Written by xiezheyuan