この数式を展開して得られる二次式が順列生成イジングモデルとなります．このイジングモデルのサイズはn³-n²であり、要求分解能は2n-4となります。

図1:　巡回セールスマン問題の解[3,1,0,2,4]とそれをone-hot表現で表す大きさ５✕5の行列

　今回開発したdual-matrix domain-wall法によるイジングモデルでは、サイズを6n²-12n+4，要求分解能を2に、それぞれ大幅に削減しました。これにより、順列型組合わせ最適化問題をイジングモデルに変換し、量子アニーラーで解を求める場合に、より大きな問題を高い精度で取り扱うことができるようになります。

この式はs_i,j，a_i,j，b_i,jを変数とする二次式に展開することができ、得られた二次式が順列生成イジングモデルとなります。4つの総和項が上の(1)から(4)の条件にそれぞれ対応します。4つの条件を満たすとき、4つの総和項のそれぞれが最小値2n，2n，0，0を取ります。したがって、このイジングモデルが最小値4nをとるとき、行列Sにone-hot表現による順列が格納されています．図2はスピン変数の行列S,A,Bが格納する値の例と、その場合の行列∆S，∆A，∆Bの各要素の値を示しています。(1)～(4)の条件をすべて満たしていることが確認できます。このイジングモデルの二次項の個数は6ⁿ²-8n、要求分解能は2となり、従来のone-hot法による順列生成イジングモデルより大幅な削減となっています。

図2：Dual-matrix domain-wall法によるイジングモデルの最適解の例．この最適解は順列[3,1,0,4,2]とその逆順列[2,1,4,0,3]を表している。

　順列型組合せ最適化問題には、例えば、巡回セールスマン問題、グラフ同型問題、二次割当問題などがあります。順列生成イジングモデルはこれらの問題に適用可能です。例えば、無向グラフの巡回セールスマン問題(＊10)では図3の重み付き無向グラフが与えられたときに、すべての頂点を一度ずつ訪問する最短巡回を求めます。図3では、巡回順を表す順列[0,1,4,8,…,3]が最適解になります。図4のグラフは、このような接続関係が平面の重み付き無向グラフについて、巡回セールスマン問題を解くイジングモデルを生成したときのサイズを表しています。グラフの頂点数が多くなると、従来のone-hot法ではサイズが急激に増加します。all-different domain-wall法を用いることにより、サイズが約半分に節約できています。一方、dual-matrix domain-wall法は、頂点数の増加にともなうサイズの増加がこれら2つの方法に比べ極めて緩やかであり、無向グラフが頂点数300のときは、従来のone-hot法に比べて25分の1以下に削減できています。

図4：無向グラフの巡回セールスマン問題を解くイジングモデルの二次項の個数

　なお、dual-matrix domain-wall法の技術内容は特許出願済みです。また、論文プレプリントでは、部分順列生成(＊11)へのdual-matrix domain-wall法の拡張や、二次割当問題、部分グラフ同型問題、マッチング問題などへ適用した場合の二次項の個数、最新の量子アニーラーD-Wave Advantageのペガサスグラフに順列生成イジングモデルを埋め込んだ場合に使用するセル数の評価を行っています。

図5：イジングモデルのグラフ表現

＊2 量子アニーラー
　量子効果に基づいて、特定のグラフ形状をもつイジングモデルの解を求める量子コンピューターを｢量子アニーラー｣と呼びます。D-Wave Systems社の開発した量子アニーラーD-Wave 2000Qは、2048個のセルがキメラグラフの形状で接続されています。各セルをイジングモデルの変数とその一次項に対応付け、セル間の接続を二次項に対応付けるよう設定すると、そのキメラグラフの形状をもつイジングモデルの解を量子焼きなまし(アニーリング)により求めることができます。図6は128セルのキメラグラフを示しており、これを16倍の大きさ(縦横それぞれ４倍)に拡張したものが2048セルのキメラグラフです。

図6：128セルのキメラグラフ

＊3 イジングモデルの量子アニーラーへの埋め込み
　量子アニーラーD-Wave 2000Qでキメラグラフとは異なる形状のイジングモデルを扱う場合は、事前に埋め込みを行います。この埋め込みでは、イジングモデルの変数をキメラグラフのセルに対応付けます。対応付けは、イジングモデルの接続関係がキメラグラフ上で維持されるように定める必要があります。
　図7は、17変数イジングモデルと33変数イジングモデルの128セルのキメラグラフへの量子アニーラーへの埋め込みをそれぞれ示しています。17変数イジングモデルは、変数のすべての組合せに二次項を持つため、密なイジングモデルとなっています。一方、33変数イジングモデルの各変数は最大で4つの二次項を構成しており、疎なイジングモデルです。変数とセルの対応は接続関係を維持するように定める必要があるため、イジングモデルの各変数がキメラグラフの複数のセルに割り当てられています。疎なイジングモデルの方が変数の個数は多いにもかかわらず、サイズ(二次項の個数)が小さいため、使用するセルの数はかなり少なくなっています。したがって、変数の個数よりもサイズを削減することが、埋め込みに必要となる量子アニーラーのセルの個数の節約につながります。

図7：密な17変数イジングモデル(二次項の個数：136)と疎な33変数イジングモデル(二次項の個数：51)の128頂点キメラグラフ(右側)への埋め込み。例えば、疎なイジングモデルの変数31は、右側のキメラグラフでは2つのセルに割り当てられています．イジングモデルでは変数31は変数14,15,30,32と接続していますが、対応付けられたセルにおいてもこれらの接続が維持されています。

＊4 量子アニーラーの計算精度とイジングモデルの要求分解能
　整数係数に限定したイジングモデルにおいて、係数の最大絶対値を要求分解能と呼びます。この値が大きいほど、量子アニーラーに高い計算精度を要求することになります。それは次のような係数圧縮が行われるためです。
　実際の量子アニーラーで扱うイジングモデルの係数は例えば範囲が-1.0～+1.0というように範囲指定された実数値を設定できます。しかし、イジングモデルが量子アニーラーに要求する計算精度を正確に評価するために、あえてイジングモデルの各項の係数を整数に限定することが行われます。そのようなイジングモデルを量子アニーラーで扱う際には、整数係数の最大絶対値(要求分解能)で各係数を除する前処理を行います。これにより、係数を-1.0～+1.0の範囲に圧縮して、量子アニーラーで扱えるようにします。例えば、上記のイジングモデルの場合、要求分解能である3で除して次のようになります：

同じ値で除しているので、最適解が変わることはありません。量子アニーラーの計算精度は限られているため、除する整数が大きくなると、絶対値の小さい係数や差の小さい2つの係数が現れることがあります。そのような小さい係数や係数の差を正しく取り扱うことができなくなるため、計算精度が限られた量子アニーラーでは最適解を得るのが困難になります。以上より、イジングモデルの係数の最大絶対値、つまり、要求分解能は小さい方が望ましいと言えます。現在の量子アニーラーの計算精度は5～6ビット程度と言われており、要求分解能は2⁵=32程度が限界です。

＊5 巡回セールスマン問題
　複数の都市が入力として与えられて、それらをすべて一度ずつ巡る訪問順で総巡回路長が最短のものを求める組合せ最適化問題。訪問順は都市の順列なので、順列型組合せ最適化問題です。

＊6 スピン値によるone-hot表現
　0からn-1の整数を表すのに、n個のスピン値を用い、k番目が+1，他が-1のとき、整数kを表すとみなします。ここで、先頭は0番目であり、例えば、 [-1,-1,-1,+1,-1]は3を表します。n個のスピン変数s₀,s₁,…,s_n-1がone-hot表現であるスピン値を格納するためのイジングモデルは、合計が(-1)×(n-1)+(+1)×1=-(n-2)となることを利用して次の二次式で得られます。

n個のスピン変数がone-hot表現のとき、この式の評価値は最小値0となります．このイジングモデルのサイズは1/2 n²-1/2 nであり、要求分解能はn-2となります．また、図１の行列のように、n×nの行列の各行がone-hot表現で整数を表し、それらの整数が異なる場合、得られる整数列は順列となります．これを順列のone-hot表現と呼びます．

＊7 スピン値によるdomain-wall表現
　0からn-1の整数を表すのに、n-1個のスピン値を用い、先頭からk番目までが+1、残りが-1のとき、整数kを表すとみなします。ここで、先頭は0番目であり、例えば、 [+1,+1,+1,-1]は3を表します。n-1個のスピン変数s₀,s₁,…,s_n-2がdomain-wall表現であるスピン値を格納するためのイジングモデルは、+1と-1が隣り合っているところが一箇所であることを利用して次の二次式で得られます。
 

ここで、s_-1=+1，s_n-1=-1として左辺を展開しています。n-1個のスピン変数がdomain-wall表現のとき、s_i-1-s_i=2となるのは一箇所で、それ以外は0となります。よって、この式の評価値はdomain-wall表現のときに、最小値4となります。このイジングモデルのサイズはn-2であり、要求分解能は1となり、one-hot表現に比べて大幅に小さくなります。

＊8 QUBOモデル
  QUBO(Quadratic Unconstrained Binary Optimization)は、 0または1の値をとる複数のバイナリ変数の二次式です。イジングモデルとは変数がスピン値かバイナリ値を取る点のみ異なります。QUBOモデルとイジングモデルは、そのグラフ形状(接続関係)を変えることなく等価に変換可能であり、また、変換によりサイズは変わりません。

＊9 順列と逆順列
　n要素の2つの順列P=[p(0),p(1),…,p(n-1)]とQ=[q(0),q(1),…,q(n-1)]が全てのiについてq(p(i))=i成り立つ時、つまり、逆関数の関係になっている時、PとQは互いに逆順列です。

＊10 無向グラフの巡回セールスマン問題を解くためのイジングモデル
　頂点数nの無向グラフの頂点集合を0,1,…,n-1とし，辺集合をE，辺(u,v)∈Eの長さをw_u,vとします。また、辺の長さの最大値より十分大きい値Mを決めます。one-hot表現による順列を格納するスピン変数の行列S=(s_i,j)としたとき、次の式により、順列が定める訪問順による巡回路長を評価することができます。

ただし、i'=(i+1)mod nとします。i番目に訪問する頂点がuであり、その次のi'番目に訪問する頂点がvである場合にs_i,u=s_i',v=+1となります。よって、(s_i,u+1)(s_i',v+1)=4が成り立ち、この式の評価値は、one-hot行列が表す訪問巡の巡回路長の4倍から4Mをマイナスした値となります。式の中で-Mは、辺がない場合の係数を0として、辺がある場合にその長さに応じて負の値とするためです。これにより辺のない訪問巡を選ばなくなります。この式を展開して得られるイジングモデルと順列を表す行列S=(s_i,j)を生成するイジングモデルをあわせることにより、最短巡回路を求めるためのイジングモデルを得ることができます。

＊11 部分順列
　n個の要素からm個(m＜n)の要素を重複なく選んで並べる順列。例えば、n=5，m=3のとき、0,1,2,3,4から3個選んで並べた[3,0,1]は部分順列です。