JP2022508009A - データ処理装置及びデータ処理方法 - Google Patents

Abstract

所与の問題を解くよう構成された、簡易なハードウエアを用いたデータ処理装置であって、所定の時間発展的処理によって状態データの更新を繰り返す状態データ処理部と、現在の状態データのコスト関数の値を演算するコスト評価部と、現在の状態データの振幅の均一性に関するエラー値を演算するエラー演算部と、を有し、状態データ処理部は、状態データについて時間発展的処理を行い、コスト関数とエラー演算部によって演算された、現在の状態データの振幅の均一性に関するエラー値とに基づいて、現在の状態データを更新する、データ処理装置である。【選択図】図１

Description

本発明は、組み合わせ最適化問題の解を求めるデータ処理装置に関する。

組み合わせ最適化問題を解決するために、古典的なデジタルコンピュータでは、汎用中央処理装置（ＣＰＵ）で動くようにコンパイルされたアルゴリズムが用いられてきた。長年にわたり、デジタルハードウエアの急速な小型化（例えば、ＣＰＵ中のトランジスタの数に対する規模の小型化）は可能であり続けた。しかし、デジタル部品の小型化は限界に近づいてきており、詰め込まれたトランジスタが、十分なエネルギー効率を保ちつつ、安定して機能するのは難しくなってきている。理論上は、このような古典的なコンピュータで用いられるコンピュータ処理は、一般にノイマン型の構成として知られている。

このような古典的なコンピュータによる計算がチューリングマシンを使って行われるものとすると、これらの古典的なコンピュータは、リソースの限界を無視さえすれば、チャーチ＝チューリングの定理として証明されるように、「あらゆるものに使える」計算能力を有する。上述のようにＣＰＵが実行できる論理演算の数にも限界があるが、それに加えて、ＣＰＵ装置とメモリとの間の帯域幅が有限であることも、古典的なコンピュータにおける情報の流れを制限する障害となっていることは良く知られている。

古典的なコンピュータにおける制限を回避するため、「非古典的な」データプロセッサが考えられてきている。このような新しいコンピュータは必ずしも何にでも使える計算能力は有さないが、特定の計算タスクに最適化されている。情報の処理を、抽象化の高いレベルではなく、物理層に近い、抽象化の低いレベルで直接実行することで、パフォーマンスを向上できる。「ソフト」データプロセッサを使えば、特定の目的に合うように動的に内部構造を再構成できるので、抽象化のレベルが最も低いレベルでの計算処理も可能となる。これによって、メモリと処理動作とが協働し、大量に並列的な計算処理が可能となる。

さらに、古典的なコンピュータで使われているバイナリ状態の信号ではなく、データプロセッサのアナログ信号を直接使用することも提案されてきている。アナログ状態の信号は、例えば閾値未満の領域（subthreshold regime）で動作する電子部品を使用するか、あるいは非線形光学現象を利用することで、電子工学領域において物理的に実装できる。このようなアナログコンピュータは、理論的にはデジタルコンピュータでシミュレートでき、またこのアナログコンピュータならば、ある種の特定の問題、例えば微分方程式のシミュレーションに係る問題に対しては、十分高速に処理することができる。このような装置を開発することへの潜在的な動機の一つは、計算に使用するハードウエアの物理的な装置で、０と１だけではなく、もっと多くの情報を符号化したいという点にある。このように、高い論理レベルではなく、より低い物理レベルで直接計算処理をすることによって、リソースを増やすことが可能となる。

具体的に、アナログ・ホップフィールド・ニューラル・ネットワーク（米国特許第４６６０１６６号）のような近年提案されたアナログコンピュータや、コヒーレント・イジングマシン（例えば、米国特許第９４１１０２６号）のような光学アナログコンピュータにより、組み合わせ最適化問題に対する近似的な解を求めることができる。

現在、高速で計算処理を行うため、このようなマシンはアナログハードウエアよって実現された並列計算を利用している（米国特許第４６６０１６６号）。このような装置とニューラルネットワークとの概念的な近似に注目するのも興味深い。計算論的神経科学の分野での最近の進歩は、脳内で起きているアナログ処理に触発された、新規なアナログ計算スキームの開発に応用することができる。

ＧＰＵ（画像処理装置）、ＴＰＵ（テンソル処理装置）、ＦＰＧＡ等の、非従来型の神経系を模した（neuro-inspired）データプロセッサは、分類の分野で成功裏に応用されており、ディープラーニングネットワークにおける最近のトレンドに例示されるように、従来のハードウエアを使用した最先端の手法をしのぐことも可能である。

しかし、上述したアナログ・ニューラルネットワークには２つの制限がある。1つは、コスト関数（又は、目的関数）をシステムのエネルギー関数（又は、通常、接続が対称的であるときに定義されるリアプノフ関数）に写像することで、組み合わせ最適化問題に対する合理的な近似解を得ることは可能だが、どのような組み合わせ最適化問題に対しても最適解を得られる保証はない。当然ながら、非凸問題の場合には、このようなシステムはエネルギー関数の局所的最小値（local minima）に捕捉されることがある。アナログ・ニューラルネットワークは一般に散逸系であり、コヒーレント・イジングマシンの構造においては、最小損失となる解のみが安定的であり、そうでない場合は不安定であるというケースで、システムのゲインを最小値に設定することで、解の質を向上する方法などが提案されてきている。しかし、アナログ変数の振幅が一般的にヘテロ（つまりすべてが等しくない）ので、エネルギー関数による目的関数のマッピングが誤った状態となり、最小損失となる領域で操作しても最適解への収束を保証できない。

２つ目の制限として、制約つき最適化問題において、通常、ペナルティ項をコスト関数に加えることでソフト制約（soft constraints）に変換される制約が、適切に考慮されないということがある。ソフト制約（つまりペナルティ項）を用いることで、収束問題に帰することは知られているが、ペナルティ項はグローバルなコスト関数の和を求める際に、互いに干渉してしまうことも知られている。さらに、ペナルティ項は一般に大きく異なったスケールを持つため、注意深く選んだ定数係数を使って訂正する必要がある。

これらの問題を解決するために種々の方法が提案されている。離散状態の確率論的探索である、一般的なシミュレーテッドアニーリング法では、最適解への収束は、システムの「温度」を、適切な割合で徐々に低下させることで保証される。この収束は解析的に証明されたものであり、現実には、与えられた組み合わせ最適化問題を効率的に解くための最適なスケジューリングを見いだすことは困難である（L. Ingber, Mathematical and computer modelling, 18, (29), 1993）。

ハードウェアの量子動力学的性質を利用することを目的とした量子アナログデバイスまで拡大された物理層に近いコンピュータ概念もある。初期の強い横磁場（transverse field）を徐々に減衰させる量子アニーリングデバイスを構成する試みもある。断熱的方法の限界として、このようなシステムは、基底状態にとどまり横磁場が消失したときに組み合わせ最適化問題の最適解の状態に到達する。例えばＤ－ウェーブ社による（米国特許第６８０３５９９号のような）、このような装置の現実的な実装は、装置内の量子状態を破壊する環境との相互作用を受け得るほか、コンポーネントを特殊な接続トポロジ（キメラグラフ）で接続する必要があり、組み合わせ最適化問題を、このトポロジに埋め込むことが要求され、リソースのコストが増大する。このような限界のため、一般的な種類の問題の解を得るためには、シミュレート・アニーリングのような従来の方法の類いに比べ、上記量子デバイスは十分な計算機的優位性を確立しているとは広く認められていないのが現状である。

米国特許第４６６０１６６号明細書 米国特許第９４１１０２６号明細書 米国特許第４６６０１６６号明細書 米国特許第６８０３５９９号明細書

上述のように、関連技術においては、デジタル及びアナログコンピュータデバイスの双方で、簡易なハードウェア構成で、特定の問題を解決する際の速度の改善が困難になっている。そこで、簡易なハードウェア構成で、迅速に特定の問題を解決するコンピュータ装置が要望されている。

上記従来例の問題点を解決するための本発明の一態様は、所与の課題を解決するためのデータ処理装置であって、所定の時間発展的処理によって状態データの更新を繰り返す状態データ処理部と、現在の状態データに対するコスト関数の値を求めるコスト評価部と、現在の状態データの振幅の均一性に関するエラー値を計算するエラー演算部、とを有し、前記状態データ処理部が、前記状態データについて所定の時間発展的処理を行い、コスト関数と、エラー演算部によって計算されたエラー値と、に基づいて現在の状態データを更新する、データ処理装置。

本発明によると、簡易なハードウェア構成で、迅速に特定の問題を解決できるコンピュータ装置が提供される。もっとも、各請求項記載の発明が有利な作用効果をすべて奏するものとして限定的に解釈されるべきではない。

また、ここでは以下の説明の理解のため、構成や効果の概要を簡潔に示した。追加的な構成及び効果が以下に、各請求項記載の構成に関連して説明される。

実施形態と効果の理解のため、以下の図面を参照しつつ説明する。

本発明の実施の形態に係るデータ処理装置の構成例を表す概略図である。 本発明の実施の形態に係るデータ処理装置の機能的構成の例を表す概略図である。 本発明の実施の形態に係るデータ処理装置のエラー演算部の構成例を表す概略図である。 本発明の実施の形態に係るデータ処理装置の状態データ処理部の例を表す説明図である。 本発明の実施の形態に係るデータ処理装置の機能的構成を表す概略図である。 本発明の実施の形態に係るデータ処理装置における処理の例を表す概略フローチャート図である。

本発明の好ましい実施態様に係るデータ処理装置１は、図１に示すように、プロセッサ１１、及び、入出力装置１３を含む。

本実施の形態のデータ処理装置１は、イジングマシンやホップフィールドニューラルネットワークと同様に、後に説明するように、組み合わせ最適化問題の二値変数をアナログ変数に変換する状態エンコード部を備える。それに加え、データ処理装置１は、もう一つのサブシステムであるエラーエンコード部を有し、このエラーエンコード部は、データ処理装置１の定常状態（steady-states）と、組み合わせ最適化問題のより低いコスト値の構成との間の変換（マッピング）を補正する。状態エンコード部は、エラーエンコード部に対して非対称に接続される。

プロセッサ１１はロジックゲートやメモリブロックを備えたＦＰＧＡでもよい。本実施態様によれば、プロセッサ１１は所定の時間発展的処理によって状態データの更新を繰り返し、状態データに対するコスト関数の値を求め、状態データの振幅の均一性に関するエラー値を計算する。状態データの更新の際には、プロセッサ１１はエラー値及びコスト関数の値も利用する。プロセッサ１１の詳細な処理については後述する。

プロセッサ１１のメモリブロックは、プロセッサ１１のロジックゲートの処理で利用される、状態データなどのデータを保持する。

入出力（Ｉ／Ｏ）装置１３は、キーボード、マウス等の入力装置を含む。また、この入出力（Ｉ／Ｏ）装置１３は、プロセッサ１１からの指示に応じて、状態データや、コスト関数の値などの情報を出力するディスプレイも含んでもよい。

以下に、プロセッサ１１での処理を記載する。以下の例では、本発明の一実施態様に係るデータ処理装置１によって解決される課題は、組み合わせ最適化問題であるものとする。組み合わせ最適化問題においては、コスト関数が規定され、組み合わせ最適化問題のコスト関数が極小化されることが、組み合わせ最適化問題が解かれたことを意味する。

ここではコスト関数はV⁽⁰⁾(σ)で表されるとし、このうち、V⁽⁰⁾(σ)の値はどのベクトルσに対しても実数となり、ベクトルσはσ={σ_i}_i (i = 1,2,…,N)、σ_i =±1であるとする。コスト関数V⁽⁰⁾(σ)は、パラメータ集合{M_0k}_k (k = 1,2,…)で規定され、M_0kはベクトルまたは行列であり、より一般的にはテンソルである。ブーリアン型変数の数（問題の大きさ）はＮで表される。

制約つき最適化問題の場合には、許容される解は、すべての配置からなる空間に対して、Ｓで表される部分集合を構成する。所与の組み合わせ最適化問題の制約に対して、部分集合Ｓは下記のような等式制約及び／又は不等式制約によって規定される。

（１）式で、k, k’ = 1,2,…,Kである。

行列及びベクトルA_k、C_k,、b_k、及びd_k’は、{N_ki}_i又は{M_ki}_iで表される別のパラメータ集合で規定される。ここで、kは制約のインデックスであり、k =1,2,…,K, i=1,2,…である。

制約は２つのカテゴリーに分けられる。制約の第一の集合は、タイプＩのソフト制約と呼ばれるもので、ペナルティ項をコスト関数に加算し、その制約によって規定される有効な部分空間にシステムを射影することによって充足される制約である。全コスト関数V^*は、これらの制約を考慮したものであり、下記のようになる。

（２）式において、U^(K)は、制約kによって与えられるペナルティ項である。q_kは正定数パラメータであり、k=1, 2, …, K_Iであり、K_Iは当該部分集合内の制約の総数を表す。ベクトルσが制約kを充足すると、ペナルティ項U^(K)(σ)は極小となる。ペナルティ項U^(K)(σ)は、パラメータ{M_ki}_iと、有効部分空間への射影Pとに依存する関数であり、ここでは、提案されたシステムに対して入力されるものである。

制約の第二の集合は、タイプＩＩのソフト制約と呼ばれるもので、エラー検出訂正フィードバックループを使って導入される制約である。パラメータ{N_ki}_i（k=1,…, K_II，ただしK_I + K_II = K）に依存するペナルティ項V^(K)も規定される必要がある。このペナルティ項V^(K)はエラー訂正に使用される。更に、g_kは制約が充足されないときに正の値をとる関数であり、エラー検出に使用される関数であるとする。制約が充足される場合、この関数g_kの値は負となるものとする。

これらのU^(k), V^(k), P, 及び g_kは、解かれるべき組み合わせ問題とその制約に合わせて選択される。言い換えれば、U^(k), V^(k), P, 及び g_kはデータ処理装置１のユーザによって決定される。

図２は、プロセッサ１１の機能的な構造の一例を示す。図２にあるように、一つの実施例では、プロセッサ１１は、機能的に、状態データプロセッサ２１、コスト評価部２２、エラー演算部２３、変調部２４、及び出力部２５を含む。

状態データプロセッサ２１は、所定の時間発展的処理によって状態データの更新を繰り返す。ここで状態データは固定小数点変数ｘ_i（ｉ＝１，２，…Ｎ）の集合であり、例えばアナログ状態を符号化して得られるものである。

本実施形態において、この状態データプロセッサ２１は所定の時間発展的処理によって状態データの更新を繰り返し、更新した状態データを有効な部分空間に射影し、射影した更新済状態データを新しい状態データとしてメモリブロックに保持する。

具体的には、状態データプロセッサ２１は処理部２１０と、ゲイン消失シミュレータ２１１と、射影部２１２とを有する。各処理部２１０は、対応する状態データｘ_i（ｉ＝１，２，…Ｎ）をそれぞれ処理する。ゲイン消失シミュレータ２１１は孤立した（非結合な）ユニットである。このゲイン消失シミュレータ計算部２１１は状態データｘ_iと、線形ゲインｐとの入力を受けて、ポテンシャルＶ_bを最急降下法により演算する。この演算は、次の常微分方程式で簡単に示すことができ、この常微分方程式は、Ｉ_ｉ＝０のときに最急降下法に相当するものとなる。
∂_ｔｘ_ｉ ＝ ｆ（ｘ_ｉ） + Ｉ_ｉ (3)
式（３）において、ｆ（ｘ_ｉ）＝－∂Ｖ_ｂ／∂ｘ_ｉ であり、V_ｂは孤立した（非結合な）ユニットのエネルギー関数またはリアプノフ関数であり、例えば、次のようなポテンシャル関数である：
Ｖ_ｂ ＝－（－１＋ｐ）ｘ_ｉ ²／２＋ｘ_ｉ ^４／４
なお、Ｉ_ｉについては後述する。

エネルギー関数Ｖ_ｂは、典型的には双安定ポテンシャル（単安定または双安定ポテンシャルの典型（archetype））であり、線形ゲインｐの値によって単安定（ｐ＜１のとき）と双安定（ｐ＞１のとき）とのいずれにもなり得る。Ｖ_ｂが双安定の場合、状態データｘ_ｉはＩ_ｉ＝０のとき、ポテンシャルＶ_ｂが最小となるところでバイナリー状態に収束する。なお、ここでＩ_ｉは、ｉ番目の処理部２１０に外部から入力されるアナログ信号を表す。外部から入力されるアナログ信号の例については、後に述べる。

これらより（３）式は下記のように書き換えられる。
∂_ｔｘ_ｉ ＝ （－１＋ｐ）ｘ_ｉ－ｘ_ｉ ³ + Ｉ_ｉ (4)
（４）式において、－ｘ_ｉ，ｐｘ_ｉ及び－ｘ_ｉ ³は、それぞれ、状態ｘの損失、線形ゲイン、及び飽和に対応する項である。

この力学系は、常微分方程式（ODE）を使用する連続時間ドメインにおけるものとなっているが、本実施の形態のシステムは、離散時間ドメインで運用することも可能である。連続時間から離散時間への変換は、上述の常微分方程式（ODE）を単にオイラー近似することなどによって行うことができる。

状態データプロセッサ２１の各処理部２１０間の結合は、下記で表されるインジェクション項Ｉ_ｉを使って実装される。

ここで、

ここで、ｋ＝１，２，…Ｋ_II。

式（６）及び（７）において、ベクトルｅ^(k)（要素ｅ_ｉ ^(k)，ｉ＝１，２，…）はエラー信号であり、Ｖ^＊はタイプＩのソフト制約を考慮したペナルティ項を含むコスト関数であり、Ｖ^(k)はタイプＩＩのｋ番目のソフト制約に関するペナルティ項である。

さらに、ε_kは、正の実数パラメータ値である。制約の種類については、後に説明する。

Ｉ_ｉを入力すると、ポテンシャルＶ^(k)（ｘ）の最急降下が行われる。しかし、勾配∂Ｖ^(k)（ｘ）／∂ｘ_ｉはｅ_ｉ ^(k)によって変調（modulate）される、すなわち、勾配ベクトルは状態空間を用いて規定され、エラー信号によってそのスケールが変更される。各エラー信号ｅ_ｉ ^(k)は、加えられる制約に応じて、空間ベクトルｘを異なるスケールに変更する。

そして、全インデックスｋについて、タイプＩＩのソフト制約を考慮した、勾配の総和を求める。これにより、複数の制約が総和の中で競合することになる。ある制約に関して、スケールが変更された勾配ベクトルの振幅が、他のものよりかなり大きい場合、その制約が最終的に勝ち残る。

射影部２１２は、状態データを、所定の部分空間に射影する。具体的に本実施の形態では、状態データベクトルｘが、タイプＩの制約に対応して予め定められた射影演算子Ｐを用いて、状態データベクトルｘを更新する繰り返し処理ごとに、有効な部分空間に射影される。

この射影は、ホップフィールドネットワークにおけるアイヤー（Aiyer）の方法と同様のものである。

射影演算子の効果は、次のように、ベクトルｘの時間発展のオイラーステップを考慮して説明される：

ここで、ｘ′_ｉ（ｔ）は、射影演算子Ｐを用いた、ｘ_ｉ（ｔ）の有効な部分空間への射影を表す：

仮に、タイプＩのソフト制約が存在しない場合、射影Ｐは恒等演算子であり、ベクトルｘ′はベクトルｘに等しくなる。

本実施の形態のシステムにおいて、時間発展は、射影Ｐを考慮に入れるため、連続時間ドメインにおいて代数的な微分方程式系によっても表現される。

コスト評価部２２は、現在の状態データに対するコスト関数の値を求める。

エラー演算部２３は、少なくとも現在の状態データの振幅の均一性（homogeneity）に関するエラー値を演算する。本実施の形態においては、このエラー信号は、
（１）状態エンコードユニットにおける振幅の非均一性を補正し、
（２）適切な制約のマッピングを許容する。

各エラーエンコードユニットは、通常、状態エンコードユニットのサブセットにのみ接続される。振幅の非均一性の補正は、アナログ空間に関する最適化問題の等式制約として解釈可能であることに注意されたい。

換言すれば、バイナリー変数の組み合わせ最適化は、すべてのアナログ状態の振幅が等しいという制約を含めた、アナログ変数の最適化問題の一種であると言える。

図３に示すように、エラー演算部２３は、エラー演算サブユニット２３１と、複数のサブユニット２３２とを含む。

エラー演算サブユニット２３１は、振幅の非均一性に係るエラーｅ⁽⁰⁾を演算する。このエラー演算サブユニット２３１は、時間発展プロセッサ２３１１とアップデータ２３１２とを含む。時間発展プロセッサ２３１１は、ａ＞０なるターゲット振幅と、エラー信号の変化率β_０との入力を得て、エラー信号ｅ⁽⁰⁾を演算する。ここで、ターゲット振幅とエラー信号の変化率β_０とは初期化の際に、ユーザが指定する。またエラー信号ｅ⁽⁰⁾は、ｅ^(k)の一つであり、ｋ＝０の場合に相当し、コスト関数Ｖ^＊の最小化に関わる。これらのエラー信号は、状態データベクトルｘの振幅の非均一性を補正する。

具体的に、時間発展プロセッサ２３１１は、エラー信号ｅ⁽⁰⁾を：

のように演算する。

アップデータ２３１２は、現在のｅ_ｉ ⁽⁰⁾に、∂_ｔｅ_ｉ ⁽⁰⁾ｄｔを加えて、更新後のｅ_ｉ ⁽⁰⁾を求めることで更新し、この更新後のｅ_ｉ ⁽⁰⁾を、メモリブロックに、次の繰り返し処理のエラー信号として保持する。

制約のエラーを演算するサブユニット２３２のそれぞれもまた、時間発展プロセッサ２３２１と、アップデータ２３２２とを備える。時間発展プロセッサ２３２１は、ターゲット関数ｇ_ｉと、このエラー信号の変化率β_ｉとの入力を受けて、エラー信号ｅ⁽ⁱ⁾を出力する。ここでターゲット関数ｇ_ｉと、エラー信号の変化率β_ｉとは、初期化の際に、ユーザが指定する。また、エラー信号ｅ⁽ⁱ⁾は、ｅ^(k)の一つで、ｋ＝ｉの場合に相当しており、コスト関数Ｖ^＊の最小化に関わる。これらのエラー信号は、状態データベクトルｘに対して問題の制約を与える。

具体的に、時間発展プロセッサ２３２１は、エラー信号ｅ⁽ⁱ⁾を：

として演算する。ここで、ｘは、状態データのベクトルであり、ｇ^(k)（ｅ_ｉ ^(k)，ｘ）は、解こうとする問題の制約に関係して定められる。

仮に、解こうとする問題に制約が存在しない場合、サブユニット２３２は必ずしも必要ではない。

アップデータ２３２２は、現在のｅ_ｉ ^(k)に、∂_ｔｅ_ｉ ^(k)ｄｔを加えて、更新後のｅ_ｉ ^(k)：

を求め、この更新後のｅ_ｉ ^(k)を、メモリブロックに、次の繰り返し処理のエラー信号として保持する。

エラー演算部２３は、現在のエラーの値ｅ_ｉ ^(k)（ｉ＝０，１，２，…）を、状態データプロセッサ２１に対して出力する。

変調部２４は、線形ゲインｐ、ターゲット振幅ａ、エラー値の変化率β_ｋなどのパラメータ値を、現在の状態データに基づいて演算する。この変調部２４がエラー演算部２３に対してターゲット振幅等のパラメータを出力することとすれば、エラー演算部２３は、変調部２４から得られる当該パラメータを、ユーザが指定した値に代えて利用できる。

この例では、変調部２４は、アナログ状態ｘを、演算の各ステップにおいて、コンピュータで演算可能なブーリアン変数σの構成に、σ＝Ｃ［ｘ］というように変換する。また、変調部２４は、この変数を、コスト関数Ｖ⁽⁰⁾の現在の値を演算するために利用する。

本実施の形態のある例では、変調部２４は、コスト関数Ｖ⁽⁰⁾の現在の値を：

として求める。ここでｈは、行列Ｍ_0iを、σの左から乗じる（いわゆるカップリング演算）などして得られる内部場演算（internal field calculation）を表すものであり、ｈ＝Ｍ_0iσと表される。なお、ドットは、ドット積を表す。

また変調部２４は、線形ゲインｐやターゲット振幅ａを：

と変調する。ここで、

である。また、Ｖ⁽⁰⁾（ｔ）は、状態ｘ（ｔ）に関連するコスト関数の値であり、Ｖopt⁽⁰⁾は、ターゲットエネルギーである、この例では、Ｖopt⁽⁰⁾は、繰り返し演算の間に見いだされたエネルギーの最低値に設定すればよい。つまり、

としてよい。あるいは、このＶopt⁽⁰⁾は、知られているのであれば、コスト関数Ｖ⁽⁰⁾の最低値に設定されてもよい。

関数φは、シグモイド（例えばハイパボリック・タンジェント）関数であり、
δ＞０，ρ₁，ρ₂
は、予め定められた定数パラメータであり、ρ₁，ρ₂はいずれも、正の値でも、負の値であってもよい。

仮に、ターゲット振幅ａや線形ゲインｐが

に従って変調されるのであれば、ρ₁＞０または、ρ₂＞０である。

一方、仮にターゲット振幅ａや線形ゲインｐが

に反して変調されるのであれば、ρ₁＜０または、ρ₂＜０である。

他のパラメータ、例えば変化率β_ｋも変調され得ることに注意されたい。

ここで説明した処理手順の効率は、α，πなどのパラメータの選択に依存するものである。

カップリング行列のスペクトル分解（最大の固有値）を用いるなどして事前の調整なしに、パラメータを選択できることが理解される。

出力部２５は、現在の状態データｘを出力する。この出力部２５は、状態データとともに、現在の状態データに対するコスト関数を併せて出力してもよい。

本実施の形態のある例では、データプロセッサは、上述のエラー補正処理を含む。エラー検出は、例えば、エラー信号などの外部アナログ動的変数を考慮して行われる。エラー補正変数は、システムによる目的関数のマッピングの誤りによって生じる振幅の非均一性を補正するために用いられる。

さらに、もう一組の変数と有効な部分空間への射影が、最適化問題に与えられる制約のために考慮される。エラー制御は、非対称なエラー補正とエラー検出のフィードバックループを用いる。

さらに、エラー信号の変動は、一般に、アナログ状態を順次エンコードして得られる現在のブーリアン変数に依存し、論理レベルでエラーを検出するために用いられる。しかしながら、エラー信号それ自体は、アナログ値であり、現在の状態エンコード変数をアナログ的に変調する。

本実施の形態のデータ処理装置１は、上述の構成を備えており、次のように動作する。

なお、本実施の形態のデータ処理装置１は、図６に例示するように、状態データｘ_ｉ（ｉ＝１，２，…Ｎ）を所定の方法で初期化する（Ｓ１）。一例として、データ処理装置１は、状態データを、ランダムなブーリアン値を生成することによって初期化する。

データ処理装置１は、状態データｘｉとコスト関数のパラメータのセット｛Ｍ_0k｝_k（ｋ＝１，２，…）とに基づいてパラメータの変調を演算する（Ｓ２）。

このステップＳ２では、データ処理装置１は、タイプＩの制約が存在する場合、変調された線形ゲインｐ、変調されたターゲット振幅ａ、コスト関数Ｖ⁽⁰⁾の現在の値：

または

を求める。

ある例では、Ｍ₀₁やＭ₀₂は行列であり、Ｍ₀₃はベクトルである。Ｍ₀₁，Ｍ₀₂，及びＭ₀₃はコスト関数Ｖ⁽⁰⁾を規定する。同様に、仮にタイプＩの制約が存在するのであれば、Ｍ₁₁やＭ₁₂は行列であり、Ｍ₁₃はベクトルである。ここでｋは、Ｕ^(k)で定義される制約ｋ（またはｋ＝０のときには、目的関数Ｖ⁽⁰⁾に対応する）を表す。Ｍ₁₁，Ｍ₁₂，Ｍ₁₃は、Ｕ⁽¹⁾で定義される第一の制約（ｋ＝１）を表す。β_kは、エラー信号の変化率を表す。

データ処理装置１は、ステップＳ２において、状態データｘ_ｉに係る所定のブーリアン変数σを演算する。このステップＳ２における詳細な演算の内容は、変調部２４の処理として既に説明した通りである。

データ処理装置１は、振幅の非均一性、及び制約に係るエラー変数を更新する（Ｓ３）。このエラー変数は、変調されたターゲット振幅、ステップＳ１で取得したエラー変数の変化率、及び、目的関数ｇを用いて：

及び、

として演算する。

ここで、関数ｇは、解こうとする問題の制約により規定される。

既に説明したように、状態データｘ_ｉのプライムサイン（′）は、それが状態データｘ_ｉの有効な部分空間への射影であることを示す。

データ処理装置１は、また、タイプＩの制約が存在する場合には、カップリング項、すなわち、例えば、

や、

を演算する。また、データ処理装置１は、

を演算する。ここで、Ｎ_kiは、タイプＩＩの制約により規定される。さらに、データ処理装置１は、ゲイン及び飽和：

を求める（Ｓ４）。

ステップＳ３，Ｓ４の処理は、並列して行われてもよいし、いずれかを先にした順序で順番に行われてもよい。

そしてデータ処理装置１は、インジェクション項を演算する（Ｓ５）。この演算は、インジェクションを総和する：
すなわち、

として、

を求めることで行われる。そしてデータ処理装置１は、状態データを更新する（Ｓ６）：

データ処理装置１は、状態データｘ_ｉに対する射影を求め（Ｓ７）、現在のエラー（現在の損失）を求める（Ｓ８）。

データ処理装置１は、そして、予め定めた条件が満足されるか否かを調べて、繰り返し処理を終了するか否かを判断する（Ｓ９）。ここで、上記予め定めた条件の一例は、時間の制限であり得る。この例では、繰り返し処理により消費した時間が、予め定められた制限時間を経過したか否かを判断することとなる。データ処理装置１は、繰り返し処理により消費した時間が予め定めた制限時間を経過していなければ（Ｓ９：Ｎｏ）、ステップＳ１に戻って処理を続ける。一方、データ処理装置１は、繰り返し処理により消費した時間が予め定めた制限時間を経過していれば（Ｓ９：Ｙｅｓ）、演算の結果（状態データｘ_ｉや、それに対応する所定のブーリアン変数σ）を、見いだされた最適な結果として出力して（Ｓ１０）、処理を終了する。

本実施の形態のデータプロセッサの処理を示すため、以下に、具体的な問題の例を挙げて説明する。

［（１）最大カット問題の例］
本実施の形態のデータ処理装置１を最大カット問題に適用した例を示す。この例の最大カット問題のためのデータ処理装置１は、重み｛ｗ_ij｝_ijが規定されたグラフのカットであって、重みの総和が最大となるような切断（カット）を見いだすよう設定される。ここでｉ，ｊは、それぞれＮ以下の自然数：１，２，…Ｎの一つであるとする。

最大カット問題の与えられる解は、カットの後に得られる２つの集合に対応する頂点ｉの分割として表され得る。頂点ｉが一方の集合または他方の集合のいずれに属するかがブーリアン値σ_ｉとしてエンコードされる。この問題の最適解、すなわち最大カットは、コスト関数

を最小化するものとなる。ここで、

であり、

は、イジング・ハミルトニアン：

であり、σは、要素がσｉ（ｉ＝１，２，…）であるようなベクトルである。

最大カット問題は、二次の非制約二値組み合わせ最適化問題（quadratic unconstrained binary combinatorial optimization problem）、あるいはイジング問題と呼ばれる問題である。コスト関数のパラメータが行列で構成されること、そしてコスト関数の総和Ｖ^＊が単一の行列：

で表されることに注意すると、目的関数は、簡潔に、

として与えられる。

この最大カット問題の例では、制約は存在しない。つまり、この例におけるデータ処理装置１は、それぞれ状態データｘ_ｉ（ｉ＝１，２，…Ｎ）に対応するＮ個の状態データプロセッサ２１のみを備える。この例のデータ処理装置１は、また、サイズＮの問題を解く際の振幅の非均一性を補正するためのエラーデータｅ_ｉ ⁽⁰⁾に対応するＮ個のエラー演算部２３のみを要する。

この例では、有効な部分空間は、すべての構成となるので、射影演算子Ｐによる射影は、恒等的なもの：

となる。

さらに、妥当な解への変換は、

として定義する。ここで、ｓｉｇｎ（ｘ）は、ｘ＞０のときｓｉｇｎ（ｘ）＝１、そうでないときにはｓｉｇｎ（ｘ）＝－１となるものとする。

そしてこの実施の形態の例の制御部１１では、状態データプロセッサ２１が所定の方法、例えばランダムなブーリアン値をｘ_ｉに設定するなどの方法で状態データｘ_ｉ（ｉ＝１，２，…Ｎ）を初期化する。

状態データが設定されると、コスト評価部２２が現在の状態データのコスト関数を評価する。この例では、コスト関数は、

と設定される。ここでハミルトニアンは、イジングハミルトニアン：

であり、ここで｛ｗ_ij｝_ijは、繰り返し処理の開始時にランダムに初期化される重みである。

一方で、エラー演算部２３は、現在の状態データの振幅の均一性に係るエラー値を演算する。

エラー演算部２３の時間発展プロセッサ２３１１は、ターゲット振幅ａ（ａ＞０）と、エラー値の変化率βとの入力を受け入れる。これらのターゲット振幅ａ及びエラー値の変化率βは、当初（ｔ＝０の時点）では、予め定めた値、例えばａ＝１．０，β（ｔ＝０）＝０．０などといった値に初期化される。

時間発展プロセッサ２３１１は、これらの値を用いて、エラー値ｅ_ｉ ⁽⁰⁾の時間発展

を演算する。ここで繰り返し処理の初回においては、エラー値ｅ_ｉ ⁽⁰⁾は、ランダムな値に設定される。

アップデータ２３１２は、現在のエラー値ｅ_ｉ ⁽⁰⁾を、∂_ｔｅ_ｉ ⁽⁰⁾ｄｔを加算することにより更新して、更新後のエラー値ｅ_ｉ ⁽⁰⁾を求め、この更新後のエラー値ｅ_ｉ ⁽⁰⁾を、状態データプロセッサ２１に出力する。

そして状態データプロセッサ２１のゲイン消失シミュレータ２１１が、状態データｘ_ｉ、線形ゲインｐ、及びエラー値ｅ_ｉ ⁽⁰⁾を得て、状態の時間発展：

を求める。

この式において、時間依存性は陽に示されていないが、状態データｘ_ｉやエラー値ｅ_ｉ ⁽⁰⁾などの値が、時間に依存して変化する。

状態データプロセッサ２１は、対応する状態データｘ_ｉと、状態データの時間発展：

とを加算して、状態データを更新する。

状態データプロセッサ２１は、更新された状態データｘ_ｉ（ｔ＋ｄｔ）を、現在の状態データとして記憶する。この問題においては、タイプＩのソフト制約は存在しないので、射影Ｐは恒等演算子であり、従ってベクトルｘ′は、ベクトルｘに等しくなり：

更新された状態データは、次の繰り返し処理における現在の状態データとしてそのまま保持される。

状態データを更新した次の繰り返し処理の前に、変調部２４は、線形ゲインｐ、ターゲット振幅ａ、及びエラー値の変化率βなどのパラメータ値を、現在の状態データに基づいて演算する。

変調部２４は、状態データｘを、可能なブーリアン変数σに変換する。この例では、状態データｘは既にブーリアン値であるので、変調部２４は、コスト関数Ｖ⁽⁰⁾の現在の値を演算し、線形ゲインｐと、ターゲット振幅ａとを、式（１２）、（１３）によって表されるように変調する：

また、変化率β₀も、同様に変調される。

上記の式において、

であり、関数φは、シグモイド関数である。また、

は、予め定められた定数パラメータである。

パラメータ値α，πは、ユーザによって予め設定される。

ここで、Ｖ⁽⁰⁾（ｔ）は、コスト評価部２２によって評価された状態データｘ（ｔ）に係るコスト関数の値であり、Ｖ_opt ⁽⁰⁾は、ターゲットエネルギーである。この例では、Ｖ_opt ⁽⁰⁾は、繰り返し処理の間に見いだされたエネルギーの最低値に設定される。すなわち：

このＶ_opt ⁽⁰⁾を得るため、変調部２４は、現在のＶ_opt ⁽⁰⁾を記憶し、現在のＶ_opt ⁽⁰⁾が記憶しているＶ_opt ⁽⁰⁾よりも低い値であるときに、Ｖ_opt ⁽⁰⁾を更新して記憶する。初期状態では、変調部２４は、Ｖ_opt ⁽⁰⁾を記憶していないので、変調部２４は、演算されたＶ⁽⁰⁾（ｔ＝０）を記憶する。

変調部２４は、更新された線形ゲインｐ、ターゲット振幅ａ、及び変化率β₀を得て、これらを出力する。

制御部１１は、現在の更新された状態データｘと、コスト関数とを出力し、次の繰り返し処理に移行して処理を続ける。

コスト評価部２２は、現在の（更新された）状態データに係るコスト関数を評価し、エラー演算部２３は、現在の状態データに係る振幅の均一性に関するエラー値を演算する。

制御部１１は、繰り返し処理を行い、エラー値の時間発展：

と、状態データの時間発展：

と、を演算し、エラー値を：

として更新し、また、状態データを：

と、更新する。

制御部１１は、更新された状態データｘ_ｉ（ｔ＋ｄｔ）を、現在の状態データとして記憶する。

そして制御部１１は、線形ゲインｐ、ターゲット振幅ａ、エラー値の変化率β_kなどのパラメータ値を、現在の状態データに基づいて演算する。

制御部１１は、現在の状態データｘと、コスト関数とを出力し、このコスト関数が、所定の条件、例えばコスト関数が予め定めたしきい値を下回るなどの条件を満足するか、ユーザが処理を停止するまで繰り返し処理を行う。

最大カット問題の場合、データ処理装置１の動作は、次のように要約できる：

これらの式は、

の勾配により与えられる。

さらに、ターゲット振幅ａは、最適解：

への収束を保証するため、次のようにして選択される。ただし、上述の式で、

またここで、

は、現在の状態（現在の状態データ）のエネルギーであり、

は、ターゲットエネルギーである。実際には、ターゲットエネルギーは、それまでに見いだされたエネルギーの最低値：

に設定されるか、または、既知であれば基底エネルギーの値に設定される。

さらに、関数ｆは、シグモイド関数であり、ρ，δの双方は正の０でない定数パラメータであって、ユーザによって予め定められる。

解を得るまでの時間をさらに短縮するため、パラメータβが時間依存であることを考慮する。パラメータβは処理中、変化率λで線形的に増加し、正の値τで表される時間の間、エネルギーが減少しないのであれば、０にリセットされる。

である場合、ただし、ｔ_ｃは、

が最低値である時点または、βがリセットされた時点を表すものとするとき、ならば、β（ｔ）の変化は、

ただし、

とする。一方で、

ならば、βは０に設定され、ｔ_ｃはｔに設定される。

また、線形ゲインであるパラメータｐは、状態に依存して：

として定められる。

［（２）二次割当問題（Quadratic Assignment Problem；ＱＡＰ）］
二次割当問題（ＱＡＰ）は、生産品を輸送するトータルのコストが最小となるように、ｎ個の工場をｎ箇所の異なる場所に割り当てる問題である。ここでコストは、工場間の距離とフローとの積の総和に依存するものとする。この問題を目的関数としてエンコードするため、ブーリアン値ｓ_ｉｕが次のように設定される。すなわち、工場ｕが場所ｉに配置されるときｓ_ｉｕ＝１とし、そうでなければｓ_ｉｕ＝０とする。次に、この問題を、コスト関数Ｖ⁽⁰⁾：

を、最小化する割り当てを見いだす問題として定式化する。ここで行列｛ｂ_ｕｖ｝_ｕｖ（（ｕ，ｖ）番目の要素が、ｂ_ｕｖとなる）は、工場ｕ，ｖ間のフローを表し、行列｛ａ_ｉｊ｝_ｉｊ（（ｉ，ｊ）番目の要素が、ａ_ｉｊとなる）が、場所ｉ，ｊ間の距離を表す。

なお、巡回セールスマン問題は、このＱＡＰの特別な場合に相当する。

変数の変化σ_i＝２（ｓ_ｉ－０．５）を考慮すると、この目的関数は、次のコスト関数にマッピングできる：

ここでＶ⁽⁰⁾は、最小とするべきコスト関数であり、Ｕ⁽¹⁾，Ｕ⁽²⁾は、それぞれ第一の制約（一つの場所あたり工場一つとの制約）と、第二の制約（一つの工場あたり一つの場所）とに係るタイプＩのソフト制約を表す。さらにｑは、ユーザによって定められる正のパラメータである。ｋ＝１，２，３…に対応する各コスト関数は：

と表すことができる。イジングカップリング及びゼーマン項は、

である。

ここで

は、Ｎ×Ｎの行列であり、その要素がすべて「１」であるような行列である。また、太字の

は、要素数Ｎのベクトルであり、その要素がすべて「１」であるようなベクトルを表す。

Ａ，Ｂは、それぞれの要素がａ_ｉｊ，ｂ_ｉｊであるような行列であり、ＩはＮ×Ｎの単位行列（identity matrix）である。

イジング・カップリングが本実施の形態のこの問題におけるコスト関数であり、コスト関数のパラメータは、行列のテンソル積に依存することに注意されたい。

コスト関数を特徴づけるパラメータは：

と表される。

さらに、第一、第二の制約を特徴づけるパラメータは：

及び、

となる。

ＱＡＰでは、すべての制約がタイプＩとなる。従って制御部１１は、サイズＮの問題の解を求める場合、状態データｘ_ｉに対応するＮ個の状態データプロセッサ２１と、振幅の非均一性を補正するためのｅ⁽⁰⁾ _ｉに対応するＮ個のエラー演算部２３とを含んで構成される。そして、コスト関数Ｖ^＊は、式（２１）で与えられたものとなる。この問題に関する有効な部分空間は：

と定義される。

ここで、Ｘ_ａは、Ｎ次元の実空間であり、Ｃ_ａは、

となるように定義された定数である。

有効な部分空間は、従って、確率行列｛ｘ_ｉｕ｝_ｉｕの集合となる。

有効な部分空間への射影演算子Ｐは、この行列の固有値分解

を考慮して定義される。

可能な解（acceptable solution）への変換は、置換行列を、各状態データ行列｛ｘ_ｉｕ｝_ｉｕに関連付けることで得られる。

ここでの例では、本実施の形態の制御部１１において、状態データプロセッサ２１は、状態データｘ_ｉ（ｉ＝１，２，…Ｎ）を、所定の方法、例えば状態データを、許される状態の一つとして（すなわち状態データが有効な部分空間にあるように）生成することによって、状態データｘ_ｉを初期化する。

ここで状態データｘｉは、ブーリアン変数ｓ_ｐｑであって、工場ｑが場所ｐに配されるときにｓ_ｐｑ＝１となり、そうでないときにはｓ_ｐｑ＝０となるようなブーリアン変数ｓ_ｐｑである。ここでｉ＝Ｎ_site（ｐ－１）＋ｑであり、Ｎ_siteは、場所の数を表す。

状態データプロセッサ２１の射影部２１２は、現在の状態データを、有効な部分空間に射影する。具体的に、本実施の形態のここでの例では、状態データベクトルｘは、タイプＩの制約に従い、予め定められた射影演算子Ｐを用いて、有効な部分空間に射影され、ｘ′となる。射影演算子Ｐは、行列の固有値分解

を考慮して定義される。ここで、

は、式（２２）において定義されたものである。

状態データが設定されると、コスト評価部２２は、現在の状態データのコスト関数を演算する。この例では、コスト関数Ｖ^＊は式（２１）により評価される。

一方で、エラー演算部２３は、現在の状態データの振幅の均一性に関するエラー値を演算する。

エラー演算部２３の時間発展プロセッサ２３１１は、当初（ｔ＝０）所定の値、例えばａ＝１．０、β（ｔ＝０）＝０．０に初期化されるターゲット振幅ａ（ａ＞０）と、エラー値の変化量βとを用いた演算を行う。

これらの値を利用して、時間発展プロセッサ２３１１は、エラー値の勾配ｅ_ｉ ⁽⁰⁾を

と求める。ここで、初回の繰り返し演算におけるエラー値ｅ_ｉ ⁽⁰⁾（ｔ＝０）は、所定の値に設定される。

アップデータ２３１２は、現在のエラー値ｅ_ｉ ⁽⁰⁾（ｔ）を、∂_ｔｅ_ｉ ⁽⁰⁾（ｔ）ｄｔを加算して更新し、次の繰り返し処理で利用するエラー値である、更新されたエラー値ｅ_ｉ ⁽⁰⁾（ｔ＋ｄｔ）を得るとともに、当該更新されたエラー値ｅ_ｉ ⁽⁰⁾（ｔ＋ｄｔ）を、状態データプロセッサ２１に出力する。すなわちアップデータ２３１２は、：

を演算して出力する。

次に、状態データプロセッサ２１のゲイン消失シミュレータ２１１が、状態データｘ′_ｉ（状態データを有効な部分空間に射影したもの）、線形ゲインｐ、及びエラー値ｅ_ｉ ⁽⁰⁾を得て、状態の時間発展：

を演算する。

ここでｈ⁽⁰⁾ _ｉ（ｔ）は、ベクトルｈ⁽⁰⁾（ｔ）のｉ番目の要素であり、

と定められる。ただし、

であり、ｑとηとはいずれも正の定数パラメータ：ｑ＞０，η＞０である。ｈ⁽⁰⁾ _ｉ（ｔ）の値は、Ｍ_０ｉ，Ｍ_１ｉ，Ｍ_２ｉ及び状態データｘ′_ｉ（有効な部分空間に射影された状態データ）を用いて演算される。

項

は、簡単な場合、ｘ′_ｉの平均絶対振幅

つまり、

に等しく、簡単な場合でなければ、項

は、

に等しいものとされる。ここで、ζは、ζ＞１を満足する値として経験的に設定されるパラメータである。

状態データプロセッサ２１は、状態データに、対応する状態データｘ_ｉと最急降下勾配とを加算して更新する：

状態データプロセッサ２１は、更新した状態データｘ_ｉ（ｔ＋ｄｔ）を、現在の状態データとして記憶する。ここでの問題を解くに際しては、タイプＩのソフト制約が存在し、状態データプロセッサ２１の射影部２１２は更新された状態データｘを、射影演算子Ｐを用いて射影する：

この射影された状態データｘ′_ｉ（ｔ＋ｄｔ）は、次の繰り返し処理で利用される、現在の状態データｘ（ｔ）として記憶される。

更新された状態データを用いた次の繰り返し処理の前に、変調部２４は、線形ゲインｐや、ターゲット振幅ａ、エラー値の変化率β_０などのパラメータ値を現在の状態データに基づいて演算する。

変調部２４は、状態データｘを、σ_pq＝２（ｓ_pq－０．５）というような式によって、可能なブーリアン値σ、及び、二次割り当て問題の制約において可能な設定、例えば

といったものに変換する。

変調部２４は、コスト関数Ｖ⁽⁰⁾の項の一つの現在の値を演算して、線形ゲインｐや、ターゲット振幅ａを、式（１２），（１３）により表されるように変調する：

そして、変化率β₀についても、同様に変調を行う。

これらの式において、

であり、関数φはシグモイド関数であり、

は、予め定められた定数パラメータである。

パラメータ値α，πは、予めユーザによって定められる。

ここで、Ｖ⁽⁰⁾（ｔ）は、状態データｘ（ｔ）に関連するコスト関数の項の一つの値であり、コスト評価部２２によって演算される。また、Ｖ_opt ⁽⁰⁾は、ターゲットエネルギーである。ここでの例では、Ｖ_opt ⁽⁰⁾は、繰り返し演算中に見いだされたエネルギーの最低値に設定される：

Ｖ_opt ⁽⁰⁾を求めるために、変調部２４は、現在のＶ_opt ⁽⁰⁾を記憶しており、現在のＶ⁽⁰⁾（ｔ）が記憶しているＶ_opt ⁽⁰⁾よりも小さい値であるときに、Ｖ_opt ⁽⁰⁾を更新する。初期の状態では、変調部２４は、Ｖ_opt ⁽⁰⁾を記憶していないので、変調部２４は、演算されたＶ⁽⁰⁾（ｔ＝０）をそのままＶ_opt ⁽⁰⁾として記憶する。

変調部２４は、更新された線形ゲインｐとターゲット振幅ａ、及び、変化率β₀を求めて、これらを出力する。

制御部１１は、現在の更新された状態データｘと、コスト関数とを出力する。そして次の繰り返し処理を実行する。

ＱＡＰの場合、上述のように、データ処理装置１の動作は、

and

と、まとめられる。

この式において、ｉ，ａは、ｉ番目の工場、ａ番目の場所に対応する。そしてここでｈ⁽⁰⁾ _ｉａは、（２４）式で定義される。

ｈ⁽⁰⁾の演算は、行列とベクトルとの乗算、例えば

を含むが、行列のサイズがＮ^２×Ｎ^２であるので、これらを直接演算することは現実的でない。さらに、この積の演算は、行列自体がテンソル積の演算の結果であることを考慮して行われるのである。

例えば、

なので、積の演算は：

のように行われる。ここで［Ｑ］は、要素がＱのそれぞれであるようなベクトルを意味する。

［（３）リード最適化（lead optimization）問題を用いた構造設計］
リード最適化問題は、元素の配置（geometry）と構成元素が既知である場合に、化合物候補の構造を見いだす問題である。つまり、この組み合わせ最適化問題の目的は、所与のたんぱく質に関して、相互作用エネルギーを最小化するよう、元素を既知の位置に配置することである。

仮に、空間的な構造が、例えば、－X1－X2－X3－X4－X5－X6－という形式であり、σ_inが、=CH－, =CH2, =O, －NH2, －OH, －CH3のうちから選択されたもの（ここで、記号「＝」または「－」は、それぞれ化学結合（chemical bond）を表す）であるとするとき、可能な解の一つは、 －CH=CH－CH=CH－CH=CH－ という鎖状の構造となる。

候補となる構造は、次の２つの条件を満足しなければならない：
（１）隣接する元素との間の化学結合の整合性は満足されなければならない、また、
（２）一つの位置あたり、元素が一つだけ割り当てられる必要がある。

以下では、簡単にするため、この問題の対象が、対象となるたんぱく質の相互作用エネルギー（interaction energy）を考慮することなく、これらの２つの条件を満足する候補を見いだすことであると限定する。

提案される構造は、制約付き組み合わせ最適化問題を解くために使われる。第一に、状態エンコード変数σ_in＝±１は、次のことを考慮してエンコードされたものとする。すなわち、元素ｎが位置ｉに配置される場合にσ_in＝１、そうでない場合、σ_in＝－１とする。既に述べた２つの制約が次の式：

で与えられるコスト関数Ｖを用いたイジング問題に変換される（H. Sakaguchi, et. al.,"Boltzmann Sampling by Degenerate Optical Parametric Oscillator Network for Structure-Based Virtual Screening", Entropy, 18, 365(2016)による）。ここで、Ｖ⁽¹⁾及びＶ⁽²⁾は、それぞれ、化学結合の整合性に係る第一の制約と、単一性（unicity）に係る第二の制約とに関連するソフト制約のコスト関数である。

具体的には、例えば、

と表される。

これらの式において、

は、変数

の変化を考慮して、ユーザにより定義されたパラメータである。なお、ｓ_iは、０または１であり得る。

条件を満足する構造を見つける問題は、上記のコスト関数Ｖ^＊を最小化することに帰着される。

現実的には、ソフト制約を表すコスト関数Ｖ⁽¹⁾及びＶ⁽²⁾は、

というようにも表記できる。

ここで、位置ｉ，ｊが隣り合っており、かつ、位置ｉに配された元素ｎの化学結合が、位置ｊに配された元素ｍの化学結合と互換性がないならば、Ｊ_injm＝１とする。

一方で、位置ｉ，ｊが隣り合っていないか、または、位置ｉに配された元素ｎの化学結合が、位置ｊに配された元素ｍの化学結合と互換性があるならば、Ｊ_injm＝０とする。

Ｃ₁及びＣ₂は、ｓ_inとは独立した定数値であり、組み合わせ最適化問題に関わらない。

データ処理装置１の動作は、既に述べた最大カット問題やＱＡＰの例におけるのと同等であるので、冗長な説明は省略する。

本実施形態のこの例では、データ処理装置１は、次の式で表現される動作を行うように設定される：

（３４）式は、化学結合の整合性を表す。また（３５）式は単一性を表している。

ここで、

であり、位置ｉに配された元素ｎが整合しないならば、Ｅin＞０とする。

β_１及びβ_２は、それぞれ、ｅ_in ⁽¹⁾とｅ_i ⁽²⁾とのエラー値の変化率である。

β_１とβ_２とは、互いに異なる値となり得る。

本実施の形態のもう一つの例を以下に説明する。

状態データプロセッサ２１の状態データは、量子動力学（quantum dynamics）を用いて表現できる。この態様では、各状態データｘ_iに対応する状態データプロセッサ２１の各部は、状態データを密度行列ρ_iとして保持し、独立した各部の動きは、量子マスター方程式（quantum master equation）によって規定される。最も一般的な場合では、この密度行列はより小さい密度行列のテンソル積として表現することはできない場合、一つの密度行列ρ_i1i2…は、複数のユニットｉ_１，ｉ_２…の状態で記述できる。

状態データは、３つの異なる方法：
（１）量子レベルにある密度行列によって
（２）古典的レベルにあるアナログ変数によって、
（３）論理レベルにあるブーリアン状態によって、
エンコードされ得る。

量子的表現から古典的表現への変換は、量子測定（quantum measurement）によって行われる。アナログからブーリアン値への変換は、アナログ－デジタル変換器によって行われる。

具体的に、ここでの状態データプロセッサ２１は、例えばＵＳ２０１７／００２４６５８の公開公報に開示されたＤＯＰＯ（Degenerate Optical Parametric Oscillator）などのイジングモデル量子計算デバイスを用いて実現される。

図４には、P. L. McMahon et al., “A fully-programmable 100-spin coherent Ising machine with all-to-all connections”, Science 354, 614 (2016)による、縮退型光パラメトリック発振器に基づくコヒーレントイジングマシン（ＣＩＭ）が示される。

上述のように、ＤＯＰＯシステムの詳細は、ＵＳ２０１７／００２４６５８の公開公報に示されているので、ここではＤＯＰＯの概括的な構成のみを示す。本実施の形態のここでの例では、状態データプロセッサ２１は、図４に例示するようなＤＯＰＯシステムを用いて実現される。この例に係る状態データプロセッサ２１は、ポンプパルス発生器（pump pulse generator：ＰＰＧ）４１と、第二高調波発生デバイス（ＳＨＧ）４２と、周期的極性導波デバイス（Periodically-poled Wave Guide Device：ＰＰＷＧ）４３と、方向性結合器（Directional Coupler）４４と、ＡＤコンバータ４５と、ＦＰＧＡデバイス４６と、ＤＡコンバータ４７と、変調器４８と、リング・キャビティ４９とを含んで構成される。

このシステムにおいて、複数のポンプパルス光波（擬スピンパルス）が、ＰＰＧ４１とＳＨＧ４２とにより生成され、ＰＰＷＧ４３を介して、リング・キャビティ４９に時分割的に導入される。

ここで、複数のポンプパルス光波は、それぞれが疑似的に、複数のイジングモデルスピンに対応しており、相互に同じ発振周波数（identical oscillation frequency）を有している。

また、隣接する光波間の時間間隔Ｔは、Ｌをリング・キャビティ４９の長さ、ｃをリング・キャビティ４９内での光速、Ｎを、Ｎ＞０なる自然数としてＬ／ｃ／Ｎに設定される。

リング・キャビティ（リング共振器）４９を伝播する光波の一部は、第１の方向性結合器４４-1を介して、ＡＤコンバータ４５に導かれる。この後、第２の方向性結合器４４-2までリング・キャビティ４９を伝播し続ける光波の他の部分は、「ターゲット波」と呼ばれる。

ＡＤコンバータ４５は、光波の強度をデジタル値に変換し、当該デジタル値を、ＦＰＧＡデバイス４６に出力する。

換言すれば、第１の方向性結合器４４-1とＡＤコンバータ４５とはリング・キャビティ４９内の複数の擬スピンパルスが周回して伝播するごとに、当該複数の擬スピンパルスの位相及び振幅を仮計測する。

ＤＯＰＯシステムにおけるＦＰＧＡデバイス４６は、ＡＤコンバータ４５の出力（前回の状態データｘ_ｉ（ｔ）を表す）を、コスト評価部２２と、エラー演算部２３と、変調部２４との出力を用いて演算した、

に加算するように構成されてもよい。

ＦＰＧＡデバイス４６は、この加算の結果をＤＡコンバータ４７に出力する。このＤＡコンバータ４７の出力は、入力パルスの変調に利用される。

変調器４８は、他のポンプパルス光波を生成し、当該光波の振幅と位相とを、ＤＡコンバータ４７の出力に基づいて変調する。ＤＡコンバータ４７の出力は、ＦＰＧＡデバイス４６の出力に対応するアナログ値である。一例では、この変調器４８は、

が正であれば位相をπ／２だけ遅延させ、

の絶対値に比例した振幅変調を行う。またこの変換器４８は、

が負であれば、位相をπ／２だけ進めて、

の絶対値に比例した振幅変調を行う。

変調器４８により変調された光波は、第２の方向性結合器４４-2を介してリング・キャビティ４９内に導入される。第２の方向性結合器４４-2は、ターゲット波が第２の方向性結合器４４-2に到達したタイミングで、当該変調された光波をリング・キャビティ４９に導入する点に注意されたい。変調された光波とターゲット波から光波が合成され、更新された擬スピンパルスとなる。

上述のように、本実施の形態のこの例では、ＦＰＧＡデバイス４６が

で表される値を出力している。

合成された光波は、リング・キャビティ４９に沿って導かれており、ＦＰＧＡデバイス４６は状態データｘｉ（ｔ）の時間発展を表す値を繰り返し出力する。

つまり、この例のＤＯＰＯシステムは、ゲイン消失シミュレータ２１１として機能しているものである。ＤＯＰＯシステムのＦＰＧＡデバイス４６の出力は、ＤＡコンバータ４７に出力されるだけでなく、状態データプロセッサ２１の他の部分、例えば射影部２１２などにも出力される。

ＤＯＰＯの動作は広く知られているため、ここでのさらなる詳しい説明は省略する。

本発明の実施の形態のさらに他の態様では、データ処理装置１は、ＣＰＵ等のデジタルプロセッサを用いて構成される。この態様では、状態データプロセッサ２１と、コスト評価部２２と、エラー演算部２３と、変調部２４は、ＣＰＵによって実行されるソフトウエアプログラムとして実装される。このプログラムは、ＣＰＵに接続されたメモリデバイス等の記憶手段にインストールされ得る。このメモリデバイス等の記憶手段はＣＰＵが利用する状態データや、エラー値等のデータを保持する。この実施態様は、ディスプレイデバイスや、キーボードなどの入力デバイスをさらに含む一般的なコンピュータデバイスによって実現されてもよい。

本発明の実施の形態のさらに別の態様では、データ処理装置１は、ＦＰＧＡ、ＧＰＵ等のデジタルプロセッサを用いて構成され得る。この態様では、状態データプロセッサ２１、コスト評価部２２、エラー演算部２３、及び変調部２４は、統合プロセスの後で得られる論理ゲートの態様で実装設計されて実現される。この実施形態は、ディスプレイデバイスや、キーボードなどの入力デバイスなどを備えた一般的なコンピュータデバイスに接続され得る。

図５は、本発明の実施形態の一態様に係る機能的構成を表す概略図である。図５に示されるように、本発明のある実施態様によれば、データ処理装置１は、状態データノード３１と、一次エラーノード３２と、高次エラーノード３３とを含む。

状態データノード３１は、ブーリアン値σ_ｉ、アナログ値ｘ_ｉ、及び量子的表現（密度行列）ρ_ｉのうちから選択されたもので表される状態データを保持する。

状態データノードは、互いに接続され、ある状態データノードが保持する値は、コスト関数Ｖ^＊を介して他の状態データノードが保持する値に影響を与える。

一次エラーノード３２の各々は、対応する状態データノード３１に接続されている。一次エラーノード３２は、対応する状態データノード３１が保持する状態データを、当該状態データの振幅の不均一性を減少させるように補正する。

高次エラーノード３３は、少なくとも一つの状態データノード３１に接続される。すべての状態データノード３１が高次エラーノード３３に接続されている必要はない。換言すると、高次エラーノード３３と状態エラーノードとの接続関係は「非対称」である。高次エラーノード３３と状態データノード３１との間の接続は、解かれるべき問題によって規定される。高次エラーノード３３は、接続されている一つまたは複数の状態データノード３１が保持する状態データを、解かれるべき問題に関する整合性を維持させるために変化させる。

上述の例は、例示に過ぎず、当業者はこれらの例示した実施形態に基づいて、新規事項や効果を追加することなく種々の変形例があることを容易に理解できるものである。すなわち、それらの変形例は本発明の実施の形態として実質的に開示されている。

Claims (14)

1. 所与の問題を解くよう構成されたデータ処理装置であって、
   所定の時間発展的処理によって状態データの更新を繰り返す状態データ処理部と、
   現在の状態データのコスト関数の値を演算するコスト評価部と、
   現在の状態データの振幅の均一性に関するエラー値を演算するエラー演算部と、を有し、
   前記状態データは、変数の集合であり、前記状態データ処理部は、状態データについて時間発展的処理を行い、前記コスト関数と前記エラー演算部によって演算された、現在の状態データの振幅の均一性に関する前記エラー値とに基づいて、現在の状態データを更新する、データ処理装置。
2. 請求項１に記載のデータ処理装置であって、
   前記状態データ処理部は、状態データについて時間発展的処理を行い、前記コスト関数の勾配と、前記エラー値とに基づいて現在の状態データを更新するデータ処理装置。
3. 請求項２に記載のデータ処理装置であって、
   前記時間発展的処理は、現在の状態データに基づくゲイン値によって単安定から双安定へ変化するアーキタイプポテンシャルに基づいて現在の状態データを更新するデータ処理装置。
4. 請求項１に記載のデータ処理装置であって、
   前記エラー演算部は、少なくとも一つの、解決するべき問題の制約に関わるエラー値をさらに演算するデータ処理装置。
5. 請求項１に記載のデータ処理装置であって、
   前記状態データ処理部は、さらに、前記状態データを所定の部分空間に射影する処理を行うデータ処理装置。
6. 請求項１に記載のデータ処理装置であって、さらに、
   前記現在の状態データに基づいてターゲット振幅を演算する変調部を有し、
   前記エラー演算部が、当該ターゲット振幅を利用して、現在の状態データの振幅の均一性に関する前記エラー値を演算するデータ処理装置。
7. 請求項６に記載のデータ処理装置であって、
   前記状態データ処理部は、さらに、前記状態データを、前記ターゲット振幅で規定された所定の部分空間に射影する処理を行うデータ処理装置。
8. 請求項１に記載のデータ処理装置であって、さらに、
   前記現在の状態データに基づいて、ゲイン値を演算する変調部を有し、
   前記状態データ処理部は、前記現在の状態データに基づいて演算されたゲイン値を利用してアーキタイプ単安定／双安定ポテンシャルも更新するデータ処理装置。
9. 請求項１に記載のデータ処理装置であって、さらに、
   エラー値の変化率を動的に定める変調部を有し、
   前記エラー演算部は、エラー値についての時間発展的処理を行い、前記エラー値の変化率と、前記エラー値と、前記状態データと、前記コスト関数とに基づいて現在のエラー値を更新するデータ処理装置。
10. 請求項８に記載のデータ処理装置であって、
    前記ゲイン値は、前記時間発展的処理の開始前に、前記コスト関数のパラメータを用いて演算された定数に依存するデータ処理装置。
11. 請求項１に記載のデータ処理装置であって、
    前記状態データ処理部が、イジングモデル量子計算デバイスを用いて、コスト関数に基づき、現在の状態データを更新して、状態データの時間発展的処理を実行し、
    前記イジングモデル量子計算デバイスは、複数のイジングモデルスピンを模した複数の擬スピンパルスであって、相互に同じ発振周波数を有する擬スピンパルスをパラメトリックに発振するパラメトリック発振器と、
    前記複数の擬スピンパルスを、その内部で周回させて伝播するリング共振器と、
    前記複数の擬スピンパルスが前記リング共振器内を周回して到来する度に、前記複数の擬スピンパルスの位相と振幅とを仮に計測し、前記複数の擬スピンパルスの擬スピンを仮計測する仮スピン計測手段と、
    前記コスト評価部の出力と前記エラー値とに従い、状態データの時間発展を得るべく前記計測された擬スピンパルスに合成されるべきデータを演算して、当該状態データの時間発展の結果を出力するＦＰＧＡデバイスと、を有するデータ処理装置。
12. 請求項１に記載のデータ処理装置であって、
    前記状態データ処理部は、状態データとエラー値の演算を行うＧＰＵまたはＦＰＧＡ、あるいは、アナログ電子コンピューティングデバイスを用いて、前記状態データの時間発展的処理を実行し、前記コスト関数に基づいてデータを更新するデータ処理装置。
13. 状態データ処理部と、コスト評価部と、エラー演算部とを備えて、所与の問題を解くデータ処理装置を用いたデータ処理方法であって、
    前記状態データ処理部が、所定の時間発展的処理によって状態データの更新を繰り返し、
    前記コスト評価部が、現在の状態データのコスト関数の値を演算し、
    前記エラー演算部が、現在の状態データの振幅の均一性に関するエラー値を演算し、
    前記状態データは、変数の集合であり、前記状態データ処理部が、状態データについて時間発展的処理を行い、前記コスト関数と前記エラー演算部によって演算された、現在の状態データの振幅の均一性に関する前記エラー値とに基づいて、現在の状態データを更新する、データ処理方法。
14. 請求項１３に記載のデータ処理方法であって、
    前記データ処理装置は、さらに、変調部を備え、
    当該変調部が、振幅の均一性と制約とに係るエラー変数を更新し、前記制約と前記現在の状態データとの間のカップリング項を演算し、ゲインと飽和とを演算し、
    前記状態データ処理部が、状態データの更新のためのインジェクション項を演算し、当該状態データ処理部が所定の時間発展的処理により状態データを繰り返し更新している間に、当該状態データを有効な部分空間に射影し、
    エラー演算部が、前記現在の状態データの振幅の均一性と制約とに関わるエラー値を演算し、
    前記繰り返し処理が所定の条件が満足されるまで繰り返されるデータ処理方法。

Images