QUBO++入門：replace関数（巡回セールスマン問題）

replace関数による二値変数への値の代入

QUBO++では，replace関数を使用してQUBO式の二値変数に値を代入することができます．以下に具体例を示します．次の4変数からなるQUBO式は，変数の合計が2のとき最小値0を取ります．

$$f(x0,x1,x2,x3)=(x0+x1+x2+x3−2)2 \begin{align*} f(x_0,x_1,x_2,x_3) &= (x_0+x_1+x_2+x_3-2)^2 \end{align*}$$

この式に対して，$x2=0x_2 = 0$ および $x3=1x_3 = 1$ を代入すると，以下のように簡約されます．

$$f(x0,x1,0,1)=(x0+x1−1)2 \begin{align*} f(x_0,x_1,0,1) &= (x_0+x_1-1)^2 \end{align*}$$

QUBO++では，このような代入は代入を定義した初期化リストを引数とするreplace関数で簡単に実現できます．

#include "qbpp.hpp"

int main() {
  auto x = qbpp::var("x", 4);
  auto f = qbpp::sqr(qbpp::sum(x) - 2);
  f.simplify_as_binary();
  std::cout << "f = " << f << std::endl;
  f.replace({{x[2], 0}, {x[3], 1}});
  f.simplify_as_binary();
  std::cout << "f = " << f << std::endl;
}

このプログラムを実行すると，以下の出力が得られます．

f = 4 -3*x[0] -3*x[1] -3*x[2] -3*x[3] +2*x[0]*x[1] +2*x[0]*x[2] +2*x[0]*x[3] +2*x[1]*x[2] +2*x[1]*x[3] +2*x[2]*x[3]
f = 1 -x[0] -x[1] +2*x[0]*x[1]

多数の変数に対する値の代入

多数の変数や代入する値が動的に決定される場合には，qbpp::MapListオブジェクトを使用すると便利です．このオブジェクトに変数と代入する値を順次追加し，replace関数でまとめて代入を実行します．

次の例では，サイズ100の変数ベクトルxに対して，$x1x_1$から$x99x_{99}$までを0に代入しています．

#include "qbpp.hpp"

int main() {
  auto x = qbpp::var("x", 100);
  auto f = qbpp::sqr(qbpp::sum(x) - 1);
  qbpp::MapList assign;
  for (size_t i = 1; i < 100; i++) {
    assign.push_back({x[i], 0});
  }
  f.replace(assign);
  f.simplify_as_binary();
  std::cout << "f = " << f << std::endl;
}

このプログラムを実行すると，以下のような出力が得られます．

f = 1 -x[0]

巡回セールスマンの開始頂点の固定

巡回セールスマン問題では，開始ノードを固定しても最適な巡回路を得ることができます．例えば，ノード0を開始ノードとして固定した場合でも，最適巡回路が計算可能です．QUBO++では，replace関数を用いて開始ノードを固定した巡回セールスマン問題を簡単に表現できます．

$nn$ノードの巡回セールスマン問題を解くQUBO式では，大きさ$n×nn \times n$の二値変数行列$xx$を用いて順列を表します．開始ノードを0に固定するため，次の制約を与えます．

• $x0,0=1x_{0,0} = 1$
• $xi,0=0x_{i,0} = 0$ （$1≤i≤n−11 \leq i \leq n-1$）
• $x0,j=0x_{0,j} = 0$ （$1≤j≤n−11 \leq j \leq n-1$）
  これにより，開始ノード順列の先頭をノード0に固定できます．

$n=4n=4$の場合，二値変数行列$xx$の固定後の形は次のようになります．
$$
\left(
\begin{array}{cccc}
1 & 0 & 0 & 0 \
0 & x_{1,1} & x_{1,2} & x_{1,3} \
0 & x_{1,2} & x_{2,2} & x_{2,3} \
0 & x_{1,3} & x_{2,3} & x_{3,3} \
\end{array}
\right)
$$
固定された要素以外の部分（未固定部分）は，残りの順列を表します．

以上を踏まえ，「QUBO++入門：巡回セールスマン問題」で示したプログラムを，開始ノードをノード0に固定するよう変更したプログラムを以下に示します．

#include "qbpp.hpp"
#include "qbpp_easy_solver.hpp"

int distance(const std::pair<int, int> &p1, const std::pair<int, int> &p2) {
  return static_cast<int>(
      std::round(std::sqrt((p1.first - p2.first) * (p1.first - p2.first) +
                           (p1.second - p2.second) * (p1.second - p2.second))));
}

int main() {
  std::vector<std::pair<int, int>> nodes = {
      {0, 0}, {10, 0}, {20, 0}, {0, 10}, {10, 10}, {20, 10}, {0, 20}, {10, 20}};

  auto x = qbpp::var("x", nodes.size(), nodes.size());
  auto permutation = qbpp::sum(qbpp::sum(x) == 1) +
                     qbpp::sum(qbpp::sum(qbpp::transpose(x)) == 1);

  auto dist = qbpp::expr();
  for (size_t i = 0; i < nodes.size(); i++) {
    for (size_t j = 0; j < nodes.size(); j++) {
      for (size_t k = 0; k < nodes.size(); k++) {
        dist += distance(nodes[j], nodes[k]) * x[i][j] *
                x[(i + 1) % nodes.size()][k];
      }
    }
  }

  auto f = simplify_as_binary(dist + permutation * 100);

  qbpp::MapList fix0 = {{x[0][0], 1}};
  for (uint32_t i = 1; i < nodes.size(); i++) {
    fix0.push_back({x[0][i], 0});
    fix0.push_back({x[i][0], 0});
  }

  auto g = qbpp::replace(f, fix0).simplify_as_binary();

  auto model = qbpp::QuadModel(g);
  std::cout << "var_count = " << model.var_count() << std::endl;
  std::cout << "linear_term_count = " << model.term_count(1) << std::endl;
  std::cout << "quadratic_term_count = " << model.term_count(2) << std::endl;

  auto solver = qbpp::easy_solver::EasySolver(model);
  solver.set_time_limit(1);
  auto sol = solver.search();

  auto tsp_sol = qbpp::Sol(f);
  tsp_sol.set(sol);
  tsp_sol.set(fix0);

  std::cout << tsp_sol << std::endl;
  auto tour = qbpp::onehot_to_int(tsp_sol.get(x));
  std::cout << str(tour, "tour") << std::endl;
}

• 開始ノードを固定しないQUBO式の構築：distとpermutationを用いて，巡回セールスマン問題のQUBO式fを構築しています．
• 開始ノードを固定する制約の追加：開始ノードをノード0に固定するため，二値変数への値の代入リストをqbpp::MapListオブジェクトfix0として構築しています．
• 開始ノード固定後のQUBO式の生成：fix0をreplace関数でfに適用し，制約が追加されたQUBO式をgとしています．
• QUBO式の解探索：EasySolverを用いてgを解き，結果をsolとして取得しています．
• 解の整理：tsp_solをfの解を格納するオブジェクトとして生成し，solとfix0をセットすることで，xのすべての二値変数の値を含むようにしています．
• 巡回路の出力：tsp_solからxの値を取得し，これを巡回路tourとして表示しています．

このプログラムを実行すると，次のような出力が得られます．

var_count = 49
linear_term_count = 49
quadratic_term_count = 1092
80:{{x[0][0],1},{x[0][1],0},{x[0][2],0},{x[0][3],0},{x[0][4],0},{x[0][5],0},{x[0][6],0},{x[0][7],0},{x[1][0],0},{x[1][1],1},{x[1][2],0},{x[1][3],0},{x[1][4],0},{x[1][5],0},{x[1][6],0},{x[1][7],0},{x[2][0],0},{x[2][1],0},{x[2][2],1},{x[2][3],0},{x[2][4],0},{x[2][5],0},{x[2][6],0},{x[2][7],0},{x[3][0],0},{x[3][1],0},{x[3][2],0},{x[3][3],0},{x[3][4],0},{x[3][5],1},{x[3][6],0},{x[3][7],0},{x[4][0],0},{x[4][1],0},{x[4][2],0},{x[4][3],0},{x[4][4],1},{x[4][5],0},{x[4][6],0},{x[4][7],0},{x[5][0],0},{x[5][1],0},{x[5][2],0},{x[5][3],0},{x[5][4],0},{x[5][5],0},{x[5][6],0},{x[5][7],1},{x[6][0],0},{x[6][1],0},{x[6][2],0},{x[6][3],0},{x[6][4],0},{x[6][5],0},{x[6][6],1},{x[6][7],0},{x[7][0],0},{x[7][1],0},{x[7][2],0},{x[7][3],1},{x[7][4],0},{x[7][5],0},{x[7][6],0},{x[7][7],0}}
{tour[0],0},{tour[1],1},{tour[2],2},{tour[3],5},{tour[4],4},{tour[5],7},{tour[6],6},{tour[7],3}

巡回路として，ノード0を開始ノードとする最適な順序0,1,2,5,4,7,6,3が得られています．開始ノードを固定したことで，QUBO式gの変数数は49，線形項は49，二次項は1092となり，開始ノードを固定しない場合と比べてQUBO式の規模が大幅に小さくなっています．

まとめ

• replace関数を使用することで，初期化リストを通じてQUBO式の変数に値を代入できます．
• 変数が多い場合や，代入する変数や値を動的に設定したい場合には，初期化リストの代わりにqbpp::MapListオブジェクトを利用し，replace関数に渡します．
• 巡回セールスマン問題を解くQUBO式では，開始ノードを固定することで，式の規模を縮小し，計算効率を向上させることが可能です．

QUBO++情報

• QUBO＋＋ツールとABS2 QUBO Solverの詳細情報：QUBO++に関する公式サイトです．ツールの特徴や使用方法，ABS2 QUBO Solverとの連携について詳しく説明されています．
• 【研究成果】組合せ最適化問題を解くためのプログラミング用ツールQUBO＋＋を無償公開します：広島大学からの公式リリース情報です．QUBO++の開発背景や研究成果について記載されています．
• 利用条件：QUBO++は非商用利用や評価研究目的に限り無償で利用可能です．商用利用や定常的な業務での使用を希望する場合は，ライセンス契約が必要となります．