Day 081-Q5 — 拡張ラグランジュ法・ADMM（凸最適化の主双対反復）

問題

グラフ $G = (V, E)$（$|V| = N$, $|E| = M$）が与えられ，各頂点 $i$ に重み $c_i > 0$ が付く。最大重み独立集合の LP 緩和を ADMM で解き，LP 最適値と整数解（丸め）の目的関数値を求めよ。

入出力例

4 4
3 5 2 4
1 2
2 3
3 4
4 1

LP: 9.00
Integer: 9

概念図: ADMM の 3 ステップ反復

ヒント

import numpy as np
rho = 2.0
A = ...  # (M, N) 辺-頂点接続行列
M_mat = A.T @ A + rho * np.eye(N)
M_inv = np.linalg.inv(M_mat)
x = np.zeros(N); z = np.ones(M); u_dual = np.zeros(M)

for _ in range(2000):
    rhs = c + rho * A.T @ (z - u_dual)
    x = np.clip(M_inv @ rhs, 0.0, 1.0)
    z = np.minimum(A @ x + u_dual, 1.0)
    u_dual = u_dual + A @ x - z

模範解答

import sys
import numpy as np

def solve():
    data = sys.stdin.read().split()
    idx = 0
    N = int(data[idx]); idx+=1
    M = int(data[idx]); idx+=1
    c = np.array([float(data[idx+i]) for i in range(N)]); idx+=N
    edges = []
    adj = [[] for _ in range(N)]
    for _ in range(M):
        u, v = int(data[idx])-1, int(data[idx+1])-1; idx+=2
        edges.append((u, v)); adj[u].append(v); adj[v].append(u)

    if M == 0:
        lp_val = float(sum(c))
        print(f"LP: {lp_val:.2f}")
        print(f"Integer: {int(lp_val)}")
        return

    A = np.zeros((M, N))
    for k, (u, v) in enumerate(edges):
        A[k, u] = 1.0; A[k, v] = 1.0

    rho = 2.0
    M_mat = A.T @ A + rho * np.eye(N)
    M_inv = np.linalg.inv(M_mat)

    x = np.zeros(N); z = np.ones(M); u_dual = np.zeros(M)

    for it in range(2000):
        x_prev = x.copy()
        rhs = c + rho * A.T @ (z - u_dual)
        x = np.clip(M_inv @ rhs, 0.0, 1.0)
        z = np.minimum(A @ x + u_dual, 1.0)
        u_dual = u_dual + A @ x - z
        if np.max(np.abs(x - x_prev)) < 1e-8: break

    lp_val = float(c @ x)

    # 整数解: x 値降順に貪欲選択
    sel = np.zeros(N, dtype=bool)
    for i in np.argsort(-x):
        if all(not sel[nb] for nb in adj[i]):
            sel[i] = True

    int_val = int(c @ sel)
    print(f"LP: {lp_val:.2f}")
    print(f"Integer: {int_val}")

solve()

Step-by-Step 解説

Step 1: ADMM の基本形

問題: $\min f(x) + g(z)$, s.t. $Ax + Bz = c$

Augmented Lagrangian: $L_\rho = f(x) + g(z) + y^\top(Ax+Bz-c) + \frac{\rho}{2}\|Ax+Bz-c\|^2$

3 ステップの交互更新で $L_\rho$ を最小化。凸問題に対して収束が保証される。

Step 2: LP 緩和への適用

辺制約 $x_i + x_j \le 1$ を補助変数 $z$ に持たせる。x-update は $[0,1]$ ボックス制約への射影（クリップ），z-update は $z \le 1$ への射影（element-wise min）。

Step 3: $M^{-1}$ の前計算

x-update の線形方程式 $(A^\top A + \rho I)x = r$ は毎回同じ左辺 → $M^{-1} = (A^\top A + \rho I)^{-1}$ を事前計算することで各 iter が $O(N^2)$ に。

Step 4: 整数解の丸め戦略

LP 最適値 $x^*$ から整数解を構築：$x_i \ge 0.5$ の頂点を降順に選び，隣接頂点が選ばれていなければ採用（貪欲）。

よくあるミス

次のステップ

• 発展問題: SDP（半正定値計画）+ ADMM による MAX-CUT 近似解法
• 参考: CVXPY（Python の凸最適化ライブラリ）で同問題を解いて ADMM 解と比較
• $\rho$ の適応的調整（Primal/Dual Residual バランシング）

自己評価