JP7341804B2 - 情報処理装置および情報処理方法 - Google Patents

Description

本発明は、情報処理技術、特に最適解の探索を行うアルゴリズムとしてアニーリングを採用する技術に関する。

最適解の探索を行う汎用アルゴリズムとしてアニーリングが知られている。アニーリングマシンはアニーリングを実行し近似最適解を出力する専用の装置である（例えば特許文献１～５、非特許文献２参照）。アニーリングマシンでは汎用的に問題を受け付け可能な計算モデルとして、イジングモデルを用いている。アニーリングマシンはイジングモデルのパラメータを入力として持つ。そのためアニーリングマシンのユーザは解きたい問題をイジングモデルに変換する必要がある。

図１Ａは、イジングモデルのエネルギーランドスケープの概念図を示す。横軸にスピンの状態（スピン配列）、縦軸にそのスピン配列におけるエネルギー関数をプロットしたものである。イジングモデルは、与えられたスピン配列、相互作用係数、および、外部磁場係数で構成される。イジングモデルのエネルギー関数Ｈ（σ_ｉ）（一般にハミルトニアンと呼ばれる）は図１Ａに示すとおりである。σ_ｉ、σ_ｊはそれぞれｉ番目とｊ番目のスピンの値で、通常上向きまたは下向きの２値をとる。Ｊ_ｉｊはｉ番目とｊ番目のスピンの間の相互作用係数、ｈ_ｉはｉ番目のスピンに対する外部磁場係数を表わす。Ｅ、Ｖは制約条件を示す。

確率的な遷移では、スピンの状態は現在の状態σの近傍のある状態σ’へ確率的に遷移を繰り返す。状態σから状態σ’へ遷移する確率を、遷移確率Ｐ（σ，σ’）と呼ぶ。遷移確率としてメトロポリス法やギブス法などが知られている（特許文献２、３、５参照）。一般に温度Ｔと呼ばれるパラメータにより、遷移確率を調整する。

近傍の状態σ’を生成する方法として、現在の状態σから一つのスピンの値を更新することが一般的である。更新するスピンを一つずつ順番に変えることで、スピン全体の探索を行う。温度Ｔを大きな値から徐々に減少させつつ遷移を実行するとき、エネルギーは局所解を経由して、エネルギーが最も低い状態（基底状態）に漸近的に収束する。これにより最小化問題の最適解または近似解を求める手法は、シミュレーティッド・アニーリング法（ＳＡ法）として知られる。

図１Ｂに、イジングモデルにおけるスピンの更新の概念を示す。イジングモデルとして問題をモデル化する場合、疎グラフと完全グラフの区別がある。図左側に示す疎グラフは、一つのノード（スピン）が限定された一部のノードとのみ接続されているモデルである。図中の上または下を向いた矢印がスピンの向きを示している。疎グラフでは同時に複数スピンを更新することができる（並列更新）。しかし、疎グラフでは結合されるノードの数に制限があるために、処理性能を向上させることが難しい。

一方、図右側に示す完全グラフは、一つのノードが他の全てのノードと接続されているモデルである。完全グラフでは、ノードは互いに相互作用（図中ｊの文字と添え字で示す）を及ぼしあっているので、複数スピンを同時に更新することができない。そこで、ＳＡ法では例えばノードを２つに分け、ひとつずつ更新する方式をとる（図右側中段。「単一更新（非並列）」）。しかし、これでは処理時間がかかるため、図１Ｂ中の式で示されるように、スピン更新は並列ではないが、瞬間磁場（mean field、cavity field、局所場などとも呼ばれる）を並列更新して、高速化する手法が提案されている（特許文献２～５、非特許文献２参照）。なお、図中ｈ_ｉとオーバーラインで示す瞬間磁場を、明細書等ではＩｈ_ｉで示す。

国際公開WO2015－132883 特開2019－16077号公報 特開2019－16129号公報 特開2019－71113号公報 特開2019－87273号公報

"Sampling from a Gibbs Measure with Pair Interaction by Means of PCA", Paolo Dai Pra et. Al., J Stat Phys (2012) 149:722－737 DOI 10.1007/s10955－012－0612－9 "An accelerator architecture for combinatorial optimization problems", Sanroku Tsukamoto et. Al., Fujitsu Sci. Tech. J., Vol. 53, No. 5 (September 2017)

図２に、発明者らが特許文献２～５、非特許文献２を検討して整理した、従来の最適解求める手順を示した。前提として、既に最適化問題は完全グラフの形でモデル化されているものとする。

処理Ｓ２０１はデータの初期化である。初期アニーリングステップと最大ステップ(ｔ_max)を設定し、スピンの更新を何回繰り返すかを設定する。逆温度β（＝１／Ｔ）を初期値β₀に設定する。モデルが含むＮ個のスピン（ノード）σのそれぞれに対してランダムに初期状態を設定する。スピンの値は２値であればよいが、ここではプラス１とマイナス１とする。各ノードに対して外部磁場係数ｈ_ｉを定め、各ノード間の相互作用係数ｊ_ｉｊを定める。瞬間磁場Ｉｈ_ｉ＝ｈ_ｉ＋Σ_{ｊ＝１～Ｎ}σ_ｊｊ_ｉｊを求める。スピンの反転を示す反転フラグＦと反転したスピンを指定する反転インデクスｊを初期値０に設定する。これらのデータは、モデルに基づいて定め、アニーリングマシンに設定される。アニーリングマシンの具体的な構成については、特許文献１～５、非特許文献２に記載がある。次に、処理Ｓ２０２～Ｓ２０４はスピンの更新フローを示している。

処理Ｓ２０２では、瞬間磁場Ｉｈ_ｉを更新する。反転フラグＦが１であれば、瞬間磁場Ｉｈ_ｉを各スピンで並列に計算して更新する。この方式では、スピンは一つずつ更新され、スピンの変化分だけを瞬間磁場に反映すればよいので、Ｉｈ_ｉ←Ｉｈ_ｉ＋２σ_ｊｊ_ｊｉとなる。

処理Ｓ２０３では、それぞれのスピンの更新確率ｐ_ｉを並列に計算して更新する。
ｐ_ｉ←sigmoid（－2σ_ｉＩｈ_ｉβ）
sigmoidはシグモイド関数である。

処理Ｓ２０４では、スピンの更新確率ｐ_ｉによって、各スピンの更新可否を決定する。更新可能なスピンがあれば、更新可能なスピンからランダムに一つσ_ｊを選択し、スピンを更新し(σ_ｊ←－σ_ｊ)、反転フラグを変更する(Ｆ←１)。

処理Ｓ２０５では、アニーリングステップを１進める（ｔ←ｔ＋１）。また、逆温度を更新する（β＝β（ｔ））。アニーリングでは、温度の逆数である逆温度は次第に大きくなっていく。

処理Ｓ２０２～Ｓ２０５は、判定処理Ｓ２０６でアニーリングステップｔがｔ_maxになるまで繰り返される。アニーリングステップｔがｔ_maxになると、スピン配列を読み出して、（疑似的に）最適化されたスピン配列として出力する（Ｓ２０７）。

上記の方式は、瞬間磁場の更新とスピンの更新確率の計算について並列的に処理ができ、高速化を図ることができる。しかし、スピンの更新については並列に行わず、スピン１つずつの順次処理となる。

そこで、アニーリングマシンにおいてスピンの更新を並列に行い得る方式の検討を行った。

本発明の好ましい一側面は、アニーリング制御部と、スピン相互作用メモリと、乱数生成部と、スピン状態更新部を備え、イジングモデルを用いて解を求める情報処理装置である。アニーリング制御部は、アニーリングのステップを制御するとともに、温度のパラメータと自己作用のパラメータを制御する。スピン相互作用メモリは、スピンの相互作用係数を記憶する。乱数生成部は、所定の乱数を生成する。スピン状態更新部は、複数のスピンの値を記憶するスピンバッファと、複数のスピンの瞬間磁場を計算する瞬間磁場計算部と、複数のスピンの更新確率を計算する確率計算部と、更新確率と乱数に基づいてスピンの値の更新を行うスピン状態決定部と、を備える。

本発明の好ましい他の一側面は、アニーリング制御部と、スピン相互作用メモリと、乱数生成部と、スピン状態更新部を備える情報処理装置により、イジングモデルを用いて解を求める情報処理方法である。アニーリング制御部は、アニーリングのステップを制御するとともに、温度のパラメータと自己作用のパラメータを更新する。スピン相互作用メモリは、スピンの相互作用係数を記憶する。乱数生成部は、所定の乱数を生成する。スピン状態更新部は、アニーリングのステップ毎に、相互作用係数を用いて複数のスピンの瞬間磁場を並列的に計算し、瞬間磁場を用いて複数のスピンの更新確率を並列的に計算し、更新確率と乱数を用いて複数のスピンの値を並列的に更新し、瞬間磁場の計算に用いる相互作用係数は、直前のアニーリングのステップで更新されたスピンに係る相互作用係数である。

アニーリングマシンにおいてスピンの更新を並列に行い得る。

アニーリングマシンの原理を示す概念図。 スピンの更新手法を説明する概念図。 スピンの更新を１つずつ行うアニーリングマシンの処理を示すフロー図。 実施例のアニーリングマシンのスピンの更新原理を示す概念図。 実施例１のアニーリングマシンの処理を示すフロー図。 実施例２のアニーリングマシンの処理を示すフロー図。 アニーリングマシンの構成例を示すブロック図。 アニーリング制御部の機能を説明するブロック図。 ステップ更新／アニーリング終了判定部の詳細ブロック図。 温度更新部と自己作用更新部の詳細ブロック図。 アニーリング制御部とＰＣＡステップ部の詳細を示すブロック図。 瞬間磁場計算部の機能ブロック図。 確率計算部の機能ブロック図。 スピン状態決定部の機能ブロック図。 スピン状態更新のタイムチャート。 反転スピンバッファおよび反転スピン選択回路のブロック図。 反転スピン選択回路の機能を説明するためのタイムチャート。 シグモイド関数によるスピン反転確率ｐを線形近似する説明図。 線形近似のグラフのスケーリング変換の説明図。 スケーリング変換した線形近似のグラフの符号変換の説明図。 スピン反転確率ｘとξの関係を示す説明図。 スピン反転確率ｐを使ったスピンの更新の処理フロー図。 スピン反転確率ｘを使ったスピンの更新の処理フロー図。 スピン反転確率ｐを用いたスピン状態決定部の機能ブロック図。 スピン反転確率ｘを用いたスピン状態決定部の機能ブロック図。 実施例の乱数の生成フロー図。 実施例の乱数生成の処理を説明する概念図。 実施例の乱数生成部のブロック図。 実施例３のシステムの概念図。 アニーリング制御部とＰＣＡステップ部の詳細を示すブロック図。 アニーリング制御部の機能を説明するブロック図。 瞬間磁場計算部を説明するブロック図。 マルチアニーラチップの概要を示す概念図。 アニーラチップのスピン相互作用メモリに格納される相互作用係数の行列を示す概念図。 マルチアニーラチップ間での変化スピン情報の共有シーケンスを説明する概念図。 マルチアニーラチップのアニーリング制御部とＰＣＡステップ部の詳細を示すブロック図。 マクロスピンの概念を示す説明図。 実施例４のアニーリングマシンの処理を示すフロー図。

実施の形態について、図面を用いて詳細に説明する。ただし、本発明は以下に示す実施の形態の記載内容に限定して解釈されるものではない。本発明の思想ないし趣旨から逸脱しない範囲で、その具体的構成を変更し得ることは当業者であれば容易に理解される。

以下に説明する発明の構成において、同一部分又は同様な機能を有する部分には同一の符号を異なる図面間で共通して用い、重複する説明は省略することがある。

同一あるいは同様な機能を有する要素が複数ある場合には、同一の符号に異なる添字を付して説明する場合がある。ただし、複数の要素を区別する必要がない場合には、添字を省略して説明する場合がある。

本明細書等における「第１」、「第２」、「第３」などの表記は、構成要素を識別するために付するものであり、必ずしも、数、順序、もしくはその内容を限定するものではない。また、構成要素の識別のための番号は文脈毎に用いられ、一つの文脈で用いた番号が、他の文脈で必ずしも同一の構成を示すとは限らない。また、ある番号で識別された構成要素が、他の番号で識別された構成要素の機能を兼ねることを妨げるものではない。

図面等において示す各構成の位置、大きさ、形状、範囲などは、発明の理解を容易にするため、実際の位置、大きさ、形状、範囲などを表していない場合がある。このため、本発明は、必ずしも、図面等に開示された位置、大きさ、形状、範囲などに限定されない。

本明細書で引用した刊行物、特許および特許出願は、そのまま本明細書の説明の一部を構成する。

本明細書において単数形で表される構成要素は、特段文脈で明らかに示されない限り、複数形を含むものとする。

図３は、以下で説明する実施例の原理を説明する概念図である。以下の実施例では、最適化問題は図３(A)のように完全グラフの形でモデル化される。そして、スピンの状態は図３(B)のように、各スピンを１対のコピーで置き換えることで並列的に更新される。スピンを並列更新するために、本実施例では新しいパラメータとして自己作用ｑＴ（ｔ）を導入する。

図４は、実施例１のアニーリングマシンがスピンの更新を並列に行う処理を示すフロー図である。前提として、既に最適化問題は完全グラフの形でモデル化されているものとする。

処理Ｓ４０１はデータの初期化である。初期アニーリングステップｔを０に設定する。最大ステップ(ｔ_max)を設定し、スピンの更新を何回繰り返すかを設定する。温度Ｔあるいは逆温度β（＝１／Ｔ）を、初期値Ｔ₀あるいはβ₀に設定する。モデルが含むＮ個のスピン（ノード）σのそれぞれに対してランダムに初期状態を設定する。スピンの値は２値であればよいが、ここではプラス１とマイナス１とする。各ノードに対して外部磁場係数ｈ_ｉを定め、各ノード間の相互作用係数ｊ_ｉｊを定める。本実施例ではパラメータとして自己作用ｑを導入し、自己作用を初期値ｑ ← ｑ_０に設定する。これらのデータは、モデルに基づいて定め、アニーリングマシンに設定される。アニーリングマシンの具体的な構成については、特許文献１～５、非特許文献２に記載がある。次に、処理Ｓ４０２～Ｓ４０４はスピンの更新フローを示している。

処理Ｓ４０２では、瞬間磁場を計算する。瞬間磁場Ｉｈ_ｉを各スピンで並列に計算する。この方式では、スピンは並列に更新されるので、すべてのスピンの変化を反映し、Ｉｈ_ｉ ← ｈ_ｉ＋Σ_{ｊ＝１～Ｎ}σ_ｊｊ_ｉｊとなる。

処理Ｓ４０３では、それぞれのスピンの更新確率ｐ_ｉを並列に計算して更新する。
ｐ_ｉ ← ρ（－βσ_ｉＩｈ_ｉ－ｑ）
ここでρはシグモイド関数を用いるのが一般的だが、他の関数でも可能である。

処理Ｓ４０４では、スピンの更新確率ｐ_ｉによって、各スピンの更新可否を決定し、並列的に更新する（σ_ｉ ← －σ_ｉ）。

処理Ｓ４０５では、アニーリングステップを１進める（ｔ ← ｔ＋１）。また、温度Ｔまたは逆温度βを更新する（（Ｔ ← Ｔ（ｔ））あるいは（β ← β（ｔ）））。また、自己作用ｑを更新する（ｑ ← ｑ（ｔ））。ここで、Ｔ（ｔ）、β（ｔ）、ｑ（ｔ）は、ユーザによって定義された関数とする。

Ｔ（ｔ）、β（ｔ）、ｑ（ｔ）の設定は、ユーザが任意に選択できる。一般に、温度Ｔはアニーリングステップが進むにつれて次第に小さく、逆数である逆温度βは次第に大きくなっていく。また、ｑ（ｔ）は原理的に、図３(B)に示した左側スピン群のｉ番目のスピンと右側スピン群のｉ番目のスピンとが、基底状態において同値となるように設定すればよい。すなわちアニーリングの最終ステップにおいて、ｑ（ｔ）が十分大きくなるような関数を選ぶ。一般に、最後のアニーリングのステップにおいて、ｑ（ｔ）は最大値を取る。

処理Ｓ４０２～Ｓ４０５は、判定処理Ｓ４０６でアニーリングステップｔがｔ_maxになるまで繰り返される。アニーリングステップｔがｔ_maxになると、スピン配列σ_ｉを読み出して、（疑似的に）最適化されたスピン配列として出力する（Ｓ４０７）。

スピンの更新に係る処理Ｓ４０２～Ｓ４０４の原理は、ＰＣＡ（Probabilistic Cellular Automaton）と呼ばれるアルゴリズムを応用したものである（非特許文献１参照）。非特許文献１では、パラメータとして自己作用ｑを導入し、パラメータｑは更新の１ステップにおけるスピン反転の平均個数を制御するとしている。ただし、非特許文献１は、ギブス測度（Gibbs Measure）に関する研究論文で、最適化問題やアニーリングに関するものではない。

イジングモデルにおいては通常、２つの結合スピンを同時に更新することはできない。これは、それらのエネルギー変化の推定値に誤差が生じるためである。よって、図３(A)に示した完全グラフを用いる場合、スピンの同時（並列）更新ができない。ＰＣＡは図３(B)に示したように、各スピンについてコピーを導入し、一対のスピン配列にするものと解釈できる。図３(B)右側のスピンの一つは左側のスピンの全てと接続されている。左側のスピンの一つは右側のスピンの全てと接続されている。一方、右側のスピン同士は接続されていない。また、左側のスピン同士は接続されていない。このため、右側のスピンは同時に更新することができる。また、左側のスピンは同時に更新することができる。

しかしながら、図３(B)の右側のスピン配列と左側のスピン配列は、異なる最適スピン配置をとることができるので、通常は図３(B)に示すように対応する左右のスピンは異なる。このため、パラメータｑを導入し、ｑを十に大きくすることにより、右側のスピン配列と左側のスピン配列が同一になるように制御する。この条件が満たされると、図３(A)に示した完全グラフを、図３(B)に示した一対のスピン配列に置き換えることが可能となり、図３(B)の一対のスピン配列を用いて並列更新を実行することができる。

実施例１では、基本的な概念を示した。この方式では、瞬間磁場の計算、スピンの更新確率の計算、スピンの更新が全て並列に実行可能である。ただし、瞬間磁場の計算は全てのノードの状態を反映しなければならない。また、シグモイド関数を用いる場合にはそれを計算する必要がある。次の実施例２では、さらに計算効率を考慮した具体的なインプリメンテーションについて説明する。

図５は、実施例２のアニーリングマシンがスピンの更新を並列に行う処理を示すフロー図である。前提として、既に最適化問題は完全グラフの形でモデル化されているものとする。

実施例２では、スピンの更新確率をｘ_ｉ＝σ_ｉＩｈ_ｉ＋ｑＴで計算する。再度図３(B)を用いて計算の例を説明する。図３(B)の左側の現在の状態σに基づいて、図３(B)の右側の更新されたスピンの状態σ′を計算するため、更新確率ｘ_ｉ＝σ_ｉＩｈ_ｉ＋ｑＴを計算する必要がある。図３(B)は、ｉ＝１のスピン（一番上のスピン）σ′_１を得るためのｘ_１の計算の一例を示している。いま単純化するために、５つのスピンで考え、ｈ_ｉ＝０，Ｊ_ｉｉ＝０とする。
Ｉｈ_ｉ＝ｈ_ｉ＋Σ_{ｊ＝１～５}σ_ｊｊ_ｉｊ
であるから、図３(B)のｘ_１を計算するため、
ｘ_１
＝σ_１Ｉｈ_１＋ｑＴ
＝ｑＴ＋σ_１Ｉｈ_１
＝ｑＴ＋σ_１（ｈ_１＋Σ_{ｊ＝１～５}σ_ｊｊ_１ｊ）
＝ｑＴ＋σ_１（σ_２Ｊ_２１＋σ_３ｊ_３１＋σ_４ｊ_４１＋σ_５ｊ_５１）
となる。ここで、スピンの値がσ_１＝σ_２＝σ_３＝σ_５＝１、σ_４＝－1とすると、
ｘ_１＝ｑＴ＋Ｊ_２１＋ｊ_３１－ｊ_４１＋ｊ_５１
となる。

＜１．処理の全体像＞
処理Ｓ５０１はデータの初期化である。アニーリングステップｔを初期値０に設定する。最大ステップ(ｔ_max)を設定し、スピンの更新を何回繰り返すかを設定する。温度２Ｔあるいは逆温度β（＝１／２Ｔ）を、初期値２Ｔ₀あるいは１／２Ｔ₀に設定する。モデルが含むＮ個のスピン（ノード）σのそれぞれに対してランダムに初期状態を設定する。スピンの値は２値であればよいが、ここではプラス１とマイナス１とする。各ノードに対して外部磁場係数ｈ_ｉを定め、各ノード間の相互作用係数ｊ_ｉｊを定める。図２と同様に瞬間磁場Ｉｈ_ｉ＝ｈ_ｉ＋Σ_{ｊ＝１～Ｎ}σ_ｊｊ_ｉｊを用意する。本実施例ではパラメータとして自己作用ｑＴを導入し、自己作用を初期値ｑＴ ← ｑＴ_０に設定する。各ノードに反転スピンフラグＦ_ｊを設定し、初期値Ｆ_ｊ＝０を設定する。これらのデータは、モデルに基づいて定め、アニーリングマシンに設定される。アニーリングマシンの具体的な構成については、特許文献１～５、非特許文献２に記載がある。次に、処理Ｓ５０２～Ｓ５０４はスピンの更新フローを示している。

実施例１との差異としては、温度がＴの２倍となっており、また自己作用ｑとして温度との積ｑＴを用いる。理由については後述する。

処理Ｓ５０２では、瞬間磁場を更新する。瞬間磁場Ｉｈ_ｉを各スピンで並列に計算して更新する。本方式では各スピンは並列的に更新されているが、瞬間磁場の計算では、反転スピンフラグＦ_ｊ＝１であるノード、すなわちスピンが更新されたノードの影響だけを反映する。すなわちＩｈ_ｉ ← Ｉｈ_ｉ＋２ｊ_ｉｊσ_ｊとなる。

処理Ｓ５０３では、それぞれのスピンの更新確率ｘ_ｉを並列に計算して更新する。
ｘ_ｉ＝ σ_ｉＩｈ_ｉ＋ｑＴ
ここでｘ_ｉにはシグモイド関数を線形近似した関数を用いる。

処理Ｓ５０４では、スピンの更新確率ｘ_ｉによって、各スピンの更新可否を決定して並列的に更新する。更新するスピンではσ_ｉ ← －σ_ｉとする。また更新されたスピン集合の反転スピンフラグを更新しＦ_ｉ ← １とする。
処理Ｓ５０５では、アニーリングステップを１進める（ｔ ← ｔ＋１）。また、温度または逆温度を更新する（（２Ｔ ← ２Ｔ（ｔ））あるいは（β／２ ← β／２（ｔ）））。温度２Ｔはステップが進むにつれて次第に小さく、逆数である逆温度β／２は次第に大きくなっていく。また、自己作用ｑＴを更新する（ｑＴ← ｑＴ（ｔ））。ここで、２Ｔ（ｔ）、β／２（ｔ）、ｑＴ（ｔ）は、ユーザによって定義された関数とする。自己作用ｑＴ（ｔ）の設定の条件は、実施例１と同様である。

処理Ｓ５０２～Ｓ５０５は、判定処理Ｓ５０６でアニーリングステップｔがｔ_maxになるまで繰り返される。アニーリングステップｔがｔ_maxになると、スピン配列σ_ｉを読み出して、（疑似的に）最適化されたスピン配列として出力する（Ｓ５０７）。

＜２．装置構成の全体像＞
図６は、図５の処理を実行するための情報処理装置の構成を示すブロック図である。情報処理装置は、制御装置６１０とアニーリングマシン６２０からなる。

＜２.１．制御装置＞
制御装置６１０は例えばパーソナルコンピュータ（ＰＣ）であり、ＰＣＡアニーラ制御部６１１と相互作用・瞬間磁場準備部６１２を含む。制御装置６１０は、アニーリングマシン６２０の上位装置として、これを制御する。また制御装置６１０は、初期化の処理Ｓ５０１を実行する。

ＰＣＡアニーラ制御部６１１は、アニーリングマシン６２０との間のデータの入出力や、アニーリングステップおよび各種パラメータの設定等、アニーリングマシン６２０の全体的な制御を行う。相互作用・瞬間磁場準備部６１２は、最適化問題をイジングモデル化し、相互作用係数ｊ_ｉｊ、外部磁場係数ｈ_ｉ、瞬間磁場Ｉｈ_ｉを計算して準備する。ＰＣＡアニーラ制御部６１１は、これらのデータをアニーリングマシン６２０に入力して設定する。

制御装置６１０は、通常のＰＣと同様に、入力装置、出力装置、処理装置、記憶装置を備える。本実施例では計算や制御等の機能は、記憶装置に格納されたプログラムが処理装置によって実行されることで、定められた処理を他のハードウェアと協働して実現される。計算機などが実行するプログラム、その機能、あるいはその機能を実現する手段を、「機能」、「手段」、「部」、「ユニット」、「モジュール」等と呼ぶ場合がある。

以上の構成は、単体のコンピュータで構成してもよいし、あるいは、入力装置、出力装置、処理装置、記憶装置の任意の部分が、ネットワークで接続された他のコンピュータで構成されてもよい。本実施例中、ソフトウエアで構成した機能と同等の機能は、ＦＰＧＡ（Field Programmable Gate Array）、ＡＳＩＣ（Application Specific Integrated Circuit）などのハードウェアでも実現できる。

＜２.２．アニーリングマシン＞
アニーリングマシン６２０は、Ｉ／Ｏインターフェース６２１とＰＣＡアニーラ６２２を含む。ＰＣＡアニーラ６２２は、メモリアクセスインターフェース６２３、アニーリング制御部７００、ＰＣＡステップ部８００を含む。アニーリングマシン６２０は、アニーリングを実行するアニーリング専用の回路である。

Ｉ／Ｏインターフェース６２１は、アニーリングマシン６２０と制御装置６１０との間のデータの送受信を司る。メモリアクセスインターフェース６２３は、アニーリングマシン６２０が含む各種のメモリに対して、データの記録および読み出しを行う。

＜３．ＰＣＡアニーラ＞
ＰＣＡアニーラ６２２は、Ｉ／Ｏインターフェース６２１を介して、制御装置６１０とデータを交換する。メモリアクセスインターフェース６２３は、Ｉ／Ｏインターフェース６２１から受け取ったデータをＰＣＡアニーラ６２２内のメモリに格納する。また、メモリからデータをＩ／Ｏインターフェース６２１経由で制御装置６１０へ読み出す。メモリへのデータの書き込みや読み出しは、例えばＳＲＡＭ（Static Random Access Memory）の技術を適用して構成することができる。このようなハードウェアについては、例えば特許文献１に記載がある。

アニーリング制御部７００は、ステップ更新／アニーリング終了判定部７１０、温度更新部７２０、自己作用更新部７３０、アニーリング制御レジスタ７４０を含む。アニーリング制御部７００は、アニーリングプロセスに伴う、スピンの更新のための各種パラメータの初期値設定や更新を行う。アニーリング制御部７００は、例えばマイクロコンピュータで構成することができる。アニーリング制御レジスタ７４０は、現在の温度２Ｔや自己作用ｑＴの値とそれらの変更パラメータ、現在のステップ数ｔと最大ステップ数ｔ_maxを格納している。温度更新部７２０や自己作用更新部７３０は、アニーリング制御レジスタ７４０の数値を用いて、温度や自己作用を更新する。

ＰＣＡステップ部８００は、スピン相互作用メモリ８１０、スピン状態更新部８２０、乱数生成部８３０を含む。スピン相互作用メモリ８１０は、相互作用係数ｊ_ｉｊを格納する。スピン状態更新部８２０は、Ｎ個のスピン更新ユニット８２１を備え、各スピン更新ユニット８２１は、瞬間磁場レジスタ８２２を備える。ＰＣＡステップ部８００は、Ｎ個のスピン更新ユニット８２１のスピンの値を並列的に更新する。また、Ｎ個の瞬間磁場レジスタ８２２の瞬間磁場の値を並列的に更新する。乱数生成部８３０は乱数を生成する。

＜３．アニーリング制御部＞
＜３．１．全体構成＞
図７Ａは、アニーリング制御部７００の機能を説明する図である。アニーリング制御部７００へは、スピン状態更新部８２０からの信号が入力される。スピン状態更新部８２０からスピン状態の更新完了を示す信号（ステップ終了信号）を得るとステップ数がインクリメントされ、ステップ更新／アニーリング終了判定部７１０、温度更新部７２０、自己作用更新部７３０により、アニーリング制御レジスタ７４０の内容が更新される。

ステップ更新／アニーリング終了判定部７１０で、アニーリングステップが最大値ｔ_maxになったかを判定して、最大値であれば、反転スピンバッファ８２６内のＮ個のスピン状態を出力する。Ｎ個のスピン状態は、メモリアクセスインターフェース６２３、Ｉ／Ｏインターフェース６２１経由で、最適化されたスピン配列として制御装置６１０に出力される。アニーリングステップが最大値でなければ、アニーリングステップｔ ← ｔ＋１とし、温度更新部７２０で温度２Ｔを、自己作用更新部７３０で自己作用ｑＴを更新し、アニーリング制御レジスタ７４０の値を更新する。更新された温度２Ｔや自己作用ｑＴは、スピン状態更新部８２０に出力される。

＜３．２．ステップ更新／アニーリング終了判定部＞
図７Ｂは、ステップ更新／アニーリング終了判定部７１０の詳細ブロック図である。ステップ更新／アニーリング終了判定部７１０は、スピン状態更新部８２０からステップ終了信号を得ると、アニーリング制御レジスタ７４０から現在のアニーリングステップ数ｔと、最大ステップ数ｔ_maxを読み込む。現在のアニーリングステップ数ｔに１加算する（ｔ ← ｔ＋１）。比較器７１２では、ｔ＋１＜ｔ_maxでなければアニーリングステップを終了し、ｔ＋１＜ｔ_maxであれば、アニーリング制御レジスタ７４０のアニーリングステップ数ｔを更新する（ｔ ← ｔ＋１）。

＜３．３．温度更新部および自己作用更新部＞
図７Ｃは、温度更新部７２０と自己作用更新部７３０の詳細ブロック図である。両者はパラメータの更新部として同様の構成をとるため、パラメータＸ更新部として説明する。図中（ｘ，ｘ_ｍｕｌｔ）はそれぞれ（２Ｔ，２Ｔ_ｍｕｌｔ）または（β／２，β／２_ｍｕｌｔ）、および（ｑＴ，ｑＴ_ｍｕｌｔ）、を表す。

スピン状態更新部８２０からステップ終了信号を得ると、アニーリング制御レジスタ７４０から現在のパラメータｘとパラメータの更新パラメータｘ_ｍｕｌｔを読み込む。パラメータｘとパラメータの更新パラメータｘ_ｍｕｌｔを積算し、ｘ ← ｘ＊ｘ_ｍｕｌｔを新たなパラメータとしてアニーリング制御レジスタ７４０に書き込む。

＜４．ＰＣＡステップ部＞
＜４．１．全体構成＞
図８Ａは、アニーリング制御部７００とＰＣＡステップ部８００の詳細を示すブロック図である。
スピン相互作用メモリ８１０には、Ｎ個のスピンに対応するノードのそれぞれに、相互作用係数ｊ_ｉｊが記憶されている。具体的なハードウェア構成は、例えばＳＲＡＭの技術を適用した特許文献１に詳しい。

アニーリング制御部７００のアニーリング制御レジスタ７４０から乱数生成部８３０へ温度２Ｔが入力され、－２Ｔ～２Ｔの範囲の値を持つ乱数が生成される。それぞれ異なる乱数が、スピン状態決定部８２５に入力される。

スピンの値は、スピン状態更新部８２０により更新される。スピン状態更新部８２０は、瞬間磁場計算部８２３、確率計算部８２４、スピン状態決定部８２５、反転スピンバッファ８２６、反転スピン選択回路８２７を含む。スピン状態更新部８２０へは、アニーリング制御部７００のアニーリング制御レジスタ７４０から自己作用ｑＴが入力される。また、乱数生成部８３０から乱数が入力される。また、スピン相互作用メモリから相互作用係数ｊ_ｉｊが入力される。これらに基づいて、各スピンの状態を更新する。

＜４．２．瞬間磁場計算部＞
図８Ｂは、瞬間磁場計算部８２３の説明のためのブロック図である。瞬間磁場計算部８２３は、入力として、スピン相互作用メモリ８１０の相互作用係数の値ｊ_ｉｊと、反転スピン選択回路のスピンの値σ_ｉを読み込む。相互作用係数の値ｊ_ｉｊは、反転スピン選択回路８２７で指定されたアドレスに対応しており、前のアニーリングステップで更新された１または複数のスピンσ_ｊに対応している。

プラス１あるいはマイナス１を持つスピンの値σ_ｉを符号拡張する。「＜＜１」は、ビットを１桁左にシフトすることを意味し、２を掛けることと同じである。これにより２σ_ｊｊ_ｉｊを計算する。瞬間磁場レジスタ８２２の現在の瞬間磁場Ｉｈ_ｉの値に２σ_ｊｊ_ｉｊを加算して、新たな瞬間磁場の値を瞬間磁場レジスタ８２２にＩｈ_ｉ ← Ｉｈ_ｉ＋２σ_ｊｊ_ｉｊとして書き込む。瞬間磁場Ｉｈ_ｉの値は、確率計算部８２４に入力される。

図５で説明したように、Ｉｈ_ｉ＋２σ_ｊｊ_ｉｊの計算において、ｊ∈{ｊ＝1,…,N｜F_j=1}である。反転スピン選択回路８２７は、前のアニーリングステップで更新されたスピンの反転フラグＦ_ｊの値を“１”に設定して（F_j=1）、更新されたスピンを指定している。よって、瞬間磁場計算部８２３は、アニーリングステップでスピンが更新（変更）されたスピンに基づいて、瞬間磁場を計算する。

＜４.３．確率計算部＞
図９Ａは、反転確率ｘ_ｉを計算するための、確率計算部８２４の詳細を説明する図である。確率計算部８２４は、入力として反転スピン選択回路のスピンの値σ_ｉ（瞬間磁場計算部８２３から引き継げばよい）、瞬間磁場計算部８２３から得た瞬間磁場Ｉｈ_ｉ、アニーリング制御部７００のアニーリング制御レジスタ７４０から得た自己作用ｑＴを得る。

プラス１あるいはマイナス１を持つスピンの値σ_ｉを符号拡張し、瞬間磁場Ｉｈ_ｉσ_ｉを計算し、自己作用ｑＴに加算する。これにより、反転確率ｘ_ｉ＝σ_ｉＩｈ_ｉ＋ｑＴを得る。

＜４.４．スピン状態決定部＞
図９Ｂは、スピンの更新を決定するための、スピン状態決定部８２５の詳細を説明する図である。スピン状態決定部８２５は、入力として確率計算部８２４から得た反転確率ｘ_ｉ、乱数生成部８３０から得た疑似乱数、反転スピン選択回路のスピンの値σ_ｉ（確率計算部８２４から引き継げばよい）を得る。乱数生成部８３０は、Ｎ個のスピン状態決定部８２５毎に異なる乱数を与える。乱数生成部８３０における乱数の生成は温度２Ｔで制御される。乱数がとる値は－２から２Ｔの範囲である。

スピン状態決定部８２５は、反転確率ｘ_ｉと乱数を比較し、反転確率ｘ_ｉが乱数より小さければ、スピンの値σ_ｉの符号を反転する。スピンを反転した場合、反転スピンバッファ８２６の反転フラグＦ_ｉを設定する。すでに述べたように、Ｆ_ｉはスピン状態決定部８２５でスピンσ_ｉが反転されたとき“１”が、スピン状態決定部８２５でスピンσ_ｉが反転されなかったとき“０”が、設定される。

＜４.５．スピン更新のタイムチャート＞
図９Ｃは、ｉ＝１のスピン状態更新のタイムチャートである。瞬間磁場計算部８２３では、前のアニーリングステップで更新されたスピンによる瞬間磁場の変化分を加算して、瞬間磁場Ｉｈ_ｉを更新（Ｉｈ_１ ← Ｉｈ_１＋２σ_ｊｊ_１ｊ）し、瞬間磁場レジスタ８２２に格納する。瞬間磁場を求めたら、アニーリング制御部７００から自己作用ｑＴを、乱数生成部８３０から乱数randを得る。確率計算部８２４が、自己作用ｑＴ、瞬間磁場Ｉｈ_１、およびスピンσ_１から反転確率ｘ_１を計算する。スピン状態決定部８２５が、反転確率ｘ_１、乱数rand、およびスピンσ_ｉからスピンの更新有無を決定し、反転スピンバッファ８２６の反転フラグＦ_１を設定しスピンσ_１を更新する。

＜４．６．反転スピンバッファおよび反転スピン選択回路＞
図１０Ａは、反転スピンバッファ８２６および反転スピン選択回路８２７の詳細を説明する図である。Ｎ個の各スピンに対応するスピン状態決定部８２５で生成された反転フラグＦ_ｉとスピンの値σ_ｉは（ｉはスピン（ノード）の番号を示す）、反転スピンバッファ８２６に並列に書き込まれる。反転フラグＦ_ｉは例えば、“１”または“０”を示す任意のビットの符号でよい。上述のように、反転フラグＦ_ｉはスピンσ_ｉが更新されたかどうかを示す。

反転スピン選択回路８２７のセレクタ１００１を介して、スピンの値σ_ｉを対応する瞬間磁場計算部８２３へ送る。また、反転フラグＦ_ｉの値が“１”であるスピンのアドレスをスピン相互作用メモリ３１０へ送る。スピン相互作用メモリ３１０は、アドレスで指定された相互作用ｊ_ｉｊを瞬間磁場計算部８２３へ送る。

プライオリティ回路１００２は、すべてのデータの送付が終了したら、アニーリング制御部７００へステップ終了信号を送信する。アニーリング制御部７００は、ステップ終了信号を受けると、次のステップを開始する。

＜４．７．反転フラグの例＞
図１０Ｂは、反転スピン選択回路８２７の機能を説明するための図である。ここでは、σ_１，σ_２，σ_３，σ_４のスピンが４個のみの系で、反転フラグの機能を説明する。また、σ_１，σ_２，σ_４はスピン状態決定部８２５でスピンが更新されており、σ_３はスピンが更新されていないものとする。すなわち、スピン状態決定部８２５は、反転フラグをＦ_１＝Ｆ_２＝Ｆ_４＝１，Ｆ_３＝０にして、反転スピンバッファ８２６に設定している。

反転スピン選択回路８２７が処理を開始するサイクル１で、反転フラグは“Ｆ_４Ｆ_３Ｆ_２Ｆ_１”＝“１０１１”に設定されている。プライオリティ回路１００２は、反転フラグ“１０１１”により、Ｆ_４＝１であり、σ_４のスピンが更新されていることを知る。よって、セレクタ１００１は、反転スピンバッファの更新されたスピンσ_４を瞬間磁場計算部８２３に送り、対応するアドレス“１１”をスピン相互作用メモリ８１０に送る。スピン相互作用メモリ８１０は、アドレス“１１”に対応するｊ_４１、ｊ_４２、ｊ_４３を瞬間磁場計算部８２３に送る。次にＦ_４＝０に設定し、反転フラグを“００１１”に変更する。

サイクル２で、反転フラグは“Ｆ_４Ｆ_３Ｆ_２Ｆ_１”＝“００１１”に変更されている。プライオリティ回路１００２は、反転フラグ“００１１”により、Ｆ_２＝１であり、σ_２のスピンが更新されていることを知る。よって、セレクタ１００１は、反転スピンバッファの更新されたスピンσ_２を瞬間磁場計算部８２３に送り、対応するアドレス“０１”をスピン相互作用メモリ８１０に送る。次にＦ_２＝０とし、反転フラグを“０００１”に変更する。

以下のサイクルも同様であり、最終的に反転フラグがオールゼロになれば、プライオリティ回路１００２は、処理が終了したことを判定してステップ終了フラグを生成する。そして、ステップ終了信号をアニーリング制御部７００に送る。

＜５．スピンの反転確率の計算＞
実施例１では、スピン更新確率ｐとしてシグモイド関数を用いるのが一般的であると述べた。しかし、シグモイド関数を用いた場合、計算の負荷が大きい。そこで、線形近似を用いて、回路規模や計算量を減少させることを検討した。

図１１Ａに、シグモイド関数を線形近似する概念を説明する。ここで、スピン更新確率ｐ_ｉを以下のシグモイド関数で規定する。
ｐ_ｉ＝sigmoid（－ｘ_ｉ／Ｔ）
ｘ_ｉ＝σ_ｉＩｈ_ｉ＋ｑＴ
もし、ｐ_ｉ＞rand’’∈（０，１）であれば、スピンを反転させ、σ_ｉ ← －σ_ｉとする。

ここでｐ_ｉを以下のように線形近似する。
ｘ_ｉ＜－２Ｔのとき、ｐ_ｉ≒１
－２Ｔ≦ｘ_ｉ≦２Ｔのとき、ｐ_ｉ≒１／２－ｘ_ｉ／４Ｔ
２Ｔ＜ｘ_ｉのとき、ｐ_ｉ≒０
もし、ｐ_ｉ＞rand’’であれば、スピンを反転させ、σ_ｉ ← －σ_ｉとする。

以上で、図１１Ａのグラフにおいて、点線（１）で示すシグモイド関数を、実線（２）で示すように線形近似することができる。

図１１Ｂは、図１１Ａの線形近似のグラフをスケーリング変換したものである。図１１Ａでは乱数rand’’のとる値は０から１の間である。これを、－２Ｔから２Ｔの間の値をとる乱数rand’にスケーリング変換する。この変換は以下のように実行できる。
rand’＝（rand’’－1／２）×４Ｔ
また、スピン更新確率ηは以下になる。
η_ｉ＝（ｐ_ｉ－1／２）×４Ｔ
線形近似により、
ｘ_ｉ＜－２Ｔのとき、η_ｉ≒２Ｔ
－２Ｔ≦ｘ_ｉ≦２Ｔのとき、η_ｉ≒－ｘ_ｉ
２Ｔ＜ｘ_ｉのとき、η_ｉ≒－２Ｔ
もし、η_ｉ＞rand’であれば、スピンを反転させ、σ_ｉ ← －σ_ｉとする。

図１１Ｃは、図１１Ｂのスケーリング変換したグラフを符号変換したものである。符号変換により、
rand＝－rand’
とする。
また、スピン更新確率ξは以下になる。
ξ_ｉ＝－η_ｉ
線形近似により、
ｘ_ｉ＜－２Ｔのとき、ξ_ｉ≒－２Ｔ
－２Ｔ≦ｘ_ｉ≦２Ｔのとき、ξ_ｉ≒ｘ_ｉ
２Ｔ＜ｘ_ｉのとき、ξ_ｉ≒２Ｔ
もし、ξ_ｉ＜randであれば、スピンを反転させ、σ_ｉ ← －σ_ｉとする。

図１１Ｄは、ξ_ｉとｘ_ｉを重ねて描いたものである。上に述べたように、－２Ｔ≦ｘ_ｉ≦２Ｔのとき、ξ_ｉ≒ｘ_ｉであるから、ξ_ｉの代わりにｘ_ｉを使っても、乱数の比較結果が変わらないことがわかる。そこで、もし、ｘ_ｉ＜randであれば、スピンを反転させ、σ_ｉ ← －σ_ｉとすることができる。

＜６．スピンの更新の処理フロー＞
図１２Ａは、図１１Ａのスピン反転確率ｐを使ったスピンの更新の処理フローを示している。反転確率を計算してｐ_ｉを得る（Ｓ１２０１ａ）。ｐ_ｉを０から１の間の値をとる乱数rand’’と比較する（Ｓ１２０２ａ）。乱数rand’’よりｐ_ｉが大きければ、スピンの状態を反転する（Ｓ１２０３ａ）。乱数rand’’よりｐ_ｉが大きくなければ、スピンの状態を維持する（Ｓ１２０４ａ）。

図１２Ｂは、図１１Ｄのスピン反転確率ｘを使ったスピンの更新の処理フローを示している。反転確率を計算してｘ_ｉを得る（Ｓ１２０１ｂ）。ｘ_ｉを－２Ｔから２Ｔの間の値をとる乱数randと比較する（Ｓ１２０２ｂ）。乱数randよりｘ_ｉが小さければ、スピンの状態を反転する（Ｓ１２０３ｂ）。乱数randよりｘ_ｉが小さくなければ、スピンの状態を維持する（Ｓ１２０４ｂ）。

＜７．スピンの更新の処理回路＞
図１２Ｃは、図１２Ａのシグモイド関数ｐを使ったスピン状態決定部８２５ａの構成を示している。入力は、乱数生成部８３０ａから得る０から１の間の値を取る乱数rand’’と、確率計算部８２４ａから得るｘ_ｉと、温度更新部７２０ａから得る逆温度β（ｔ）＝１／Ｔ（ｔ）である。ｐ_ｉ＝sigmoid（－ｘ_ｉ／Ｔ）を得るために、ｘ_ｉの符号を反転して積算し、シグモイド関数ｐ_ｉを計算する必要がある。

図１２Ｄは、図１２Ｂの線形近似とスケール変換した関数ｘを使ったスピン状態決定部８２５ｂの構成を示している。入力は、乱数生成部８３０ｂから得る－２Ｔから２Ｔの間の値を取る乱数randと、確率計算部８２４ｂから得るｘ_ｉである。この方式では、乱数randとｘ_ｉを比較するだけでよい。符号拡張および積算器は不要であり、シグモイド関数ｐ_ｉを計算する必要もない。

＜８．乱数の生成方法＞
図１３Ａは、乱数生成部８３０における、－２Ｔから２Ｔの間の値をとる乱数の生成フローである。

図１３Ｂは、乱数生成部８３０における、乱数生成の処理を説明する概念図である。乱数生成部８３０は、各スピン状態決定部８２５に異なる乱数を与える。ここでは、そのうちの一つの乱数の生成について説明するが、同様の構成がスピンの数Ｎだけ並列に動作すればよい。

まず公知の方法で所定ビット（たとえば３２ビット）の乱数を生成する（Ｓ１３０１）。図１３Ｂ(A)に生成した３２ビットの乱数（以下の説明で「乱数(A)」という）を示す。

次に、乱数生成部８３０は２Ｔマスクを作成する（Ｓ１３０２）。２Ｔマスクは、２Ｔの１の存在する最上位ビット以下のビットを１に変更したものである。たとえば、２Ｔが“01010”の場合、２Ｔマスクは“01111”である。

そして、２Ｔマスクと乱数(A)の論理積（AND）をとる（Ｓ１３０３）。図１３Ｂ(B)に、２Ｔマスクと生成した乱数のANDの結果（以下の説明で「AND結果(B)」という）を示す。２ＴマスクとANDをとることにより、乱数(A)は最大で２Ｔマスクの最上位ビットと同じ桁までが切り出される。しかし、AND結果(B)は、まだ２Ｔより大きい可能性がある。ここで、AND結果(B)の０から２Ｔまでの値をＰ１、２Ｔから２^ceil(log(2T)までの値をＰ２とする。ceilは、引数として与えた数以上の最小の整数を返す関数である。

次に、AND結果(B)と２Ｔを比較する（Ｓ１３０４）。AND結果が２Ｔより小さくない場合には、AND結果(B)と２Ｔマスクの排他的論理和（XOR）をとる。図１３Ｂ(C)に、AND結果(B)と２ＴマスクのXORの結果（以下の説明で「XOR結果(C)」という）を示す。２ＴマスクとXORをとることにより、AND結果(B)の最上位ビットが”1”だった場合には、”0”に変更されるので、XOR結果(C)は２Ｔより小さくなる。XOR結果(C)を符号拡張して、－２Ｔから２Ｔの間の値(D)を得る。図１３Ｂ(D)に、符号拡張した結果(D)を示す。

AND結果(B)が２Ｔより小さい場合には、AND結果(B)をそのまま符号拡張して－２Ｔから２Ｔの乱数を得る。

図１３Ｃは、上記の処理を実現する乱数生成部８３０の構成例である。疑似乱数列生成法の１つであるXOR-shift８３１を用いて乱数(A)を生成する（Ｓ１３０１）。任意に定めた２Ｔを入力する。マスク作成部８３２は、２Ｔに基づいて２Ｔマスクを作成する（Ｓ１３０２）。AND回路８３３で、２Ｔマスクと乱数(A)の論理積（AND）をとる（Ｓ１３０３）。

比較器８３５で、AND結果(B)と２Ｔを比較する（Ｓ１３０４）。比較結果はセレクタ８３６に入力される。セレクタ８３６は、AND結果(B)が２Ｔより小さい場合には、AND結果(B)を符号拡張器８３７に入力する。セレクタ８３６は、AND結果(B)が２Ｔより小さくない場合には、XOR回路８３４で生成された、AND結果(B)と２Ｔマスクの排他的論理和（XOR）を符号拡張器８３７に入力する。符号拡張器８３７は－２Ｔから２Ｔの乱数を生成する。

以上の説明したように、本実施例では並列処理により計算速度が高速化される。また、並列処理の個々の計算量も圧縮することにより、回路構成が簡単になり消費電力を低減するとともにさらに処理速度が高速化される。

実施例２では、図１０Ｂで説明したように、更新されたスピンの処理は、セレクタ１００１で１つずつシーケンシャルに実行されている。すなわち反転スピン選択回路８２７による、スピン状態の瞬間磁場計算部８２３への送信、スピン相互作用メモリ８１０の読み出しの処理はスピンの数に比例したサイクル数を必要とする。実施例３では、この部分の並列化により処理時間を圧縮する例を示す。

図１４は、実施例３のシステムの概念図である。スピン相互作用メモリ８１０ａ，８１０ｂとスピン状態更新部８２０ａ，８２０ｂのように、並列動作可能な複数の回路に分割する。スピン相互作用メモリ８１０ａは、ｉ＝１～Ｎ／２の番号のスピンの相互作用を格納し、スピン相互作用メモリ８１０ｂはＮ／２＋１～Ｎの番号のスピンの相互作用を格納する。スピン状態更新部８２０ａはｉ＝１～Ｎ／２の番号のスピンの更新を行い、スピン状態更新部８２０ｂはＮ／２＋１～Ｎの番号のスピンの更新を行う。

図１５は、アニーリング制御部７００－２とＰＣＡステップ部８００－２の詳細を示すブロック図である。基本的に図８Ａと同様だが、添字ａとｂで示す並列動作可能な２ブロックに分割されている。

図１６Ａは、アニーリング制御部７００－２の機能を説明する図である。基本的に図７Ａと同様だが、並列動作可能な２ブロックを調整する機能を含む。すなわち、一対のスピン状態更新部１（８２０ａ）とスピン状態更新部２（８２０ｂ）からのステップ終了信号の論理積を取り、両方の信号が得られてから、すなわち両方のアニーリングステップが終了してから、次のアニーリングステップに移るように構成されている。

図１６Ｂは、瞬間磁場計算部８２３－２の説明のためのブロック図である。並列的にスピン相互作用計算を行う。瞬間磁場計算部８２３－２は、入力として、反転スピン選択回路８２７ａ，８２７ｂでアドレスを指定されたスピン相互作用メモリ８１０ａ，８１０ｂの相互作用係数の値ｊ_ｋｉ,ｊ_ｍｉと、反転スピン選択回路８２７ａ，８２７ｂのスピンの値σ_ｋ、σ_ｍを並列に読み込む。プラス１あるいはマイナス１を持つスピンの値σ_ｋ、σ_ｍを符号拡張する。「＜＜１」は、ビットを１桁左にシフトすることを意味し、２を掛けることと同じである。これにより２σ_ｋｊ_ｉｋと２σ_ｍｊ_ｉｍを計算してこれらを加算する。瞬間磁場レジスタ８２２の現在の瞬間磁場Ｉｈ_ｉの値に２σ_ｋｊ_ｉｋ＋２σ_ｍｊ_ｉｍを加算して、新たな瞬間磁場の値を瞬間磁場レジスタ８２２にＩｈ_ｉ ← Ｉｈ_ｉ＋２σ_ｋｊ_ｉｋ＋２σ_ｍｊ_ｉｍとして書き込む。瞬間磁場Ｉｈ_ｉの値は、確率計算部８２４－２に入力される。このような構成により、２つのブロックで更新されたスピンの情報を瞬間磁場計算部８２３－２で得ることができる。

多数のスピンを含む大規模なイジングモデルの基底状態を探索可能なアニーリングマシンを構築するためには、単位素子をスピン数に応じた数だけ半導体チップに搭載する必要がある。このような半導体装置は、チップサイズが大きく、製造コストも高くなる。従って、このような半導体装置を実現するに際しては、ある程度の数の単位素子が搭載された半導体チップを複数接続して構築することが望ましい。実施例４は、実施例１で説明したアニーリングマシンを複数接続したものである。ここで、アニーリングマシンは半導体製造技術で作成された半導体チップであり、アニーラチップと呼ぶことにする。

図１７Ａは、実施例４のマルチアニーラチップの概要を示す概念図である。アニーラチップ１７００は、アニーリングマシン６２０の構成を基本的に含む半導体チップである。ここでは、Ｌ個のアニーラチップを接続した例を示している。各アニーラチップは、Ｎ個のスピンをＬ分割して処理を担当する。例えば、アニーラチップ１７００ａはσ_１～σ_Ｎ／Ｌのスピンを担当しアニーラチップ１７００ｂはσ_{Ｎ／Ｌ＋１}～σ_２Ｎ／Ｌのスピンを担当する。各マルチアニーラチップは、データを循環して転送可能な構成になっている。

図１７Ｂは、例としてＬ＝４の場合の各アニーラチップのスピン相互作用メモリ８１０に格納される相互作用係数の行列を示す概念図である。

マルチアニーラチップの構成では、チップ相互でスピン状態σ_ｉと更新されたスピンの情報を共有する必要がある。このため、これらの情報は、アニーラチップ１７００間で送信データとして送られる。送信データの内容は、更新されたスピンの値σ_ｉと相互作用係数ｊ_ｉｊである。

図１７Ｃは、Ｌ＝４の場合の各アニーラチップにおける処理の流れを説明する図である。図５のアニーリングステップの、最初の更新確率計算Ｓ５０３とスピンの更新Ｓ５０４で、更新されるスピンとその値が各アニーラチップで計算される（最初の段階では瞬間磁場は初期値である）。瞬間磁場の更新処理Ｓ５０２のためには、更新されたスピンの相互作用係数とスピンの値を他のアニーラチップに知らせなければならない。

図１７Ｃでは、ステップ１で、各アニーラチップ０～３では、自分のチップで更新されたスピンに関する情報（変化スピン０～３）を得る。続くステップ２～４において、変化スピン０～３を順次隣接するアニーラチップに転送する。ステップ４で、すべての更新されたスピンに関する情報が、各アニーラチップに行き渡るので、瞬間磁場の更新処理Ｓ５０２が可能となる。

図１８は、マルチアニーラチップ構成における、アニーリング制御部７００とＰＣＡステップ部８００の詳細を示すブロック図である。図８Ａとほぼ同様な構成であるが、ＰＣＡステップ部８００は、反転スピンキュー１８０１を備える。反転スピンキュー１８０１は、図１７Ｃのステップ１では、自チップの変化スピンの情報を隣接チップに転送する。ステップ２以降ではステップごと、隣接するアニーラチップから送られてくる変化スピンの情報を記憶し、さらに隣接するアニーラ―チップに転送する。アニーラチップがＮ個あれば、Ｎステップ後に瞬間磁場の更新処理Ｓ５０２が可能となる。

この実施例では、スピンに代えてマクロスピンを用いた例を示す。マクロスピンは、個々のスピンの値に代えて、複数のスピンの値の多数派の値を適用するものである。

図１９Ａに、マクロスピンの概念を示す。Ｎ個のスピンがＭ個のグループに分割されている。ｉ番目のスピンの値がσ_ｉで表される。各グループに属する複数のスピンの値の多数決によって、ｍ番目のグループのマクロスピンσ’_ｍが決定される。マクロスピンσ’_ｍは、ｍ番目のグループに属するスピン全ての値として適用される。

例えば、あるグループに属するスピンの数が６だとする。６個のスピンの値が1,1,－1,1,－1,1であるとすると、”1”が４個、”－1”が２個なので、このグループのマクロスピンは”1”ということになる。当該グループのスピンの値は全て”1”として扱われる。マクロスピンを用いることにより、計算量を圧縮することができる。

図１９Ｂは、実施例５のアニーリングマシンがマクロスピンを用いて行う処理を示すフロー図である。前提として、既に最適化問題は完全グラフの形でモデル化されているものとする。基本的な処理は図５に示した実施例と同様なので、特に差異について説明する。

初期化処理Ｓ１９０１は、図５の初期化処理Ｓ５０１と同様である。ただし、初期化処理Ｓ１９０１では、初期マクロスピンサイズＭ＝Ｍ_ｉｎｉｔ（２≦Ｍ_ｉｎｉｔ≦Ｎ）を設定する。

マクロスピンの状態決定処理Ｓ１９０２では、上述の多数決によりマクロスピンの値を決定する。例えば、グループｍ内のスピンの値σ_ｉの合計が０より大きければ、マクロスピンσ’_ｍ＝１、そうでなければσ’_ｍ＝－1とする。

スピンの更新処理Ｓ１９０３では、マクロスピンを用いてそれぞれのスピンの瞬間磁場Ｉｈと更新確率ｘを計算しスピンを更新する。スピンの代わりにマクロスピンを用いること以外は、図５の実施例と同様でよい。

処理Ｓ１９０４では、アニーリングパラメータを更新する。更新の内容は図５の実施例と同様でよい。マクロスピンサイズＭも更新する。マクロスピンサイズＭは、例えば徐々に大きくしてＮに近づけていく。アニーリングステップの最後ではＭ＝Ｎとして、マクロスピンは用いなくすることが望ましい。

判定処理Ｓ１９０５でアニーリングステップｔがｔ_maxになるまで繰り返される。アニーリングステップｔがｔ_maxになると、スピン配列σ_ｉを読み出して、（疑似的に）最適化されたスピン配列として出力する（Ｓ１９０６）。

本発明は上記した実施形態に限定されるものではなく、様々な変形例が含まれる。例えば、ある実施例の構成の一部を他の実施例の構成に置き換えることが可能であり、また、ある実施例の構成に他の実施例の構成を加えることが可能である。また、各実施例の構成の一部について、他の実施例の構成の追加・削除・置換をすることが可能である。

制御装置６１０、アニーリングマシン６２０、ＰＣＡアニーラ制御部６１１、相互作用・瞬間磁場準備部６１２、アニーリング制御部７００、ステップ更新／アニーリング終了判定部７１０、温度更新部７２０、自己作用更新部７３０、アニーリング制御レジスタ７４０、ＰＣＡステップ部８００、スピン相互作用メモリ８１０、スピン状態更新部８２０、乱数生成部８３０

Claims (15)

1. アニーリング制御部と、スピン相互作用メモリと、乱数生成部と、スピン状態更新部を備え、イジングモデルを用いて解を求める情報処理装置であって、
   前記アニーリング制御部は、アニーリングのステップを制御するとともに、温度のパラメータと自己作用のパラメータを制御し、
   前記スピン相互作用メモリは、スピンの相互作用係数を記憶し、
   前記乱数生成部は、所定の乱数を生成し、
   前記スピン状態更新部は、
   複数のスピンの値を記憶するスピンバッファと、
   複数のスピンの瞬間磁場を計算する瞬間磁場計算部と、
   複数のスピンの更新確率を計算する確率計算部と、
   前記更新確率と前記乱数に基づいてスピンの値の更新を行うスピン状態決定部と、
   を備える、情報処理装置。
2. 前記スピン状態更新部は、
   複数のスピンの瞬間磁場を並列的に計算する複数の瞬間磁場計算部と、
   複数のスピンの更新確率を並列的に計算する複数の確率計算部と、
   前記更新確率と前記乱数に基づいて複数のスピンの値の更新を並列的に行う複数のスピン状態決定部と、
   を備える、
   請求項１記載の情報処理装置。
3. アニーリングのステップの開始前において、前記瞬間磁場は、前記スピンの値と、外部磁場係数と、前記相互作用係数とに基づいて、初期設定され、
   前記瞬間磁場計算部は、
   前記アニーリングのステップの進行に伴って前記瞬間磁場の値を更新する、
   請求項２記載の情報処理装置。
4. 前記瞬間磁場計算部は、
   瞬間磁場の値を格納する瞬間磁場レジスタを備え、
   前記アニーリングのステップの進行に伴って前記瞬間磁場の値を更新する際に、直前のアニーリングのステップで変更されたスピンに関し、スピンの値と当該スピンに関する相互作用係数と、前記瞬間磁場レジスタに格納されている更新前の瞬間磁場の値を用いて、当該スピンの更新後の瞬間磁場を計算する、
   請求項３記載の情報処理装置。
5. 反転スピンバッファと反転スピン選択回路を備え、
   前記スピン状態決定部は、
   前記アニーリングのステップの進行に伴って、スピンの値の更新を行う際に、更新したスピンの値と更新したスピンを特定する反転フラグを前記反転スピンバッファに格納し、
   前記反転スピン選択回路は、
   前記更新したスピンの値を前記瞬間磁場計算部に送り、前記反転フラグに基づいて更新したスピンを特定する情報を前記スピン相互作用メモリに送り、
   前記スピン相互作用メモリは、
   前記更新したスピンを特定する情報に基づいて、当該スピンに関連する相互作用係数を、前記瞬間磁場計算部に送る、
   請求項４記載の情報処理装置。
6. 前記スピン状態決定部は、
   所定の関数と前記所定の乱数とを比較して、比較結果に基づいて前記スピンを反転させて更新するか、あるいは、スピンの状態を保持するかを決定する、
   請求項２記載の情報処理装置。
7. 前記所定の関数は、シグモイド関数を線形近似した関数である、
   請求項６記載の情報処理装置。
8. 前記アニーリング制御部は、
   前記温度のパラメータとして２Ｔの値を制御し、前記自己作用のパラメータとしてｑＴの値を制御し、
   前記乱数生成部は、
   －２Ｔから２Ｔの間の値を持つ乱数を生成し、
   前記所定の関数は、スピンの値と、瞬間磁場と、ｑＴを引数とする、
   請求項６記載の情報処理装置。
9. 前記自己作用のパラメータは、最後のアニーリングのステップで最大値を取る、
   請求項１記載の情報処理装置。
10. 前記反転スピン選択回路は、
    所定のアニーリングのステップにおいて、全ての前記更新したスピンの値を前記瞬間磁場計算部に送り、前記反転フラグに基づいて更新した全てのスピンを特定する情報を前記スピン相互作用メモリに送ると、ステップ終了信号を前記アニーリング制御部に送り、
    前記アニーリング制御部は、
    前記ステップ終了信号を受けると、アニーリングのステップを１進めるとともに、温度のパラメータと自己作用のパラメータを更新する、
    請求項５記載の情報処理装置。
11. 前記スピン相互作用メモリとして、第１のスピン相互作用メモリと第２のスピン相互作用メモリを備え、
    前記スピン状態更新部として、第１のスピン状態更新部と第２のスピン状態更新部を備え、
    前記第１のスピン状態更新部は、第１の瞬間磁場計算部と第１の反転スピンバッファと第１の反転スピン選択回路を備え、
    前記第２のスピン状態更新部は、第２の瞬間磁場計算部と第２の反転スピンバッファと第２の反転スピン選択回路を備え、
    前記第１の瞬間磁場計算部は、前記第１のスピン相互作用メモリと前記第２のスピン相互作用メモリから相互作用係数を受信し、前記第１の反転スピンバッファと前記第２の反転スピンバッファからスピンの値を受信し、
    前記第２の瞬間磁場計算部は、前記第１のスピン相互作用メモリと前記第２のスピン相互作用メモリから相互作用係数を受信し、前記第１の反転スピンバッファと前記第２の反転スピンバッファからスピンの値を受信し、
    前記アニーリング制御部は、
    前記第１の反転スピン選択回路からの第１のステップ終了信号と、前記第２の反転スピン選択回路からの第２のステップ終了信号との両方を受け取ると、アニーリングのステップを１進めるとともに、温度のパラメータと自己作用のパラメータを更新する、
    請求項１０記載の情報処理装置。
12. 請求項１記載の情報処理装置を実装した半導体チップを複数備え、複数の前記半導体チップの間で、スピン相互作用メモリの相互作用係数とスピンの値を共有する、
    情報処理装置。
13. 複数のスピンをグルーピングし、１つのグループの中のスピンは、当該グループのスピンの値の多数決で定まる値に統一する、
    請求項１記載の情報処理装置。
14. アニーリング制御部と、スピン相互作用メモリと、乱数生成部と、スピン状態更新部を備える情報処理装置により、イジングモデルを用いて解を求める情報処理方法であって、
    前記アニーリング制御部は、アニーリングのステップを制御するとともに、温度のパラメータと自己作用のパラメータを更新し、
    前記スピン相互作用メモリは、スピンの相互作用係数を記憶し、
    前記乱数生成部は、所定の乱数を生成し、
    前記スピン状態更新部は、前記アニーリングのステップ毎に、
    前記相互作用係数を用いて複数のスピンの瞬間磁場を並列的に計算し、
    前記瞬間磁場を用いて複数のスピンの更新確率を並列的に計算し、
    前記更新確率と前記乱数を用いて複数のスピンの値を並列的に更新し、
    前記瞬間磁場の計算に用いる前記相互作用係数は、直前のアニーリングのステップで更新されたスピンに係る相互作用係数である、
    情報処理方法。
15. 前記温度のパラメータを２Ｔとし、
    前記自己作用のパラメータをｑＴとし、
    前記所定の乱数は－２Ｔから２Ｔの間の値を取り、
    前記更新確率は前記瞬間磁場の値あるいは前記瞬間磁場の符号を反転した値に前記ｑＴを加えた値である、
    請求項１４記載の情報処理方法。

Images