範囲制約と整数線形計画問題の解法

範囲制約の多項式による定式化

\(f\) を2値変数の多項式とします。
範囲制約は、\(l \leq f \leq u~(l < u)\) の形で表されます。ここでの目的は、この範囲制約が満たされる場合に限り最小値 0 をとる多項式を設計することです。

基本的な考え方は、補助整数変数 \(a\) を導入し、\(a\) が \(l \leq a \leq u\) の範囲の値をとるようにすることです。次の式を考えます。

\[ g=(f − a)² \]

この式 \(g\) は、\(f = a\) のときにちょうど最小値 0 をとります。\(a\) は \(l \leq a \leq u\) の任意の整数値をとることができるため、この式 \(g\) は \(f\) 自身が同じ範囲内の整数値をとる場合に限り 0 を達成します。

この補助変数の手法を用いて、Hi-QUBO では範囲制約を実装しています。
\(f\) が線形式である場合、\(g\) は QUBO 式になります。
\(f\) が3次以上の式である場合、\(g\) は HUBO 式になります。

NOTE

Hi-QUBO は内部的に軽量な改良を用いており、範囲制約をわずかに少ない2値変数で符号化できるようにしています。

PyQBPP は内部的に軽量な改良を用いており、範囲制約をわずかに少ない2値変数で符号化できるようにしています。

詳細は 比較演算子 で説明されています。

整数線形計画問題の解法

整数線形計画問題のインスタンスは、目的関数と複数の線形制約から構成されます。
例えば、次の整数線形計画問題は、2つの変数、1つの目的関数、および2つの制約を持ちます。

\[\begin{split} \mathrm{Maximize:}&~5x+4y,\\ \mathrm{Subject}~\mathrm{to:}&~2x+3y \leq 24,\\ &~7x+5y \leq 54. \end{split}\]

この問題の最適解は、
\(x = 1, y = 2\) であり、目的関数値は 40 です。

以下の Hi-QUBO プログラムは、Easy Solver を用いてこの最適解を求めます。

range-constraints-and-ilp-program1.cpp
#include <qbpp/qbpp.hpp>
#include <qbpp/easy_solver.hpp>

int main() {
  auto x = 0 <= qbpp::var_int("x") <= 10;
  auto y = 0 <= qbpp::var_int("y") <= 10;
  auto f = 5 * x + 4 * y;
  auto c1 = 0 <= 2 * x + 3 * y <= 24;
  auto c2 = 0 <= 7 * x + 5 * y <= 54;
  auto g = -f + 100 * (c1 + c2);
  g.simplify_as_binary();
  auto solver = qbpp::EasySolver(g);
  auto sol = solver.search({{"time_limit", 1.0}});
  std::cout << "x = " << sol(x) << ", y = " << sol(y) << std::endl;
  std::cout << "f = " << sol(f) << std::endl;
  std::cout << "c1 = " << sol(c1) << ", c2 = " << sol(c2) << std::endl;
  std::cout << "c1.body(sol) = " << c1.body(sol) << ", c2.body(sol) = " << c2.body(sol) << std::endl;
}

この Hi-QUBO プログラムでは、

  • f は目的関数を表します。

  • c1 および c2 は範囲制約を表します。

  • g はそれらを1つの最適化式にまとめたものです。

目的は最大化であるため、目的関数は -f として符号を反転させています。また、制約 c1 および c2 には重み 100 のペナルティを与え、制約が高い優先度で満たされるようにしています。

g に対して Easy Solver のインスタンスを作成し、時間制限 1.0 秒で探索を実行します。最適解 sol を取得した後、プログラムは x、y、f、c1、c2、*c1、および *c2 の値を出力します。

このプログラムの出力は次のとおりです。

x = 4, y = 5
f = 40
c1 = 0, c2 = 0
*c1 = 23, *c2 = 53

ここで、

  • c1 は制約 0 ≤ 2x + 3y ≤ 24 を表す式です。

  • *c1 は線形式 2x + 3y を表します。

以上より、ソルバーが正しく最適解を見つけていることが確認できます。

以下の PyQBPP プログラムは、Easy Solver を用いてこの最適解を求めます。

range-constraints-and-ilp-program1.py
import pyqbpp as qbpp

x = qbpp.var("x", between=(0, 10))
y = qbpp.var("y", between=(0, 10))
f = 5 * x + 4 * y
c1 = qbpp.constrain(2 * x + 3 * y, between=(0, 24))
c2 = qbpp.constrain(7 * x + 5 * y, between=(0, 54))
g = -f + 100 * (c1 + c2)
g.simplify_as_binary()

solver = qbpp.EasySolver(g)
sol = solver.search(time_limit=1.0)

print(f"x = {sol(x)}, y = {sol(y)}")
print(f"f = {sol(f)}")
print(f"c1 = {sol(c1)}, c2 = {sol(c2)}")
print(f"2x+3y = {sol(c1.body)}, 7x+5y = {sol(c2.body)}")

この PyQBPP プログラムでは、

  • f は目的関数を表します。

  • c1 および c2 は範囲制約を表します。

  • g はそれらを1つの最適化式にまとめたものです。

目的は最大化であるため、目的関数は -f として符号を反転させています。また、制約 c1 および c2 には重み 100 のペナルティを与え、制約が高い優先度で満たされるようにしています。

g に対して Easy Solver のインスタンスを作成し、時間制限 1.0 秒で探索を実行します。最適解 sol を取得した後、プログラムは x、y、f、c1、c2、および制約本体の式の値を出力します。

このプログラムの出力は次のとおりです。

x = 4, y = 5
f = 40
c1 = 0, c2 = 0
*c1 = 23, *c2 = 53

ここで、

  • c1 は制約 0 ≤ 2x + 3y ≤ 24 を表す式です。

  • c1.body は線形式 2x + 3y を表します。

以上より、ソルバーが正しく最適解を見つけていることが確認できます。