0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

この記事について

整数計画で困ることのひとつが「非線形なものを扱いたい」という場面です。区分関数($\max, \min$)、非凸な非線形関数、2次の積項…。これらを線形のMIPで表現できれば、既存のソルバー(Gurobi、CPLEX、HiGHS等)の強力なエンジンに乗せられます。

この記事では、区分線形近似2値変数の積の線形化という2つの代表的なテクニックを紹介します。


非凸な非線形関数の区分線形近似

非凸な非線形関数 $f(x)$ を最小化したい場合、そのままだと整数計画では扱えません。そこで 区分線形関数 $g(x)$ で近似します。

20260608_piecewise_linear_approx.png

ブレークポイント $a_1, a_2, \ldots, a_m$ を選び、その点での値 $f(a_i)$ を結んだ折れ線で近似します。ブレークポイントが多いほど近似精度は上がりますが、変数も増えるので計算が重くなります。

線分上の点の表現

区分線形関数の任意の点は「いずれかの線分上にある」ので、その線分の2端点の重み付き平均で書けます。

(x, g(x)) = t_i (a_i, f(a_i)) + t_{i+1} (a_{i+1}, f(a_{i+1}))
t_i + t_{i+1} = 1, \quad t_i, t_{i+1} \geq 0

これを全体に拡張すると、ブレークポイント $a_1, \ldots, a_m$ の重みベクトル $(t_1, \ldots, t_m)$ で:

(x, g(x)) = \sum_{i=1}^{m} t_i (a_i, f(a_i)), \quad \sum_{i=1}^{m} t_i = 1, \quad t_i \geq 0

$t_i$ のうち高々2つの隣り合うものだけが正でなければなりません。これがないと、隣同士でない2点の重み付き平均になってしまい、関数 $f$ と値がずれます。

「隣り合う2つだけが正」を整数計画で表現

各線分 $i$ に対応する2値変数 $z_i \in {0, 1}$ を導入し、「点が $i$ 番目の線分上にある」ことを表します。

\begin{aligned}
t_1 &\leq z_1 \\
t_i &\leq z_{i-1} + z_i, \quad i = 2, \ldots, m-1 \\
t_m &\leq z_{m-1} \\
\sum_{i=1}^{m-1} z_i &= 1 \\
z_i &\in \{0, 1\}
\end{aligned}
  • $z_i = 1$ なら点は $i$ 番目の線分上にある($t_i$ と $t_{i+1}$ のみが正になり得る)
  • $\sum z_i = 1$ により、ちょうど1つの線分が選ばれる

これは SOS2(Special Ordered Set type 2) という制約タイプとして Gurobi・CPLEX に標準搭載されています。直接 SOS2 制約として渡せば、より効率的に処理してくれます。

Pythonで区分線形近似を実装

import numpy as np
import pulp

# 近似したい非線形関数(例:sin に少し小細工を加えた非凸関数)
def f(x):
    return np.sin(x) + 0.1 * x

# ブレークポイント
xs = np.linspace(0, 10, 11)   # 11点
fs = [f(x) for x in xs]
m = len(xs)

prob = pulp.LpProblem("piecewise", pulp.LpMinimize)

# 区分線形近似に使う変数
t = [pulp.LpVariable(f"t{i}", lowBound=0, upBound=1) for i in range(m)]
z = [pulp.LpVariable(f"z{i}", cat='Binary') for i in range(m-1)]

# x と g(x) を t で表現
x_var = pulp.lpSum(xs[i] * t[i] for i in range(m))
g_x   = pulp.lpSum(fs[i] * t[i] for i in range(m))

# 目的関数:g(x) を最小化
prob += g_x

# SOS2 制約(隣り合う2つの t のみが正)
prob += t[0] <= z[0]
for i in range(1, m-1):
    prob += t[i] <= z[i-1] + z[i]
prob += t[m-1] <= z[m-2]

prob += pulp.lpSum(t) == 1
prob += pulp.lpSum(z) == 1

prob.solve(pulp.PULP_CBC_CMD(msg=0))

x_opt = pulp.value(x_var)
g_opt = pulp.value(g_x)
print(f"近似解 x = {x_opt:.4f}")
print(f"近似 g(x) = {g_opt:.4f}")
print(f"真の f(x) = {f(x_opt):.4f}")

実行結果例:

近似解 x = 5.0000
近似 g(x) = -0.4589
真の f(x) = -0.4589

これで非凸関数の大域最適解に近い点を見つけられます。


2値変数の積を線形化する

2次計画問題(QP)や QUBO(Quadratic Unconstrained Binary Optimization)を MIP に変換するときに重要な技法です。

問題

\text{minimize} \quad \sum_{i=1}^{n} \sum_{j=1}^{n} q_{ij} x_i x_j
\text{subject to} \quad x_i \in \{0, 1\}

$x_i x_j$ の積項が非線形です。これを線形化したいです。

古典的な線形化

新しい2値変数 $y_{ij} \in {0, 1}$ を導入して $y_{ij} = x_i x_j$ を表します。$x_i, x_j \in {0,1}$ の組み合わせは4通り:

$x_i$ $x_j$ $x_i x_j$
0 0 0
1 0 0
0 1 0
1 1 1

これを満たすには3本の線形制約で十分です。

\begin{aligned}
y_{ij} &\leq x_i \\
y_{ij} &\leq x_j \\
y_{ij} &\geq x_i + x_j - 1 \\
y_{ij} &\in \{0, 1\}
\end{aligned}

確認

$x_i$ $x_j$ 上限 $\min(x_i, x_j)$ 下限 $x_i + x_j - 1$ $y_{ij}$
0 0 0 -1 0
1 0 0 0 0
0 1 0 0 0
1 1 1 1 1

線形化後の MIP

\begin{aligned}
\text{minimize} \quad & \sum_{i=1}^{n} \sum_{j=1}^{n} q_{ij} y_{ij} \\
\text{subject to} \quad & y_{ij} \leq x_i \\
& y_{ij} \leq x_j \\
& y_{ij} \geq x_i + x_j - 1 \\
& x_i \in \{0, 1\} \\
& y_{ij} \in \{0, 1\}
\end{aligned}

これで標準的なMIPソルバーで解けます。

Pythonでの実装:QUBO → MIP変換

import numpy as np
import pulp

# QUBO 行列の例(対称・非対角に交互作用あり)
np.random.seed(42)
n = 5
Q = np.random.randint(-5, 5, (n, n))
Q = (Q + Q.T) / 2   # 対称化

prob = pulp.LpProblem("qubo_to_mip", pulp.LpMinimize)

x = [pulp.LpVariable(f"x{i}", cat='Binary') for i in range(n)]
y = {(i, j): pulp.LpVariable(f"y{i}_{j}", cat='Binary')
     for i in range(n) for j in range(n)}

# 線形化された目的関数
prob += pulp.lpSum(Q[i][j] * y[(i, j)] for i in range(n) for j in range(n))

# 線形化制約
for i in range(n):
    for j in range(n):
        prob += y[(i, j)] <= x[i]
        prob += y[(i, j)] <= x[j]
        prob += y[(i, j)] >= x[i] + x[j] - 1

prob.solve(pulp.PULP_CBC_CMD(msg=0))

x_opt = [int(pulp.value(xi)) for xi in x]
print(f"x* = {x_opt}")
print(f"目的関数値 = {pulp.value(prob.objective):.2f}")

# 検証:直接 xQx を計算
x_np = np.array(x_opt)
direct = x_np @ Q @ x_np
print(f"直接計算 xQx = {direct:.2f}")

実行結果例:

x* = [0, 1, 0, 1, 1]
目的関数値 = -12.00
直接計算 xQx = -12.00

線形化した目的関数値と直接計算した $\mathbf{x}^\top Q \mathbf{x}$ が一致しています。これで非線形の2次問題を線形のMIPで解けました。


まとめ

  • 区分線形近似:非凸な非線形関数 $f(x)$ をブレークポイントで折れ線近似
  • 「隣り合う2つの重みのみが正」を SOS2 制約 で表現
  • 2値変数の積の線形化:$y_{ij} = x_i x_j$ を3本の線形制約で表現
  • これによりQUBO(2次の0-1最適化)が MIP として解ける
0
0
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?