斯坦纳树(Steiner Tree)问题是图论中经典的NP-hard组合优化问题:给定无向带权图和终端节点集合,寻找连接所有终端节点的最小权重子图。与最小生成树(MST)不同,斯坦纳树允许引入额外的Steiner点来降低总代价。该问题在VLSI布线、网络设计、生物信息学等领域有广泛应用。本文用Java实现斯坦纳树的精确DP算法与启发式近似解法。
核心原理
问题定义
给定无向图G=(V, E, w)和终端节点集合T⊆V,斯坦纳树问题要求找到一棵包含T中所有节点的最小权重树。当|T|=2时退化为最短路径;当T=V时退化为最小生成树。
关键性质
- Steiner点度数:最优解中每个Steiner点的度数至少为3。
- 边数上界:n个终端的斯坦纳树最多包含2n-2个节点和2n-3条边。
- NP-hard性:一般图上的斯坦纳树问题是NP-hard,但小规模实例可用动态规划精确求解。
Dreyfus-Wagner算法
对于终端数k较小的场景,使用状态压缩动态规划:
– 状态dp[mask][v]表示以mask描述的终端集合连通、且树以v为根的最小代价
– 转移:将mask拆分为两个非空子集sub和mask\sub,dp[mask][v] = min(dp[sub][v] + dp[mask\sub][v])
– 对每个状态运行Dijkstra在图上松弛边权
近似算法(Kou-Markowsky-Berman)
- 计算所有终端点对之间的最短路径
- 构造完全图( metric closure),边权为最短路径距离
- 求该完全图的最小生成树
- 将MST边展开为原图中的最短路径
- 在展开子图上求MST,即为2-近似解
完整Java实现
import java.util.*;
/**
* 斯坦纳树(Steiner Tree)实现
* 包含:Dreyfus-Wagner精确DP、KMB 2-近似算法、最短路径计算
*/
public class SteinerTree {
private final int n; // 节点数
private final List<Edge>[] graph; // 邻接表
private final int INF = Integer.MAX_VALUE / 4;
public SteinerTree(int n) {
this.n = n;
this.graph = new ArrayList[n];
for (int i = 0; i < n; i++) graph[i] = new ArrayList<>();
}
/**
* 添加无向边
*/
public void addEdge(int u, int v, int w) {
graph[u].add(new Edge(v, w));
graph[v].add(new Edge(u, w));
}
/**
* 边类
*/
static class Edge {
int to, weight;
Edge(int to, int weight) { this.to = to; this.weight = weight; }
}
/**
* Dijkstra单源最短路径
*/
public int[] dijkstra(int src) {
int[] dist = new int[n];
Arrays.fill(dist, INF);
dist[src] = 0;
PriorityQueue<int[]> pq = new PriorityQueue<>(Comparator.comparingInt(a -> a[1]));
pq.offer(new int[]{src, 0});
while (!pq.isEmpty()) {
int[] cur = pq.poll();
int u = cur[0], d = cur[1];
if (d > dist[u]) continue;
for (Edge e : graph[u]) {
int v = e.to;
if (dist[u] + e.weight < dist[v]) {
dist[v] = dist[u] + e.weight;
pq.offer(new int[]{v, dist[v]});
}
}
}
return dist;
}
/**
* Floyd-Warshall全源最短路径
*/
public int[][] allPairsShortestPath() {
int[][] dist = new int[n][n];
for (int i = 0; i < n; i++) {
Arrays.fill(dist[i], INF);
dist[i][i] = 0;
}
for (int u = 0; u < n; u++) {
for (Edge e : graph[u]) {
dist[u][e.to] = Math.min(dist[u][e.to], e.weight);
}
}
for (int k = 0; k < n; k++) {
for (int i = 0; i < n; i++) {
for (int j = 0; j < n; j++) {
if (dist[i][k] + dist[k][j] < dist[i][j]) {
dist[i][j] = dist[i][k] + dist[k][j];
}
}
}
}
return dist;
}
/**
* Dreyfus-Wagner精确算法:状态压缩DP求解斯坦纳树
* 时间复杂度: O(3^k * n + 2^k * n * log n)
* 适用于终端数k <= 10-12的情况
*/
public Result dreyfusWagner(int[] terminals) {
int k = terminals.length;
int fullMask = (1 << k) - 1;
// dp[mask][v]: 连接mask中所有终端、且以v为根的最小代价
int[][] dp = new int[1 << k][n];
for (int[] row : dp) Arrays.fill(row, INF);
// 初始化:单个终端到自身的代价为0
for (int i = 0; i < k; i++) {
dp[1 << i][terminals[i]] = 0;
}
// 按mask大小递增处理
for (int mask = 1; mask <= fullMask; mask++) {
// 子集划分转移:mask = sub | (mask ^ sub)
for (int sub = (mask - 1) & mask; sub > 0; sub = (sub - 1) & mask) {
int other = mask ^ sub;
if (other == 0) continue;
// 确保不重复划分(sub < other)
if (sub > other) continue;
for (int v = 0; v < n; v++) {
if (dp[sub][v] + dp[other][v] < dp[mask][v]) {
dp[mask][v] = dp[sub][v] + dp[other][v];
}
}
}
// 对当前mask运行Dijkstra-like松弛(SPFA风格)
// 找到所有dp[mask][v] < INF的节点,在图上传播
boolean[] inQueue = new boolean[n];
Queue<Integer> queue = new LinkedList<>();
for (int v = 0; v < n; v++) {
if (dp[mask][v] < INF) {
queue.offer(v);
inQueue[v] = true;
}
}
while (!queue.isEmpty()) {
int u = queue.poll();
inQueue[u] = false;
for (Edge e : graph[u]) {
int v = e.to;
if (dp[mask][u] + e.weight < dp[mask][v]) {
dp[mask][v] = dp[mask][u] + e.weight;
if (!inQueue[v]) {
queue.offer(v);
inQueue[v] = true;
}
}
}
}
}
// 找到最优根节点
int minCost = INF, bestRoot = -1;
for (int v = 0; v < n; v++) {
if (dp[fullMask][v] < minCost) {
minCost = dp[fullMask][v];
bestRoot = v;
}
}
// 回溯重建树边
Set<int[]> treeEdges = new HashSet<>();
if (bestRoot != -1) {
reconstructTree(dp, fullMask, bestRoot, terminals, treeEdges, new boolean[1 << k][n]);
}
return new Result(minCost, bestRoot, treeEdges);
}
/**
* 递归回溯重建斯坦纳树的边
*/
private void reconstructTree(int[][] dp, int mask, int root, int[] terminals,
Set<int[]> edges, boolean[][] visited) {
if (visited[mask][root]) return;
visited[mask][root] = true;
int k = terminals.length;
// 检查是否是子集合并点
boolean foundSplit = false;
for (int sub = (mask - 1) & mask; sub > 0; sub = (sub - 1) & mask) {
int other = mask ^ sub;
if (other == 0 || sub > other) continue;
if (dp[sub][root] + dp[other][root] == dp[mask][root]) {
reconstructTree(dp, sub, root, terminals, edges, visited);
reconstructTree(dp, other, root, terminals, edges, visited);
foundSplit = true;
break;
}
}
if (!foundSplit) {
// 是边上移:找到前驱节点
for (Edge e : graph[root]) {
int prev = e.to;
if (dp[mask][prev] + e.weight == dp[mask][root]) {
int a = Math.min(root, prev);
int b = Math.max(root, prev);
edges.add(new int[]{a, b, e.weight});
reconstructTree(dp, mask, prev, terminals, edges, visited);
break;
}
}
}
}
/**
* Kou-Markowsky-Berman (KMB) 2-近似算法
* 时间复杂度: O(k * n^2 + k^2 * log k)
* 适用于大规模图,保证结果不超过最优解的2倍
*/
public Result kmbApproximation(int[] terminals) {
int k = terminals.length;
// 步骤1:计算所有终端之间的最短路径
int[][] allDist = allPairsShortestPath();
// 步骤2:构造metric closure(完全图,边权为最短路径距离)
// 步骤3:求完全图的最小生成树
List<int[]> metricEdges = new ArrayList<>();
for (int i = 0; i < k; i++) {
for (int j = i + 1; j < k; j++) {
metricEdges.add(new int[]{i, j, allDist[terminals[i]][terminals[j]]});
}
}
metricEdges.sort(Comparator.comparingInt(e -> e[2]));
// Kruskal求MST
int[] parent = new int[k];
for (int i = 0; i < k; i++) parent[i] = i;
List<int[]> mstEdges = new ArrayList<>();
for (int[] e : metricEdges) {
int ri = find(parent, e[0]);
int rj = find(parent, e[1]);
if (ri != rj) {
parent[ri] = rj;
mstEdges.add(e);
}
}
// 步骤4:展开MST边为原图中的最短路径,收集所有节点和边
Set<Integer> steinerNodes = new HashSet<>();
Set<String> steinerEdges = new HashSet<>();
for (int[] me : mstEdges) {
int u = terminals[me[0]];
int v = terminals[me[1]];
// 展开u到v的最短路径
List<Integer> path = getShortestPath(u, v, allDist);
for (int node : path) steinerNodes.add(node);
for (int i = 0; i < path.size() - 1; i++) {
int a = path.get(i), b = path.get(i + 1);
int w = getEdgeWeight(a, b);
steinerEdges.add(a + "#" + b + "#" + w);
}
}
// 步骤5:在展开子图上求MST
// 构建展开子图的邻接表
Map<Integer, List<Edge>> subGraph = new HashMap<>();
for (String es : steinerEdges) {
String[] parts = es.split("#");
int a = Integer.parseInt(parts[0]);
int b = Integer.parseInt(parts[1]);
int w = Integer.parseInt(parts[2]);
subGraph.computeIfAbsent(a, x -> new ArrayList<>()).add(new Edge(b, w));
subGraph.computeIfAbsent(b, x -> new ArrayList<>()).add(new Edge(a, w));
}
// Prim算法求MST
int totalCost = 0;
Set<Integer> mstNodes = new HashSet<>();
Set<String> finalEdges = new HashSet<>();
PriorityQueue<int[]> pq = new PriorityQueue<>(Comparator.comparingInt(e -> e[1]));
int start = terminals[0];
mstNodes.add(start);
for (Edge e : subGraph.getOrDefault(start, Collections.emptyList())) {
pq.offer(new int[]{e.to, e.weight, start});
}
while (!pq.isEmpty() && mstNodes.size() < steinerNodes.size()) {
int[] cur = pq.poll();
int v = cur[0], w = cur[1], u = cur[2];
if (mstNodes.contains(v)) continue;
mstNodes.add(v);
totalCost += w;
int a = Math.min(u, v), b = Math.max(u, v);
finalEdges.add(a + "#" + b + "#" + w);
for (Edge e : subGraph.getOrDefault(v, Collections.emptyList())) {
if (!mstNodes.contains(e.to)) {
pq.offer(new int[]{e.to, e.weight, v});
}
}
}
// 转换边格式
Set<int[]> edgeSet = new HashSet<>();
for (String es : finalEdges) {
String[] parts = es.split("#");
edgeSet.add(new int[]{Integer.parseInt(parts[0]),
Integer.parseInt(parts[1]),
Integer.parseInt(parts[2])});
}
return new Result(totalCost, start, edgeSet);
}
/**
* 获取a到b的最短路径(基于Floyd结果回溯)
*/
private List<Integer> getShortestPath(int a, int b, int[][] dist) {
// 简化为Dijkstra路径回溯
int[] prev = new int[n];
Arrays.fill(prev, -1);
int[] d = new int[n];
Arrays.fill(d, INF);
d[a] = 0;
PriorityQueue<int[]> pq = new PriorityQueue<>(Comparator.comparingInt(x -> x[1]));
pq.offer(new int[]{a, 0});
while (!pq.isEmpty()) {
int[] cur = pq.poll();
int u = cur[0];
if (u == b) break;
if (cur[1] > d[u]) continue;
for (Edge e : graph[u]) {
if (d[u] + e.weight < d[e.to]) {
d[e.to] = d[u] + e.weight;
prev[e.to] = u;
pq.offer(new int[]{e.to, d[e.to]});
}
}
}
List<Integer> path = new ArrayList<>();
int cur = b;
while (cur != -1) {
path.add(cur);
cur = prev[cur];
}
Collections.reverse(path);
return path;
}
private int getEdgeWeight(int u, int v) {
for (Edge e : graph[u]) {
if (e.to == v) return e.weight;
}
return INF;
}
private int find(int[] parent, int x) {
if (parent[x] != x) parent[x] = find(parent, parent[x]);
return parent[x];
}
/**
* 最小生成树(Prim算法)——用于对比
*/
public int mstWeight() {
boolean[] used = new boolean[n];
int[] minEdge = new int[n];
Arrays.fill(minEdge, INF);
minEdge[0] = 0;
int total = 0;
for (int i = 0; i < n; i++) {
int v = -1;
for (int j = 0; j < n; j++) {
if (!used[j] && (v == -1 || minEdge[j] < minEdge[v])) {
v = j;
}
}
used[v] = true;
total += minEdge[v];
for (Edge e : graph[v]) {
if (!used[e.to] && e.weight < minEdge[e.to]) {
minEdge[e.to] = e.weight;
}
}
}
return total;
}
/**
* 结果封装
*/
static class Result {
int cost; // 总代价
int root; // 根节点(参考)
Set<int[]> edges; // 树的边集
Result(int cost, int root, Set<int[]> edges) {
this.cost = cost;
this.root = root;
this.edges = edges;
}
@Override
public String toString() {
StringBuilder sb = new StringBuilder();
sb.append("Steiner Tree Cost: ").append(cost).append("\n");
sb.append("Edges (").append(edges.size()).append("):\n");
for (int[] e : edges) {
sb.append(" ").append(e[0]).append(" - ").append(e[1])
.append(" (w=").append(e[2]).append(")\n");
}
return sb.toString();
}
}
/**
* 主函数:演示斯坦纳树的精确解与近似解
*/
public static void main(String[] args) {
// 构造测试图:网格状拓扑
SteinerTree st = new SteinerTree(12);
// 构建图(0-1-2-3 横向,0-4-8 纵向等)
int[][] edges = {
{0,1,2}, {1,2,3}, {2,3,2},
{0,4,4}, {1,5,1}, {2,6,3}, {3,7,2},
{4,5,2}, {5,6,1}, {6,7,3},
{4,8,2}, {5,9,4}, {6,10,1}, {7,11,2},
{8,9,3}, {9,10,2}, {10,11,1}
};
for (int[] e : edges) st.addEdge(e[0], e[1], e[2]);
// 终端节点:0, 3, 8, 11(四个角落)
int[] terminals = {0, 3, 8, 11};
System.out.println("=== 斯坦纳树问题演示 ===");
System.out.println("图节点数: 12, 终端节点: " + Arrays.toString(terminals));
System.out.println("全图MST权重: " + st.mstWeight());
// 精确算法
System.out.println("\n--- Dreyfus-Wagner 精确算法 ---");
long start1 = System.currentTimeMillis();
Result exact = st.dreyfusWagner(terminals);
long time1 = System.currentTimeMillis() - start1;
System.out.println("耗时: " + time1 + "ms");
System.out.println(exact);
// 2-近似算法
System.out.println("\n--- KMB 2-近似算法 ---");
long start2 = System.currentTimeMillis();
Result approx = st.kmbApproximation(terminals);
long time2 = System.currentTimeMillis() - start2;
System.out.println("耗时: " + time2 + "ms");
System.out.println(approx);
// 对比
System.out.println("\n=== 结果对比 ===");
System.out.printf("精确解代价: %d%n", exact.cost);
System.out.printf("近似解代价: %d%n", approx.cost);
System.out.printf("近似比: %.2f%n", (double) approx.cost / exact.cost);
// 测试更多终端的情况
System.out.println("\n=== 更多终端测试 (6个终端) ===");
int[] moreTerminals = {0, 2, 5, 7, 9, 11};
Result exact2 = st.dreyfusWagner(moreTerminals);
Result approx2 = st.kmbApproximation(moreTerminals);
System.out.printf("6终端精确解: %d, 近似解: %d, 近似比: %.2f%n",
exact2.cost, approx2.cost, (double) approx2.cost / exact2.cost);
}
}
关键设计解析
状态压缩DP的子集枚举
Dreyfus-Wagner算法的核心在于for (int sub = (mask-1) & mask; sub > 0; sub = (sub-1) & mask)这一枚举技巧。它按降序枚举mask的所有非空真子集,每个mask被划分为两个互补子集sub和mask^sub。通过sub > other的剪枝避免重复划分,将时间复杂度控制在O(3^k)。
SPFA式图松弛
在DP状态转移后,需要通过图上的边松弛来传播代价改善。不同于对每次状态运行完整Dijkstra,这里使用SPFA(Shortest Path Faster Algorithm)风格:将所有已更新节点入队,逐层松弛邻接边。这是因为状态间的依赖具有局部性,SPFA在小规模图上效率更高。
Metric Closure的构造
KMB算法的精髓在于将原图”度量闭包化”:终端间的最短路径距离满足三角不等式,构成一个度量空间。在这个度量空间上的MST,展开回原始图后代价不会超过最优斯坦纳树的2倍。这一证明基于树度量到图度量的Lipschitz嵌入理论。
复杂度分析
| 算法 | 时间复杂度 | 空间复杂度 | 适用场景 |
|---|---|---|---|
| Dreyfus-Wagner | O(3^k · n + 2^k · m) | O(2^k · n) | k ≤ 12的精确解 |
| KMB近似 | O(k·n² + k²·log k) | O(n²) | 大规模图的快速近似 |
| 全对最短路径 | O(n³)或O(k·(n+m)log n) | O(n²) | KMB前置步骤 |
其中k为终端数,n为节点数,m为边数。
扩展方向
- 终端数大于12的精确算法:结合整数线性规划(ILP)求解器处理中等规模实例。
- Steiner森林:多个不相交终端组分别连通,推广为Steiner森林问题。
- 有向斯坦纳树:处理有向图上的连通问题,近似比更高。
- ** prize-collecting变体**:允许不连接某些终端,但需支付惩罚代价。
- GPU并行:使用CUDA加速全对最短路径和DP状态计算。