Hi-QUBO

部分グラフ同型問題

二つの無向グラフ $G_H=(V_H,E_H)$(ホストグラフ)と $G_G=(V_G,E_G)$(ゲストグラフ)が与えられたとき、部分グラフ同型問題とは、$G_H$ に $G_G$ と同型な部分グラフが含まれるかどうかを問う問題である。

より形式的に言えば、目的は単射写像 $\sigma:V_G\rightarrow V_H$ を求めることであり、 そのとき、任意の辺 $(u,v)\in E_G$ に対して、組 $(\sigma(u),\sigma(v))$ もホストグラフの辺となる、すなわち $(\sigma(u),\sigma(v))\in E_H$ となる。

例えば、以下のホストグラフとゲストグラフを考えてみよう:

Host Graph
ホストグラフ $G_H=(V_H,E_H)$ の例(10ノード)

Guest Graph
6つの頂点を持つゲストグラフ $G_G=(V_G,E_G)$ の例

一つの解 $\sigma$ は次の通り:

ノード $i$ 0 1 2 3 4 5
$\sigma(i)$ 1 4 6 7 9 8

この解法は以下のように可視化される:

The solution of the subgraph isomorphism problem
部分グラフ同型問題の解法

部分グラフ同型問題のQUBO定式化

ゲストグラフ $G_G=(V_G,E_G)$ は $m$ 個の頂点 $0, 1, \ldots, m-1$ を持ち、 ホストグラフ $G_H=(V_H,E_H)$ は $n$ 個の頂点 $0, 1, \ldots, n-1$ を持つと仮定する。 $m\times n$ の二値行列 $X=(x_{i,j})$ ($0\leq i\leq m-1, 0\leq j\leq n-1$) を導入する。これは $mn$ 個の二値変数を持つ。 この行列は単射写像 $\sigma:V_G\rightarrow V_H$ を表し、 $x_{i,j}=1$ となるのは $\sigma(i)=j$ である場合に限る。

例えば、部分グラフ同型問題の解は次の $6\times 10$ バイナリ行列で表現できる:

  0 1 2 3 4 5 6 7 8 9
0 0 1 0 0 0 0 0 0 0 0
1 0 0 0 0 1 0 0 0 0 0
2 0 0 0 0 0 0 1 0 0 0
3 0 0 0 0 0 0 0 1 0 0
4 0 0 0 0 0 0 0 0 0 1
5 0 0 0 0 0 0 0 0 1 0

$X$ が単射写像を表すため、以下の制約を満たさなければならない:

これらを以下のQUBO++形式の制約に統合できる。全制約が満たされた時に最小値を取る:

\[\begin{aligned} \text{制約} &= \sum_{i=0}^{m-1}\Bigl(\sum_{j=0}^{n-1}x_{i,j} = 1\Bigr)+\sum_{j=0}^{m-1}\Bigl(0\leq \sum_{i=0}^{n-1}x_{i,j} \leq 1\Bigr) \end{aligned}\]

QUBO形式では、同じ制約を次のように表現できる:

\[\begin{aligned} \text{制約条件} &= \sum_{i=0}^{m-1}\Bigl(\sum_{j=0}^{n-1}x_{i,j} - 1\Bigr)^2+\sum_{j=0}^{m-1}\sum_{i=0}^{n-1}x_{i,j}\Bigl(\sum_{i=0}^{n-1}x_{i,j}-1\Bigr) \end{aligned}\]

次に、目的関数をホストエッジにマッピングされるゲストエッジの数として定義する:

\[\begin{aligned} \text{目的関数} &= \sum_{(u_G,v_G)\in E_G}\sum_{(u_H,v_H)\in E_H} (x_{u_G,u_H}x_{v_G,v_H}+x_{u_G,v_H}x_{v_G,u_H}) \end{aligned}\]

ここで、無向ゲストエッジ $(u_G,v_G)\in E_G$ は、ホストエッジ $(u_H,v_H)\in E_H$ に2つの対称的な方法で対応させることができる:

したがって、二次項 $x_{u_G,u_H}x_{v_G,v_H}$ と $x_{u_G,v_H}x_{v_G,u_H}$ の両方を包含する。

最後に、目的関数と制約を単一のQUBO式に統合する:

\[\begin{aligned} f &= -\text{目的関数} + mn\times \text{制約条件} \end{aligned}\]

ペナルティ係数 $mn$ は、目的関数の改善よりも制約条件の充足を優先するように選択される。 制約項がゼロで目的関数がゲストエッジの数に等しいときに、$f$ の最良値が達成される。

部分グラフ同型問題のQUBO++プログラム

上記のQUBO定式化に基づき、以下のQUBO++プログラムはゲストグラフの頂点数$M=6$、ホストグラフの頂点数$N=10$の場合のサブグラフ同型問題を解く:

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

int main() {
  const size_t N = 10;
  std::vector<std::pair<size_t, size_t>> host = {
      {0, 1}, {0, 2}, {1, 3}, {1, 4}, {1, 6}, {2, 5}, {3, 7}, {4, 6},
      {4, 7}, {5, 6}, {5, 8}, {6, 8}, {6, 7}, {7, 9}, {8, 9}};

  const size_t M = 6;
  std::vector<std::pair<size_t, size_t>> guest = {
      {0, 1}, {0, 2}, {1, 2}, {1, 3}, {2, 3}, {2, 5}, {3, 4}, {4, 5}};

  auto x = qbpp::var("x", M, N);

  auto host_assigned = qbpp::vector_sum(qbpp::transpose(x));

  auto constraint =
      qbpp::sum(qbpp::vector_sum(x) == 1) + qbpp::sum(0 <= host_assigned <= 1);

  auto objective = qbpp::toExpr(0);

  for (const auto& e_g : guest) {
    for (const auto& e_h : host) {
      objective += x[e_g.first][e_h.first] * x[e_g.second][e_h.second] +
                   x[e_g.first][e_h.second] * x[e_g.second][e_h.first];
    }
  }

  auto f = -objective + constraint * (M * N);

  f.simplify_as_binary();

  auto solver = qbpp::easy_solver::EasySolver(f);
  solver.target_energy(-static_cast<int>(guest.size()));
  auto sol = solver.search();

  std::cout << "sol(x) = " << sol(x) << std::endl;

  std::cout << "sol(objective) = " << sol(objective) << std::endl;
  std::cout << "sol(constraint) = " << sol(constraint) << std::endl;

  auto guest_to_host = qbpp::onehot_to_int(sol(x));
  std::cout << "guest_to_host = " << guest_to_host << std::endl;

  auto host_to_guest = qbpp::onehot_to_int(sol(qbpp::transpose(x)));
  std::cout << "host_to_guest = " << host_to_guest << std::endl;

  qbpp::graph::GraphDrawer guest_graph;
  for (size_t i = 0; i < M; ++i) {
    guest_graph.add_node(qbpp::graph::Node(i));
  }
  for (const auto& e : guest) {
    guest_graph.add_edge(qbpp::graph::Edge(e.first, e.second));
  }
  guest_graph.write("guest_graph.svg");

  qbpp::graph::GraphDrawer host_graph;
  for (size_t i = 0; i < N; ++i) {
    host_graph.add_node(qbpp::graph::Node(i));
  }
  for (const auto& e : host) {
    host_graph.add_edge(qbpp::graph::Edge(e.first, e.second));
  }
  host_graph.write("host_graph.svg");

  qbpp::graph::GraphDrawer graph;
  for (size_t i = 0; i < N; ++i) {
    graph.add_node(qbpp::graph::Node(i).color(sol(host_assigned[i])));
  }

  std::vector<std::vector<bool>> guest_adj(N, std::vector<bool>(N, false));
  for (auto [u, v] : guest) {
    guest_adj[u][v] = guest_adj[v][u] = true;
  }
  for (const auto& e_h : host) {
    auto u = host_to_guest[e_h.first];
    auto v = host_to_guest[e_h.second];
    if (u != -1 && v != -1 &&
        guest_adj[static_cast<size_t>(u)][static_cast<size_t>(v)]) {
      graph.add_edge(
          qbpp::graph::Edge(e_h.first, e_h.second).color(1).penwidth(2.0f));
    } else {
      graph.add_edge(qbpp::graph::Edge(e_h.first, e_h.second));
    }
  }

  graph.write("subgraph_isomorphism.svg");
}

ゲストグラフとホストグラフは、それぞれ辺リスト guest と host として与えられる。 $M\times N$ の二値行列 xを定義し、式 constraint, objective, および f を構築する。

対象エネルギーは fが作成され、目標エネルギーは $−∣E_G|$(ゲストエッジの負の数)に設定される。これは、 -objective すべてのゲストエッジがホストエッジにマッピングされた場合の最適値である。 得られた解は solに格納される。 の値は x, objective, および constraintsol の値がそれぞれ出力される。

関数 qbpp::onehot_to_int()を使用すると、プログラムはゲストノードからホストノードへのマッピング(guest_to_host)およびホストノードからゲストノードへのマッピング(host_to_guest).

ゲストとホストのグラフは guest_graph.svghost_graph.svgにそれぞれ保存される。 最後に、解は subgraph_isomorphism.svgで可視化され、マッピングによって選択されたホストノードと、ゲストエッジに対応するホストエッジが強調表示されます。

このプログラムは以下の出力を生成する:

sol(x) = {{0,1,0,0,0,0,0,0,0,0},{0,0,0,0,1,0,0,0,0,0},{0,0,0,0,0,0,1,0,0,0},{0,0,0,0,0,0,0,1,0,0},{0,0,0,0,0,0,0,0,0,1},{0,0,0,0,0,0,0,0,1,0}}
sol(objective) = 8
sol(constraint) = 0
guest_to_host = {1,4,6,7,9,8}
host_to_guest = {-1,0,-1,-1,1,-1,2,3,5,4}

目的関数値はゲスト辺の数($|E_G|=8$)に等しく、全ての制約が満たされる(constraint = 0)。 したがって、このプログラムは有効な部分グラフ同型に対応する最適解を見出している。 なお、ホストノードがゲストノードからマッピングされていない場合、host_to_guestのエントリは -1 である場合、対応するホストノードはどのゲストノードからもマッピングされていないことを示す。


最終更新日: 2026.01.13