(57)【要約】【課題】トンネル効果を用いて連続変数の最適化問題を高精度に解けるようにした連続最適化問題の大域的探索装置及びプログラムを提供する。【解決手段】コンピュータは、評価関数Ｖ（）の微小変化に沿う勾配法により連続変数ｘを更新し（Ｓ４）、ボルツマン分布に従って調和振動子の固有状態を選択し（Ｓ５）、選択された固有状態の存在確率を用いて当該固有状態の値を連続的なノイズとして連続変数ｘに加算し（Ｓ９ａ）、ノイズが加算された連続変数ｘを用いて勾配法による更新を繰り返す（Ｓ４〜Ｓ１０）。【選択図】図１０
2019159782.tif 000002
【特許請求の範囲】
【請求項１】
連続変数を用いて生成された評価関数が最小値又は最大値となる条件を満たす最適解を探索する連続最適化問題の大域的探索装置（１）であって、
前記評価関数の微小変化に沿う勾配法により前記連続変数を更新する更新部（６）と、
ボルツマン分布に従って調和振動子の固有状態を選択する選択部（７）と、
前記選択された前記固有状態の存在確率を用いて前記固有状態の値を連続的なノイズとして前記連続変数に加算する加算部（９）と、を備え、
前記更新部は、前記加算部によりノイズが加算された前記連続変数を用いて勾配法による更新を繰り返す連続最適化問題の大域的探索装置。
【請求項２】
連続変数を用いて生成された評価関数が最小値又は最大値となる条件を満たす最適解を探索する連続最適化問題の大域的探索装置であって、
前記評価関数の微小変化に沿う勾配法により前記連続変数を更新する更新部（６）と、
ボルツマン分布に従って調和振動子の固有状態を選択し、当該選択された前記固有状態の存在確率がピークとなる条件を満たす値を離散ノイズとしてランダムに選択する選択部（７）と、
前記離散ノイズを加算する前後のエネルギー差を計算し、前記評価関数に依存して予め設定される温度に依存した確率で受理するか否かを判定する判定部（８）と、
前記受理されないときには前記離散ノイズを０とし前記受理されたときには前記選択部により選択された前記離散ノイズを前記連続変数に加算する加算部（９）と、を備え、
前記更新部は、前記加算部により前記離散ノイズが加算された前記連続変数を用いて勾配法により更新を繰り返す連続最適化問題の大域的探索装置。
【請求項３】
初期状態では前記温度を所定の変数値とし、前記更新部が前記連続変数の更新を繰り返すときに、前記温度を徐々に減少させ（Ｓ１０ａ）、前記加算部は前記連続変数に前記離散ノイズと共に前記温度に依存するガウスノイズを加算する（Ｓ９ａ）請求項２記載の連続最適化問題の大域的探索装置。
【請求項４】
前記選択部が、前記調和振動子の固有状態を選択するときには、一定範囲の励起状態の中から選択する請求項１から３の何れか一項に記載の連続最適化問題の大域的探索装置。
【請求項５】
前記選択部が、前記調和振動子の固有状態を選択するときには、エネルギーが最低の基底状態から有限個の励起状態を選定しこの中から選択する請求項４記載の連続最適化問題の大域的探索装置。
【請求項６】
前記選択部は、前記調和振動子の固有状態として常に第一励起状態を選択する請求項１から５の何れか一項に記載の連続最適化問題の大域的探索装置。
【請求項７】
連続変数を用いて生成された評価関数が最小値又は最大値となる条件を満たす最適解を探索するプログラムであって、
連続最適化問題の大域的探索装置（１）に、
前記評価関数の微小変化に沿う勾配法により前記連続変数を更新する手順（Ｓ４）と、
ボルツマン分布に従って調和振動子の固有状態を選択する手順（Ｓ５）と、
前記選択された前記固有状態の存在確率を用いて前記固有状態の値を連続的なノイズとして前記連続変数に加算する手順（Ｓ９ａ）と、実行させると共に、
前記ノイズが加算された前記連続変数を用いて前記勾配法による更新を繰り返すように実行させるプログラム。
【請求項８】
連続変数を用いて生成された評価関数が最小値又は最大値となる条件を満たす最適解を探索するプログラムであって、
連続最適化問題の大域的探索装置（１）に、
前記評価関数の微小変化に沿う勾配法により前記連続変数を更新する手順（Ｓ４）と、
ボルツマン分布に従って調和振動子の固有状態を選択し、当該選択された前記固有状態の存在確率がピークとなる条件を満たす値を離散ノイズとしてランダムに選択する手順（Ｓ６）と、
前記離散ノイズを加算する前後のエネルギー差を計算し、前記評価関数に依存して設定される温度に依存した確率で受理するか否かを判定する手順（Ｓ８）と、
前記受理されないときには前記離散ノイズを０とし前記受理されたときには前記選択された前記離散ノイズを前記連続変数に加算する手順（Ｓ９）と、を実行させると共に、
前記離散ノイズが加算された前記連続変数を用いて前記勾配法による更新を繰り返すように実行させるプログラム。
【請求項９】
初期状態では前記温度を所定の変数値とし、前記連続変数を更新を繰り返すときに前記温度を徐々に減少させ（Ｓ１０ａ）、前記離散ノイズと共に前記温度に依存するガウスノイズを加算する手順（Ｓ９ａ）、を備える請求項８記載のプログラム。
【請求項１０】
前記調和振動子の固有状態を選択するときには、一定範囲の励起状態の中から選択する請求項７から９の何れか一項に記載のプログラム。
【請求項１１】
前記調和振動子の固有状態を選択するときには、エネルギーが最低の基底状態から有限個の励起状態を選定しこの中から選択する請求項１０記載のプログラム。
【請求項１２】
前記調和振動子の固有状態として常に第一励起状態を選択する請求項７から１１の何れか一項に記載のプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、連続最適化問題の大域的探索装置及びプログラムに関する。
【背景技術】
【０００２】
局所的な最小点を複数個備えた誤差関数を適用し、量子トンネル効果を利用して大域的な最小点を検索する試みがなされている（例えば、特許文献１参照）。この特許文献１記載の技術によれば、力学システムのプランク定数に相当する部分に置換されるべき設計定数と、力学システムの時間発展を規定する時間微分方程式における該力学システムの散面の変化率を規定する摩擦係数とを具備し、量子力学的トンネル効果を使うことにより真の最小値を計算するようにしている。
【０００３】
しかしながら、この特許文献１記載の技術によれば、量子力学的トンネル効果を用いることについて言及されているものの、具体的に、どのようにして変数を更新して局所解を脱出するのかについて記載も示唆もない。例えば、最適化変数を物理系の座標、評価関数を物理系のポテンシャル関数と見做したとき、どのような量子揺らぎによりトンネル効果を発生させるか記載も示唆もなされておらず、また、力学システムの時間発展を規定する時間微分方程式の具体的な形も示唆されていない。すなわち、特許文献１記載の内容は、トンネル効果を利用することを提起していることに留まるものであり、当該トンネル効果を利用した最適化方法について何ら示されているものではない。
【先行技術文献】
【特許文献】
【０００４】
【特許文献１】特開２００６−５９２３７号公報
【発明の概要】
【発明が解決しようとする課題】
【０００５】
本発明の目的は、トンネル効果を用いて連続変数の最適化問題を高精度に解けるようにした連続最適化問題の大域的探索装置及びプログラムを提供することにある。
【課題を解決するための手段】
【０００６】
請求項１又は７に記載した発明は、連続変数を用いて生成された評価関数が最小値又は最大値となる条件を満たす最適解を探索する連続最適化問題の大域的探索装置を対象としている。この請求項１又は７に記載した発明によれば、評価関数の微小変化に沿う勾配法により前記連続変数を更新し、ボルツマン分布に従って調和振動子の固有状態を選択し、選択された固有状態の存在確率を用いて当該固有状態の値を連続的なノイズとして連続変数に加算し、ノイズが加算された連続変数を用いて勾配法による更新を繰り返している。このため、連続的なノイズを連続変数に加算することで、トンネル効果を用いて局所解を脱出できるようになり、連続変数の最適化問題を高精度に解くことができる。
【０００７】
また、請求項２又は８に記載した発明によれば、評価関数の微小変化に沿う勾配法により連続変数を更新し、ボルツマン分布に従って調和振動子の固有状態を選択し、当該選択された固有状態の存在確率がピークとなる条件を満たす値を離散ノイズとしてランダムに選択し、離散ノイズを加算する前後のエネルギー差を計算し、評価関数に依存して予め設定される温度に依存した確率で受理するか否かを判定し、受理されないときには離散ノイズを０とし受理されたときには選択された離散ノイズを連続変数に加算し、離散ノイズが加算された連続変数を用いて勾配法により更新を繰り返している。このため、離散ノイズを連続変数に加算することで、トンネル効果を用いて局所解を脱出できるようになり、連続変数の最適化問題を高精度に解くことができる。
【図面の簡単な説明】
【０００８】
【図１Ａ】第１実施形態を示す電気的構成図
【図１Ｂ】機能的構成図
【図２】評価関数の例
【図３】量子力学的な調和振動子の固有値と固有状態を表す図
【図４】固有状態のピーク位置
【図５】最適解の導出処理内容を示すフローチャート
【図６】勾配法の処理イメージを表す説明図
【図７】トンネル効果による局所解の脱出イメージを表す説明図
【図８】第２実施形態を示す最適解の導出処理内容を示すフローチャート
【図９】シミュレーテッドアニーリング法による処理イメージを表す説明図
【図１０】第３実施形態における最適解の導出処理内容を示すフローチャート
【発明を実施するための形態】
【０００９】
以下、本発明の連続最適化問題の大域的探索装置、及びプログラムの幾つかの実施形態について図面を参照して説明する。以下の実施形態中では、各実施形態間において、同一機能又は類似機能を備えた部分に同一符号又は類似符号（例えば添え字「ａ」）を付して説明を行い、同一又は類似機能を備えた連携動作説明を必要に応じて省略する。
【００１０】
（第１実施形態）
図１Ａから図７は、第１実施形態の説明図を示している。図１Ａに示す装置１は、量子力学的な性質を利用して最適化問題の最適化処理のシミュレーションを実行する連続最適化問題の大域的探索装置として構成される。
【００１１】
この装置１は、ＣＰＵ２と、ＲＯＭ、ＲＡＭ等のメモリ３と、入出力インタフェース４とをバス接続した汎用のコンピュータ５を用いて構成される。このコンピュータ５は、ＣＰＵ２によってメモリ３に記憶された変換プログラムを実行し、各種手順を実行することで大域的探索処理を実行する。メモリ３は、非遷移的実体的記録媒体として用いられる。
【００１２】
コンピュータ５が実行する大域的探索処理は、１以上のＮ次元を備えたユークリッド空間からなる探索空間を想定し、この探索空間の中で、複数の要求や制約によって生成された評価関数Ｖ（）が最小値となる条件を満たす連続変数ｘ、すなわち最適解（図２のＡ３）を求める処理である。図１Ｂに示すように、コンピュータ５は、その実現する機能として、更新部６、選択部７、判定部８、及び加算部９としての各種機能を備える。
【００１３】
評価関数Ｖ（）は、例えば図２に示すように、複数の要求や制約によって生成されており、１以上のＮ個の連続変数ｘ、すなわちパラメータを用いて生成された数式による関数を示すものであり、例えば任意の多項式、有理関数、無理関数、指数関数、対数関数やその加減乗除等による組み合わせなどを挙げることができる。
【００１４】
図２に示すように、評価関数Ｖ（）は、連続変数ｘに応じて変化する関数であり、数々の極小値を含む関数である。この条件下において、コンピュータ５は、評価関数Ｖ（）の極小値の中でもその最低値を満たす連続変数ｘの最適解Ａ３として求めることになるが、評価関数Ｖ（）が極小値となる条件を満たす連続変数ｘの局所解Ａ１、Ａ２、Ａ４が多数存在する。このため、コンピュータ５がこの問題を解いても局所解Ａ１、Ａ２、Ａ４に陥りやすい。このため、本実施形態においては、コンピュータ５は、量子力学的なトンネル効果を用いて局所解Ａ１、Ａ２、Ａ４を脱出して最適解Ａ３を求めるようにしている。
【００１５】
＜量子揺らぎの概念の導入＞
評価関数Ｖ（）の評価値Ｖ（ｘ）に量子力学的なトンネル効果を生じさせることで局所解（例えばＡ４）を脱出するため、本実施形態において、量子揺らぎの概念を導入する。本実施形態では、量子アニーリングのハミルトニアンＨ＾（ｍ）を、下記の（１）式に示すように与える。ｍは質量を示す。
【数１】
2019159782.tif 000003
【００１６】
この（１）式では、ｘを連続変数とすると共に、評価関数Ｖ（）をポテンシャルとした量子アニーリングを導入している。また（１）式の右辺第２項は、運動量ｐの演算子ｐ＾による量子揺らぎの導入項を示している。この（１）式において、初期状態では質量ｍを十分小さい値に設定することで、量子揺らぎの導入項、すなわち（１）式の右辺第２項の影響を強くすることが望ましい。そして、探索処理を進めるに連れて質量ｍを大きくすることで（１）式の右辺第１項の評価関数Ｖ（）の影響を強めると共に、右辺第２項の量子揺らぎの導入項の影響を弱めるようにすると良い。すると、探索当初は量子揺らぎの影響を受けて例えば大域的に連続変数ｘが移動するようになり、探索処理を進めるに連れて評価関数Ｖ（）の影響を大きく受け例えば局所的に最適解Ａ３を導出できるようになる。
【００１７】
＜連続変数ｘの更新処理＞
連続変数ｘの更新処理を量子系に適用する場合には、時間に依存するシュレーディンガー方程式により記述できるが、シュレーディンガー方程式を解くのは膨大な計算が必要なため非現実的である。このため、量子アニーリングの性能を評価するときには、シュレーディンガー方程式を直接解く方法を採用することは稀であり、実際には例えば温度Ｔを極低温条件とした平衡状態を求めることが望ましい。連続変数ｘの更新処理においては、平衡状態に収束するように実行する。モンテカルロ法を用いて平衡状態に収束するように計算処理を行うことで、シュレーディンガー方程式を解くよりはるかに少ない計算量で最適化変数を更新できるようになる。
【００１８】
＜分配関数と量子揺らぎの解釈＞
（１）式のハミルトニアンＨ＾（ｍ）を用いると、（２）式のように分配関数を表すことができる。
【数２】
2019159782.tif 000004
【００１９】
この（２）式において、βは熱ノイズ（＝１／Ｔ）を示す。また、Ｔｒはトレースを示しており行列の対角和を表している。そして、この（２）式を変数分離すると、（３）式のように分配関数を表すことができる。この（３）式では、ｋを無限大にしてその指数関数ｅｘｐの中身の値の極限値を取得することで｜｜ｚ−ｘ｜｜＾２のＬ２ノルムにより等式制約化していることになる。
【数３】
2019159782.tif 000005
【００２０】
この（３）式の分配関数は、評価関数Ｖ（ｚ）と、ｚを中心とする量子力学的な調和振動子との和と解釈できる。このことから、連続変数ｘを更新処理するときに、評価関数Ｖ（）の微小変化に沿って連続変数ｘを更新する勾配法を適用しつつ、量子力学的な調和振動子によるノイズ成分を加算することで、量子トンネル効果により連続変数ｘが局所解Ａ１、Ａ２、Ａ４を脱出し、最適解Ａ３に至らせることが可能になる。
【００２１】
＜量子力学的な調和振動子の説明＞
量子力学的な調和振動子の固有値と固有状態を図３に示す。各固有状態の第ｎ励起状態の曲線は、各状態の存在確率Ｐｃを表している。基底状態を第ゼロ励起状態と定義すればｎ≧０を満たす。また、基底状態〜第三励起状態において、存在確率Ｐｃがピーク条件を満たす調和振動子のｚ−ｘの位置を図４に示す。
【００２２】
すなわち図４に示すように、基底状態であればピーク条件を満たす位置は０である。また、第一励起状態において存在確率Ｐｃがピーク条件を満たす位置は（４）式である。
【数４】
2019159782.tif 000006
ここで、ｍは質量、ｋはばね定数を示す。第二励起状態において存在確率Ｐｃがピーク条件を満たす位置は、下記の（５）式である。
【数５】
2019159782.tif 000007
また、第三励起状態において存在確率Ｐｃがピーク条件を満たす位置は、下記の（６）式である。
【数６】
2019159782.tif 000008
【００２３】
＜調和振動子の固有状態を選択＞
このような調和振動子の固有状態を考慮した上で、温度Ｔ（＝１／β）のボルツマン分布に従って固有状態を所定の確率で選択すると良い。このボルツマン分布によれば、第ｎ励起状態（ｎ≧０）の選択確率Ｐosc（ｎ）は、下記の（７−１）式により表すことができる。ここで、Ｚoscは（７−２）式により表すことができる。
【数７】
2019159782.tif 000009
【００２４】
調和振動子の固有状態は、理論上、無限個存在することになるが、全固有状態を考慮すると、必要な精度に対して計算量が大幅に増えてしまうため、必要な精度を考慮して一定範囲の励起状態の中から選択するすることが望ましい。また、さらにはエネルギーの最低の基底状態から有限個Ｎoscの励起状態を選定し、この中から選択することが望ましい。
【００２５】
＜調和振動子に基づく離散ノイズΔquantumの加算方法＞
調和振動子によるノイズは、（７−１）式、及び（７−２）式のボルツマン分布により固有状態を選択した後、選択した固有状態の存在確率Ｐｃがピークとなる条件を満たす値を離散ノイズΔquantumとして連続変数ｘに加算することが望ましい。図３に示すように、存在確率Ｐｃがピークとなる条件を満たす値以外にも高い確率条件を満たす値も存在するが、ピークとなる条件だけ考慮することで計算量を少なくできるためである。しかも、離散ノイズΔquantumを加えることで、トンネル効果により局所解Ａ１、Ａ２、Ａ４を容易に脱出できるようになる。
【００２６】
＜最適解Ａ３の導出方法＞
以下では、このような技術的意義の下で、コンピュータ５が実際に最適解Ａ３を導出するための実際の方法について説明する。図５は最適解Ａ３の導出処理内容をフローチャートにより概略的に示している。
【００２７】
コンピュータ５は、図５のＳ１において温度Ｔとばね定数ｋとを定数として初期設定し、Ｓ２において質量ｍを変数として初期設定する。温度Ｔとばね定数ｋは、評価関数Ｖ（）に依存して定まるパラメータであるため、例えばシミュレーションを用いて予め定数として算出することが望ましい。また初期状態では、質量ｍを予め小さい所定の変数値に設定することが望ましい。
さらにコンピュータ５は、Ｓ３において連続変数ｘの初期値を例えばランダムに設定する。そしてコンピュータ５は、評価関数Ｖ（）に連続変数ｘの初期値を代入して評価値Ｖ（ｘ）を算出し、その後、Ｓ４において勾配法を用いて連続変数ｘを更新する。勾配法では、下記の（８）式に示すように、評価関数Ｖ（）の微小変化に沿って連続変数ｘを更新することが望ましい。
【数８】
2019159782.tif 000010
【００２８】
ここで、ηは勾配法で用いられる所定の係数を示しており、ｘは更新前の連続変数を示し、ｘ＾＊は勾配法による更新後の連続変数を示している。図６は、勾配法による連続変数ｘの更新イメージを示している。この図６に示すように、評価関数Ｖ（）の勾配に沿って低下する方向に連続変数ｘを更新することになる。この後、コンピュータ５は、Ｓ５においてボルツマン分布に従って調和振動子の固有状態として第ｎ励起状態を選択する。このとき、前述の（７−１）式、（７−２）式のボルツマン分布に従って第ｎ励起状態を選択する。
【００２９】
前述したように、調和振動子の固有状態は、理論上、無限個存在することになるが、全固有状態を考慮すると、必要な精度に対して計算量が大幅に増えてしまうため、必要な精度を考慮して、一定範囲の第ｎ励起状態の中から選択するすることが望ましい。また、さらには、エネルギーの最低の基底状態から有限個Ｎoscの励起状態を選定し、この中から選択することが望ましい。すると計算量を削減できる。
例えば、コンピュータ５が、Ｓ５において第一励起状態、すなわちｎ＝１を選択したときには、Ｓ６において当該第一励起状態の（４）式が示す２つのピークとなる条件を満たす値のうち何れかの値をランダムに選択し、離散ノイズΔquantumとする。このとき、コンピュータ５は、選択すべき複数のピークを互いに同一確率、この場合５０％の確率で選択し、選択された値を離散ノイズΔquantumとする。その後、コンピュータ５は、Ｓ７において連続変数ｘに離散ノイズΔquantumを加算する前後のエネルギー変化ΔＶを下記の（９）式のように計算する。
【数９】
2019159782.tif 000011
【００３０】
そしてコンピュータ５は、このエネルギー変化ΔＶについて、評価関数Ｖ（）に依存して設定される温度Ｔに依存した確率で受理判定を行うと良い。この受理判定方法は、メトロポリス法を用いても良いし熱浴法を用いても良い。例えば、メトロポリス法を用いる場合には、コンピュータ５は、例えばΔＶ≦０であるときに１００％受理し、ΔＶ＞０であるときに温度Ｔに依存した例えばｅｘｐ（−ΔＶ／Ｔ）の確率で受理し、その他の場合、破棄する。コンピュータ５が、この内容を受理した場合にはＳ８でＹＥＳと判定し、連続変数ｘに離散ノイズΔquantumを加算し更新する。
【００３１】
そしてコンピュータ５は、Ｓ１０において質量ｍを大きくする。質量ｍが大きくなると、（１）式はその右辺第１項の評価関数Ｖ（）の影響が強くなり、同時に右辺第２項の量子揺らぎの導入項の影響が弱くなる。
【００３２】
この後、コンピュータ５は、これらのＳ４〜Ｓ１０の処理を繰り返すが、質量ｍを増加させながらこれらのＳ４〜Ｓ１０の処理を繰り返すため、（１）式の右辺第１項に相当する評価関数Ｖ（）の影響を徐々に強くしながら、（１）式の右辺第２項に示す量子揺らぎの導入項の影響を徐々に弱めることができる。
【００３３】
その後、コンピュータ５は、Ｓ１１において終了条件を満たしたときに最適化したと見做してＳ１２において解を出力して処理を終了する。Ｓ１１の終了条件としては、Ｓ１０にて徐々に増加している質量ｍが上限値に達することを条件としても良いし、処理を開始してから所定時間経過したことを条件としても良いし、Ｓ４〜Ｓ１０の処理を所定回数以上繰り返したことを条件としても良いし、又は、Ｓ７で算出したエネルギー変化ΔＶが所定値以下となる条件を満たすことを条件としても良い。すなわちＳ１１の終了条件は様々な条件を適用できる。
【００３４】
＜技術的イメージの説明＞
コンピュータ５が、Ｓ４において勾配法により連続変数ｘを更新すると、図６に技術的イメージを示したように、評価関数Ｖ（）が低下する方向にだけ連続変数ｘを更新することになる。このため、一旦、図６に示す局所解Ａ４に嵌ると当該局所解Ａ４から抜け出すことができない。しかしながら、コンピュータ５が、本実施形態のＳ５〜Ｓ１０の処理を実行し、Ｓ８において受理判定されれば、連続変数ｘに離散ノイズΔquantumを加算したエネルギー変化ΔＶに基づくトンネル効果を発生させることができ、図７にそのイメージを示すように、トンネル効果により局所解Ａ４を脱出でき、さらに勾配法を繰り返すことで最適解Ａ３に導かれるようになる。特に量子揺らぎによるトンネル効果を模擬することで、鋭く深い局所解Ａ４に嵌った場合においてもこの局所解Ａ４を効率的に抜け出すことができる。
【００３５】
＜本実施形態のまとめ、効果＞
以上説明したように、本実施形態によれば、コンピュータ５は、評価関数Ｖ（）の微小変化に沿う勾配法により連続変数ｘを更新し、ボルツマン分布に従って調和振動子の固有状態を選択し、当該選択された第ｎ励起状態の存在確率Ｐｃがピークとなる条件を満たす値を離散ノイズΔquantumとしてランダムに選択し、離散ノイズΔquantumを加算する前後のエネルギー差を計算し、評価関数Ｖ（）に依存して予め設定される温度Ｔに依存した確率で受理するか否かを判定し、受理されないときには離散ノイズΔquantumを０とし受理されたときには選択された離散ノイズΔquantumを連続変数ｘに加算し、離散ノイズΔquantumが加算された連続変数ｘを用いて勾配法により更新を繰り返すようにしている。このため、トンネル効果を用いて局所解Ａ１、Ａ２、Ａ４を脱出して最適解Ａ３を導出できるようになり、連続変数ｘの最適化問題を高精度に解けるようになる。
【００３６】
（変形例）
前述では、コンピュータ５が、Ｓ５において（７−１）式、（７−２）式のボルツマン分布に従って固有状態を選択する形態を示したが、この確率的な選択処理に代えて、調和振動子の固有状態として常に第一励起状態を選択するようにしても良い。このとき、固有状態を選択するための計算量を削減しながら、離散ノイズΔquantumのトンネル効果を用いて局所解Ａ１、Ａ２、Ａ４を脱出でき、連続変数ｘの最適化問題を高精度に解ける。
【００３７】
（第２実施形態）
図８は、第２実施形態の追加説明図を示している。第２実施形態が、第１実施形態と異なるところは、シミュレーテッド・アニーリング法を適用したところにある。また、温度Ｔを変数としつつ、離散ノイズΔquantumに併せてガウスノイズΔthermalを加算したところにある。第１実施形態と同一部分には同一符号を付して説明を省略し、以下、異なる部分について説明する。
【００３８】
図８は、最適解Ａ３の導出処理内容をフローチャートにより示している。コンピュータ５は、図８のＳ１ａに示すようにばね定数ｋを定数として設定し、Ｓ２ａに示すように、質量ｍと温度Ｔを変数として初期設定する。本実施形態では、ばね定数ｋが、評価関数Ｖ（）に依存して定まるパラメータであるため、例えばシミュレーションを用いて予め定数として算出することが望ましい。
また初期状態では、質量ｍを予め小さい所定の変数値に設定し、温度Ｔを予め高い所定値に設定すると良い。その後、コンピュータ５は、Ｓ３において連続変数ｘの初期値を例えばランダムに設定する。そしてコンピュータ５は、評価関数Ｖ（）に連続変数ｘの初期値を代入して評価値Ｖ（ｘ）を算出し、その後、Ｓ４において勾配法を用いて連続変数ｘを更新する。勾配法は第１実施形態で説明した方法と同様であるため説明を省略する。本実施形態では、コンピュータ５は、Ｓ４ａにおいて更新した連続変数ｘにガウスノイズΔthermalを加算する。ここで、このガウスノイズΔthermalは、下記の（１０）式のように表すことができる。
【数１０】
2019159782.tif 000012
この（１０）式において、Ｔは温度、ηは勾配法の係数、Ｎ（０，１）は平均０，分散１のガウス分布を示している。
【００３９】
その後、コンピュータ５は、Ｓ５においてボルツマン分布に従って調和振動子の固有状態を所定の確率で選択する。このときコンピュータ５は、例えば前述の（７−１）式、（７−２）式に示すボルツマン分布に従って固有状態を選択すると良い。コンピュータ５が、Ｓ５において例えば第一励起状態を選択したときには、Ｓ６において当該第一励起状態の（４）式が示す２つのピークのうち何れかのピークをランダムに選択する。このときコンピュータ５は、選択すべき複数のピークを互いに同一確率、この場合、５０％の確率で選択し、選択された値を離散ノイズΔquantumとする。
【００４０】
その後、コンピュータ５は、Ｓ７において連続変数ｘに離散ノイズΔquantumを加算する前後のエネルギー変化ΔＶを（９）式のように計算し、前述実施形態と同様にＳ８において受理判定する。すなわち、勾配法により更新された直後の連続変数ｘをｘ＾＊とすれば、離散ノイズΔquantumを加算する前後のエネルギー変化ΔＶは例えば下記の（１１）式のように計算される。
【数１１】
2019159782.tif 000013
【００４１】
その後、コンピュータ５は、このエネルギー変化ΔＶについて、温度Ｔに依存した確率で受理判定を行う。この受理判定方法は、メトロポリス法を用いても良いし熱浴法を用いても良い。例えばメトロポリス法を用いる場合には、コンピュータ５は、ΔＶ≦０であるときに１００％受理し、ΔＶ＞０であるときにｅｘｐ（−ΔＶ／Ｔ）の確率で受理し、その他の場合、破棄する。コンピュータ５が、この内容を受理した場合にはＳ８でＹＥＳと判定し、Ｓ９において連続変数ｘ＾＊＋Δthermalに離散ノイズΔquantumを加算、更新する。
【００４２】
そしてコンピュータ５は、Ｓ１０ａにおいて質量ｍを大きくしつつ温度Ｔを減少させる。第１実施形態でも説明したように、質量ｍが大きくなると、（１）式はその右辺第１項の評価関数Ｖ（）の影響が強くなり、同時に右辺第２項の量子揺らぎの導入項の影響が弱くなる。また、温度Ｔが減少すると、（１０）式に示すガウスノイズΔthermalの影響も弱くなる。
【００４３】
この後、コンピュータ５は、これらのＳ４〜Ｓ１０ａの処理を繰り返すが、質量ｍを増加させると共に温度Ｔを減少させながら、これらのＳ４〜Ｓ１０の処理を繰り返すため、（１）式の右辺第１項に相当する評価関数Ｖ（）の影響を徐々に強くしながら、（１）式の右辺第２項に示す量子揺らぎの導入項の影響を徐々に弱めることができ、更にガウスノイズΔthermalの影響も徐々に弱めることができる。
【００４４】
コンピュータ５は、これらのＳ４〜Ｓ１０ａの処理を繰り返し、Ｓ１１において終了条件を満たしたときに最適化したと見做してＳ１２において解を出力して処理を終了する。Ｓ１１の終了条件としては、第１実施形態と同一条件を用いれば良いため、説明を省略する。
【００４５】
＜技術的イメージの説明＞
コンピュータ５が、Ｓ４において勾配法により連続変数ｘを更新すると、図６にイメージを示したように、評価関数Ｖ（）が低下する方向にだけ連続変数ｘを更新することになる。例えば、図９に示すように、評価関数Ｖ（）が比較的緩やかに変化する場合を想定したとしても、一旦、局所解Ａ５に嵌ると当該局所解Ａ５から抜け出すことができない。しかしながら、コンピュータ５が、連続変数ｘにガウスノイズΔthermalを加入したシミュレーテッドアニーリング法を用いることで、例えば図９に示すように、評価関数Ｖ（）が比較的緩やかに変化する連続変数ｘの領域においても、評価関数Ｖ（）が緩やかに上昇する方向へ更新させることができ、評価関数Ｖ（）の極値の山を昇ることができ、局所解Ａ５から脱出できるようになる。この結果、緩やかで幅が広い谷に対してもガウスノイズΔthermalを加算することで効率的に局所解Ａ５を脱出できる。
【００４６】
また本実施形態では、離散ノイズΔquantumと共にガウスノイズΔthermalを導入しているため、鋭くて高い谷と緩やかで幅が広い谷が混在する評価関数Ｖ()においても高精度な探索が可能となる。
【００４７】
以上説明したように、本実施形態によれば、コンピュータ５が、連続変数ｘの更新を繰り返すときに温度Ｔを徐々に低下させると共に、連続変数ｘに離散ノイズΔquantumと共に温度Ｔに依存するガウスノイズΔthermalを加算するようにしているため、評価関数Ｖ（）の極値の山を昇ることで局所解Ａ５から脱出でき、しかも鋭くて高い谷と緩やかで幅が広い谷が混在する評価関数Ｖ（）でも高精度に探索できるようになる。
【００４８】
（第３実施形態）
図１０は、第３実施形態の追加説明図を示している。第３実施形態が、第１実施形態と異なるところは、第ｎ励起状態の値を連続的なノイズとして連続変数ｘに加算したところにある。第１実施形態と同一部分には同一符号を付して説明を省略し、以下、異なる部分について説明する。
【００４９】
図１０は、最適解Ａ３の導出処理内容をフローチャートにより示している。コンピュータ５は、第１実施形態に示したように、図１０のＳ１〜Ｓ５の処理を実行する。ここでコンピュータ５は、Ｓ５においてボルツマン分布に従って調和振動子の固有状態を所定の確率で選択する。このとき、前述の（７−１）式、（７−２）式に示すボルツマン分布に従って固有状態を選択する。その後、コンピュータ５は、選択された固有状態の存在確率Ｐｃを用いて調和振動子の第ｎ励起状態の値を連続的なノイズとして連続変数ｘに加算する（Ｓ９ａ）。ここでは受理／破棄の判定を行うことなくノイズを加算するため判定処理を削減できる。
【００５０】
そしてコンピュータ５は、Ｓ１０において質量ｍを大きくする。質量ｍが大きくなると、（１）式はその右辺第１項の評価関数Ｖ（）の影響が強くなり、同時に右辺第２項の量子揺らぎの導入項の影響が弱くなる。この後、コンピュータ５は、これらのＳ４〜Ｓ１０の処理を繰り返すが、質量ｍを増加させながらこれらのＳ４〜Ｓ１０の処理を繰り返すため、（１）式の右辺第１項に相当する評価関数Ｖ（）の影響を徐々に強くしながら、（１）式の右辺第２項に示す量子揺らぎの導入項の影響を徐々に弱めることができる。
【００５１】
その後、コンピュータ５は、Ｓ１１において終了条件を満たしたときに最適化したと見做してＳ１２において解を出力して処理を終了する。Ｓ１１の終了条件としては第１実施形態と同一条件を用いれば良いため説明を省略する。
【００５２】
以上説明したように、本実施形態によれば、評価関数Ｖ（）の微小変化に沿う勾配法により連続変数ｘを更新し、ボルツマン分布に従って調和振動子の固有状態を選択し、選択された固有状態の存在確率Ｐｃを用いて当該第ｎ励起状態の値を連続的なノイズとして連続変数ｘに加算し、ノイズが加算された連続変数ｘを用いて勾配法による更新を繰り返すようにしている。このような処理を行ったとしても第１実施形態と同様の作用効果を奏し、トンネル効果を用いて高精度に最適解Ａ３を導出できるようになる。
【００５３】
（他の実施形態）
本開示は、前述実施形態に限定されるものではなく、例えば、以下に示す変形又は拡張が可能である。
評価関数Ｖ（）の最小値を最適解Ａ３として探索する形態を示したが、最大値を最適解Ａ３として探索する形態に適用しても良い。
前述実施形態の構成、処理内容を組み合わせて構成することもできる。また、特許請求の範囲に記載した括弧内の符号は、本発明の一つの態様として前述する実施形態に記載の具体的手段との対応関係を示すものであって、本発明の技術的範囲を限定するものではない。前述実施形態の一部を、課題を解決できる限りにおいて省略した態様も実施形態と見做すことが可能である。また、特許請求の範囲に記載した文言によって特定される発明の本質を逸脱しない限度において、考え得るあらゆる態様も実施形態と見做すことが可能である。
【００５４】
また本発明は、前述した実施形態に準拠して記述したが、本発明は当該実施形態や構造に限定されるものではないと理解される。本発明は、様々な変形例や均等範囲内の変形をも包含する。加えて、様々な組み合わせや形態、さらには、それらに一要素、それ以上、あるいはそれ以下、を含む他の組み合わせや形態をも、本開示の範畴や思想範囲に入るものである。
【符号の説明】
【００５５】
図面中、１は装置（連続最適化問題の大域的探索装置）、５はコンピュータ、６は更新部、７は選択部、８は判定部、９は加算部、を示す。
【図１Ａ】
2019159782.tif 000014
【図１Ｂ】
2019159782.tif 000015
【図２】
2019159782.tif 000016
【図３】
2019159782.tif 000017
【図４】
2019159782.tif 000018
【図５】
2019159782.tif 000019
【図６】
2019159782.tif 000020
【図７】
2019159782.tif 000021
【図８】
2019159782.tif 000022
【図９】
2019159782.tif 000023
【図１０】
2019159782.tif 000024