(57)【要約】      （修正有）【課題】組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算するにあたり、光学回路を安定化する必要をなくし、機械学習のボルツマンサンプリングを行なうために、機械学習を連続値スピンモデルへマッピングし、連続値スピンモデルの平衡状態を計算するにあたり、光学回路を安定化する必要をなくす。【解決手段】結合振動子系のダイナミクスを記述する微分方程式を、光学回路で実装するのではなく、電子回路で実装し、上記の微分方程式のイテレーションの反復後に、上記の微分方程式の定常状態における各振動子毎の位相を計算することにより、結合振動子系の相互結合に対応するイジングモデル又は連続値スピンモデルの基底状態又は平衡状態を計算する。【選択図】図３
2019008502.tif 000002
【特許請求の範囲】
【請求項１】
結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーション部と、
前記イテレーション部により反復されたイテレーションの後に、前記微分方程式の定常状態における各振動子毎の位相を計算することにより、前記結合振動子系の相互結合に対応するスピン系の基底状態又は平衡状態を計算するスピン系計算部と、
を備え、
前記イテレーション部は、
前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算部と、
前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算部と、
前記微分方程式の各イテレーションにおいて、前記利得損失計算部及び前記相互結合計算部により各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算部と、
を備えることを特徴とする結合振動子系の計算装置。
【請求項２】
前記微分方程式の各イテレーションにおいて、各振動子毎の雑音に起因する各振動子毎の振幅及び位相の微小時間変化を計算する雑音計算部、
をさらに備え、
前記振幅位相計算部は、前記微分方程式の各イテレーションにおいて、前記利得損失計算部、前記相互結合計算部及び前記雑音計算部により各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する
ことを特徴とする、請求項１に記載の結合振動子系の計算装置。
【請求項３】
結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーション部、
を備え、
前記イテレーション部は、
前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算部と、
前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算部と、
前記微分方程式の各イテレーションにおいて、各振動子毎の雑音に起因する各振動子毎の振幅及び位相の微小時間変化を計算する雑音計算部と、
前記微分方程式の各イテレーションにおいて、前記利得損失計算部、前記相互結合計算部及び前記雑音計算部により各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算部と、
を備えることを特徴とする結合振動子系の計算装置。
【請求項４】
前記雑音計算部は、前記微分方程式の初期状態のみにおいて、各振動子毎の雑音に起因する各振動子毎の振幅及び位相の微小時間変化を計算する
ことを特徴とする、請求項２又は３に記載の結合振動子系の計算装置。
【請求項５】
前記イテレーション部は、前回のイテレーションにおける各振動子毎の微小時間後の振幅及び位相が、前記微分方程式の定常状態に所定の精度で接近したときに、次回のイテレーションにおける各振動子毎の微小時間後の振幅及び位相の計算を停止する
ことを特徴とする、請求項１から４のいずれかに記載の結合振動子系の計算装置。
【請求項６】
前記イテレーション部は、前記微分方程式の各イテレーションにおいて、各振動子毎の並列処理により、各振動子毎の微小時間後の振幅及び位相を計算する
ことを特徴とする、請求項１から５のいずれかに記載の結合振動子系の計算装置。
【請求項７】
前記イテレーション部は、複数のイテレーション部を備え、
前記複数のイテレーション部の各々は、前記微分方程式の各イテレーションにおいて、前記結合振動子系のうちの各々に特定のグループに属する振動子の微小時間後の振幅及び位相のデータを、前記複数のイテレーション部のその他との間で共有する
ことを特徴とする、請求項６に記載の結合振動子系の計算装置。
【請求項８】
前回のイテレーションにおける前記データの共有処理が、次回のイテレーションにおける計算処理を停止／中断させないように、前記複数のイテレーション部の各々は、前記結合振動子系のうちの各々に特定のグループに属する振動子を割り当てられる
ことを特徴とする、請求項７に記載の結合振動子系の計算装置。
【請求項９】
前記振幅位相計算部は、前記微分方程式の各イテレーションにおいて、所定の精度で量子化したうえで、各振動子毎の微小時間後の振幅及び位相を計算する
ことを特徴とする、請求項１から８のいずれかに記載の結合振動子系の計算装置。
【請求項１０】
前記振幅位相計算部により計算される各振動子毎の微小時間後の振幅についての量子化カットオフは、前記振幅位相計算部により計算される各振動子毎の微小時間後の振幅についての量子化の前記所定の精度を満足するように設定されており、かつ、前記微分方程式の定常状態における各振動子毎の振幅と比べて大きく設定されている
ことを特徴とする、請求項９に記載の結合振動子系の計算装置。
【請求項１１】
前記相互結合計算部により計算される各振動子毎の振幅の微小時間変化についての量子化ビット数は、前記振幅位相計算部により計算される各振動子毎の微小時間後の振幅についての量子化ビット数と比べて、少なくともｌｏｇ２Ｍ（ただし、Ｍは、前記結合振動子系の全振動子数である。）だけ大きく設定されている
ことを特徴とする、請求項９又は１０に記載の結合振動子系の計算装置。
【請求項１２】
結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーションステップと、
前記イテレーションステップにより反復されたイテレーションの後に、前記微分方程式の定常状態における各振動子毎の位相を計算することにより、前記結合振動子系の相互結合に対応するスピン系の基底状態又は平衡状態を計算するスピン系計算ステップと、
をコンピュータに実行させ、
前記イテレーションステップは、
前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算ステップと、
前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算ステップと、
前記微分方程式の各イテレーションにおいて、前記利得損失計算ステップ及び前記相互結合計算ステップにより各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算ステップと、
を備えることを特徴とする結合振動子系の計算プログラム。
【請求項１３】
結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーションステップと、
前記イテレーションステップにより反復されたイテレーションの後に、前記微分方程式の定常状態における各振動子毎の位相を計算することにより、前記結合振動子系の相互結合に対応するスピン系の基底状態又は平衡状態を計算するスピン系計算ステップと、
を備え、
前記イテレーションステップは、
前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算ステップと、
前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算ステップと、
前記微分方程式の各イテレーションにおいて、前記利得損失計算ステップ及び前記相互結合計算ステップにより各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算ステップと、
を備えることを特徴とする結合振動子系の計算方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本開示は、スピン系の基底状態又は平衡状態を計算することにより、そのスピン系へマッピングされた組み合わせ最適化問題の最適解を求めて、そのスピン系へマッピングされた機械学習のボルツマンサンプリングを行なう技術に関する。
【背景技術】
【０００２】
組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算する技術が知られている。機械学習のボルツマンサンプリングを行なうために、機械学習を連続値スピンモデルへマッピングし、連続値スピンモデルの平衡状態を計算する技術も知られている。
【０００３】
特許文献１、２では、複数のレーザパルスで複数のイジングスピンを実装する。ここで、複数のレーザパルスの間の相互注入を行ない、複数のイジングスピンの間の相互結合を実装する。そして、複数のレーザパルス全体の閾値利得が最低となる発振モードを立ち上げ、複数のイジングスピン全体のエネルギーが最低となる基底状態を実装する。
【先行技術文献】
【特許文献】
【０００４】
【特許文献１】特表２０１６−５２８６１１号公報
【特許文献２】国際公開第２０１５／１５６１２６号
【発明の概要】
【発明が解決しようとする課題】
【０００５】
特許文献１のイジングモデルの計算装置Ｃ１の構成を図１に示す。レーザ発振部１１は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・を生成する。リング共振器１２は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・を周回させる。
【０００６】
パルス遅延線１４は、１パルス間隔分の遅延量を有する。相互結合実装部１３は、パルス遅延線１４を伝搬して相互結合係数Ｊ１２、Ｊ２３、Ｊ３４、Ｊ４１を反映した振幅及び位相を有するレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４を、パルス遅延線１４を伝搬せずリング共振器１２を伝搬するレーザパルスＰ２、Ｐ３、Ｐ４、Ｐ１にそれぞれ相互注入する。
【０００７】
パルス遅延線１６は、２パルス間隔分の遅延量を有する。相互結合実装部１５は、パルス遅延線１６を伝搬して相互結合係数Ｊ１３、Ｊ２４、Ｊ３１、Ｊ４２を反映した振幅及び位相を有するレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４を、パルス遅延線１６を伝搬せずリング共振器１２を伝搬するレーザパルスＰ３、Ｐ４、Ｐ１、Ｐ２にそれぞれ相互注入する。
【０００８】
パルス遅延線１８は、３パルス間隔分の遅延量を有する。相互結合実装部１７は、パルス遅延線１８を伝搬して相互結合係数Ｊ１４、Ｊ２１、Ｊ３２、Ｊ４３を反映した振幅及び位相を有するレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４を、パルス遅延線１８を伝搬せずリング共振器１２を伝搬するレーザパルスＰ４、Ｐ１、Ｐ２、Ｐ３にそれぞれ相互注入する。
【０００９】
スピン測定部１９は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・が定常状態に到達した後に、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・の位相を測定することにより、複数のレーザパルスの位相に対応する複数のイジングスピンの方向を測定する。
【００１０】
しかし、特許文献１では、リング共振器１２の光路長を安定化しなければ、複数のレーザパルスを安定化することができない。そして、パルス遅延線１４、１６、１８の光路長を安定化しなければ、複数のレーザパルスの間の相互注入を正確に実行することができず、複数のイジングスピンの間の相互結合を正確に実装することができない。
【００１１】
特許文献２のイジングモデルの計算装置Ｃ２の構成を図２に示す。レーザ発振部２１は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・を生成する。リング共振器２２は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・を周回させる。
【００１２】
スピン測定部２３は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・の位相を測定する。相互結合実装部２４は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・の間の相互注入を行ない、複数のイジングスピンの間の相互結合を実装する。
【００１３】
相互結合実装部２４は、相互結合係数ΣＪ１ｋσｋ（ｋ＝２、３、４）を電子回路上で計算し、相互結合係数ΣＪ１ｋσｋ（ｋ＝２、３、４）を反映した振幅及び位相を有するレーザパルスを生成し、生成したレーザパルスをレーザパルスＰ１に注入する。
【００１４】
相互結合実装部２４は、相互結合係数ΣＪ２ｋσｋ（ｋ＝１、３、４）を電子回路上で計算し、相互結合係数ΣＪ２ｋσｋ（ｋ＝１、３、４）を反映した振幅及び位相を有するレーザパルスを生成し、生成したレーザパルスをレーザパルスＰ２に注入する。
【００１５】
相互結合実装部２４は、相互結合係数ΣＪ３ｋσｋ（ｋ＝１、２、４）を電子回路上で計算し、相互結合係数ΣＪ３ｋσｋ（ｋ＝１、２、４）を反映した振幅及び位相を有するレーザパルスを生成し、生成したレーザパルスをレーザパルスＰ３に注入する。
【００１６】
相互結合実装部２４は、相互結合係数ΣＪ４ｋσｋ（ｋ＝１、２、３）を電子回路上で計算し、相互結合係数ΣＪ４ｋσｋ（ｋ＝１、２、３）を反映した振幅及び位相を有するレーザパルスを生成し、生成したレーザパルスをレーザパルスＰ４に注入する。
【００１７】
スピン測定部２３は、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・が定常状態に到達した後に、複数のレーザパルスＰ１、Ｐ２、Ｐ３、Ｐ４、・・・の位相を測定することにより、複数のレーザパルスの位相に対応する複数のイジングスピンの方向を測定する。
【００１８】
つまり、特許文献２では、相互結合実装部２４が有する電子回路を用いて、複数のレーザパルスの間の相互注入を正確に実行することができ、複数のイジングスピンの間の相互結合を正確に実装することができる。しかし、リング共振器２２の光路長を安定化しなければ、複数のレーザパルスを安定化することができない。
【００１９】
そこで、前記課題を解決するために、本開示は、組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算するにあたり、光学回路を安定化する必要をなくすことを目的とする。
【００２０】
そして、前記目的に加えて、本開示は、機械学習のボルツマンサンプリングを行なうために、機械学習を連続値スピンモデルへマッピングし、連続値スピンモデルの平衡状態を計算するにあたり、光学回路を安定化する必要をなくすことも目的とする。
【課題を解決するための手段】
【００２１】
前記課題を解決するために、結合振動子系のダイナミクスを記述する微分方程式を、光学回路で実装するのではなく、電子回路で実装することとした。そのうえで、上記の微分方程式のイテレーションの反復後に、上記の微分方程式の定常状態における各振動子毎の位相を計算することにより、結合振動子系の相互結合に対応するイジングモデル又は連続値スピンモデルの基底状態又は平衡状態を計算することとした。
【００２２】
具体的には、本開示は、結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーション部と、前記イテレーション部により反復されたイテレーションの後に、前記微分方程式の定常状態における各振動子毎の位相を計算することにより、前記結合振動子系の相互結合に対応するスピン系の基底状態又は平衡状態を計算するスピン系計算部と、を備え、前記イテレーション部は、前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算部と、前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算部と、前記微分方程式の各イテレーションにおいて、前記利得損失計算部及び前記相互結合計算部により各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算部と、を備えることを特徴とする結合振動子系の計算装置である。
【００２３】
また、本開示は、結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーションステップと、前記イテレーションステップにより反復されたイテレーションの後に、前記微分方程式の定常状態における各振動子毎の位相を計算することにより、前記結合振動子系の相互結合に対応するスピン系の基底状態又は平衡状態を計算するスピン系計算ステップと、をコンピュータに実行させ、前記イテレーションステップは、前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算ステップと、前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算ステップと、前記微分方程式の各イテレーションにおいて、前記利得損失計算ステップ及び前記相互結合計算ステップにより各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算ステップと、を備えることを特徴とする結合振動子系の計算プログラムである。
【００２４】
また、本開示は、結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーションステップと、前記イテレーションステップにより反復されたイテレーションの後に、前記微分方程式の定常状態における各振動子毎の位相を計算することにより、前記結合振動子系の相互結合に対応するスピン系の基底状態又は平衡状態を計算するスピン系計算ステップと、を備え、前記イテレーションステップは、前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算ステップと、前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算ステップと、前記微分方程式の各イテレーションにおいて、前記利得損失計算ステップ及び前記相互結合計算ステップにより各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算ステップと、を備えることを特徴とする結合振動子系の計算方法である。
【００２５】
この構成によれば、組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。
【００２６】
そして、機械学習のボルツマンサンプリングを行なうために、機械学習を連続値スピンモデルへマッピングし、連続値スピンモデルの平衡状態を計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。
【００２７】
また、本開示は、前記微分方程式の各イテレーションにおいて、各振動子毎の雑音に起因する各振動子毎の振幅及び位相の微小時間変化を計算する雑音計算部、をさらに備え、前記振幅位相計算部は、前記微分方程式の各イテレーションにおいて、前記利得損失計算部、前記相互結合計算部及び前記雑音計算部により各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算することを特徴とする結合振動子系の計算装置である。
【００２８】
この構成によれば、結合振動子系のダイナミクスを記述する微分方程式の定常状態として、複数の結合振動子の振幅が０である状態を回避することができる。
【００２９】
また、本開示は、結合振動子系のダイナミクスを記述する微分方程式の各イテレーションにおいて、各振動子毎の微小時間後の振幅及び位相を計算するイテレーション部、を備え、前記イテレーション部は、前記微分方程式の各イテレーションにおいて、各振動子毎の利得及び損失に起因する各振動子毎の振幅及び位相の微小時間変化を計算する利得損失計算部と、前記微分方程式の各イテレーションにおいて、各振動子の間の相互結合に起因する各振動子毎の振幅及び位相の微小時間変化を計算する相互結合計算部と、前記微分方程式の各イテレーションにおいて、各振動子毎の雑音に起因する各振動子毎の振幅及び位相の微小時間変化を計算する雑音計算部と、前記微分方程式の各イテレーションにおいて、前記利得損失計算部、前記相互結合計算部及び前記雑音計算部により各々計算された各振動子毎の振幅及び位相の微小時間変化に基づいて、各振動子毎の微小時間後の振幅及び位相を計算する振幅位相計算部と、を備えることを特徴とする結合振動子系の計算装置である。
【００３０】
この構成によれば、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。
【００３１】
また、本開示は、前記雑音計算部は、前記微分方程式の初期状態のみにおいて、各振動子毎の雑音に起因する各振動子毎の振幅及び位相の微小時間変化を計算することを特徴とする結合振動子系の計算装置である。
【００３２】
この構成によれば、複数の結合振動子の振幅が０である初期状態から脱出することができた後には、雑音計算部の計算負担をなくすことができる。
【００３３】
また、本開示は、前記イテレーション部は、前回のイテレーションにおける各振動子毎の微小時間後の振幅及び位相が、前記微分方程式の定常状態に所定の精度で接近したときに、次回のイテレーションにおける各振動子毎の微小時間後の振幅及び位相の計算を停止することを特徴とする結合振動子系の計算装置である。
【００３４】
この構成によれば、複数の結合振動子の振幅が有限値である定常状態に接近することができた後には、イテレーション部の計算負担をなくすことができる。
【００３５】
また、本開示は、前記イテレーション部は、前記微分方程式の各イテレーションにおいて、各振動子毎の並列処理により、各振動子毎の微小時間後の振幅及び位相を計算することを特徴とする結合振動子系の計算装置である。
【００３６】
結合振動子系の全振動子数をＭとすると、各振動子毎の並列処理を行わないときには、２体相互結合の計算量は、Ｍ２のオーダーに増えるが、各振動子毎の並列処理を行なうときには、２体相互結合の計算量は、ｌｏｇ２Ｍのオーダーに減らせる。
【００３７】
また、本開示は、前記イテレーション部は、複数のイテレーション部を備え、前記複数のイテレーション部の各々は、前記微分方程式の各イテレーションにおいて、前記結合振動子系のうちの各々に特定のグループに属する振動子の微小時間後の振幅及び位相のデータを、前記複数のイテレーション部のその他との間で共有することを特徴とする結合振動子系の計算装置である。
【００３８】
この構成によれば、結合振動子系の全振動子数が膨大な数になるときでも、各々のイテレーション部を各々のプロセッサに実装することができる。
【００３９】
また、本開示は、前回のイテレーションにおける前記データの共有処理が、次回のイテレーションにおける計算処理を停止／中断させないように、前記複数のイテレーション部の各々は、前記結合振動子系のうちの各々に特定のグループに属する振動子を割り当てられることを特徴とする結合振動子系の計算装置である。
【００４０】
この構成によれば、各々のイテレーション部を各々のプロセッサに実装するときでも、各々のイテレーション部の間のデータの共有処理を円滑に実行することができる。
【００４１】
また、本開示は、前記振幅位相計算部は、前記微分方程式の各イテレーションにおいて、所定の精度で量子化したうえで、各振動子毎の微小時間後の振幅及び位相を計算することを特徴とする結合振動子系の計算装置である。
【００４２】
この構成によれば、結合振動子系の計算装置の計算負担を減らすことができるとともに、結合振動子系の計算装置のメモリ容量を減らすことができる。
【００４３】
また、本開示は、前記振幅位相計算部により計算される各振動子毎の微小時間後の振幅についての量子化カットオフは、前記振幅位相計算部により計算される各振動子毎の微小時間後の振幅についての量子化の前記所定の精度を満足するように設定されており、かつ、前記微分方程式の定常状態における各振動子毎の振幅と比べて大きく設定されていることを特徴とする結合振動子系の計算装置である。
【００４４】
この構成によれば、各振動子毎の振幅の分解能を低くし過ぎず、かつ、各振動子毎の振幅の制限を厳しくし過ぎず、結合振動子系の計算装置の計算負担を減らすことができるとともに、結合振動子系の計算装置のメモリ容量を減らすことができる。
【００４５】
また、本開示は、前記相互結合計算部により計算される各振動子毎の振幅の微小時間変化についての量子化ビット数は、前記振幅位相計算部により計算される各振動子毎の微小時間後の振幅についての量子化ビット数と比べて、少なくともｌｏｇ２Ｍ（ただし、Ｍは、前記結合振動子系の全振動子数である。）だけ大きく設定されていることを特徴とする結合振動子系の計算装置である。
【００４６】
この構成によれば、各振動子毎の振幅の分解能と比べて、各振動子毎の相互結合の大きさの分解能を低くし過ぎず、結合振動子系の計算装置の計算負担を減らすことができるとともに、結合振動子系の計算装置のメモリ容量を減らすことができる。
【発明の効果】
【００４７】
このように、本開示によれば、組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算するにあたり、光学回路を安定化する必要をなくすことができる。
【００４８】
そして、本開示によれば、機械学習のボルツマンサンプリングを行なうために、機械学習を連続値スピンモデルへマッピングし、連続値スピンモデルの平衡状態を計算するにあたり、光学回路を安定化する必要をなくすことができる。
【図面の簡単な説明】
【００４９】
【図１】特許文献１のイジングモデルの計算装置の構成を示す図である。
【図２】特許文献２のイジングモデルの計算装置の構成を示す図である。
【図３】第１実施形態の結合振動子系の計算装置の構成を示す図である。
【図４】第１実施形態の利得損失計算部の処理内容を示す図である。
【図５】第１実施形態の結合振動子系の計算装置の計算結果を示す図である。
【図６】第１実施形態のイテレーション部の第１構成を示す図である。
【図７】第１実施形態のイテレーション部の第２構成を示す図である。
【図８】第１実施形態のイテレーション部の第３構成を示す図である。
【図９】第１実施形態の振幅位相計算部の量子化を示す図である。
【図１０】第１実施形態の結合振動子系の計算装置の計算結果を示す図である。
【図１１】第１実施形態の結合振動子系の計算装置の計算結果を示す図である。
【図１２】第１実施形態の結合振動子系の計算装置の計算結果を示す図である。
【図１３】第２実施形態の結合振動子系の計算装置の構成を示す図である。
【図１４】第３実施形態の結合振動子系の計算装置の構成を示す図である。
【発明を実施するための形態】
【００５０】
添付の図面を参照して本開示の実施形態を説明する。以下に説明する実施形態は本開示の実施の例であり、本開示は以下の実施形態に制限されるものではない。
【００５１】
（第１実施形態：計算原理）
第１実施形態では、組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算する。または、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算する。
【００５２】
そこで、複数の振動子で複数のイジングスピンを実装する。ここで、複数の振動子の間の相互結合を行なうことにより、複数のイジングスピンの間の相互結合を実装する。そして、複数の振動子全体の閾値利得が最低となる発振モードを立ち上げることにより、複数のイジングスピン全体のエネルギーが最低となる基底状態を実装する。
【００５３】
まず、複数の振動子全体の閾値利得が最低となる発振モードを立ち上げることにより、複数のイジングスピン全体のエネルギーが最低となる基底状態を実装することができる理由について説明する。イジングモデルのハミルトニアンは、数１で表される。ここで、Ｊｉｊは、スピンｉ、ｊの間の相互結合係数である。そして、σｉ、σｊ（＝＋１又は−１）は、それぞれ、スピンｉ、ｊの位相（０ｒａｄ又はπｒａｄ）を表わす。
【数１】
2019008502.tif 000003
【００５４】
各振動子ｉにおいて、Ｉ成分強度ｃｉ及びＱ成分強度ｓｉについて、レート方程式は、ファンデルポール方程式に対応して、数２、３のようになる。
【数２】
2019008502.tif 000004
【数３】
2019008502.tif 000005
【００５５】
ｔは、無次元時間であり、ｔ＝γＳτ／２である。τは、実時間である。γＳは、振動子の強度減衰率である。ｃｉ及びｓｉは、それぞれ、規格化後のＩ成分及びＱ成分の強度であり、ｃｉ＝Ｃｉ／Ａｓ及びｓｉ＝Ｓｉ／Ａｓである。Ｃｉ及びＳｉは、それぞれ、規格化前のＩ成分及びＱ成分の強度である。規格化因子Ａｓは、ｐ（後述する規格化後のポンプレート）＝２における振動子の強度であり、Ａｓ＝√（γＳγｐ／２κ２）である。γｐは、ポンピングの強度減衰率である。κは、縮退パラメトリックゲインである。ｐは、規格化後のポンプレートであり、ｐ＝Ｆｐ／Ｆｔｈである。Ｆｐは、規格化前のポンプレートである。規格化因子Ｆｔｈは、閾値ポンプレートであり、Ｆｔｈ＝γＳ√（γｐ）／４κである。
【００５６】
数２の−ｃｉ及び数３の−ｓｉは、損失に関わる項である。数２の＋ｐｃｉ及び数３の−ｐｓｉは、線形利得に関わる項である。数２の−（ｃｉ２＋ｓｉ２）ｃｉ及び数３の−（ｃｉ２＋ｓｉ２）ｓｉは、飽和利得に関わる項である。これらの項は、２体の相互結合等による摂動項を含まない、ファンデルポール方程式を構成する。
【００５７】
数２、３のξｉｊ（Ｊｉｊに比例）が関わる項は、２体の相互結合に関わる項であり、ファンデルポール方程式に対する、２体の相互結合による摂動項である。数２、３のｄＷ（確率微分）が関わる項は、振動子強度の増幅段階等での雑音に関わる項であり、ファンデルポール方程式に対する、振動子強度の増幅段階等での雑音による摂動項である。
【００５８】
定常状態において、数２、３は、それぞれ、数４、５のようになる。
【数４】
2019008502.tif 000006
【数５】
2019008502.tif 000007
【００５９】
数２のｐ−（ｃｉ２＋ｓｉ２）は、各振動子ｉについての飽和利得である。ここで、定常状態において、複数の振動子全体についての飽和利得は、複数の振動子全体についての強度減衰率に等しく、Ｉ成分ｃｉは、有限値であるが、Ｑ成分ｓｉは、０である。よって、複数の振動子全体についての強度減衰率Γは、数６のようになる。
【数６】
2019008502.tif 000008
【００６０】
ここで、数６の最右辺の第１項は、数４の左辺の第３項を摂動項としたときにおける、摂動の０次の寄与を示す。そして、数６の最右辺の第２項は、数４の左辺の第３項を摂動項としたときにおける、摂動の１次の寄与を示す。さらに、σｉ＝ｓｇｎ（ｃｉ）〜ｓｇｎ（ｃｉ（０））（ｃｉ（０）は摂動の０次の寄与）を用いている。
【００６１】
ポンプレートｐを漸増制御するときには、複数の振動子全体として、「最初に」最小の強度減衰率Γを実現する発振位相状態｛σｉ｝が選択される。つまり、複数の振動子全体として、１個の特定の発振モードが選択される。そして、発振モードの間の競合に起因して、１個の特定の発振モードは、他の発振モードを抑制する。つまり、複数の振動子全体として、数６のΓは最小化される。一方で、複数の振動子全体として、数６のＭは一定である。よって、複数の振動子全体として、数６の−Σξｉｊσｉσｊは最小化される。つまり、数１のイジングハミルトニアンを最小化する基底状態が実現されたことになる。
【００６２】
（第１実施形態：計算方法）
次に、第１実施形態の結合振動子系の計算装置の構成を図３に示す。結合振動子系の計算装置Ｃ３は、イテレーション部３１及びスピン系計算部３２から構成される。結合振動子系の計算装置Ｃ３に対して、結合振動子系の計算プログラムをインストールすることにより、イテレーション部３１及びスピン系計算部３２を実装することができる。
【００６３】
イテレーション部３１は、結合振動子系のダイナミクスを記述する微分方程式（数２）の各イテレーションにおいて、各振動子ｉ毎の微小時間後の振幅及び位相（数２のｃｉの絶対値及び符号）を計算する。よって、イテレーション部３１は、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算することができる。
【００６４】
スピン系計算部３２は、イテレーション部３１により反復されたイテレーションの後に、微分方程式（数２）の定常状態における各振動子ｉ毎の位相（数２のｃｉの符号）を計算することにより、結合振動子系の相互結合（数２のξｉｊ）に対応するイジングモデル（数１）の基底状態を計算する。よって、スピン系計算部３２は、イジングモデルの基底状態を計算することにより、組み合わせ最適化問題の最適解を求めることができる。
【００６５】
ここで、結合振動子系の全振動子数をＭとすると、数２の微分方程式のうち、右辺第１項の利得損失項及び右辺第３項の雑音項の計算量は、Ｍに比例するが、右辺第２項の相互結合項の計算量は、Ｍ２に比例する。そこで、数２の微分方程式のうち、右辺第２項の相互結合項の計算時間が、右辺第１項の利得損失項及び右辺第３項の雑音項の計算時間と比べて、十分に長いときには、図３に示したように各イテレーションを行なえばよい。
【００６６】
つまり、ある回のイテレーションの後でのｃｉを用いて計算した右辺第１項の利得損失項及び右辺第３項の雑音項と、その回のイテレーションよりも数回（例えば、１回のみ。）遅れた回のイテレーションの後でのｃｉを用いて計算した右辺第２項の相互結合項と、を加算することにより、その回のイテレーションでの左辺第１項の微小時間変化を計算するのである。よって、右辺第２項の相互結合項を余裕を持って計算しながら、右辺第１項の利得損失項及び右辺第３項の雑音項を高速に計算することができる。
【００６７】
結合振動子系の計算装置Ｃ３のうち、イテレーション部３１は、利得損失計算部３１１、雑音計算部３１２、相互結合計算部３１３及び振幅位相計算部３１４から構成される。
【００６８】
利得損失計算部３１１は、数２の微分方程式の各イテレーションにおいて、各振動子ｉ毎の利得及び損失に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、利得損失計算部３１１は、数２の微分方程式の各イテレーションにおいて、第ｎ回目のイテレーションの後でのｃｉ（ｎ）を入力され、ｃｉ＝Ａｔａｎｈ（ｃｉ／Ｂ）を出力する。ここで、第１回目のイテレーションの前でのｃｉ（０）は、０に初期設定される。そして、Ａ及びＢは、数２の右辺第１項のポンプレートｐに応じて設定される。
【００６９】
第１実施形態の利得損失計算部の処理内容を図４に示す。図４の横軸は、利得損失計算部３１１の入力ｃｉｎを示し、図４の縦軸は、利得損失計算部３１１の出力ｃｏｕｔを示す。２つの曲線のうちの一方の曲線は、数２の右辺第１項を計算したときにおける、利得損失計算部３１１の入力ｃｉｎと出力ｃｏｕｔとの間の関係を示す。２つの曲線のうちの他方の曲線は、数２の右辺第１項をｔａｎｈで近似したときにおける、利得損失計算部３１１の入力ｃｉｎと出力ｃｏｕｔとの間の関係を示す。２つの曲線のうち、上記の他方の曲線は、上記の一方の曲線に対して、正当な近似となっていることが分かる。
【００７０】
そこで、利得損失計算部３１１は、図４に示した入力ｃｉｎと出力ｃｏｕｔとの間の関係をルックアップテーブルとして記憶していれば、数２の右辺第１項を効率よく計算することができる。或いは、利得損失計算部３１１は、数２の右辺第１項をルンゲクッタ法又はオイラー丸山法で計算すれば、数２の右辺第１項を精度高く計算することができる。
【００７１】
雑音計算部３１２は、数２の微分方程式の各イテレーションにおいて、各振動子ｉ毎の雑音に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、雑音計算部３１２は、数２の微分方程式の各イテレーションにおいて、利得損失計算部３１１からの出力ｃｉを入力され、ｃｉ＝ｃｉ＋Ｎｏｉｓｅを出力する。ここで、Ｎｏｉｓｅは、数２の右辺第３項に応じて設定される。そして、Ｎｏｉｓｅは、数２の微分方程式の各イテレーションのうち、いずれの段階において付加されてもよい。
【００７２】
相互結合計算部３１３は、数２の微分方程式の各イテレーションにおいて、各振動子ｉの間の相互結合に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、相互結合計算部３１３は、数２の微分方程式の各イテレーションにおいて、雑音計算部３１２からの出力ｃｉを入力され、Σξｉｊｃｊを出力する。ここで、相互結合計算部３１３でのΣξｉｊｃｊの計算時間は、利得損失計算部３１１でのｃｉ＝Ａｔａｎｈ（ｃｉ／Ｂ）の計算時間及び雑音計算部３１２でのｃｉ＝ｃｉ＋Ｎｏｉｓｅの計算時間の合計時間と比べて、イテレーション数回分（例えば、１回分。）だけ長くなっている。
【００７３】
ただし、相互結合計算部３１３は、ξｉｊを表わす行列のうち０となる要素を省略し、ξｉｊを表わす行列を圧縮すれば、数２の右辺第２項を効率よく計算することができる。
【００７４】
振幅位相計算部３１４は、数２の微分方程式の各イテレーションにおいて、利得損失計算部３１１、雑音計算部３１２及び相互結合計算部３１３により各々計算された各振動子ｉ毎の振幅及び位相の微小時間変化に基づいて、各振動子ｉ毎の微小時間後の振幅及び位相を計算する。具体的には、振幅位相計算部３１４は、数２の微分方程式の各イテレーションにおいて、雑音計算部３１２からの出力ｃｉを入力され、相互結合計算部３１３からの出力Σξｉｊｃｊを入力され、雑音計算部３１２からの出力ｃｉ及び相互結合計算部３１３からの出力Σξｉｊｃｊを加算し、第（ｎ＋１）回目のイテレーションの後でのｃｉ（ｎ＋１）を出力する。ここで、雑音計算部３１２からの出力ｃｉは、第ｎ回目のイテレーションの後でのｃｉ（ｎ）を用いて計算されたものであるが、相互結合計算部３１３からの出力Σξｉｊｃｊは、第（ｎ−ｍ）回目のイテレーションの後でのｃｉ（ｎ−ｍ）を用いて計算されたものである。ただし、ｍは、イテレーションｍ回分（例えば、１回分。）を表わす。
【００７５】
スピン系計算部３２は、イテレーション部３１により反復されたイテレーションの後に、数２の微分方程式の定常状態における各振動子ｉ毎の位相を計算することにより、結合振動子系の相互結合に対応するイジングモデルの基底状態を計算する。
【００７６】
このように、組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。そして、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。
【００７７】
（第１実施形態：計算結果）
次に、第１実施形態の結合振動子系の計算装置の計算結果を図５に示す。図５の左欄は、完全グラフＫ４（全振動子数Ｍ＝４）の計算結果を示し、図５の右欄は、ランダムグラフＧ１（全振動子数Ｍ＝８００）の計算結果を示す。図５の横軸は、イテレーション回数を示し、図５の縦軸は、振幅位相計算部３１４の出力ｃｉ（ｎ）を示す。
【００７８】
図５の左欄の完全グラフＫ４では、第１回目のイテレーションの前でのｃｉ（０）を０に初期設定し、数２の右辺第１項のポンプレートをｐ＝２．０に設定し、利得損失計算部３１１からの出力をｃｏｕｔ＝２．２ｔａｎｈ（ｃｉｎ／１．８）に設定した。すると、イテレーション回数〜３０回目において、ｃｉ（ｎ）が０から立ち上がり、イテレーション回数〜５０回目において、ｃｉ（ｎ）が定常状態に落ち着き、完全グラフＫ４の最大のカット数が得られ、完全グラフＫ４に対応するイジングモデルの基底状態が得られた。
【００７９】
図５の右欄のランダムグラフＧ１では、第１回目のイテレーションの前でのｃｉ（０）を０に初期設定し、数２の右辺第１項のポンプレートをｐ＝１．１に設定し、利得損失計算部３１１からの出力をｃｏｕｔ＝２．３５ｔａｎｈ（ｃｉｎ／２．４６）に設定した。すると、イテレーション回数〜１０回目において、ｃｉ（ｎ）が０から立ち上がり、イテレーション回数〜１００回目において、ｃｉ（ｎ）が定常状態に落ち着き、ランダムグラフＧ１の最大のカット数が得られ、ランダムグラフＧ１に対応するイジングモデルの基底状態が得られた。
【００８０】
ここで、雑音計算部３１２は、数２の微分方程式の初期状態のみにおいて、各振動子ｉ毎の雑音に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算してもよい。具体的には、雑音計算部３１２は、図５の左欄の完全グラフＫ４では、イテレーション回数〜３０回目において、計算を停止してもよく、図５の右欄のランダムグラフＧ１では、イテレーション回数〜１０回目において、計算を停止してもよい。
【００８１】
よって、数２の微分方程式の定常状態として、複数の結合振動子の振幅が０である状態を回避することができる。そして、複数の結合振動子の振幅が０である初期状態から脱出することができた後には、雑音計算部３１２の計算負担をなくすことができる。
【００８２】
さらに、イテレーション部３１は、前回のイテレーションにおける各振動子ｉ毎の微小時間後の振幅及び位相が、数２の微分方程式の定常状態に所定の精度で接近したときに、次回のイテレーションにおける各振動子ｉ毎の微小時間後の振幅及び位相の計算を停止してもよい。具体的には、イテレーション部３１は、図５の左欄の完全グラフＫ４では、イテレーション回数〜５０回目において、計算を停止してもよく、図５の右欄のランダムグラフＧ１では、イテレーション回数〜１００回目において、計算を停止してもよい。
【００８３】
より一般的には、イテレーション部３１は、所定の割合（例えば、全振動子の９０％。）の振動子について、ｃｉ（ｎ）が所定の閾値（例えば、飽和振幅の９０％。）を超えたときに、計算を停止してもよい。よって、複数の結合振動子の振幅が有限値である定常状態に接近することができた後には、イテレーション部３１の計算負担をなくすことができる。
【００８４】
（第１実施形態：並列計算）
次に、第１実施形態のイテレーション部の第１構成を図６に示す。イテレーション部３１は、複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍ及びＦＩＦＯ３３から構成される。複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍ及びＦＩＦＯ３３は、１プロセッサ又は１ＦＰＧＡに収容される。
【００８５】
複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍは、数２の微分方程式の各イテレーションにおいて、各振動子１、・・・、ｍ、・・・、Ｍ毎の並列処理により、各振動子１、・・・、ｍ、・・・、Ｍ毎の微小時間後の振幅及び位相を計算する。
【００８６】
図６の以下の説明では、イテレーション部３１−ｍについて説明する。イテレーション部３１−ｍは、利得損失計算部３１１−ｍ、雑音計算部３１２−ｍ、相互結合計算部３１３−ｍ及び振幅位相計算部３１４−ｍから構成される。
【００８７】
利得損失計算部３１１−ｍは、ＦＩＦＯ３３（詳細については後述する。）に格納されたｃ１（ｎ）、・・・、ｃｍ（ｎ）、・・・、ｃＭ（ｎ）をこの順序で入力され、ｃ１＝Ａｔａｎｈ（ｃ１／Ｂ）、・・・、ｃｍ＝Ａｔａｎｈ（ｃｍ／Ｂ）、・・・、ｃＭ＝Ａｔａｎｈ（ｃＭ／Ｂ）をこの順序で出力する。ここで、利得損失計算部３１１−ｍは、複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍの間で共有されていてもよい。
【００８８】
雑音計算部３１２−ｍは、利得損失計算部３１１−ｍからの出力ｃ１、・・・、ｃｍ、・・・、ｃＭをこの順序で入力され、ｃ１＝ｃ１＋Ｎｏｉｓｅ、・・・、ｃｍ＝ｃｍ＋Ｎｏｉｓｅ、・・・、ｃＭ＝ｃＭ＋Ｎｏｉｓｅをこの順序で出力する。ここで、雑音計算部３１２−ｍは、複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍの間で共有されていてもよい。
【００８９】
相互結合計算部３１３−ｍは、雑音計算部３１２−ｍからの出力ｃ１、・・・、ｃｍ、・・・、ｃＭをこの順序で入力され、Σξｍｊｃｊを出力する。ここで、相互結合計算部３１３−ｍは、ＦＩＦＯ３１５−ｍ、乗算器３１６−ｍ及び加算器３１７−ｍから構成される。
【００９０】
乗算器３１６−ｍは、雑音計算部３１２−ｍからの出力ｃ１、・・・、ｃｍ、・・・、ｃＭをこの順序で入力され、ＦＩＦＯ３１５−ｍに格納されたξｍ１、・・・、ξｍｍ（＝０）、・・・、ξｍＭをこの順序で入力され、ξｍ１ｃ１、・・・、ξｍｍｃｍ（＝０）、・・・、ξｍＭｃＭをこの順序で出力する。ここで、乗算器３１６−ｍでのξｍ１ｃ１、・・・、ξｍｍｃｍ（＝０）、・・・、ξｍＭｃＭの計算時間は、利得損失計算部３１１−ｍでのｃｍ＝Ａｔａｎｈ（ｃｍ／Ｂ）の計算時間及び雑音計算部３１２−ｍでのｃｍ＝ｃｍ＋Ｎｏｉｓｅの計算時間の合計時間と比べて、イテレーション数回分（例えば、１回分。）だけ長くなっている。
【００９１】
加算器３１７−ｍは、乗算器３１６−ｍからの出力ξｍ１ｃ１、・・・、ξｍｍｃｍ（＝０）、・・・、ξｍＭｃＭをこの順序で入力され、Σξｍｊｃｊを出力する。
【００９２】
振幅位相計算部３１４−ｍは、雑音計算部３１２−ｍからの出力ｃｍを入力され、加算器３１７−ｍからの出力Σξｍｊｃｊを入力され、雑音計算部３１２−ｍからの出力ｃｍ及び加算器３１７−ｍからの出力Σξｍｊｃｊを加算し、ｃｍ（ｎ＋１）を出力する。ここで、雑音計算部３１２−ｍからの出力ｃｍは、第ｎ回目のイテレーションの後でのｃｍ（ｎ）を用いて計算されたものであるが、加算器３１７−ｍからの出力Σξｍｊｃｊは、第（ｎ−ｍ）回目のイテレーションの後でのｃｍ（ｎ−ｍ）を用いて計算されたものである。ただし、ｍは、イテレーションｍ回分（例えば、１回分。）を表わす。
【００９３】
最初に、イテレーション部３１−１は、ｃ１（ｎ＋１）をＦＩＦＯ３３に格納する。・・・次に、イテレーション部３１−ｍは、ｃｍ（ｎ＋１）をＦＩＦＯ３３に格納する。・・・最後に、イテレーション部３１−Ｍは、ｃＭ（ｎ＋１）をＦＩＦＯ３３に格納する。
【００９４】
ここで、結合振動子系の全振動子数をＭとすると、各振動子毎の並列処理を行わないときには、２体相互結合の計算量は、Ｍ２のオーダーに増える。しかし、各振動子毎の並列処理を行なうときには、２体相互結合の計算量は、ｌｏｇ２Ｍのオーダーに減らすことができる。
【００９５】
第１実施形態のイテレーション部の第２構成を図７に示す。イテレーション部３１は、複数のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩから構成される。各々のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩは、１プロセッサ又は１ＦＰＧＡに収容される。
【００９６】
図６に示した第１構成に加えて、各々のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩは、数２の微分方程式の各イテレーションにおいて、結合振動子系のうちの各々に特定のグループに属する振動子の微小時間後の振幅及び位相のデータを、その他のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩとの間で共有する。
【００９７】
イテレーション部３１−Ｉは、ＦＩＦＯ３３−Ｉ、イテレーション部３４−Ｉ及びＦＩＦＯ３５−Ｉから構成される。イテレーション部３４−Ｉは、ｃ１、・・・、ｃ３Ｍのうち、ｃ１、ｃ２、・・・、ｃＭについて、図６に示した複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍと同様の処理を実行する。
【００９８】
イテレーション部３４−Ｉは、ＦＩＦＯ３３−Ｉ（詳細については後述する。）に格納されたｃ１（ｎ）、・・・、ｃＭ（ｎ）、ｃＭ＋１（ｎ）、・・・、ｃ２Ｍ（ｎ）、ｃ２Ｍ＋１（ｎ）、・・・、ｃ３Ｍ（ｎ）をこの順序で入力され、ｃ１（ｎ＋１）、ｃ２（ｎ＋１）、・・・、ｃＭ（ｎ＋１）をこの順序で出力し、ＦＩＦＯ３３−Ｉに格納する前にＦＩＦＯ３５−Ｉに格納する。
【００９９】
イテレーション部３１−ＩＩは、ＦＩＦＯ３３−ＩＩ、イテレーション部３４−ＩＩ及びＦＩＦＯ３５−ＩＩから構成される。イテレーション部３４−ＩＩは、ｃ１、・・・、ｃ３Ｍのうち、ｃＭ＋１、ｃＭ＋２、・・・、ｃ２Ｍについて、図６に示した複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍと同様の処理を実行する。
【０１００】
イテレーション部３４−ＩＩは、ＦＩＦＯ３３−ＩＩ（詳細については後述する。）に格納されたｃ１（ｎ）、・・・、ｃＭ（ｎ）、ｃＭ＋１（ｎ）、・・・、ｃ２Ｍ（ｎ）、ｃ２Ｍ＋１（ｎ）、・・・、ｃ３Ｍ（ｎ）をこの順序で入力され、ｃＭ＋１（ｎ＋１）、ｃＭ＋２（ｎ＋１）、・・・、ｃ２Ｍ（ｎ＋１）をこの順序で出力し、ＦＩＦＯ３３−ＩＩに格納する前にＦＩＦＯ３５−ＩＩに格納する。
【０１０１】
イテレーション部３１−ＩＩＩは、ＦＩＦＯ３３−ＩＩＩ、イテレーション部３４−ＩＩＩ及びＦＩＦＯ３５−ＩＩＩから構成される。イテレーション部３４−ＩＩＩは、ｃ１、・・・、ｃ３Ｍのうち、ｃ２Ｍ＋１、ｃ２Ｍ＋２、・・・、ｃ３Ｍについて、図６に示した複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍと同様の処理を実行する。
【０１０２】
イテレーション部３４−ＩＩＩは、ＦＩＦＯ３３−ＩＩＩ（詳細については後述する。）に格納されたｃ１（ｎ）、・・・、ｃＭ（ｎ）、ｃＭ＋１（ｎ）、・・・、ｃ２Ｍ（ｎ）、ｃ２Ｍ＋１（ｎ）、・・・、ｃ３Ｍ（ｎ）をこの順序で入力され、ｃ２Ｍ＋１（ｎ＋１）、ｃ２Ｍ＋２（ｎ＋１）、・・・、ｃ３Ｍ（ｎ＋１）をこの順序で出力し、ＦＩＦＯ３３−ＩＩＩに格納する前にＦＩＦＯ３５−ＩＩＩに格納する。
【０１０３】
最初に、ｃ１（ｎ＋１）、ｃ２（ｎ＋１）、・・・、ｃＭ（ｎ＋１）が、ＦＩＦＯ３５−ＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。次に、ｃＭ＋１（ｎ＋１）、ｃＭ＋２（ｎ＋１）、・・・、ｃ２Ｍ（ｎ＋１）が、ＦＩＦＯ３５−ＩＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。最後に、ｃ２Ｍ＋１（ｎ＋１）、ｃ２Ｍ＋２（ｎ＋１）、・・・、ｃ３Ｍ（ｎ＋１）が、ＦＩＦＯ３５−ＩＩＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。ここで、ＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩは、単一のＦＩＦＯであってもよい。
【０１０４】
このように、結合振動子系の全振動子数が膨大な数になるときでも、各々のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩを各々のプロセッサに実装することができる。
【０１０５】
第１実施形態のイテレーション部の第３構成を図８に示す。イテレーション部３１は、複数のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩから構成される。各々のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩは、１プロセッサ又は１ＦＰＧＡに収容される。
【０１０６】
図７に示した第２構成に加えて、前回のイテレーションにおける上記のデータ（ｃ１（ｎ＋１）、・・・ｃ３Ｍ（ｎ＋１））の共有処理が、次回のイテレーションにおける計算処理を停止／中断させないように、各々のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩは、結合振動子系のうちの各々に特定のグループに属する振動子を割り当てられる。
【０１０７】
イテレーション部３１−Ｉは、ＦＩＦＯ３３−Ｉ、イテレーション部３４−Ｉ及びＦＩＦＯ３５−Ｉから構成される。イテレーション部３４−Ｉは、ｃ１、・・・、ｃ３Ｍのうち、ｃ１、ｃ４、・・・、ｃ３Ｍ−２について、図６に示した複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍと同様の処理を実行する。
【０１０８】
イテレーション部３４−Ｉは、ＦＩＦＯ３３−Ｉ（詳細については後述する。）に格納されたｃ１（ｎ）、・・・、ｃＭ（ｎ）、ｃＭ＋１（ｎ）、・・・、ｃ２Ｍ（ｎ）、ｃ２Ｍ＋１（ｎ）、・・・、ｃ３Ｍ（ｎ）をこの順序で入力され、ｃ１（ｎ＋１）、ｃ４（ｎ＋１）、・・・、ｃ３Ｍ−２（ｎ＋１）をこの順序で出力し、ＦＩＦＯ３３−Ｉに格納する前にＦＩＦＯ３５−Ｉに格納する。
【０１０９】
イテレーション部３１−ＩＩは、ＦＩＦＯ３３−ＩＩ、イテレーション部３４−ＩＩ及びＦＩＦＯ３５−ＩＩから構成される。イテレーション部３４−ＩＩは、ｃ１、・・・、ｃ３Ｍのうち、ｃ２、ｃ５、・・・、ｃ３Ｍ−１について、図６に示した複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍと同様の処理を実行する。
【０１１０】
イテレーション部３４−ＩＩは、ＦＩＦＯ３３−ＩＩ（詳細については後述する。）に格納されたｃ１（ｎ）、・・・、ｃＭ（ｎ）、ｃＭ＋１（ｎ）、・・・、ｃ２Ｍ（ｎ）、ｃ２Ｍ＋１（ｎ）、・・・、ｃ３Ｍ（ｎ）をこの順序で入力され、ｃ２（ｎ＋１）、ｃ５（ｎ＋１）、・・・、ｃ３Ｍ−１（ｎ＋１）をこの順序で出力し、ＦＩＦＯ３３−ＩＩに格納する前にＦＩＦＯ３５−ＩＩに格納する。
【０１１１】
イテレーション部３１−ＩＩＩは、ＦＩＦＯ３３−ＩＩＩ、イテレーション部３４−ＩＩＩ及びＦＩＦＯ３５−ＩＩＩから構成される。イテレーション部３４−ＩＩＩは、ｃ１、・・・、ｃ３Ｍのうち、ｃ３、ｃ６、・・・、ｃ３Ｍについて、図６に示した複数のイテレーション部３１−１、・・・、３１−ｍ、・・・、３１−Ｍと同様の処理を実行する。
【０１１２】
イテレーション部３４−ＩＩＩは、ＦＩＦＯ３３−ＩＩＩ（詳細については後述する。）に格納されたｃ１（ｎ）、・・・、ｃＭ（ｎ）、ｃＭ＋１（ｎ）、・・・、ｃ２Ｍ（ｎ）、ｃ２Ｍ＋１（ｎ）、・・・、ｃ３Ｍ（ｎ）をこの順序で入力され、ｃ３（ｎ＋１）、ｃ６（ｎ＋１）、・・・、ｃ３Ｍ（ｎ＋１）をこの順序で出力し、ＦＩＦＯ３３−ＩＩＩに格納する前にＦＩＦＯ３５−ＩＩＩに格納する。
【０１１３】
最初に、ｃ１（ｎ＋１）が、ＦＩＦＯ３５−ＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。次に、ｃ２（ｎ＋１）が、ＦＩＦＯ３５−ＩＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。次に、ｃ３（ｎ＋１）が、ＦＩＦＯ３５−ＩＩＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。・・・次に、ｃ３Ｍ−２（ｎ＋１）が、ＦＩＦＯ３５−ＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。次に、ｃ３Ｍ−１（ｎ＋１）が、ＦＩＦＯ３５−ＩＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。最後に、ｃ３Ｍ（ｎ＋１）が、ＦＩＦＯ３５−ＩＩＩからＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩへと転送される。ここで、ＦＩＦＯ３３−Ｉ、３３−ＩＩ、３３−ＩＩＩは、単一のＦＩＦＯであってもよい。
【０１１４】
このように、各々のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩを各々のプロセッサに実装するときでも、各々のイテレーション部３１−Ｉ、３１−ＩＩ、３１−ＩＩＩの間の上記のデータ（ｃ１（ｎ＋１）、・・・ｃ３Ｍ（ｎ＋１））の共有処理を円滑に実行することができる。
【０１１５】
（第１実施形態：ビット量子化）
次に、第１実施形態の振幅位相計算部の量子化を図９に示す。図９の左欄は、振幅位相計算部３１４の量子化の第１形態を示し、図９の右欄は、振幅位相計算部３１４の量子化の第２形態を示す。図９の横軸は、振幅位相計算部３１４の量子化前の入力を示し、図９の縦軸は、振幅位相計算部３１４の量子化後の出力を示す。
【０１１６】
振幅位相計算部３１４は、数２の微分方程式の各イテレーションにおいて、所定の精度で量子化したうえで、各振動子ｉ毎の微小時間後の振幅及び位相を計算する。具体的には、振幅位相計算部３１４は、数２の微分方程式の各イテレーションにおいて、雑音計算部３１２からの出力ｃｉ及び相互結合計算部３１３からの出力Σξｉｊｃｊの加算結果を所定の精度で量子化し、第（ｎ＋１）回目のイテレーションの後でのｃｉ（ｎ＋１）を出力する。
【０１１７】
図９の左欄に示した第１形態では、振幅位相計算部３１４の出力ビット数ｘは、４である。そして、振幅位相計算部３１４の出力は、０を含まず、振幅位相計算部３１４の出力の階調数は、２ｘ＝２４＝１６である。よって、数２の微分方程式の定常状態として、複数の結合振動子の振幅が０である状態を回避することができる。
【０１１８】
図９の右欄に示した第２形態でも、振幅位相計算部３１４の出力ビット数ｘは、４である。しかし、振幅位相計算部３１４の出力は、０を含み、振幅位相計算部３１４の出力の階調数は、２ｘ−１＝２４−１＝１５である。よって、数２の微分方程式の初期状態において、結合振動子の位相が正答の位相からずれた位相になる状態を回避することができる（図５の右欄に示したｃｉ（ｎ）の正／負から負／正への変化を参照。）。
【０１１９】
このように、結合振動子系の計算装置Ｃ３の計算負担を減らすことができるとともに、結合振動子系の計算装置Ｃ３のメモリ（ＦＩＦＯ）容量を減らすことができる。
【０１２０】
振幅位相計算部３１４により計算される各振動子ｉ毎の微小時間後の振幅についての量子化カットオフは、振幅位相計算部３１４により計算される各振動子ｉ毎の微小時間後の振幅についての量子化の上記の所定の精度を満足するように設定されている。具体的には、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化カットオフは、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化の上記の所定の精度を満足する。
【０１２１】
図９の左欄及び右欄に示した第１、２形態では、振幅位相計算部３１４の量子化カットオフ±Ｎは、±２である。よって、振幅位相計算部３１４の量子化精度は、Ｎ×２／２ｘ＝２×２／２４＝０．２５程度である。もっとも、図１０〜１２に示すように、振幅位相計算部３１４の量子化精度は、０．２５程度では正答を得るには低過ぎる。
【０１２２】
振幅位相計算部３１４により計算される各振動子ｉ毎の微小時間後の振幅についての量子化カットオフは、数２の微分方程式の定常状態における各振動子ｉ毎の振幅と比べて大きく設定されている。具体的には、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化カットオフは、振幅位相計算部３１４からの定常出力ｃｉ（ｎ＋１）と比べて大きい。
【０１２３】
図９の左欄及び右欄に示した第１、２形態では、振幅位相計算部３１４の量子化カットオフ±Ｎは、±２である。そして、図５に示したように、振幅位相計算部３１４からの定常出力ｃｉ（ｎ＋１）は、±１〜２程度である。よって、振幅位相計算部３１４の量子化カットオフＮは、振幅位相計算部３１４からの定常出力ｃｉ（ｎ＋１）と比べて大きい。
【０１２４】
このように、各振動子毎の振幅の分解能を低くし過ぎず、かつ、各振動子毎の振幅の制限を厳しくし過ぎず、結合振動子系の計算装置Ｃ３の計算負担を減らすことができるとともに、結合振動子系の計算装置Ｃ３のメモリ（ＦＩＦＯ）容量を減らすことができる。
【０１２５】
相互結合計算部３１３により計算される各振動子ｉ毎の振幅の微小時間変化についての量子化ビット数は、振幅位相計算部３１４により計算される各振動子ｉ毎の微小時間後の振幅についての量子化ビット数と比べて、少なくともｌｏｇ２Ｍだけ大きく設定されている。ただし、Ｍは、結合振動子系の全振動子数である。具体的には、相互結合計算部３１３からの出力Σξｉｊｃｊについての量子化ビット数は、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化ビット数と比べて、少なくともｌｏｇ２Ｍだけ大きい。
【０１２６】
仮に、相互結合計算部３１３からの出力Σξｉｊｃｊについての量子化ビット数を、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化ビット数と比べて、同等にする。すると、相互結合計算部３１３からの出力Σξｉｊｃｊは、Ｍ個のｃｊを重み付け加算したものであるため、相互結合計算部３１３からの出力Σξｉｊｃｊについての量子化精度は、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化精度と比べて、低くなる。
【０１２７】
そこで、相互結合計算部３１３からの出力Σξｉｊｃｊについての量子化ビット数を、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化ビット数と比べて、少なくともｌｏｇ２Ｍだけ大きくする。すると、相互結合計算部３１３からの出力Σξｉｊｃｊは、Ｍ個のｃｊを重み付け加算したものであるところ、相互結合計算部３１３からの出力Σξｉｊｃｊについての量子化精度は、振幅位相計算部３１４からの出力ｃｉ（ｎ＋１）についての量子化精度と比べて、少なくとも同等になる。
【０１２８】
このように、各振動子毎の振幅の分解能と比べて、各振動子毎の相互結合の大きさの分解能を低くし過ぎず、結合振動子系の計算装置Ｃ３の計算負担を減らすことができるとともに、結合振動子系の計算装置Ｃ３のメモリ（ＦＩＦＯ）容量を減らすことができる。
【０１２９】
第１実施形態の結合振動子系の計算装置における、完全グラフＫ４（全振動子数Ｍ＝４）の計算結果を図１０に示す。１００００回の試行のうち、正答率が高かったときほど、塗りつぶしを薄くし、正答率が低かったときほど、塗りつぶしを濃くする。
【０１３０】
図１０の左欄の横軸は、振幅位相計算部３１４の出力ビット数ｘを示し、図１０の左欄の縦軸は、利得損失計算部３１１の出力ビット数ｙを示す。ただし、相互結合計算部３１３の出力ビット数ｚは、３２に固定されており、振幅位相計算部３１４の量子化カットオフＮは、２．０に固定されている。ｘ及びｙがともに小さいときや、ｘが大きくてもｙが小さいときや、ｙが大きくてもｘが小さいときは、正答率が低くなった。ｘ及びｙがほぼ等しくともに大きい（例えば、１５程度。）ときは、正答率が高くなった。
【０１３１】
図１０の中欄の横軸は、振幅位相計算部３１４の出力ビット数ｘ＝利得損失計算部３１１の出力ビット数ｙを示し、図１０の中欄の縦軸は、相互結合計算部３１３の出力ビット数ｚを示す。ただし、振幅位相計算部３１４の量子化カットオフＮは、２．０に固定されている。ｘ＝ｙが大きくなるにつれて、ｚがより大きくならなければ、相互結合が正確に計算されず正答率が低くなった。ｘ＝ｙが大きい（例えば、１５程度。）うえで、ｚがｘ＝ｙにほぼ等しいときは、相互結合がやや正確に計算され正答率がやや高くなった。ｘ＝ｙが大きい（例えば、１５程度。）うえで、ｚがｘ＝ｙよりｌｏｇ２Ｍ＝ｌｏｇ２４＝２だけ大きいときは、相互結合がより正確に計算され正答率がより高くなった。
【０１３２】
図１０の右欄の横軸は、振幅位相計算部３１４の出力ビット数ｘを示し、図１０の右欄の縦軸は、振幅位相計算部３１４の量子化カットオフｌｏｇ１０Ｎを示す。ただし、利得損失計算部３１１の出力ビット数ｙは、振幅位相計算部３１４の出力ビット数ｘに等しく、相互結合計算部３１３の出力ビット数ｚは、３２に固定されている。ｘが小さくなるにつれて、Ｎがより小さくならなければ、量子化精度が低くなり正答率が低くなった。ｘが図１０の右欄に示したいずれの値であっても、Ｎが小さすぎるときは、振動子振幅が制限され正答率が低くなった。ｘが大きい（例えば、１５程度。）うえで、Ｎが１程度であるときは、量子化精度が高くなり振動子振幅が制限されず、正答率がより高くなった。
【０１３３】
第１実施形態の結合振動子系の計算装置における、メビウスラダーＭ１６（全振動子数Ｍ＝１６）の計算結果を図１１に示す。１００００回の試行のうち、正答率が高かったときほど、塗りつぶしを薄くし、正答率が低かったときほど、塗りつぶしを濃くする。
【０１３４】
図１１の左欄の横軸は、振幅位相計算部３１４の出力ビット数ｘを示し、図１１の左欄の縦軸は、利得損失計算部３１１の出力ビット数ｙを示す。ただし、相互結合計算部３１３の出力ビット数ｚは、３２に固定されており、振幅位相計算部３１４の量子化カットオフＮは、２．０に固定されている。ｘ及びｙがともに小さいときや、ｘが大きくてもｙが小さいときや、ｙが大きくてもｘが小さいときは、正答率が低くなった。ｘ及びｙがほぼ等しくともに大きい（例えば、１５程度。）ときは、正答率が高くなった。
【０１３５】
図１１の中欄の横軸は、振幅位相計算部３１４の出力ビット数ｘ＝利得損失計算部３１１の出力ビット数ｙを示し、図１１の中欄の縦軸は、相互結合計算部３１３の出力ビット数ｚを示す。ただし、振幅位相計算部３１４の量子化カットオフＮは、２．０に固定されている。ｘ＝ｙが大きくなるにつれて、ｚがより大きくならなければ、相互結合が正確に計算されず正答率が低くなった。ｘ＝ｙが大きい（例えば、１５程度。）うえで、ｚがｘ＝ｙにほぼ等しいときは、相互結合がやや正確に計算され正答率がやや高くなった。ｘ＝ｙが大きい（例えば、１５程度。）うえで、ｚがｘ＝ｙよりｌｏｇ２Ｍ＝ｌｏｇ２１６＝４だけ大きいときは、相互結合がより正確に計算され正答率がより高くなった。
【０１３６】
図１１の右欄の横軸は、振幅位相計算部３１４の出力ビット数ｘを示し、図１１の右欄の縦軸は、振幅位相計算部３１４の量子化カットオフｌｏｇ１０Ｎを示す。ただし、利得損失計算部３１１の出力ビット数ｙは、振幅位相計算部３１４の出力ビット数ｘに等しく、相互結合計算部３１３の出力ビット数ｚは、３２に固定されている。ｘが小さくなるにつれて、Ｎがより小さくならなければ、量子化精度が低くなり正答率が低くなった。ｘが図１１の右欄に示したいずれの値であっても、Ｎが小さすぎるときは、振動子振幅が制限され正答率が低くなった。ｘが大きい（例えば、１５程度。）うえで、Ｎが１程度であるときは、量子化精度が高くなり振動子振幅が制限されず、正答率がより高くなった。
【０１３７】
第１実施形態の結合振動子系の計算装置における、ランダムグラフＧ１（全振動子数Ｍ＝８００）の計算結果を図１２に示す。１０００回の試行のうち、カット数が高かったときほど、塗りつぶしを薄くし、カット数が低かったときほど、塗りつぶしを濃くする。
【０１３８】
図１２の左欄の横軸は、振幅位相計算部３１４の出力ビット数ｘを示し、図１２の左欄の縦軸は、利得損失計算部３１１の出力ビット数ｙを示す。ただし、相互結合計算部３１３の出力ビット数ｚは、１６に固定されており、振幅位相計算部３１４の量子化カットオフＮは、１．０に固定されている。ｘ及びｙがともに小さいときや、ｘが大きくてもｙが小さいときや、ｙが大きくてもｘが小さいときは、カット数が低くなった。ｘ及びｙがほぼ等しくともに大きい（例えば、１５程度。）ときは、カット数が高くなった。
【０１３９】
図１２の中欄の横軸は、振幅位相計算部３１４の出力ビット数ｘ＝利得損失計算部３１１の出力ビット数ｙを示し、図１２の中欄の縦軸は、相互結合計算部３１３の出力ビット数ｚを示す。ただし、振幅位相計算部３１４の量子化カットオフＮは、１．０に固定されている。ｘ＝ｙが大きくなるにつれて、ｚがより大きくならなければ、相互結合が正確に計算されずカット数が低くなった。ｘ＝ｙが大きい（例えば、１５程度。）うえで、ｚがｘ＝ｙにほぼ等しいときは、相互結合がやや正確に計算されカット数がやや高くなった。ｘ＝ｙが大きい（例えば、１５程度。）うえで、ｚがｘ＝ｙよりｌｏｇ２Ｍ＝ｌｏｇ２８００〜１０だけ大きいときは、相互結合がより正確に計算されカット数がより高くなった。
【０１４０】
図１２の右欄の横軸は、振幅位相計算部３１４の出力ビット数ｘを示し、図１２の右欄の縦軸は、振幅位相計算部３１４の量子化カットオフｌｏｇ１０Ｎを示す。ただし、利得損失計算部３１１の出力ビット数ｙは、振幅位相計算部３１４の出力ビット数ｘに等しく、相互結合計算部３１３の出力ビット数ｚは、３２に固定されている。ｘが小さくなるにつれて、Ｎがより小さくならなければ、量子化精度が低くなりカット数が低くなった。ｘが図１２の右欄に示したいずれの値であっても、Ｎが小さすぎるときは、振動子振幅が制限されカット数が低くなった。ｘが大きい（例えば、１５程度。）うえで、Ｎが１程度であるときは、量子化精度が高くなり振動子振幅が制限されず、カット数がより高くなった。
【０１４１】
（第２実施形態：計算方法）
第２実施形態の結合振動子系の計算装置の構成を図１３に示す。結合振動子系の計算装置Ｃ４は、イテレーション部４１及びスピン系計算部４２から構成される。結合振動子系の計算装置Ｃ４に対して、結合振動子系の計算プログラムをインストールすることにより、イテレーション部４１及びスピン系計算部４２を実装することができる。
【０１４２】
イテレーション部４１は、結合振動子系のダイナミクスを記述する微分方程式（数２）の各イテレーションにおいて、各振動子ｉ毎の微小時間後の振幅及び位相（数２のｃｉの絶対値及び符号）を計算する。よって、イテレーション部４１は、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算することができる。
【０１４３】
スピン系計算部４２は、イテレーション部４１により反復されたイテレーションの後に、微分方程式（数２）の定常状態における各振動子ｉ毎の位相（数２のｃｉの符号）を計算することにより、結合振動子系の相互結合（数２のξｉｊ）に対応するイジングモデル（数１）の基底状態を計算する。よって、スピン系計算部４２は、イジングモデルの基底状態を計算することにより、組み合わせ最適化問題の最適解を求めることができる。
【０１４４】
ここで、結合振動子系の全振動子数をＭとすると、数２の微分方程式のうち、右辺第１項の利得損失項及び右辺第３項の雑音項の計算量は、Ｍに比例するが、右辺第２項の相互結合項の計算量は、Ｍ２に比例する。しかし、数２の微分方程式のうち、右辺第２項の相互結合項の計算時間が、右辺第１項の利得損失項及び右辺第３項の雑音項の計算時間と比べて、同程度であるときには、図１３に示したように各イテレーションを行なえばよい。
【０１４５】
つまり、ある回のイテレーションの後でのｃｉを用いて計算した、右辺第１項の利得損失項、右辺第２項の相互結合項及び右辺第３項の雑音項を加算することにより、その回のイテレーションでの左辺第１項の微小時間変化を計算するのである。よって、右辺第２項の相互結合項を図３のような遅延（例えば、イテレーション数回分。）を伴わず計算しながら、右辺第１項の利得損失項及び右辺第３項の雑音項を高速に計算することができる。
【０１４６】
結合振動子系の計算装置Ｃ４のうち、イテレーション部４１は、利得損失計算部４１１、相互結合計算部４１２、雑音計算部４１３及び振幅位相計算部４１４から構成される。
【０１４７】
利得損失計算部４１１は、数２の微分方程式の各イテレーションにおいて、各振動子ｉ毎の利得及び損失に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、利得損失計算部４１１は、数２の微分方程式の各イテレーションにおいて、第ｎ回目のイテレーションの後でのｃｉ（ｎ）を入力され、数２の右辺第１項を出力する。ここで、第１回目のイテレーションの前でのｃｉ（０）は、０に初期設定される。
【０１４８】
雑音計算部４１３は、数２の微分方程式の各イテレーションにおいて、各振動子ｉ毎の雑音に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、雑音計算部４１３は、数２の微分方程式の各イテレーションにおいて、第ｎ回目のイテレーションの後でのｃｉ（ｎ）を入力され、数２の右辺第３項を出力する。
【０１４９】
相互結合計算部４１２は、数２の微分方程式の各イテレーションにおいて、各振動子ｉの間の相互結合に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、相互結合計算部４１２は、数２の微分方程式の各イテレーションにおいて、第ｎ回目のイテレーションの後でのｃｉ（ｎ）を入力され、数２の右辺第２項を出力する。ここで、相互結合計算部４１２での数２の右辺第２項の計算時間は、利得損失計算部４１１及び雑音計算部４１３での数２の右辺第１、３項の計算時間と比べて、同程度である。
【０１５０】
振幅位相計算部４１４は、数２の微分方程式の各イテレーションにおいて、利得損失計算部４１１、雑音計算部４１３及び相互結合計算部４１２により各々計算された各振動子ｉ毎の振幅及び位相の微小時間変化に基づいて、各振動子ｉ毎の微小時間後の振幅及び位相を計算する。具体的には、振幅位相計算部４１４は、数２の微分方程式の各イテレーションにおいて、利得損失計算部４１１、雑音計算部４１３及び相互結合計算部４１２からの数２の右辺第１、３、２項を入力され、利得損失計算部４１１、雑音計算部４１３及び相互結合計算部４１２からの数２の右辺第１、３、２項を加算し、第（ｎ＋１）回目のイテレーションの後でのｃｉ（ｎ＋１）を出力する。ここで、利得損失計算部４１１、雑音計算部４１３及び相互結合計算部４１２からの数２の右辺第１、３、２項は、第ｎ回目のイテレーションの後でのｃｉ（ｎ）を用いて計算されたものである。
【０１５１】
スピン系計算部４２は、イテレーション部４１により反復されたイテレーションの後に、数２の微分方程式の定常状態における各振動子ｉ毎の位相を計算することにより、結合振動子系の相互結合に対応するイジングモデルの基底状態を計算する。
【０１５２】
このように、組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。そして、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。
【０１５３】
第２実施形態においても、第１実施形態の図６〜８と同様に、並列計算を行なうことができ、第１実施形態の図９〜１２と同様に、ビット量子化を行なうことができる。
【０１５４】
（第３実施形態：計算原理）
第３実施形態では、機械学習のボルツマンサンプリングを行なうために、機械学習をＸＹモデルへマッピングし、ＸＹモデルの平衡状態を計算する。または、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算する。
【０１５５】
そこで、複数の振動子で複数のＸＹ面内スピンを実装する。ここで、複数の振動子の間の相互結合を行なうことにより、複数のＸＹ面内スピンの間の相互結合を実装する。そして、複数の振動子全体の様々な発振モードのサンプリングを行なうことにより、複数のＸＹ面内スピン全体の様々なエネルギーのサンプリングを行なう。
【０１５６】
まず、複数の振動子全体の様々な発振モードのサンプリングを行なうことにより、複数のＸＹ面内スピン全体の様々なエネルギーのサンプリングを行なうことができる理由について説明する。ＸＹモデルのハミルトニアンは、数７で表される。ここで、Ｊｉｊは、スピンｉ、ｊの間の相互結合係数である。そして、θｉ、θｊは、それぞれ、スピンｉ、ｊのＸＹ面内での連続値をとる位相（０ｒａｄから２πｒａｄまで）を表わす。
【数７】
2019008502.tif 000009
【０１５７】
複数の振動子の時間発展を表わす複素ランジュバン方程式は、数８により表わされる。ここで、ｄＡｉは、振動子ｉの信号成分についての複素振幅変化を表わし、ｄｔは、時間変化を表わし、ｄＷｉは、振動子ｉの雑音成分についての複素ウィーナー過程を表わす。
【数８】
2019008502.tif 000010
【０１５８】
ｇ（Ａｉ）は、振動子ｉの強度増幅率であり、数９により表わされる。ここで、ｇ０は、小信号利得を表わし、ｎ０は、飽和強度を表わす。よって、（１／２）ｇ（Ａｉ）Ａｉｄｔは、振動子ｉの信号成分について、利得に起因する複素振幅変化を表わす。
【数９】
2019008502.tif 000011
【０１５９】
γｃは、振動子の強度減衰率であり、数１０により表わされる。ここで、ωは、振動子の発振周波数を表わし、Ｑは、振動子の共振Ｑ値を表わす。ただし、ω及びＱは、すべての振動子について、同一の値をとる。よって、−（１／２）γｃＡｉｄｔは、振動子ｉの信号成分について、損失に起因する複素振幅変化を表わす。
【数１０】
2019008502.tif 000012
【０１６０】
γｉｎｊは、相互結合率を表わす。よって、（γｉｎｊ／２）ＪｉｊＡｊｄｔは、振動子ｉの信号成分について、振動子ｊから振動子ｉへの相互結合に起因する複素振幅変化を表わす。そして、（γｉｎｊ／２）ΣＪｉｊＡｊｄｔは、振動子ｉの信号成分について、振動子ｉ以外のすべての振動子から振動子ｉへの相互結合に起因する複素振幅変化を表わす。
【０１６１】
Ｄは、振動子の振幅の拡散係数を表わす。よって、√（Ｄ）ｄＷｉは、振動子ｉの信号成分について、振幅拡散に起因する複素振幅変化を表わす。
【０１６２】
複数の振動子のポテンシャル関数は、数１１により表わされる。ここで、Ｈ〜（Ａ）は、ポテンシャル関数を表わし、太字のＡは、複数の振動子の信号成分について、複素振幅のベクトル表示を表わし、Ａｉ＊は、Ａｉの複素共役を表わす。
【数１１】
2019008502.tif 000013
【０１６３】
複数の振動子の時間発展を表わす複素ランジュバン方程式は、数１２のように書き換えられる。複数の振動子の信号成分について、複素振幅の定常状態の分布関数は、ポテンシャル関数を用いて、数１３により表わされ、より具体的には、数１４により表わされる。ここで、Ｐｓｔ（Ａ）は、定常状態の分布関数を表わし、Ｃは、規格化因子を表わす。
【数１２】
2019008502.tif 000014
【数１３】
2019008502.tif 000015
【数１４】
2019008502.tif 000016
【０１６４】
ここで、複数の振動子の間の相互結合の強度が低いため、すべての振動子の強度が独立して数１５により表わされる定常状態の同じ強度に安定化する、と仮定する。なお、ｎｓは、定常状態の同じ強度を表わす。すると、数１４の第１のｅｘｐの項のＤが小さいため、数１４の第１のｅｘｐの項が数１６により表わされるΠδ（｜Ａｉ｜２−ｎｓ）に近似される。なお、δは、δ関数を表わし、Ｃ’は、新たな規格化因子を表わす。
【数１５】
2019008502.tif 000017
【数１６】
2019008502.tif 000018
【０１６５】
ここで、振動子ｉの信号成分の複素振幅を、振幅√（ｎｓ）及び位相θｉを用いて、数１７のように表わす。すると、定常状態の分布関数は、数１８により表わされる。
【数１７】
2019008502.tif 000019
【数１８】
2019008502.tif 000020
【０１６６】
さて、複数の振動子の系と複数のスピンの系を対応付けるために、数１９で表わされる疑似的な逆温度β＝１／（ｋＢＴ）を定義する。ここで、ｋＢは、ボルツマン定数を表わし、Ｔは、複数のスピンの系の疑似的な温度を表わし、Ｄθ＝Ｄ／ｎｓは、振動子の位相の拡散係数を表わす。
【数１９】
2019008502.tif 000021
【０１６７】
すると、定常状態の位相の分布関数は、数２０により表わされる。ここで、太字のθは、複数の振動子の信号成分について、位相のベクトル表示を表わし、Ｈ（θ）は、数２１により表わされる、ＸＹモデルのハミルトニアンを表わす。つまり、複数の振動子の系の定常状態の分布関数は、複数のスピンの系の疑似的な熱平衡状態のボルツマン分布に対応するのである。
【数２０】
2019008502.tif 000022
【数２１】
2019008502.tif 000023
【０１６８】
（第３実施形態：計算方法）
次に、第３実施形態の結合振動子系の計算装置の構成を図１４に示す。結合振動子系の計算装置Ｃ５は、イテレーション部５１及びスピン系計算部５２から構成される。結合振動子系の計算装置Ｃ５に対して、結合振動子系の計算プログラムをインストールすることにより、イテレーション部５１及びスピン系計算部５２を実装することができる。
【０１６９】
イテレーション部５１は、結合振動子系のダイナミクスを記述する微分方程式（数８）の各イテレーションにおいて、各振動子ｉ毎の微小時間後の振幅及び位相（数８のＡｉの絶対値及び位相）を計算する。よって、イテレーション部５１は、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算することができる。
【０１７０】
スピン系計算部５２は、イテレーション部５１により反復されたイテレーションの後に、微分方程式（数８）の定常状態における各振動子ｉ毎の位相（数８のＡｉの位相）を計算することにより、結合振動子系の相互結合（数７のＪｉｊ）に対応するＸＹモデル（数７）の平衡状態を計算する。よって、スピン系計算部５２は、ＸＹモデルの平衡状態を計算することにより、機械学習のボルツマンサンプリングを行なうことができる。
【０１７１】
ここで、結合振動子系の全振動子数をＭとすると、数８の微分方程式のうち、右辺第１項の利得損失項及び右辺第３項の雑音項の計算量は、Ｍに比例するが、右辺第２項の相互結合項の計算量は、Ｍ２に比例する。
【０１７２】
しかし、数８の微分方程式のうち、右辺第２項の相互結合項の計算時間が、右辺第１項の利得損失項及び右辺第３項の雑音項の計算時間と比べて、同程度であるときには、図１４に示したように各イテレーションを行なえばよい。
【０１７３】
つまり、ある回のイテレーションの後でのＡｉを用いて計算した、右辺第１項の利得損失項、右辺第２項の相互結合項及び右辺第３項の雑音項を加算することにより、その回のイテレーションでの左辺第１項の微小時間変化を計算するのである。よって、右辺第２項の相互結合項を図３のような遅延（例えば、イテレーション数回分。）を伴わず計算しながら、右辺第１項の利得損失項及び右辺第３項の雑音項を高速に計算することができる。
【０１７４】
ただし、数８の微分方程式のうち、右辺第２項の相互結合項の計算時間が、右辺第１項の利得損失項及び右辺第３項の雑音項の計算時間と比べて、十分に長いときには、図３と同じように各イテレーションを行なえばよい。
【０１７５】
つまり、ある回のイテレーションの後でのＡｉを用いて計算した右辺第１項の利得損失項及び右辺第３項の雑音項と、その回のイテレーションよりも数回（例えば、１回のみ。）遅れた回のイテレーションの後でのＡｉを用いて計算した右辺第２項の相互結合項と、を加算することにより、その回のイテレーションでの左辺第１項の微小時間変化を計算するのである。よって、右辺第２項の相互結合項を余裕を持って計算しながら、右辺第１項の利得損失項及び右辺第３項の雑音項を高速に計算することができる。
【０１７６】
結合振動子系の計算装置Ｃ５のうち、イテレーション部５１は、利得損失計算部５１１、相互結合計算部５１２、雑音計算部５１３及び振幅位相計算部５１４から構成される。
【０１７７】
利得損失計算部５１１は、数８の微分方程式の各イテレーションにおいて、各振動子ｉ毎の利得及び損失に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、利得損失計算部５１１は、数８の微分方程式の各イテレーションにおいて、第ｎ回目のイテレーションの後でのＡｉ（ｎ）を入力され、数８の右辺第１項を出力する。ここで、第１回目のイテレーションの前でのＡｉ（０）は、０に初期設定される。
【０１７８】
雑音計算部５１３は、数８の微分方程式の各イテレーションにおいて、各振動子ｉ毎の雑音に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、雑音計算部５１３は、数８の微分方程式の各イテレーションにおいて、第ｎ回目のイテレーションの後でのＡｉ（ｎ）を入力され、数８の右辺第３項を出力する。
【０１７９】
相互結合計算部５１２は、数８の微分方程式の各イテレーションにおいて、各振動子ｉの間の相互結合に起因する各振動子ｉ毎の振幅及び位相の微小時間変化を計算する。具体的には、相互結合計算部５１２は、数８の微分方程式の各イテレーションにおいて、第ｎ回目のイテレーションの後でのＡｉ（ｎ）を入力され、数８の右辺第２項を出力する。ここで、相互結合計算部５１２での数８の右辺第２項の計算時間は、利得損失計算部５１１及び雑音計算部５１３での数８の右辺第１、３項の計算時間と比べて、同程度である。
【０１８０】
振幅位相計算部５１４は、数８の微分方程式の各イテレーションにおいて、利得損失計算部５１１、雑音計算部５１３及び相互結合計算部５１２により各々計算された各振動子ｉ毎の振幅及び位相の微小時間変化に基づいて、各振動子ｉ毎の微小時間後の振幅及び位相を計算する。具体的には、振幅位相計算部５１４は、数８の微分方程式の各イテレーションにおいて、利得損失計算部５１１、雑音計算部５１３及び相互結合計算部５１２からの数８の右辺第１、３、２項を入力され、利得損失計算部５１１、雑音計算部５１３及び相互結合計算部５１２からの数８の右辺第１、３、２項を加算し、第（ｎ＋１）回目のイテレーションの後でのＡｉ（ｎ＋１）を出力する。ここで、利得損失計算部５１１、雑音計算部５１３及び相互結合計算部５１２からの数８の右辺第１、３、２項は、第ｎ回目のイテレーションの後でのＡｉ（ｎ）を用いて計算されたものである。
【０１８１】
スピン系計算部５２は、イテレーション部５１により反復されたイテレーションの後に、数８の微分方程式の定常状態における各振動子ｉ毎の位相を計算することにより、結合振動子系の相互結合に対応するＸＹモデルの平衡状態を計算する。
【０１８２】
このように、機械学習のボルツマンサンプリングを行なうために、機械学習を連続値スピンモデルへマッピングし、連続値スピンモデルの平衡状態を計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。そして、完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算するにあたり、電子回路で実装することにより、光学回路を安定化する必要をなくすことができる。
【０１８３】
第３実施形態においても、第１実施形態の図６〜８と同様に、並列計算を行なうことができ、第１実施形態の図９〜１２と同様に、ビット量子化を行なうことができる。
【産業上の利用可能性】
【０１８４】
本開示の結合振動子系の計算装置、プログラム及び方法は、次の目的に利用可能である：（１）組み合わせ最適化問題の最適解を求めるために、組み合わせ最適化問題をイジングモデルへマッピングし、イジングモデルの基底状態を計算するにあたり、光学回路を安定化する必要をなくす、（２）機械学習のボルツマンサンプリングを行なうために、機械学習を連続値スピンモデルへマッピングし、連続値スピンモデルの平衡状態を計算するにあたり、光学回路を安定化する必要をなくす、（３）完全結合系のニューロチップ等において、結合振動子系のダイナミクスを計算するにあたり、光学回路を安定化する必要をなくす。
【符号の説明】
【０１８５】
Ｃ１、Ｃ２：イジングモデルの計算装置
Ｃ３、Ｃ４、Ｃ５：結合振動子系の計算装置
Ｐ１、Ｐ２、Ｐ３、Ｐ４：レーザパルス
１１、２１：レーザ発振部
１２、２２：リング共振器
１３、１５、１７、２４：相互結合実装部
１４、１６、１８：パルス遅延線
１９、２３：スピン測定部
３１、４１、５１：イテレーション部
３２、４２、５２：スピン系計算部
３３：ＦＩＦＯ
３１１、４１１、５１１：利得損失計算部
３１２、４１３、５１３：雑音計算部
３１３、４１２、５１２：相互結合計算部
３１４、４１４、５１４：振幅位相計算部
【図１】
2019008502.tif 000024
【図２】
2019008502.tif 000025
【図３】
2019008502.tif 000026
【図４】
2019008502.tif 000027
【図５】
2019008502.tif 000028
【図６】
2019008502.tif 000029
【図７】
2019008502.tif 000030
【図８】
2019008502.tif 000031
【図９】
2019008502.tif 000032
【図１０】
2019008502.tif 000033
【図１１】
2019008502.tif 000034
【図１２】
2019008502.tif 000035
【図１３】
2019008502.tif 000036
【図１４】
2019008502.tif 000037