特許 7273288 サンプリング装置及びサンプリング装置の制御方法

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】

(24)【登録日】2023-05-02

(45)【発行日】2023-05-15

(54)【発明の名称】サンプリング装置及びサンプリング装置の制御方法

(51)【国際特許分類】

   G06N 99/00 20190101AFI20230508BHJP

【ＦＩ】

G06N99/00 180

【請求項の数】 7

(21)【出願番号】P 2019031685

(22)【出願日】2019-02-25

(65)【公開番号】P2020135727

(43)【公開日】2020-08-31

【審査請求日】2021-11-09

(73)【特許権者】

【識別番号】000005223

【氏名又は名称】富士通株式会社

(74)【代理人】

【識別番号】110002918

【氏名又は名称】弁理士法人扶桑国際特許事務所

(72)【発明者】

【氏名】土手 暁

(72)【発明者】

【氏名】田村 泰孝

【審査官】金田 孝之

(56)【参考文献】

【文献】特開２０１８－０６３６２６（ＪＰ，Ａ）

【文献】特開平０６－１４９８６６（ＪＰ，Ａ）

【文献】特開２０１８－２０６１２７（ＪＰ，Ａ）

(58)【調査した分野】(Int.Cl.，ＤＢ名)

Ｇ０６Ｎ ９９／００

(57)【特許請求の範囲】

【請求項1】

エネルギーを表す評価関数に含まれる複数の状態変数の値をそれぞれ保持するとともに、所定試行回数毎に前記複数の状態変数の値を出力する状態保持部と、
前記複数の状態変数の値の何れかが変化することに応じて状態遷移が起こる場合、更新インデックス値に基づいて選択された重み値に基づき、エネルギーの変化値を各状態遷移に対して計算するエネルギー変化計算部と、
複数の前記エネルギーの変化値に基づいて前記各状態遷移のうちの少なくとも１つが許容されるようにオフセット値を決定し、複数の前記エネルギーの変化値のそれぞれに前記オフセット値を付加した複数の評価値を出力するとともに、前記所定試行回数毎に前記オフセット値を出力するオフセット制御部と、
温度を示す温度値と乱数値とに基づいて決まる閾値と、前記複数の評価値のそれぞれとを比較した結果に基づき、前記各状態遷移について状態遷移を許容するか否かを示す複数のフラグ値を出力する比較部と、
前記複数のフラグ値のうち、状態遷移を許容することを示すフラグ値の数を計数した計数値を、前記所定試行回数毎に出力するフラグ計数部と、
前記複数のフラグ値、または複数の前記エネルギーの変化値に基づいて、前記各状態遷移のうちの１つに対応するインデックス値を前記更新インデックス値として出力する選択部と、
前記オフセット制御部が出力した前記オフセット値と、前記フラグ計数部が出力した前記計数値に基づき、マルコフ連鎖モンテカルロ法にしたがった状態遷移における１つの状態に滞在する試行回数の期待値を計算する試行回数計算部と、
を有するサンプリング装置。

【請求項2】

前記オフセット制御部は、複数の前記エネルギーの変化値のうちの最小値が０以上である場合、前記最小値を前記オフセット値として決定し、複数の前記エネルギーの変化値のそれぞれから前記最小値を差し引くことで、前記複数の評価値を算出する、請求項１に記載のサンプリング装置。

【請求項3】

前記選択部は、前記複数のフラグ値のそれぞれに、それぞれ独立な複数の乱数値の何れかを加算した加算結果の大小関係に基づき、前記複数のフラグ値の何れかを選択する、請求項１または２に記載のサンプリング装置。

【請求項4】

前記複数のフラグ値に基づき、複数の前記エネルギーの変化値のうち、許容されない状態遷移であるとされている複数の状態遷移に対する複数の非許容エネルギー変化値を出力するエネルギー変化置き換え部と、
前記複数の非許容エネルギー変化値に基づき前記複数の状態遷移のうちの少なくとも１つが許容されるように前記オフセット値とは別のオフセット値を決定し、前記複数の非許容エネルギー変化値のそれぞれに前記別のオフセット値を付加した、前記複数の評価値とは別の複数の評価値を出力するとともに、前記所定試行回数毎に前記別のオフセット値を出力する他のオフセット制御部と、
前記閾値と、前記別の複数の評価値のそれぞれとを比較した結果に基づき、前記複数の状態遷移のそれぞれについて状態遷移を許容するか否かを示す、前記複数のフラグ値とは別の複数のフラグ値を出力する他の比較部と、
前記別の複数のフラグ値のうち、状態遷移を許容することを示すフラグ値の数を計数した、前記計数値とは別の計数値を、前記所定試行回数毎に出力する他のフラグ計数部と、
前記別の計数値、前記オフセット値及び前記別のオフセット値に基づいて、前記計数値を補正する計数値補正部と、
をさらに有する請求項１に記載のサンプリング装置。

【請求項5】

前記試行回数計算部は、補正された前記計数値と、前記オフセット値と、前記別のオフセット値と、に基づいて、前記期待値を計算する請求項４に記載のサンプリング装置。

【請求項6】

前記選択部は、複数の前記エネルギーの変化値のうち、０より小さい変化値を０に更新させた複数の更新エネルギー変化値を生成し、前記複数の更新エネルギー変化値のそれぞれに、前記温度値と互いに独立な複数の乱数値とに基づいてそれぞれ算出される複数の第２の閾値の何れかを加算した複数の第２の評価値のうちで、最小となる第２の評価値に対応する前記インデックス値を前記更新インデックス値として出力する、請求項１乃至５の何れか一項に記載のサンプリング装置。

【請求項7】

サンプリング装置の制御方法において、
前記サンプリング装置が有する状態保持部が、エネルギーを表す評価関数に含まれる複数の状態変数の値をそれぞれ保持するとともに、所定試行回数毎に前記複数の状態変数の値を出力し、
前記サンプリング装置が有するエネルギー変化計算部が、前記複数の状態変数の値の何れかが変化することに応じて状態遷移が起こる場合、更新インデックス値に基づいて選択された重み値に基づき、エネルギーの変化値を各状態遷移に対して計算し、
前記サンプリング装置が有するオフセット制御部が、複数の前記エネルギーの変化値に基づいて前記各状態遷移のうちの少なくとも１つが許容されるようにオフセット値を決定し、複数の前記エネルギーの変化値のそれぞれに前記オフセット値を付加した複数の評価値を出力するとともに、前記所定試行回数毎に前記オフセット値を出力し、
前記サンプリング装置が有する比較部が、温度を示す温度値と乱数値とに基づいて決まる閾値と、前記複数の評価値のそれぞれとを比較した結果に基づき、前記各状態遷移について状態遷移を許容するか否かを示す複数のフラグ値を出力し、
前記サンプリング装置が有するフラグ計数部が、前記複数のフラグ値のうち、状態遷移を許容することを示すフラグ値の数を計数した計数値を、前記所定試行回数毎に出力し、
前記サンプリング装置が有する選択部が、前記複数のフラグ値、または複数の前記エネルギーの変化値に基づいて、前記各状態遷移のうちの１つに対応するインデックス値を前記更新インデックス値として出力する、
前記サンプリング装置が有する試行回数計算部が、前記オフセット制御部が出力した前記オフセット値と、前記フラグ計数部が出力した前記計数値に基づき、マルコフ連鎖モンテカルロ法にしたがった状態遷移における１つの状態に滞在する試行回数の期待値を計算する、
サンプリング装置の制御方法。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、サンプリング装置及びサンプリング装置の制御方法に関する。

【背景技術】

【0002】

組合せ最適化問題を解くための手法として、組合せ最適化問題を磁性体のスピンの振る舞いを表すイジングモデルに変換し、マルコフ連鎖モンテカルロ法を用いて、イジングモデルの状態をエネルギーが低い状態に遷移させていく手法がある。以下、マルコフ連鎖モンテカルロ法を、ＭＣＭＣ（Markov-Chain Monte Carlo）法と略す。イジングモデルの状態は、複数の状態変数の値の組合せにより表現でき、状態変数の数がＮ個の場合、ｘ＝（ｘ_１，ｘ_２，…，ｘ_Ｎ）などと表せる。各状態変数の値として、０または１を用いることができる。

【0003】

イジングモデルのエネルギーを表すイジング型のエネルギー関数Ｅ（ｘ）は、たとえば、以下の式（１）で定義される。

【0004】

【数1】

【0005】

右辺の１項目は、イジングモデルの全状態変数の全組合せについて、漏れと重複なく、２つの状態変数の値（０または１）と重み値（２つの状態変数の間の相互作用の強さを表す）との積を積算したものである。ｘ_ｉは、識別情報（以下インデックス値という）がｉの状態変数、ｘ_ｊは、インデックス値＝ｊの状態変数であり、Ｗ_ｉｊは、インデックス値＝ｉ，ｊの状態変数間の相互作用の大きさを示す重み値である。

【0006】

右辺の２項目は、各インデックス値についてのバイアス係数と状態変数との積の総和を求めたものである。ｂ_ｉは、インデックス値＝ｉについてのバイアス係数を示している。
また、状態変数ｘ_ｉが変化して１－ｘ_ｉとなると、状態変数ｘ_ｉの増加分は、Δｘ_ｉ＝（１－ｘ_ｉ）－ｘ_ｉ＝１－２ｘ_ｉと表せる。スピン反転（状態変数の値の変化）に伴うエネルギー変化ΔＥ_ｉは、以下の式（２）で表される。

【0007】

【数2】

【0008】

式（２）において、状態変数ｘ_ｉが１から０に変化するとき、Δｘ_ｉは－１となり、状態変数ｘ_ｉが０から１に変化するとき、Δｘ_ｉは１となる。なお、ｈ_ｉはローカルフィールド（局所場）と呼ばれ、Δｘ_ｉに応じてローカルフィールドｈ_ｉに符号（＋１または－１）を乗じたものがエネルギー変化ΔＥ_ｉとなる。

【0009】

ある状態遷移に伴うエネルギー変化ΔＥに対するその状態遷移の許容確率として、以下の式（３）で表されるメトロポリス法またはギブス法で規定される許容確率Ａ（ΔＥ）を用いることができる。

【0010】

【数3】

【0011】

式（３）においてβは、逆温度（温度を表す温度値の逆数）である。エネルギーが増加する状態遷移についても確率的に許容される。
通常のＭＣＭＣ法では、ランダムまたはインデックス値の順（シーケンシャル）に、状態変数が選択され、その状態変数の値が変化する状態遷移に伴うエネルギー変化ΔＥに基づいて、上記の許容確率Ａ（ΔＥ）でその状態遷移を許容する。そしてその状態遷移を許容する場合には状態変数の値が更新される。このような処理が所定の試行回数、繰り返される。また、最低エネルギーとなる状態（最適解）を探索するために、温度を徐々に下げていくシミュレーテッド・アニーリング法が用いられる場合もある。

【0012】

一方、組合せ最適化問題を高速に解くためのハードウェアとして、デジタル回路を用いた最適化装置がある（たとえば、特許文献１参照）。最適化装置は、以下に示すように複数の状態遷移を同時に遷移候補として、１つの状態遷移を選択する並列探索を行う。

【0013】

最適化装置は、各状態遷移に伴うエネルギー変化と温度値に基づいて、上記の許容確率Ａ（ΔＥ）で各状態遷移を許容する。そして、最適化装置は複数の状態遷移の中から、許容された状態遷移を優先的に１つ選択し状態を更新する。最適化装置は、以上のような処理を、所定の試行回数、繰り返し行う。また、従来の最適化装置では、状態が局所解に拘束されるなどによって状態遷移が生じない場合、エネルギー変化にオフセット値を付加していた。

【先行技術文献】

【特許文献】

【0014】

【文献】特開２０１８－６３６２６号公報

【発明の概要】

【発明が解決しようとする課題】

【0015】

ところで、平衡状態では各状態の占有確率を示す確率分布はボルツマン分布となる。そのため、ＭＣＭＣ法を用いて温度一定とした状態遷移を繰り返す過程により得られた状態または状態に基づいた値をサンプルとして出力することで、ボルツマン分布にしたがうサンプルを発生するサンプラーが実現できる。発生したサンプルは、たとえば、機械学習などでの期待値の計算に用いられる。

【0016】

しかし、従来のデジタル回路を用いた最適化装置は、前述のような並列探索処理及びオフセット値の付加による計算処理の高速化を行っているため、確率分布がボルツマン分布からずれる可能性がある。従来のデジタル回路を用いた最適化装置をサンプラー（以下サンプリング装置という）として用いる場合、確率分布のボルツマン分布からのずれによって、サンプリング精度が通常のＭＣＭＣ法を用いた場合よりも悪化してしまう可能性があった。

【0017】

１つの側面では、本発明は、計算処理の高速性を損なうことなく高精度のサンプリングが可能なサンプリング装置及びサンプリング装置の制御方法を提供することを目的とする。

【課題を解決するための手段】

【0018】

１つの実施態様では、エネルギーを表す評価関数に含まれる複数の状態変数の値をそれぞれ保持するとともに、所定試行回数毎に前記複数の状態変数の値を出力する状態保持部と、前記複数の状態変数の値の何れかが変化することに応じて状態遷移が起こる場合、更新インデックス値に基づいて選択された重み値に基づき、エネルギーの変化値を各状態遷移に対して計算するエネルギー変化計算部と、複数の前記エネルギーの変化値に基づいて前記各状態遷移のうちの少なくとも１つが許容されるようにオフセット値を決定し、複数の前記エネルギーの変化値のそれぞれに前記オフセット値を付加した複数の評価値を出力するとともに、前記所定試行回数毎に前記オフセット値を出力するオフセット制御部と、温度を示す温度値と乱数値とに基づいて決まる閾値と、前記複数の評価値のそれぞれとを比較した結果に基づき、前記各状態遷移について状態遷移を許容するか否かを示す複数のフラグ値を出力する比較部と、前記複数のフラグ値のうち、状態遷移を許容することを示すフラグ値の数を計数した計数値を、前記所定試行回数毎に出力するフラグ計数部と、前記複数のフラグ値、または複数の前記エネルギーの変化値に基づいて、前記各状態遷移のうちの１つに対応するインデックス値を前記更新インデックス値として出力する選択部と、を有するサンプリング装置が提供される。

【0019】

また、１つの実施態様では、サンプリング装置の制御方法が提供される。

【発明の効果】

【0020】

１つの側面では、本発明は、計算処理の高速性を損なうことなく高精度のサンプリングが可能となる。

【図面の簡単な説明】

【0021】

【図1】第１の実施の形態のサンプリング装置の一例を示す図である。

【図2】各状態に留まる試行回数を１回とした場合の選択確率の分布例を示す図である。

【図3】状態保持部及びエネルギー変化計算部の一例を示す図である。

【図4】オフセット制御部及び比較部の一例を示す図である。

【図5】オフセット値の適用例を示す図である。

【図6】フラグ計数部の一例を示す図である。

【図7】フラグ計数部の他の例を示す図である。

【図8】選択部の一例を示す図である。

【図9】サンプリングタイミングの例を示す図である。

【図10】サンプリング動作の一例の流れを示すフローチャートである。

【図11】最適化処理動作の一例の流れを示すフローチャートである。

【図12】第２の実施の形態のサンプリング装置の一例を示す図である。

【図13】２種類のフラグ値についての計数値Ｎ_ｆ１，Ｎ_ｆ２を発生させる動作の一例の流れを示すフローチャートである。

【図14】第３の実施の形態のサンプリング装置の一例を示す図である。

【図15】閾値生成部の一例を示す図である。

【図16】温度とエネルギーとの関係についてのシミュレーション結果の例を示す図である。

【図17】試行回数と最低エネルギー状態が得られる確率との関係についてのシミュレーション結果の例を示す図である。

【図18】レプリカ交換法を利用したサンプリング装置の例を示す図である。

【図19】レプリカ交換法を利用したサンプリング装置の動作例を示すフローチャートである。

【図20】サンプリング装置の他の変形例を示す図である。

【図21】期待値計算処理を行うサンプリング装置の動作例を示すフローチャートである。

【発明を実施するための形態】

【0022】

以下、発明を実施するための形態を、図面を参照しつつ説明する。
なお、イジングモデルを用いたシミュレーテッド・アニーリングにおいては、状態遷移に伴い変化する状態変数は１つだけである。そこで、以下では各状態遷移をそれぞれ識別するインデックス値は、１つの状態変数のインデックス値と等しいものとして説明を行う。しかし、状態遷移のインデックス値と状態遷移に伴い変化する状態遷移のインデックス値が一致する形態に限定されるものではない。

【0023】

従来の最適化装置において確率分布がボルツマン分布からずれる原因は、並列探索処理及びオフセット値の付加による計算処理の高速化を行っているため、１つの状態に滞在する試行回数が、通常のＭＣＭＣ法を用いた場合よりも少なくなるためである。従来の最適化装置において、理想的なオフセット値が適用された場合、１つの状態に滞在する試行回数は、１回となる。

【0024】

ボルツマン分布にしたがったサンプルを得るために、上記のようなオフセット値を用い、サンプリングを行うたび、通常のＭＣＭＣ法であればその状態に留まったであろう試行回数の値またはその期待値により、サンプルに重み付けを行えばよい。つまり、ある状態が１回観測された（サンプリングされた）場合、その状態が通常のＭＣＭＣ法を用いた場合の試行回数だけサンプリングされたとするような重み付けを行えばよい。

【0025】

通常のＭＣＭＣ法において、ある状態からインデックス値＝ｉの状態変数の値が変化した別の状態に遷移する確率は、Γ_ｉ＝Ａ（ΔＥ_ｉ）／Ｎである。Ａ（ΔＥ_ｉ）は、インデックス値＝ｉの状態変数の値が変化する状態遷移の許容確率である。Ｎは、状態変数の総数である。試行毎に別の状態に遷移する確率は、以下の式（４）で表せる。

【0026】

【数4】

【0027】

式（４）で表される確率の逆数が、通常のＭＣＭＣ法において、１つの状態に滞在する試行回数の期待値となる。
なお、何回かの試行を行って最終的に状態変数ｘ_ｉの値の変化が生じる確率Ｐ_ｉは、以下の式（５）で表せる。

【0028】

【数5】

【0029】

以下に示す第１の実施の形態のサンプリング装置は、所定試行回数毎に、通常のＭＣＭＣ法において１つの状態に滞在する試行回数の期待値を、その状態とともに出力する装置である。

【0030】

（第１の実施の形態）
図１は、第１の実施の形態のサンプリング装置の一例を示す図である。
サンプリング装置１０は、状態保持部１１、エネルギー変化計算部１２、オフセット制御部１３、比較部１４、フラグ計数部１５、選択部１６、制御部１７、試行回数計算部１８を有する。

【0031】

状態保持部１１は、たとえば、エネルギーを表す評価関数に含まれる複数の状態変数の値をそれぞれ保持するとともに、所定試行回数毎に複数の状態変数の値を出力する。評価関数は、たとえば、式（１）に示したようなエネルギー関数Ｅ（ｘ）である。

【0032】

以下では、複数の状態変数の数（状態変数の総数）をＮ個とし、Ｎ個の状態変数を状態変数ｘ_ｉ（ｉ＝１～Ｎ）または状態変数ｘ_１～ｘ_Ｎと表記する。また、状態変数ｘ_１～ｘ_Ｎの各値の組合せを状態ｘという。なお、以下の説明では、状態保持部１１は、式（２）に示したローカルフィールドｈ_ｉ（ｉ＝１～Ｎ）についても保持するものとして説明するが、ローカルフィールドｈ_ｉは、エネルギー変化計算部１２が保持していてもよい。

【0033】

状態保持部１１は、たとえば、レジスタやＳＲＡＭ（Static Random Access Memory）などである。
エネルギー変化計算部１２は、状態変数ｘ_１～ｘ_Ｎの値の何れかが変化することに応じて状態遷移が起こる場合、更新インデックス値に基づいて選択された重み値に基づき、エネルギーの変化値（以下エネルギー変化ΔＥ_ｉという）を各状態遷移に対して計算する。

【0034】

エネルギー変化ΔＥ_ｉは、前述の式（２）のように表せ、エネルギー変化計算部１２は、Δｘ_ｉに応じてローカルフィールドｈ_ｉに符号（＋１または－１）を乗じることでエネルギー変化ΔＥ_ｉを計算できる。また、状態変数ｘ_ｊが０から１に変化したときのローカルフィールドｈ_ｉの変化分Δｈ_ｉは、＋Ｗ_ｉｊ、状態変数ｘ_ｊが１から０に変化したときの変化分Δｈ_ｉは、－Ｗ_ｉｊと表せる。したがって、ローカルフィールドｈ_ｉは行列演算により毎回計算しなおす必要はなく、状態遷移にともなって変化のあったビットによる変化分だけＷ_ｉｊを加減算すればよい。エネルギー変化計算部１２は、計算（または更新）したローカルフィールドｈ_ｉを状態保持部１１に保持させる。

【0035】

このようなエネルギー変化計算部１２は、たとえば、重み値Ｗ_ｉｊなどを記憶する記憶部（レジスタやＳＲＡＭなど）、セレクタ、乗算器、加算器などを用いて実現できる。なお、以下の説明ではエネルギー変化計算部１２は、Ｎ個のエネルギー変化ΔＥ_ｉ（ｉ＝１～Ｎ）（以下エネルギー変化ΔＥ_１～ΔＥ_Ｎと表記する場合もある）とともに、インデックス値＝ｉ（１～Ｎ）についても出力するものとする。

【0036】

オフセット制御部１３は、エネルギー変化ΔＥ_１～ΔＥ_Ｎに基づいて各状態遷移のうちの少なくとも１つが許容されるようにオフセット値Ｅ_ｏｆｆを決定し、エネルギー変化ΔＥ_１～ΔＥ_Ｎのそれぞれにオフセット値Ｅ_ｏｆｆを付加した複数の評価値を出力する。また、オフセット制御部１３は、所定試行回数毎にオフセット値Ｅ_ｏｆｆを出力する。オフセット制御部１３の例やオフセット値Ｅ_ｏｆｆの決定方法の例については後述する。

【0037】

比較部１４は、温度値Ｔと乱数値とに基づいて決まる閾値と、複数の評価値のそれぞれとを比較した結果に基づき、各状態遷移（の候補）について状態遷移を許容するか否かを示す複数のフラグ値Ｆ_１，Ｆ_２，…，Ｆ_Ｎを出力する。比較部１４は、このような処理によって、式（３）に示した許容確率Ａ（ΔＥ）で状態遷移を許容させることができる。比較部１４の例や、許容確率Ａ（ΔＥ）で状態遷移を許容させることができる理由については後述する。なお、以下の説明では、サンプリング動作は一定の温度で行われる（温度値Ｔが固定値となる）ものとするが、温度は変更されるようにしてもよい。

【0038】

フラグ計数部１５は、フラグ値Ｆ_１～Ｆ_Ｎのうち、状態遷移を許容することを示すフラグ値の数を計数し、所定試行回数毎に計数したフラグ値の数（計数値Ｎ_ｆ）を出力する。フラグ計数部１５の例については後述する。

【0039】

選択部１６はフラグ値Ｆ_１～Ｆ_Ｎの何れかを選択するとともに、選択した何れかのフラグ値に対応するインデックス値を更新インデックス値として出力する。選択部１６の例についても後述する。

【0040】

制御部１７は、サンプリング装置１０の各部を制御する。制御部１７は、たとえば、試行回数が所定試行回数に達したか否かを判定し、所定試行回数に達したときに、そのときの状態を状態保持部１１に出力させ、オフセット値Ｅ_ｏｆｆをオフセット制御部１３に出力させる。さらに、制御部１７は、フラグ計数部１５にそのときの計数値Ｎ_ｆを計数させる。また、制御部１７は、比較部１４に温度値Ｔを供給する。

【0041】

制御部１７は、たとえば、ＡＳＩＣ（Application Specific Integrated Circuit）やＦＰＧＡ（Field Programmable Gate Array）などの特定用途の電子回路にて実現できる。なお、制御部１７は、ＣＰＵ（Central Processing Unit）やＤＳＰ（Digital Signal Processor）などのプロセッサであってもよい。その場合、プロセッサは、図示しないメモリに記憶されたプログラムを実行することで、上記の処理を行う。

【0042】

試行回数計算部１８は、オフセット制御部１３が出力したオフセット値Ｅ_ｏｆｆと、フラグ計数部１５が出力した計数値Ｎ_ｆに基づき、通常のＭＣＭＣ法において、１つの状態に滞在する試行回数Ｎ_{ｔｒｉａｌ}の期待値＜Ｎ_{ｔｒｉａｌ}＞（近似値）を計算する。試行回数計算部１８は、たとえば、ＡＳＩＣやＦＰＧＡなどの特定用途の電子回路、ＣＰＵやＤＳＰなどのプロセッサにて実現できる。なお、試行回数計算部１８は、サンプリング装置１０の外部に設けられていてもよいし、制御部１７が、オフセット値Ｅ_ｏｆｆと、計数値Ｎ_ｆに基づいて、期待値＜Ｎ_{ｔｒｉａｌ}＞を計算してもよい。

【0043】

オフセット値Ｅ_ｏｆｆと計数値Ｎ_ｆを、期待値＜Ｎ_{ｔｒｉａｌ}＞を求めるための値として用いることができる理由は以下の通りである。
図１には、通常のＭＣＭＣ法における遷移確率Γ_ｉの分布の例が示されている。通常のＭＣＭＣ法では、遷移確率Γ_ｉを全状態変数について積算した値（試行毎に別の状態に遷移する確率（式（４））は、１に満たない。前述のように、式（４）で表される確率の逆数が、通常のＭＣＭＣ法において、１つの状態に滞在する試行回数Ｎ_{ｔｒｉａｌ}の期待値＜Ｎ_{ｔｒｉａｌ}＞となる。期待値＜Ｎ_{ｔｒｉａｌ}＞は、現在の状態に留まっている時間（滞在時間）の期待値に相当する。

【0044】

一方、サンプリング装置１０において、各状態に留まる試行回数を１回とした場合、各状態変数が反転する対象として選ばれる選択確率の分布は、たとえば、以下のようになる。

【0045】

図２は、各状態に留まる試行回数を１回とした場合の選択確率の分布例を示す図である。
ｐ_１～ｐ_Ｎは、状態変数ｘ_ｉ（ｉ＝１～Ｎ）の選択確率である。各状態に留まる試行回数を１回とした場合、各試行において状態変数ｘ_ｉ（ｉ＝１～Ｎ）の何れかの変化が生じるため、ｐ_１～ｐ_Ｎの積算結果は１となる。選択確率ｐ_ｉは、以下の式（６）で表せる。

【0046】

【数6】

【0047】

式（６）の右辺は、式（５）と同じになる。
サンプリング装置１０は、各状態遷移のうちの少なくとも１つが許容されるようなオフセット値Ｅ_ｏｆｆを用いることで、上記のように各状態に留まる試行回数が１回となるようにする。オフセット値Ｅ_ｏｆｆを適用したときの許容確率をＡ（ΔＥ_ｉ－Ｅ_ｏｆｆ）とすると、Ａ（ΔＥ_ｉ－Ｅ_ｏｆｆ）は、前述の式（３）により（メトロポリス法を用いた場合）、以下の式（７）のように表せる。

【0048】

【数7】

【0049】

逆温度βとオフセット値Ｅ_ｏｆｆが正の場合、許容確率Ａ（ΔＥ_ｉ－Ｅ_ｏｆｆ）は、本来の許容確率Ａ（ΔＥ_ｉ）に対して、１以上の係数ｅｘｐ（β・Ｅ_ｏｆｆ）を乗じた確率となる。このため、許容確率Ａ（ΔＥ_ｉ－Ｅ_ｏｆｆ）が用いられる場合、本来の許容確率Ａ（ΔＥ_ｉ）を用いた場合に対して、各状態遷移の許容確率の相互の比率は変わらない。

【0050】

期待値＜Ｎ_{ｔｒｉａｌ}＞は、式（７）を用いて、以下の式（８）のように表せる。

【0051】

【数8】

【0052】

すなわち、状態変数ｘ_１～ｘ_ＮについてのＡ（ΔＥ_ｉ－Ｅ_ｏｆｆ）を積算した値＜Ｎ_ｆ＞は、状態遷移を許容することを示すフラグ値の数である計数値Ｎ_ｆにより近似できる。
式（８）において、逆温度βは一定値（温度値Ｔが固定値であるため）であり、状態変数の総数＝Ｎも固定値である。そのため、サンプリングタイミングにおけるオフセット値Ｅ_ｏｆｆと計数値Ｎ_ｆが得られれば、そのサンプリングタイミングにおける期待値＜Ｎ_{ｔｒｉａｌ}＞を近似的に求めることができる。

【0053】

サンプリング装置１０は、オフセット値Ｅ_ｏｆｆを用いて各状態に留まる試行回数が１回となるようにして、所定試行回数毎、状態と、オフセット値Ｅ_ｏｆｆと、計数値Ｎ_ｆとを出力（サンプリング）する。これにより、上記のようにオフセット値Ｅ_ｏｆｆと、計数値Ｎ_ｆとから得られる期待値＜Ｎ_{ｔｒｉａｌ}＞による重み付けにより、得られた状態が通常のＭＣＭＣ法ではどのくらいの試行回数、その状態に留まるのかがわかり、通常のＭＣＭＣ法の確率過程を再現できる。

【0054】

図１には、通常のＭＣＭＣ法を用いた場合と、サンプリング装置１０を用いた場合のエネルギーＥと試行回数との関係の一例を示すグラフが示されている。図１に示すように、たとえば、あるエネルギーＥａをもつ状態ｘの滞在時間は、サンプリング装置１０を用いた場合は１試行回数分であるのに対し、通常のＭＣＭＣ法を用いた場合、期待値＜Ｎ_{ｔｒｉａｌ}＞分である。期待値＜Ｎ_{ｔｒｉａｌ}＞がわかれば、上記のように通常のＭＣＭＣ法の確率過程を再現できる。

【0055】

つまり、通常のＭＣＭＣ法にしたがったサンプリングを再現できるため、従来の最適化装置をサンプリング装置として用いた場合に比べて高精度なサンプリングが可能になり、ボルツマン分布にしたがったサンプルが得られる。

【0056】

機械学習の理論は、サンプルがボルツマン分布になっていることを根拠にしているものがほとんどであり、上記のようなサンプリング装置１０を用いることで、機械学習に適したサンプルを得ることができる。

【0057】

また、各試行で状態遷移が生じることと、１回のサンプリングで実効的に複数回サンプリングするのと同じ精度が得られることにより、従来の最適化装置をサンプリング装置として用いた場合に比べて計算処理の高速性が損なわれることもない。

【0058】

以下、サンプリング装置１０の各部の例を説明する。
（状態保持部１１及びエネルギー変化計算部１２の例）
図３は、状態保持部及びエネルギー変化計算部の一例を示す図である。

【0059】

状態保持部１１は、状態変数ｘ_１～ｘ_Ｎの値とローカルフィールドｈ_１～ｈ_Ｎとを保持するレジスタ１１ａを有する。レジスタ１１ａは、クロック信号ｃｌｋに同期して、エネルギー変化計算部１２が更新した、状態変数ｘ_１～ｘ_Ｎの何れかの値とローカルフィールドｈ_１～ｈ_Ｎとを取り込む。クロック信号ｃｌｋは、図示しないクロック信号生成回路から供給される。状態保持部１１は、たとえば、図１に示した制御部１７の制御のもと、所定試行回数毎に複数の状態変数の値を、サンプリング装置１０の外部に出力する。

【0060】

エネルギー変化計算部１２は、ΔＥ計算回路１２ａ１，１２ａ２，…，１２ａＮ、更新回路１２ｂを有する。
ΔＥ計算回路１２ａ１～１２ａＮのそれぞれは、状態変数ｘ_１～ｘ_Ｎの値の何れかとローカルフィールドｈ_１～ｈ_Ｎの何れかを状態保持部１１から読み出し、エネルギー変化ΔＥ_１～Ｅ_Ｎの何れかを計算する。たとえば、ΔＥ計算回路１２ａ１は、状態変数ｘ_１の値とローカルフィールドｈ_１とを状態保持部１１から読み出し、状態変数ｘ_１の値が変化することによるエネルギー変化ΔＥ_１を計算する。ΔＥ計算回路１２ａＮは、状態変数ｘ_Ｎの値とローカルフィールドｈ_Ｎとを状態保持部１１から読み出し、状態変数ｘ_Ｎの値が変化することによるエネルギー変化ΔＥ_Ｎを計算する。

【0061】

エネルギー変化ΔＥ_ｉ（ｉ＝１～Ｎ）は、式（２）のように表せるため、ΔＥ計算回路１２ａ１～１２ａＮのそれぞれは、たとえば、状態変数ｘ_１～ｘ_Ｎの値に応じてローカルフィールドｈ_１～ｈ_Ｎに符号（＋１または－１）を乗じる乗算器を用いて実現できる。

【0062】

なお、ΔＥ計算回路１２ａ１～１２ａＮのそれぞれは、インデックス値＝ｉ（１～Ｎ）の何れかを出力する。たとえば、ΔＥ計算回路１２ａ１は、インデックス値として１を出力し、ΔＥ計算回路１２ａＮは、インデックス値としてＮを出力する。ΔＥ計算回路１２ａ１～１２ａＮは、インデックス値を保持する保持部（たとえば、レジスタ）を有する。

【0063】

更新回路１２ｂは、選択部１６より更新インデックス値を受けると、更新インデックス値に対応した状態変数の値を変化させるとともに、その変化に基づいて、ローカルフィールドｈ_１～ｈ_Ｎを更新する。更新後の各値は、状態保持部１１に書き込まれる。

【0064】

たとえば、更新回路１２ｂは、更新インデックス値がｊの場合、状態変数ｘ_ｊを１－ｘ_ｊに変える。そして、更新回路１２ｂは、ローカルフィールドｈ_ｉ（ｉ＝１～Ｎ）を、ｈ_ｉ＋Ｗ_ｉｊ・Δｘ_ｊに更新する。

【0065】

なお、更新回路１２ｂは、上記のような更新処理を行うために、現在の状態変数ｘ_ｉ（ｉ＝１～Ｎ）の値とローカルフィールドｈ_ｉを、状態保持部１１から取得する。重み値Ｗ_ｉｊは、図示しない記憶部に記憶されており、更新回路１２ｂは、その記憶部から重み値Ｗ_ｉｊを取得する。

【0066】

このような更新回路１２ｂは、重み値Ｗ_ｉｊなどを記憶する記憶部（レジスタやＳＲＡＭなど）、セレクタ、乗算器、加算器などを用いて実現できる。
（オフセット制御部１３及び比較部１４の例）
図４は、オフセット制御部及び比較部の一例を示す図である。なお、図４では、前述のエネルギー変化ΔＥ_１～ΔＥ_Ｎがまとめて｛ΔＥ_ｉ｝と表記され、前述のフラグ値Ｆ_１～Ｆ_Ｎがまとめて｛Ｆ_ｉ｝と表記されている。

【0067】

オフセット制御部１３は、最小値検出回路１３ａ、加算回路１３ｂ、選択信号生成回路１３ｃ、セレクタ１３ｄを有する。
最小値検出回路１３ａは、エネルギー変化ΔＥ_１～ΔＥ_Ｎのうちの最小値ΔＥ_ｍｉｎを検出する。

【0068】

加算回路１３ｂは、エネルギー変化ΔＥ_１～ΔＥ_Ｎから、最小値ΔＥ_ｍｉｎをオフセット値として差し引く（エネルギー変化ΔＥ_１～ΔＥ_Ｎに－ΔＥ_ｍｉｎを加算する、ということもできる）演算回路である。加算回路１３ｂの代りに減算回路を用いることもできる。

【0069】

選択信号生成回路１３ｃは、選択信号ｓｉｇｎを生成する。選択信号生成回路１３ｃは、最小値ΔＥ_ｍｉｎが０より小さい場合、セレクタ１３ｄにエネルギー変化ΔＥ_１～ΔＥ_Ｎを選択させる選択信号ｓｉｇｎ（以下では、０とする）を生成する。選択信号生成回路１３ｃは、最小値ΔＥ_ｍｉｎが０以上の場合、セレクタ１３ｄに加算回路１３ｂが算出した値を選択させる選択信号ｓｉｇｎ（以下では、１とする）を生成する。

【0070】

セレクタ１３ｄは、選択信号ｓｉｇｎが０の場合、エネルギー変化ΔＥ_１～ΔＥ_Ｎ（｛ΔＥ_ｉ｝）を選択し、選択信号ｓｉｇｎが１の場合、加算回路１３ｂが算出した値、すなわち｛ΔＥ_ｉ｝－ΔＥ_ｍｉｎを選択して、評価値として出力する。

【0071】

図４に示すようなオフセット制御部１３を用いた場合、オフセット値Ｅ_ｏｆｆは、Ｅ_ｏｆｆ＝ｍａｘ［０，ｍｉｎ（ΔＥ_ｉ）］と表せる。つまり、エネルギー変化ΔＥ_１～ΔＥ_Ｎのうちの最小値ΔＥ_ｍｉｎが０以上の場合、Ｅ_ｏｆｆ＝ΔＥ_ｍｉｎ、それ以外の場合、Ｅ_ｏｆｆ＝０である。

【0072】

このようなオフセット値Ｅ_ｏｆｆを用いることで、状態変数ｘ_１～ｘ_Ｎのうち、少なくとも１つは状態遷移が許容される。その理由を以下に示す。
図５は、オフセット値の適用例を示す図である。図５において、横軸は状態遷移の候補となる状態変数ｘ_ｉを表し、縦軸はエネルギー変化ΔＥを表す。

【0073】

波形２０ａは、状態変数ｘ_１～ｘ_Ｎのそれぞれが変化することによるエネルギー変化｛ΔＥ_ｉ｝の例を示している。図５の例では、エネルギー変化｛ΔＥ_ｉ｝のうち、状態変数ｘ_ｍｉｎの値が変化することによるエネルギー変化ΔＥが最小値ΔＥ_ｍｉｎであり、０より大きい値である。波形２０ｂは、波形２０ａで示されるエネルギー変化｛ΔＥ_ｉ｝から最小値ΔＥ_ｍｉｎを引いた場合のエネルギー変化ΔＥを示す。

【0074】

波形２０ｂにおけるエネルギー変化ΔＥの最小値は０となるため、状態変数ｘ_ｍｉｎの値が変化する状態遷移の許容確率Ａ（ΔＥ）が１となる。つまり、状態変数ｘ_１～ｘ_Ｎのうち、少なくとも状態変数ｘ_ｍｉｎは状態遷移が許容される。

【0075】

このように、エネルギー変化ΔＥ_１～ΔＥ_Ｎのうちの最小値ΔＥ_ｍｉｎが０より大きい場合でも、Ｅ_ｏｆｆ＝ΔＥ_ｍｉｎをエネルギー変化｛ΔＥ_ｉ｝から差し引くことで、状態変数ｘ_１～ｘ_Ｎのうち、少なくとも１つの状態遷移が許容される。

【0076】

図４の説明に戻る。
比較部１４は、乱数発生回路１４ｂ１、選択法則適用部１４ｂ２、乗算回路１４ｂ３、符号反転回路１４ｂ４、比較回路１４ｂ５を有する。

【0077】

乱数発生回路１４ｂ１は、０より大きく、１以下の一様乱数である乱数値ｒを発生する。乱数発生回路１４ｂ１は、たとえば、メルセンヌツイスタ、ＬＦＳＲ（Linear Feedback Shift Register）などを用いて実現できる。

【0078】

選択法則適用部１４ｂ２は、選択法則（前述のメトロポリス法またはギブス法）に基づいた値を出力する。
前述の式（３）の許容確率Ａ（ΔＥ）で、エネルギー変化ΔＥを引き起こす状態遷移を許容することを示すフラグ（＝１）を出力する回路は、許容確率Ａ（ΔＥ）と、乱数値ｒとの比較結果に基づいた値を出力する比較器によって実現できる。

【0079】

ただ、次のような変形を行っても同じ機能が実現できる。２つの数に同じ単調増加関数を作用させても大小関係は変化しない。したがって比較器の２つの入力に同じ単調増加関数を作用させても比較器の出力は変わらない。たとえば、許容確率Ａ（ΔＥ）を、ｆ（－β・ΔＥ）と表記する。ｆ（－β・ΔＥ）に作用させる単調増加関数としてｆ（－β・ΔＥ）の逆関数ｆ^－１（－β・ΔＥ）、乱数値ｒに作用させる単調増加関数としてｆ^－１（－β・ΔＥ）の－β・ΔＥをｒとしたｆ^－１（ｒ）を用いることができる。その場合、上記の比較器と同様の機能を有する回路は、－β・ΔＥがｆ^－１（ｒ）より大きいとき１を出力する回路でよいことがわかる。さらにβ＝１／Ｔで、温度値Ｔが正であることから、その回路は、－ΔＥがＴ・ｆ^－１（ｒ）より大きいとき１を出力する回路でよい。

【0080】

選択法則適用部１４ｂ２は、入力される乱数値ｒを上記のｆ^－１（ｒ）の値に変換する変換テーブルを用いて、ｆ^－１（ｒ）の値を出力する。メトロポリス法が適用される場合、ｆ^－１（ｒ）は、ｌｏｇ（ｒ）である。変換テーブルは、たとえば、ＲＡＭ（Random Access Memory）、フラッシュメモリなどのメモリに記憶されている。

【0081】

乗算回路１４ｂ３は、制御部１７から供給される温度値Ｔと、ｆ^－１（ｒ）との積（Ｔ・ｆ^－１（ｒ））を閾値として出力する。Ｔ・ｆ^－１（ｒ）は、熱励起エネルギーに相当する。

【0082】

符号反転回路１４ｂ４は、評価値である｛ΔＥ_ｉ｝または｛Ｅ_ｉ｝－ΔＥ_ｍｉｎの符号を反転し、－｛ΔＥ_ｉ｝または－｛ΔＥ_ｉ｝＋ΔＥ_ｍｉｎを出力する。
比較回路１４ｂ５は、－｛ΔＥ_ｉ｝または－｛ΔＥ_ｉ｝＋ΔＥ_ｍｉｎと、Ｔ・ｆ^－１（ｒ）とを比較する。そして、比較回路１４ｂ５は、－｛ΔＥ_ｉ｝または－｛ΔＥ_ｉ｝＋ΔＥ_ｍｉｎがＴ・ｆ^－１（ｒ）より大きい場合、フラグ値｛Ｆ_ｉ｝として１を出力し、－｛ΔＥ_ｉ｝または－｛ΔＥ_ｉ｝＋ΔＥ_ｍｉｎがＴ・ｆ^－１（ｒ）以下の場合、フラグ値｛Ｆ_ｉ｝として０を出力する。

【0083】

（フラグ計数部１５の例）
図１に示したフラグ計数部１５の例として、以下に、２種のフラグ計数部１５ａ，１５ｂを示す。

【0084】

図６は、フラグ計数部の一例を示す図である。
フラグ計数部１５ａは、ツリー状に複数段、配置された複数の加算回路（加算回路１５ａ１，１５ａ２，１５ａ３，１５ａ４，１５ａ５など）を有する。

【0085】

初段の各加算回路は、２つのフラグ値を加算した加算結果を出力する。たとえば、加算回路１５ａ１は、フラグ値Ｆ_１，Ｆ_２を加算した加算結果を出力し、加算回路１５ａ２は、フラグ値Ｆ_３，Ｆ_４を加算した加算結果を出力する。２段目以降の加算回路は、前段の２つの加算回路の出力値を加算した加算結果を出力する。たとえば、２段目の加算回路１５ａ４は、初段の加算回路１５ａ１，１５ａ２の出力値を加算した加算結果を出力する。そして、最後段の加算回路１５ａ５が出力する加算結果が、計数値Ｎ_ｆとなる。

【0086】

初段の各加算回路は、２つの１ビットのフラグ値を加算して、２ビットの加算結果を出力し、２段目の各加算回路は、２つの加算回路のそれぞれの２ビットの出力値を加算して、３ビットの加算結果を出力する。Ｎ＝１０２４の場合、最後段（１０段目）の加算回路１５ａ５は、前段の２つの加算回路のそれぞれの１０ビットの出力値を加算して、１１ビットの計数値Ｎ_ｆを出力する。

【0087】

図７は、フラグ計数部の他の例を示す図である。
フラグ計数部１５ｂは、マイクロコントローラ１５ｂ１とプログラムメモリ１５ｂ２とを有する。

【0088】

マイクロコントローラ１５ｂ１は、プログラムメモリ１５ｂ２に記憶されたプログラムを読み出して実行することで、フラグ値Ｆ_１～Ｆ_Ｎのうち値が１である数を計数し、その計数結果である計数値Ｎ_ｆを出力する。

【0089】

サンプリング装置１０によるサンプリングは、通常１０００クロックサイクル以上の間隔（サンプリング間隔）で行われるため、比較的処理速度が遅いマイクロコントローラ１５ｂ１及び、比較的小規模のプログラムメモリ１５ｂ２を用いることができる。

【0090】

（選択部１６の例）
図８は、選択部の一例を示す図である。
選択部１６は、乱数発生回路１６ａ、加算回路１６ｂ、選択回路部１６ｃを有する。

【0091】

乱数発生回路１６ａは、たとえば、０＜ｒ_ｉ＜１の一様乱数である乱数値ｒ_ｉを、インデックス値＝ｉ（１～Ｎ）のそれぞれについて異なるシードを用いて発生する。乱数発生回路１６ａは、たとえば、メルセンヌツイスタ、ＬＦＳＲなどを用いて実現できる。

【0092】

加算回路１６ｂは、フラグ値Ｆ_１～Ｆ_Ｎのそれぞれに対して、複数の乱数値ｒ_ｉ（ｉ＝１～Ｎ）の何れかを加算した加算結果を、フラグ値Ｆ_１～Ｆ_Ｎに対応するインデックス値とともに出力する。

【0093】

選択回路部１６ｃは、ツリー状に複数段、配置された複数の選択回路（選択回路１６ｃ１，１６ｃ２，１６ｃ３，１６ｃ４，１６ｃ５など）を有する。
初段の各選択回路は、２つのフラグ値についての加算結果のうち、大きい方の加算結果とそれに対応したインデックス値を選択して出力する。たとえば、選択回路１６ｃ１は、Ｆ_１＋ｒ_１がＦ_２＋ｒ_２より大きい場合、Ｆ_１＋ｒ_１と、インデックス値＝１を出力し、Ｆ_１＋ｒ_１がＦ_２＋ｒ_２より小さい場合、Ｆ_２＋ｒ_２と、インデックス値＝２を出力する。

【0094】

２段目以降の選択回路は、前段の２つの選択回路が選択した加算結果のうち、大きい方の加算結果とそれに対応したインデックス値を選択して出力する。たとえば、２段目の選択回路１６ｃ４は、初段の選択回路１６ｃ１，１６ｃ２が選択した加算結果のうち、大きい方の加算結果とそれに対応したインデックス値を選択して出力する。最後段の選択回路１６ｃ５は、前段の２つの選択回路が選択した加算結果のうち、大きい方の加算結果に対応したインデックス値＝ｊを更新インデックス値として出力する。

【0095】

このように、選択部１６は、加算結果の大小関係に基づいて、フラグ値Ｆ_１～Ｆ_Ｎの何れかを選択することで、フラグ値Ｆ_１～Ｆ_Ｎのうち、値が１のものを、ランダムに等確率で選択することができる。

【0096】

（動作例）
以下、サンプリング装置１０の動作例を説明する。
なお、以下の例ではサンプリング動作は、温度値Ｔが一定の値で行われるものとする。

【0097】

図９は、サンプリングタイミングの例を示す図である。
図９において、期間Ｔｕは、１回の試行（状態更新処理）にかかる期間を示す。サンプリング装置１０は、所定の試行回数（バーンイン期間）の試行が行われたのちに、サンプリングを開始する。バーンイン期間は、状態が平衡状態に達するまでの期間であり、バーンイン期間に対応する試行回数が予め設定される。

【0098】

サンプリング期間では、サンプリング装置１０は、所定の試行回数（サンプリング間隔Ｔｓ）毎に、サンプリング（たとえば、状態ｘ、オフセット値Ｅ_ｏｆｆ、計数値Ｎ_ｆなどの出力）を行う。サンプリング間隔Ｔｓに対応する試行回数は予め設定される。

【0099】

図１０は、サンプリング動作の一例の流れを示すフローチャートである。
まず、初期化やパラメータ設定が行われる（ステップＳ１）。ステップＳ１の処理において、制御部１７は、たとえば、重み値、バイアス係数、状態変数ｘ_ｉ（ｉ＝１～Ｎ）の初期値に基づいて、ローカルフィールドｈ_ｉ（ｉ＝１～Ｎ）の初期値の計算を行う。重み値、バイアス係数、状態変数ｘ_ｉ（ｉ＝１～Ｎ）の初期値は、計算対象の組合せ最適化問題を変換したイジングモデルの情報として、たとえば、図示しない記憶部に予め記憶されている。状態変数ｘ_ｉの初期値やローカルフィールドｈ_ｉの初期値は、状態保持部１１に保持される。

【0100】

また、バーンイン期間に対応する試行回数Ｎ_Ｂ、サンプリング間隔Ｔｓに対応する試行回数Ｎ_ＳＩ、全サンプルが得られる試行回数Ｎ_ＳＮについても、たとえば、図示しない記憶部に予め記憶されている。そして、ステップＳ１の処理において、制御部１７にそれらのパラメータが設定される。

【0101】

その後、エネルギー変化計算部１２は、状態変数ｘ_ｉの値とローカルフィールドｈ_ｉを、状態保持部１１から読み出す（ステップＳ２）。
そして、エネルギー変化計算部１２は、前述のように、エネルギー変化ΔＥ_ｉ（ｉ＝１～Ｎ）を計算する（ステップＳ３）。

【0102】

オフセット制御部１３は、エネルギー変化計算部１２から供給されるエネルギー変化ΔＥ_ｉに対して、前述のようなオフセット値Ｅ_ｏｆｆを付加する（ステップＳ４）。
そして、前述の比較部１４と選択部１６の処理により、状態遷移を許容する状態変数を識別する更新インデックス値が出力される（ステップＳ５）。

【0103】

エネルギー変化計算部１２は、更新インデックス値により識別される状態変数ｘ_ｉの値を更新する。また、エネルギー変化計算部１２は、その更新に基づいて、ローカルフィールドｈ_ｉを、前述の処理により更新する（ステップＳ６）。

【0104】

更新後の状態変数ｘ_ｉとローカルフィールドｈ_ｉの各値は、状態保持部１１に書き込まれる（ステップＳ７）。
制御部１７は、ステップＳ２～Ｓ７の処理（試行）の回数が試行回数Ｎ_Ｂに達したか否かを判定する（ステップＳ８）。ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_Ｂに達していない場合、ステップＳ２からの処理が繰り返される。

【0105】

ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_Ｂに達した場合、制御部１７は、ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_ＳＩに達したか否かを判定する（ステップＳ９）。ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_ＳＩに達していない場合、ステップＳ２からの処理が繰り返される。

【0106】

ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_ＳＩに達した場合、制御部１７は、フラグ計数部１５に、計数値Ｎ_ｆを計数させる（ステップＳ１０）。そして、制御部１７は、フラグ計数部１５に計数値Ｎ_ｆ、状態保持部１１に状態ｘ、オフセット制御部１３にオフセット値Ｅ_ｏｆｆをそれぞれ出力させる（ステップＳ１１）。出力されたオフセット値Ｅ_ｏｆｆと、計数値Ｎ_ｆに基づいて、試行回数計算部１８は、式（８）に示される期待値＜Ｎ_{ｔｒｉａｌ}＞を計算する。

【0107】

制御部１７は、ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_ＳＮに達したか否かを判定する（ステップＳ１２）。ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_ＳＮに達していない場合、ステップＳ２からの処理が繰り返される。ステップＳ２～Ｓ７の処理の回数が試行回数Ｎ_ＳＮに達した場合、制御部１７は、サンプリング動作を終了させる。

【0108】

ところで、サンプリング装置１０は、最適化装置としても使用できる。最適化装置として使用する場合の動作例を以下に示す。
図１１は、最適化処理動作の一例の流れを示すフローチャートである。

【0109】

まず、初期化やパラメータ設定が行われる（ステップＳ２０）。ステップＳ２０の処理において、制御部１７は、サンプリング動作時と同様に、ローカルフィールドｈ_ｉ（ｉ＝１～Ｎ）の初期値の計算などを行う。

【0110】

また、ステップＳ２０の処理において、たとえば、図示しない記憶部に予め記憶された温度変更スケジュールにしたがって、制御部１７は、温度値Ｔの初期値を比較部１４に設定する。また、たとえば、図示しない記憶部には予め、温度値Ｔの更新を行う試行回数Ｎ１、上限の試行回数Ｎ２が記憶されており、ステップＳ２０の処理において、制御部１７にそれらのパラメータが設定される。

【0111】

その後のステップＳ２１～Ｓ２６の処理は、図１０に示したサンプリング動作時におけるステップＳ２～Ｓ７の処理とほぼ同様である。ただし、たとえば、制御部１７は、式（１）で表せるエネルギーＥの初期値を計算するとともに、更新インデックス値に対応するエネルギー変化（たとえば、更新インデックス値＝ｊの場合、エネルギー変化ΔＥ_ｊ）を用いて、エネルギーＥを更新する。そして、制御部１７は、これまでの最低エネルギーとそのときの状態を図示しない記憶部に保持させる。

【0112】

ステップＳ２６の後、制御部１７は、ステップＳ２１～Ｓ２６の処理（試行）の回数が試行回数Ｎ１に達したか否かを判定する（ステップＳ２７）。ステップＳ２１～Ｓ２６の処理の回数が試行回数Ｎ１に達していない場合、ステップＳ２１からの処理が繰り返される。

【0113】

ステップＳ２１～Ｓ２６の処理の回数が試行回数Ｎ１に達した場合、制御部１７は、シミュレーテッド・アニーリングを実現するために、温度値Ｔを、温度変更スケジュールにしたがって小さくするように更新する（ステップＳ２８）。

【0114】

その後、制御部１７は、ステップＳ２１～Ｓ２６の処理の回数が試行回数Ｎ２に達したか否かを判定する（ステップＳ２９）。ステップＳ２１～Ｓ２６の処理の回数が試行回数Ｎ２に達していない場合、ステップＳ２１からの処理が繰り返される。ステップＳ２１～Ｓ２６の処理の回数が試行回数Ｎ２に達した場合、制御部１７は、その時点で記憶部に保持させている最低エネルギーのときの状態を、組合せ最適化問題の解として出力し（ステップＳ３０）、最適化処理動作を終了させる。

【0115】

サンプリング装置１０を最適化装置として用いた場合、前述のようなオフセット値Ｅ_ｏｆｆを用いることで各試行において状態遷移が生じ、探索が高速化する。
（第２の実施の形態）
以下、第２の実施の形態のサンプリング装置を説明する。

【0116】

前述のように値が１であるフラグ値の数を数えることは、許容確率Ａ（ΔＥ_ｉ）（ｉ＝１～Ｎ）を積算したものを、仮数部１ビット、指数部をｅｘｐ（－β・Ｅ_ｏｆｆ）で表すことに相当する。以下に示す第２の実施の形態のサンプリング装置では、２種類のフラグ値から計数値Ｎ_ｆを求めることで、より精度のよい期待値＜Ｎ_{ｔｒｉａｌ}＞を計算するものである。

【0117】

図１２は、第２の実施の形態のサンプリング装置の一例を示す図である。図１２において、図１に示した要素と同じ要素については同じ符号が付されている。なお、図１２では、図１に示した制御部１７については図示が省略されている。また、以下では、図１においてオフセット制御部１３が出力するオフセット値Ｅ_ｏｆｆをオフセット値Ｅ_ｏｆｆ１といい、図１においてフラグ計数部１５が出力する計数値Ｎ_ｆを計数値Ｎ_ｆ１という。

【0118】

第２の実施の形態のサンプリング装置３０は、ΔＥ置き換え部３１、オフセット制御部３２、比較部３３、フラグ計数部３４、選択部３５、計数値補正部３６、試行回数計算部３７を有する。

【0119】

ΔＥ置き換え部３１は、フラグ値Ｆ_１～Ｆ_Ｎに基づき、エネルギー変化ΔＥ_１～ΔＥ_Ｎのうち、許容されない状態遷移であるとされている複数の状態遷移に対する複数の非許容エネルギー変化値（以下エネルギー変化ΔＥ_ｎａ１～ΔＥ_ｎａＭという）を出力する。ΔＥ置き換え部３１は、たとえば、エネルギー変化ΔＥ_１～ΔＥ_Ｎと、フラグ値Ｆ_１～Ｆ_Ｎを入力し、値が０のフラグ値のインデックス値に対応するインデックス値をもつエネルギー変化をエネルギー変化ΔＥ_１～ΔＥ_Ｎから選択する回路によって実現できる。

【0120】

オフセット制御部３２は、エネルギー変化ΔＥ_ｎａ１～ΔＥ_ｎａＭに基づき、フラグ値Ｆ_１～Ｆ_Ｎでは許容されない状態遷移とされている複数の状態遷移のうちの少なくとも１つが許容されるようにオフセット値Ｅ_ｏｆｆ１とは別のオフセット値Ｅ_ｏｆｆ２を決定する。そして、オフセット制御部３２は、エネルギー変化ΔＥ_ｎａ１～ΔＥ_ｎａＭのそれぞれにオフセット値Ｅ_ｏｆｆ２を付加した、オフセット制御部１３が出力する複数の評価値とは別の複数の評価値を出力する。また、オフセット制御部３２は、所定試行回数毎にオフセット値Ｅ_ｏｆｆ２を出力する。オフセット制御部３２は、オフセット制御部１３と同様の回路にて実現できる。

【0121】

比較部３３は、温度値Ｔと乱数値とに基づいて決まる閾値と、比較部３３が出力する複数の評価値のそれぞれとを比較した結果に基づき複数の状態遷移のそれぞれについて状態遷移を許容するか否かを示す、フラグ値Ｆ_１～Ｆ_Ｎとは別の複数のフラグ値を出力する。以下、比較部３３が出力する複数のフラグ値をフラグ値Ｆ_ｎａ１～Ｆ_ｎａＭという。比較部３３は、比較部１４と同様の回路にて実現できる。

【0122】

フラグ計数部３４は、比較部３３が出力するフラグ値Ｆ_ｎａ１～Ｆ_ｎａＭのうち、状態遷移を許容することを示すフラグ値（値が１）の数を計数した、計数値Ｎ_ｆ１とは別の計数値Ｎ_ｆ２を、所定試行回数毎に出力する。フラグ計数部３４は、フラグ計数部１５と同様の回路にて実現できる。

【0123】

選択部３５は、計数値Ｎ_ｆ１，Ｎ_ｆ２、オフセット値Ｅ_ｏｆｆ１，Ｅ_ｏｆｆ２に基づき、Ｎ_ｆ１：Ｎ_ｆ２ｅｘｐ［－β（Ｅ_ｏｆｆ１－Ｅ_ｏｆｆ２）］の比率で、フラグ値Ｆ_１～Ｆ_Ｎのうち値が１のグループと、フラグ値Ｆ_ｎａ１～Ｆ_ｎａＭのうち値が１のグループの何れか一方を選択する。そして、選択部３５は、選択したグループから等確率に１つのフラグ値を選択し、選択したフラグ値に対応したインデックス値を更新インデックス値として出力する。選択部３５は、たとえば、論理回路、ＡＳＩＣやＦＰＧＡなどの特定用途の電子回路、ＣＰＵやＤＳＰなどのプロセッサ、マイクロコントローラなどで実現できる。

【0124】

計数値補正部３６は、計数値Ｎ_ｆ２とオフセット値Ｅ_ｏｆｆ１，Ｅ_ｏｆｆ２に基づいて、計数値Ｎ_ｆ１を補正し、計数値Ｎ_ｆとして出力する。計数値補正部３６は、たとえば、Ｎ_ｆ＝Ｎ_ｆ１＋Ｎ_ｆ２ｅｘｐ［β（Ｅ_ｏｆｆ１－Ｅ_ｏｆｆ２）］を計算して出力する。計数値補正部３６は、たとえば、ＡＳＩＣやＦＰＧＡなどの特定用途の電子回路、ＣＰＵやＤＳＰなどのプロセッサにて実現できる。

【0125】

試行回数計算部３７は、オフセット値Ｅ_ｏｆｆ１，Ｅ_ｏｆｆ２と、計数値Ｎ_ｆに基づいて、たとえば、前述の式（８）を用いて期待値＜Ｎ_{ｔｒｉａｌ}＞を計算する。ただし、式（８）のＥ_ｏｆｆとして、（Ｅ_ｏｆｆ１－Ｅ_ｏｆｆ２）が用いられる。試行回数計算部３７は、たとえば、ＡＳＩＣやＦＰＧＡなどの特定用途の電子回路、ＣＰＵやＤＳＰなどのプロセッサにて実現できる。

【0126】

その他の構成及び動作については、第１の実施の形態のサンプリング装置１０と同じである。
上記のような第２の実施の形態のサンプリング装置３０によれば、たとえば、フラグ値Ｆ_１～Ｆ_Ｎのうち、１つだけが１であるような（１つの状態遷移だけが許容されるような）極端な場合でも、他の状態遷移についてもある割合で許容されるようになる。これにより、サンプリング装置３０を用いた場合の各状態遷移（状態変数）間の遷移確率の比を、通常のＭＣＭＣ法で得られるような各状態遷移間の遷移確率の比に近づけることができる。そして、上記のように補正された計数値Ｎ_ｆを用いて計算される期待値＜Ｎ_{ｔｒｉａｌ}＞は、通常のＭＣＭＣ法で得られる１つの状態に滞在する試行数により近いものになり、通常のＭＣＭＣ法により近い確率過程が得られる。

【0127】

なお、上記の第２の実施の形態のサンプリング装置３０の動作は、１つのオフセット制御部１３、１つの比較部１４、１つのフラグ計数部１５を用いても実行できる。その場合、以下のような動作によって、２種類のフラグ値（フラグ値Ｆ_１～Ｆ_Ｎとフラグ値Ｆ_ｎａ１～Ｆ_ｎａＭ）についての計数値Ｎ_ｆ１，Ｎ_ｆ２が得られる。

【0128】

図１３は、２種類のフラグ値についての計数値Ｎ_ｆ１，Ｎ_ｆ２を発生させる動作の一例の流れを示すフローチャートである。
以下に示す処理は、たとえば、図１に示した制御部１７（図１２では図示が省略されている）によって制御される。

【0129】

まず、制御部１７は、ｎ_ＭＡＸ＝２、ｎ＝１とする（ステップＳ４０）。そして、制御部１７は、オフセット制御部１３に、オフセット値Ｅ_ｏｆｆｎ（ｎ＝１の場合はオフセット値Ｅ_ｏｆｆ１）を発生させる（ステップＳ４１）。オフセット制御部１３は、エネルギー変化ΔＥ_ｉ（ｉ＝１～Ｎ）のそれぞれからオフセット値Ｅ_ｏｆｆｎを差し引いた評価値を計算する（ステップＳ４２）。図４のようなオフセット制御部１３が用いられる場合、ΔＥ_ｉの最小値ΔＥ_ｍｉｎが０より小さい場合、オフセット値Ｅ_ｏｆｆｎは０であり、最小値ΔＥ_ｍｉｎが０以上の場合、オフセット値Ｅ_ｏｆｆｎは最小値ΔＥ_ｍｉｎである。

【0130】

比較部１４は、温度値Ｔと乱数値とに基づいて決まる閾値と、複数の評価値のそれぞれとを比較した結果に基づき、フラグ値（ｎ＝１の場合、フラグ値Ｆ_１～Ｆ_Ｎ）を発生する（ステップＳ４３）。

【0131】

そして、フラグ計数部１５は、フラグ値が１である数を計数し（ステップＳ４４）、計数値Ｎ_ｆｎ（ｎ＝１の場合は計数値Ｎ_ｆ１）を出力する（ステップＳ４５）。
制御部１７は、ｎ≧ｎ_ＭＡＸであるか否かを判定し（ステップＳ４６）、ｎ≧ｎ_ＭＡＸの場合、計数値Ｎ_ｆ１，Ｎ_ｆ２を発生させる処理を終了する。ｎ＜ｎ_ＭＡＸの場合、制御部１７は、ｎ＝ｎ＋１とする（ステップＳ４７）。また、制御部１７は、前述のΔＥ置き換え部３１の機能により、エネルギー変化ΔＥ_ｉを、エネルギー変化ΔＥ_ｎａ１～ΔＥ_ｎａＭに置き換え（ステップＳ４８）、ステップＳ４１からの処理を繰り返す。

【0132】

なお、上記ステップＳ４４，Ｓ４５の処理以外は、各試行において行われる。
（第３の実施の形態）
図１４は、第３の実施の形態のサンプリング装置の一例を示す図である。図１４において、図１に示した要素と同じ要素については図示が省略されている。

【0133】

第３の実施の形態のサンプリング装置４０では、図１に示した選択部１６の代りに、選択部４１が設けられている。
選択部４１は、エネルギー変化ΔＥ_１～ΔＥ_Ｎのうち、０より小さいものを０に更新させた複数の更新エネルギー変化値を生成する。そして、選択部４１は、複数の更新エネルギー変化値のそれぞれに、複数の閾値の何れかを加算した複数の評価値のうちで、最小となる評価値に対応するインデックス値を更新インデックス値として出力する。複数の閾値は、温度値Ｔと互いに独立な複数の乱数値とに基づいてそれぞれ算出される。

【0134】

選択部４１は、更新エネルギー変化値算出回路４１ａ、閾値生成部４１ｂ、加算回路４１ｃ１，４１ｃ２，…，４１ｃＮ、最小値選択回路４１ｄを有する。
更新エネルギー変化値算出回路４１ａは、エネルギー変化ΔＥ_１～ΔＥ_Ｎのうち、０以上のものをそのまま出力し、０より小さいものを０に更新することで、Ｎ個の更新エネルギー変化値を生成する。更新エネルギー変化値算出回路４１ａの処理は、ｍａｘ［０，ΔＥ_ｉ］を計算することに相当する。

【0135】

閾値生成部４１ｂは、Ｎ個の閾値を生成する。閾値生成部４１ｂの例については後述する。
加算回路４１ｃ１～４１ｃＮは、更新エネルギー変化値算出回路４１ａが出力したＮ個の更新エネルギー変化値のそれぞれに、複数の閾値の何れかを加算したＮ個の評価値を出力する。

【0136】

最小値選択回路４１ｄは、加算回路４１ｃ１～４１ｃＮが出力するＮ個の評価値のうちで、最小となる評価値に対応するインデックス値を更新インデックス値として出力する。なお、インデックス値は、たとえば、エネルギー変化計算部１２から供給される。最小値選択回路４１ｄは、たとえば、ツリー状に複数段、配置された複数の選択回路を有し、各選択回路は２つの評価値のうち、小さい方を選択して出力する。最後段の選択回路から出力される評価値に対応したインデックス値が更新インデックス値となる。

【0137】

図１５は、閾値生成部の一例を示す図である。
閾値生成部４１ｂは、メルセンヌツイスタ５０ａ１，５０ａ２，…，５０ａｎ、変換部５０ｂ１，５０ｂ２，…，５０ｂｎ、レジスタ５０ｃ１１～５０ｃｎｍ、乗算部５０ｄを有する。

【0138】

メルセンヌツイスタ５０ａ１～５０ａｎは、それぞれ０＜ｒ_ｉ＜１の一様乱数である乱数値ｒ_ｉを、インデックス値＝ｉ（１～Ｎ）のそれぞれについて異なるシードを用いて発生する。乱数値ｒ_ｉは、たとえば、１６ビットの値で表される。なお、乱数値ｒ_ｉは、１クロックサイクル毎に更新される。

【0139】

変換部５０ｂ１～５０ｂｎは、変換テーブルを用いて、乱数値ｒ_ｉをｌｏｇ（－ｌｏｇｒ_ｉ）に変換する。ｌｏｇ（－ｌｏｇｒ_ｉ）は、たとえば、２７ビットの値で表される。変換テーブルは、たとえば、ＲＡＭ、フラッシュメモリなどのメモリに記憶されている。

【0140】

レジスタ５０ｃ１１～５０ｃｎｍは、図示しないクロック信号に同期して、変換部５０ｂ１～５０ｂｎが出力した値を遅延して出力する。たとえば、直列に接続されたレジスタ５０ｃ１１，５０ｃ１２，…，５０ｃ１ｍにより、変換部５０ｂ１が出力した値から、ｍ個の異なる乱数値が生成される。また、直列に接続されたレジスタ５０ｃ２１，５０ｃ２２，…，５０ｃ２ｍにより、変換部５０ｂ２が出力した値から、ｍ個の異なる乱数値が生成される。また、直列に接続されたレジスタ５０ｃｎ１，５０ｃｎ２，…，５０ｃｎｍにより、変換部５０ｂｎが出力した値から、ｍ個の異なる乱数値が生成される。

【0141】

乗算部５０ｄは、変換部５０ｂ１～５０ｂｎが出力した値と、レジスタ５０ｃ１１～５０ｃｎｍのそれぞれが出力した値に、温度値Ｔを乗じた値を閾値として出力する。乗算部５０ｄは、Ｎ個の乗算回路を有する。

【0142】

たとえば、Ｎ＝１０２４の場合、メルセンヌツイスタ５０ａ１～５０ａｎと変換部５０ｂ１～５０ｂｎは、それぞれ３２個（ｎ＝３２）、レジスタ５０ｃ１１～５０ｃｎｍの数（ｎ×ｍ）は、３２×３１である。これにより、１０２４個の閾値が生成される。

【0143】

選択部１６の代りに上記の選択部４１を用いることで、試行毎に適切な許容確率で状態遷移が発生する。その理由を以下に説明する。
前述した０＜ｒ_ｉ＜１の一様乱数である乱数値ｒ_ｉの累積分布関数Ｆ（ｒ）＝Ｐｒｏｂ（ｒ_ｉ≦ｒ）の値は、ｒ≦０の場合に０、０＜ｒ＜１の場合にｒ、ｒ≧１の場合に１となる。ここで、乱数値ｒ_ｉから発生させる正の乱数値ｙ_ｉを、ｙ_ｉ＝－ｌｏｇ（ｒ_ｉ）／Ａ_ｉ（Ａ_ｉ＞０）とする。乱数値ｙ_ｉがｙ（＞０）より大きくなる確率Ｐｒｏｂ（ｙ_ｉ≧ｙ）は、Ｐｒｏｂ（ｙ_ｉ≧ｙ）＝Ｐｒｏｂ（ｒ_ｉ≦ｅｘｐ（－Ａ_ｉｙ））＝Ｆ（ｅｘｐ（－Ａ_ｉｙ））＝ｅｘｐ（－Ａ_ｉｙ）と表せる。したがって乱数値ｒ_ｉの確率密度関数ｐ（ｙ_ｉ）は、以下の式（９）のように表せる。

【0144】

【数9】

【0145】

あるインデックス値＝ｉに対する乱数値ｙ_ｉが他の全ての乱数値ｙ_ｊ（ｊ≠ｉ）より小さい確率は、以下の式（１０）のように表せる。

【0146】

【数10】

【0147】

式（１０）の右辺は、式（５）の右辺と同じ形である。
ここで、Ａ_ｉをエネルギー変化ΔＥ_ｉに対する許容確率とし、Ａ_ｉ＝ｍｉｎ［１，ｅｘｐ（－βΔＥ_ｉ）］とすれば、式（１０）が１となるときのＡ_ｉが求めたい許容確率となる。なお、ｌｏｇ（Ａ_ｉ）／β＝－ｍａｘ［０，ΔＥ_ｉ］となる。

【0148】

あるインデックス値＝ｉに対する乱数値ｙ_ｉが他の全ての乱数値ｙ_ｊ（ｊ≠ｉ）より小さいことはｌｏｇ（ｙ_ｉ）／βがｌｏｇ（ｙ_ｊ）／β（ｊ≠ｉ）より小さいことと同じである。ここで、ｌｏｇ（ｙ_ｉ）／β＝ｌｏｇ（－ｌｏｇ（ｒ_ｉ））／β＋ｍａｘ［０，ΔＥ_ｉ］、β＝１／Ｔであるから、ｍａｘ［０，ΔＥ_ｉ］にＴ・ｌｏｇ（－ｌｏｇ（ｒ_ｉ））を加算した値の中で最小のものを選べば、試行毎に適切な許容確率で状態遷移が発生する。

【0149】

その他の構成及び動作については、第１の実施の形態のサンプリング装置１０と同じである。
図１６は、温度とエネルギーとの関係についてのシミュレーション結果の例を示す図である。

【0150】

また、図１７は、試行回数と最低エネルギー状態が得られる確率との関係についてのシミュレーション結果の例を示す図である。図１６、図１７で扱われている問題は、１６都市の巡回セールスマン問題である。

【0151】

図１６では、横軸は温度値Ｔ、縦軸は式（１）に示したエネルギーＥを表している。特性６０ａは、通常のＭＣＭＣ法を用いた場合の温度とエネルギーとの関係を示し、特性６０ｂは、第３の実施の形態のサンプリング装置４０を用いた場合の温度とエネルギーとの関係を示している。また、特性６０ｃは、デジタル回路を用いた従来の最適化装置を用いた場合の温度とエネルギーとの関係を示している。

【0152】

図１６に示すように、第３の実施の形態のサンプリング装置４０は、通常のＭＣＭＣ法を用いた場合とほぼ同様の特性が得られている。
図１７では、横軸は試行回数、縦軸は最低エネルギー状態が得られる確率を表している。特性６１ａは、通常のＭＣＭＣ法を用いた場合の試行回数と上記確率との関係を示し、特性６１ｂは、第３の実施の形態のサンプリング装置４０を用いた場合の試行回数と上記確率との関係を示している。また、特性６１ｃは、デジタル回路を用いた従来の最適化装置を用いた場合の試行回数と上記確率との関係を示している。

【0153】

図１７に示すように、第３の実施の形態のサンプリング装置４０は、従来の最適化装置や通常のＭＣＭＣ法を用いた場合よりも、高速に（早い試行回数で）最低エネルギー状態が得られることがわかる。

【0154】

以上のように、第３の実施の形態のサンプリング装置４０は、理想的に１回の試行で状態遷移が行われ、理論的にも通常のＭＣＭＣ法と同等のサンプル分布が得られるため、高速で高精度のサンプリングが可能となる。

【0155】

（変形例）
以下に示すサンプリング装置は、レプリカ交換法を利用したものである。
レプリカ交換法は複数の温度を用いたＭＣＭＣ法を同時に行い、ある試行回数毎に、それぞれ状態のエネルギーを比較し、適切な確率で２つの温度に対する状態（または温度）を交換するという操作を行う方法である。以下の例では、温度を交換するものとして説明するが、状態を交換しても同じことである。

【0156】

図１８は、レプリカ交換法を利用したサンプリング装置の例を示す図である。
図１８に示すサンプリング装置７０は、サンプリング処理部７１とレプリカ交換コントローラ７２とを有する。

【0157】

サンプリング処理部７１は、図１に示したサンプリング装置１０の各部を有しており、パイプライン処理によって、異なる温度値Ｔ（図１８の例では逆温度β_１，β_２，…，β_Ｍ）が設定されたレプリカ７１ａ１，７１ａ２，…，７１ａＭのそれぞれについて、前述の試行を行う。

【0158】

状態ｘやエネルギーＥは、レプリカ７１ａ１～７１ａＭのそれぞれについて保持され、試行が行われるたびに更新される。
レプリカ交換コントローラ７２は、所定試行回数毎にレプリカ７１ａ１～７１ａＭのエネルギーＥを比較して、その比較結果に基づいてレプリカ７１ａ１～７１ａＭのうちの２つにおいて設定されている逆温度（または温度）を、所定の交換確率ｐ_ｉｊで交換する。交換確率ｐ_ｉｊは、ｐ_ｉｊ＝ｆ（（β_ｉ－β_ｊ）（Ｅ_ｉ－Ｅ_ｊ））と表せる。β_ｉはｉ番目のレプリカに設定された逆温度、β_ｊはｊ番目のレプリカに設定された逆温度、Ｅ_ｉはｉ番目のレプリカにおけるエネルギー、Ｅ_ｊはｊ番目のレプリカにおけるエネルギーである。また、関数ｆは式（３）のβの代りにβ_ｉ－β_ｊが用いられ、式（３）のΔＥの代りにＥ_ｉ－Ｅ_ｊが用いられる。

【0159】

レプリカ交換コントローラ７２は、たとえば、ＡＳＩＣやＦＰＧＡなどの特定用途の電子回路、または、ＣＰＵやＤＳＰなどのプロセッサなどである。
なお、上記のレプリカ交換コントローラ７２の処理を、図１に示した制御部１７が行ってもよい。

【0160】

サンプリング動作は、所定の試行回数毎に行われ、たとえば、温度値Ｔが最低のレプリカについての状態ｘと計数値Ｎ_ｆとオフセット値Ｅ_ｏｆｆが出力される。
図１９は、レプリカ交換法を利用したサンプリング装置の動作例を示すフローチャートである。なお、図１０に示した第１の実施の形態のサンプリング装置１０の動作と同じ処理については、図示及び説明を省略する。

【0161】

レプリカ交換法を利用する場合、初期化（ステップＳ５０）の際に、レプリカ交換コントローラ７２に、最高逆温度β_ｍａｘ、最低逆温度β_ｍｉｎ、交換間隔を示す試行回数Ｎ_{ｉ＿ｔｏｔ}、総サンプリング回数Ｎ_{ｓ＿ｔｏｔ}が設定される。レプリカ交換コントローラ７２は、これらのパラメータを、たとえば、図示しない記憶部から読み出してもよい。また、初期化の際に、レプリカ交換コントローラ７２は、最高逆温度β_ｍａｘ、最低逆温度β_ｍｉｎに基づいて、レプリカ７１ａ１～７１ａＭに逆温度β_１～β_Ｍの初期値を設定する。

【0162】

次に、たとえば、レプリカ交換コントローラ７２は、図示しないカウンタにサンプリング回数Ｎ_ｓとして０を設定し（ステップＳ５１）、別の図示しないカウンタに試行回数Ｎ_ｉ＝０を設定する（ステップＳ５２）。その後、試行が行われ各レプリカの状態ｘとエネルギーＥが更新される（ステップＳ５３）。レプリカ交換コントローラ７２は、試行回数Ｎ_ｉを１増加させ（ステップＳ５４）、Ｎ_ｉ≧Ｎ_{ｉ＿ｔｏｔ}であるか否かを判定する（ステップＳ５５）。Ｎ_ｉ≧Ｎ_{ｉ＿ｔｏｔ}ではない場合には、ステップＳ５３からの処理が繰り返される。

【0163】

Ｎ_ｉ≧Ｎ_{ｉ＿ｔｏｔ}である場合、レプリカ交換コントローラ７２は、前述の交換確率を用いてレプリカ間の温度（逆温度β）を交換するとともに、サンプリング回数Ｎ_Ｓを１増加させる（ステップＳ５６）。そして、レプリカ交換コントローラ７２は、温度値Ｔが最低（逆温度βが最高）のレプリカについての計数値Ｎ_ｆと状態ｘとオフセット値Ｅ_ｏｆｆを出力させる（ステップＳ５７）。また、出力された計数値Ｎ_ｆとオフセット値Ｅ_ｏｆｆとに基づいて、図１に示した試行回数計算部１８は、期待値＜Ｎ_{ｔｒｉａｌ}＞を計算する。

【0164】

その後、レプリカ交換コントローラ７２は、Ｎ_ｓ≧Ｎ_{ｓ＿ｔｏｔ}であるか否かを判定する（ステップＳ５８）。Ｎ_ｓ≧Ｎ_{ｓ＿ｔｏｔ}ではない場合には、ステップＳ５２からの処理が繰り返される。

【0165】

Ｎ_ｓ≧Ｎ_{ｓ＿ｔｏｔ}である場合、レプリカ交換コントローラ７２は、サンプリング動作を終了させる。
このようなレプリカ交換法を利用したサンプリング装置７０は、局所解からの脱出が早く、温度値Ｔの最適化が容易であるという利点がある。

【0166】

次に、サンプリング装置の他の変形例を示す。
図２０は、サンプリング装置の他の変形例を示す図である。
図２０に示すサンプリング装置８０は、サンプリング処理部８１と期待値計算部８２を有する。

【0167】

サンプリング処理部８１は、図１に示したサンプリング装置１０の各部を有しており、所定の試行回数毎に期待値＜Ｎ_{ｔｒｉａｌ}＞と状態ｘを出力する。
期待値計算部８２は、期待値＜Ｎ_{ｔｒｉａｌ}＞と状態ｘに基づいて、サンプル平均により求めたい量ｆ_ｎ＝ｆ（ｘ（ｉ））の期待値を計算し、出力する。関数ｆ（ｘ（ｉ））は、たとえば、ｘ_ｉ、ｘ_ｉ×ｘ_ｊなどである。期待値は、以下の式（１１）で表せる。

【0168】

【数11】

【0169】

式（１１）において、ｎはサンプリング回数を示し、Ｍは総サンプリング回数を示す。また、τ_ｎは、ｎ回目のサンプリングにおける期待値＜Ｎ_{ｔｒｉａｌ}＞を示す。
なお、上記の期待値計算部８２の処理を、図１に示した制御部１７が行ってもよい。

【0170】

図２１は、期待値計算処理を行うサンプリング装置の動作例を示すフローチャートである。なお、図１０に示した第１の実施の形態のサンプリング装置１０の動作と同じ処理については、図示及び説明を省略する。

【0171】

期待値計算処理が行われる場合、初期化（ステップＳ６０）の際に、期待値計算部８２は、式（１１）のｆ_ｎτ_ｎの積算値を示す変数Ｓと、τ_ｎの積算値を示す変数ＴＡを０にする。また、期待値計算部８２に総サンプリング回数Ｍが設定される。さらに、サンプリング回数ｎを計数する図示しないカウンタの値が０にリセットされる。

【0172】

その後、サンプリングタイミングにおいて、期待値計算部８２は、サンプリング処理部８１が出力した状態ｘと期待値＜Ｎ_{ｔｒｉａｌ}＞とを読み込み（ステップＳ６１）、ｆ_ｎ＝ｆ（ｘ（ｉ））を計算する（ステップＳ６２）。その後、期待値計算部８２は、現在の変数Ｓにｆ_ｎτ_ｎを加算することで変数Ｓを更新し（ステップＳ６３）、現在の変数ＴＡにτ_ｎを加算することで変数ＴＡを更新する（ステップＳ６４）。そして、期待値計算部８２は、図示しないカウンタにサンプリング回数ｎを１増加させ（ステップＳ６５）、ｎ＞Ｍであるか否かを判定する（ステップＳ６６）。ｎ＞Ｍではない場合、期待値計算部８２は、ステップＳ６１からの処理を繰り返す。

【0173】

ｎ＞Ｍの場合、期待値計算部８２は、Ｓ／ＴＡを計算し、式（１１）で示される期待値として出力し（ステップＳ６７）、期待値計算処理を終える。
上記のような期待値計算処理を行う期待値計算部８２は、たとえば、積和演算回路や記憶部（レジスタなど）などを用いて実現できる。また、期待値計算部８２は、たとえば、ＡＳＩＣやＦＰＧＡなどの特定用途の電子回路、ＣＰＵやＤＳＰなどのプロセッサなどであってもよい。

【0174】

期待値計算処理は比較的演算量が大きいため、期待値計算処理をサンプリング装置８０の外部で行う場合、サンプリング装置８０と外部装置間での状態変数ｘなどの入出力によるボトルネックが生じる可能性がある。しかし、上記のようにサンプリング装置８０内で期待値計算処理を行うことで、このボトルネックを解消できる利点がある。

【0175】

なお、上記各変形例は、第１乃至第３の実施の形態のサンプリング装置１０，３０，４０のそれぞれと組合せることができる。
以上、実施の形態に基づき、本発明のサンプリング装置及びサンプリング装置の制御方法の一観点について説明してきたが、これらは一例にすぎず、上記の記載に限定されるものではない。

【符号の説明】

【0176】

１０ サンプリング装置
１１ 状態保持部
１２ エネルギー変化計算部
１３ オフセット制御部
１４ 比較部
１５ フラグ計数部
１６ 選択部
１７ 制御部
１８ 試行回数計算部

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】

【図12】

【図13】

【図14】

【図15】

【図16】

【図17】

【図18】

【図19】

【図20】

【図21】