(57)【要約】【課題】  １つの乱数生成器をＧＰＵ上で並列化し、逐次的に乱数生成処理を行った場合と同一の乱数列をより効率よく生成する。【解決手段】  １個以上のコアからなるブロックを複数備えたＧＰＵを用いて、状態ベクトルを更新する更新処理と、更新された状態ベクトルを別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う方法は、ある１つのブロックが前記更新処理を実行するステップ（Ｓ３）と、前記更新処理の結果を受けて、複数のブロックの各々が前記変換処理を並列に実行するステップ（Ｓ９）とを含む。【選択図】  図１
2012094040.tif 000002
【特許請求の範囲】
【請求項１】
１個以上のコアからなるブロックを複数備えたＧＰＵを用いて、状態ベクトルを更新する更新処理と、更新された状態ベクトルを別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う方法であって、
ある１つのブロックが前記更新処理を実行するステップと、
前記更新処理の結果を受けて、複数のブロックの各々が前記変換処理を並列に実行するステップと
を含む方法。
【請求項２】
１個以上のコアからなるブロックを複数備えたＧＰＵを用いて、状態ベクトルを更新する更新処理と、更新された状態ベクトルを整数乱数へと変換する調律処理と、変換の結果生成された整数乱数を別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う方法であって、
ある１つのブロックが前記更新処理を実行するステップと、
前記更新処理の結果を受けて、複数のブロックの各々が前記調律処理及び前記変換処理を並列に実行するステップと
を含む方法。
【請求項３】
前記変換処理あるいは前記調律処理及び前記変換処理を実行する複数のブロックが当該処理を実行する際に生成するスレッドの数が、前記更新処理を実行する１つのブロックが当該処理を実行する際に生成するスレッドの数よりも多い、請求項１または２に記載の方法。
【請求項４】
１個以上のコアからなるブロックを複数備えたＧＰＵに対して、状態ベクトルを更新する更新処理と、更新された状態ベクトルを別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を実行させるコンピュータプログラムであって、
ある１つのブロックが前記更新処理を実行するステップと、
前記更新処理の結果を受けて、複数のブロックの各々が前記変換処理を並列に実行するステップと
を含むコンピュータプログラム。
【請求項５】
１個以上のコアからなるブロックを複数備えたＧＰＵに対して、状態ベクトルを更新する更新処理と、更新された状態ベクトルを整数乱数へと変換する調律処理と、変換の結果生成された整数乱数を別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を実行させるコンピュータプログラムであって、
ある１つのブロックが前記更新処理を実行するステップと、
前記更新処理の結果を受けて、複数のブロックの各々が前記調律処理及び前記変換処理を並列に実行するステップと
を含むコンピュータプログラム。
【請求項６】
前記変換処理あるいは前記調律処理及び前記変換処理を実行する複数のブロックが当該処理を実行する際に生成するスレッドの数が、前記更新処理を実行する１つのブロックが当該処理を実行する際に生成するスレッドの数よりも多い、請求項４または５に記載のコンピュータプログラム。
【請求項７】
１個以上のコアからなるブロックを複数備えたＧＰＵを有しており、状態ベクトルを更新する更新処理と、更新された状態ベクトルを別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う計算装置であって、
ある１つのブロックが前記更新処理を実行し、
前記更新処理の結果を受けて、複数のブロックの各々が前記変換処理を並列に実行する、計算装置。
【請求項８】
１個以上のコアからなるブロックを複数備えたＧＰＵを有しており、状態ベクトルを更新する更新処理と、更新された状態ベクトルを整数乱数へと変換する調律処理と、変換の結果生成された整数乱数を別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う計算装置であって、
ある１つのブロックが前記更新処理を実行し、
前記更新処理の結果を受けて、複数のブロックの各々が前記調律処理及び前記変換処理を並列に実行する、計算装置。
【請求項９】
前記変換処理あるいは前記調律処理及び前記変換処理を実行する複数のブロックが当該処理を実行する際に生成するスレッドの数が、前記更新処理を実行する１つのブロックが当該処理を実行する際に生成するスレッドの数よりも多い、請求項７または８に記載の計算装置。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、ＧＰＵ（Graphics Processing Unit）を用いた乱数の生成処理に関する。
【背景技術】
【０００２】
ＧＰＵは、ＣＧ（Computer Graphics）などの画像処理を高速に実行することに特化したデバイスであり、通常はパーソナルコンピュータやワークステーション等の計算装置に内蔵されて画像の表示のために用いられている。このＧＰＵは、安価でありながら汎用のＣＰＵ（Central Processing Unit）を大幅に超える理論性能を持つという性質から、本来の画像処理だけではなく科学技術計算等にも利用されている（特許文献１を参照）。このようなＧＰＵを用いた計算手法は、ＧＰＵコンピューティングと呼ばれている。
【０００３】
図４はＧＰＵを備えた計算装置のブロック図である。図示されているように、計算装置は、ＣＰＵ１と、メインメモリ２と、ＧＰＵ３と、ＧＰＵメモリ４と、キーボード及びマウスを含む入力装置５と、液晶ディスプレイなどの表示装置６と、これらを接続するバス７とを備えている。ＧＰＵ３は、図示しないコアと呼ばれる計算ユニットからなるブロックを複数持っている。そして、ブロック内のコア同士は高速にデータを共有することができる。他方、あるブロック内のコアと別のブロック内のコアとによるデータの共有は、グローバルメモリであるＧＰＵメモリ４を介して行わなければならないため、比較的に時間がかかる。また、ＧＰＵメモリ４は可変のデータ記憶領域を有しており、ＧＰＵ３から直接アクセスされてデータの読み書きをすることができる。ＣＰＵ１は、ＧＰＵメモリ４に対して直接にアクセスすることはできないものの、ＧＰＵ３を介してＧＰＵメモリ４にアクセスすることができる。
【０００４】
ＧＰＵは、ＣＰＵと比べて多くのコアを有している。そのため、ＧＰＵはＣＰＵよりも優れた計算性能を持つ。その一方で、ＧＰＵ内のそれぞれのコアは同一時刻に同一の演算しか実行できないため、ＣＰＵと比べて並列化の際の制約が大きく、ＧＰＵの性能を発揮できる問題が限られる。また、ＧＰＵのアーキテクチャはＣＰＵのそれと大きく異なるため、ＣＰＵ用のプログラムコードを直接ＧＰＵに利用してもＧＰＵの性能を十分に利用できない。そこで、ＧＰＵに適したアルゴリズムや処理順序に合わせてＣＰＵ用のプログラムコードを書き換える必要がある。
【０００５】
ＧＰＵがその性能を最も発揮できる分野の一つとして確率的現象の解析等に用いられるモンテカルロ法（Monte Carlo method）が挙げられる。モンテカルロ法の実行には大量の乱数が用いられることと、乱数の質がモンテカルロ法の実行結果に影響を与えることとから、質の良い乱数列を高速に生成することのできる乱数生成器が必要となる。また、通常、乱数生成器が生成する乱数は一様乱数であるが、モンテカルロ法で用いる乱数はそれぞれ必要に応じた分布に従わなければならず、乱数生成器が出力する一様乱数を適切な分布に従う乱数へと変換しなければならない。
【０００６】
乱数生成器として高速かつ質の良い乱数を生成するメルセンヌ・ツイスター（Mersenne Twister）乱数生成器は１９９６年から１９９８年にかけて開発、発表されたものである（非特許文献１を参照）。ＣＰＵ用にＣ言語で書かれた２００２年版標準コード（以下、ｍｔ１９９３７ａｒ）が代表的な例である。加えて、ＧＰＵ用に書かれたＭＴＧＰ（Mersenne Twister for Graphic Processors）が２００９年に公開されている。詳細は非特許文献２を参照されたい。ところが、このＭＴＧＰはＣＰＵ向けのｍｔ１９９３７ａｒとは処理手順が異なる。そのため、ＭＴＧＰによって生成される乱数列は、ｍｔ１９９３７ａｒによって生成される乱数列とは通常は異なる。
【０００７】
非特許文献２のＭＴＧＰは、パラメータの異なった複数のメルセンヌ・ツイスター乱数生成器を用意し、ＧＰＵ内のそれぞれのブロックで別々の乱数生成器を実行している。つまり、１つのブロックが１つの乱数生成器に相当する。さらに、ＭＴＧＰは、ブロック内の通信の高速性を利用して、ブロック内で乱数生成器の並列化を行っている。
【０００８】
図５は上記ＭＴＧＰにより乱数を生成する際のＧＰＵ３の構成を模式的に示した図である。図示されているように、ＧＰＵ３はｋ個のブロック、すなわち第１ブロック〜第ｋブロックを有している。そして、各ブロックは図示しない複数のコアを有している。これら複数のコアは、第１スレッド〜第ｎスレッドの実行を担当する。ここで、スレッドとは、プロセス内の処理の流れのことを指す。なお、１個のコアが１個のスレッドを担当してもよいし、２個以上のスレッドを担当しても構わない。このようにして、第１ブロック〜第ｋブロックに対して乱数生成器ＭＴ１〜ＭＴｋをそれぞれ対応させている。つまり、第１ブロック〜第ｋブロックが乱数生成器ＭＴ１〜ＭＴｋにそれぞれ相当する。そして、各乱数生成器の処理は、各ブロック内のｎ個のスレッドにより並列に行われる。
【０００９】
図６は、上記乱数生成器ＭＴ１〜ＭＴｋの各々が実行する乱数生成処理のフローチャートである。メルセンヌ・ツイスター法に限らず、多くの乱数生成器における乱数生成処理は、状態ベクトルと呼ばれるベクトルを更新する更新処理と、更新された状態ベクトルを性質のよい整数乱数へと変換する調律処理と、その整数乱数を別の分布を持った乱数へと変換する変換処理とを含む。まず、ステップＳ５１にて乱数生成処理を開始すると同時に、ＧＰＵ３は乱数生成処理プログラムをＧＰＵメモリ４から読み込むなどの準備を行う。ステップＳ５２において変数Ｍの値を変数Ｎに代入する。この変数Ｍは、後述する変数Ｌとともに、メルセンヌ・ツイスター法の実行の態様によって定まるものである。次にステップＳ５３に進み、状態ベクトルの更新処理を行う。具体的には、次に示す漸化式を用いて状態ビットＲｉ（ｉは自然数）を更新する。
【数１】
2012094040.tif 000003
ただし、状態ビットＲｉはＫビットの変数であり、Ｆはビット演算及び四則演算などからなる関数であり、Ｌ及びＭはＬ＜Ｍを満たす正の整数である。これらＬ及びＭはメルセンヌ・ツイスター法の実行の態様によって定まる。
【００１０】
ステップＳ５４では調律処理を行う。具体的には、次式により、ステップＳ５３で求めた状態ビットＲｉを整数値Ｓｉへ変換する。
【数２】
2012094040.tif 000004
ただし、Ｇはビット演算からなる関数である。
【００１１】
続いて、変換処理を行うステップＳ５５は、必須のステップではなく、必要に応じて実行するステップである。具体的には、次式により、ステップＳ５４で求めた整数値Ｓｉを一様乱数Ｔｉへと変換してＧＰＵメモリ４に保存する。
【数３】
2012094040.tif 000005
ただし、Ｈは一様乱数への変換関数である。入力は複数必要な場合もある。また、ＭｅｍはＧＰＵメモリ４への保存を示している。このように、ステップＳ５５は、ステップＳ５４で求めた整数値Ｓｉを一様乱数Ｔｉへ変換したい場合に行うステップである。
【００１２】
ステップＳ５６では、必要とされる個数の乱数が生成されたか否かを各ブロックが判断する。必要とされる個数の乱数が生成されていない場合には、ステップＳ５７に進み、変数Ｎに並列実行数ｎを加えた値を変数Ｎに代入する。そして、ステップＳ５３〜Ｓ５５を実行する。このようにして、必要とされる個数の乱数が生成されると、ステップＳ５８にて乱数生成処理を終了する。
【００１３】
上述したように、乱数生成器ＭＴ１〜ＭＴｋのそれぞれがステップＳ５１〜Ｓ５８を実行する。また、図６にも示したように、ステップＳ５３〜Ｓ５７の処理はｎ個のスレッドにより並列に行われる。
【００１４】
また、ＧＰＵ開発を行っている米ＮＶＩＤＩＡ社も、乱数生成器の並列化のサンプルとしてメルセンヌ・ツイスター乱数生成器を公開している（非特許文献３を参照）。米ＮＶＩＤＩＡ社の態様では、１つのスレッドに対して１つの乱数生成器を持たせ、並列に乱数を生成させている。図７は米ＮＶＩＤＩＡ社の態様により乱数を生成する際のＧＰＵの構成を模式的に示した図である。図示しているように、第１ブロック内の第１〜第ｍのスレッドは、第１〜第ｍの乱数生成器ＭＴ１〜ＭＴｍをそれぞれ担当する。第２ブロック内の第１〜第ｍのスレッドは、第（ｍ＋１）〜第２ｍの乱数生成器ＭＴｍ＋１〜ＭＴ２ｍをそれぞれ担当する。同様にして、第ｋブロック内の第１〜第ｍのスレッドは、第（ｋ−１）ｍ＋１〜第ｋｍの乱数生成器ＭＴ（ｋ−１）ｍ＋１〜ＭＴｋｍを担当する。このようにして、１つのスレッドに対して１つの乱数生成器を持たせている。つまり、１つのスレッドが１つの乱数生成器に相当する。各スレッド内の処理は通常のメルセンヌ・ツイスター・アルゴリズムと同一である。
【先行技術文献】
【特許文献】
【００１５】
【特許文献１】特開２００９−８６２０２号公報
【非特許文献１】M. Matsumoto and T. Nishimura, "Mersenne Twister: A 623-dimensionally equidistributed uniform pseudorandom number generator", ACM Trans. on Modeling and Computer Simulation, January, 1998, Vol. 8, No. 1, pp.3-30
【非特許文献２】松本  眞、”Mersenne Twister for Graphic Processors (MTGP): a new variant of Mersenne Twister”、[online]、［平成２２年９月２５日検索］、インターネット（ＵＲＬ：http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/MTGP/index-jp.html）
【非特許文献３】米ＮＶＩＤＩＡ社、”CUDA ZONE”、[online]、［平成２２年９月２５日検索］、インターネット（http://www.nvidia.co.jp/object/cuda_get_samples_jp.html）
【発明の概要】
【発明が解決しようとする課題】
【００１６】
図５に示したように、既存のＭＴＧＰは、複数の乱数生成器を並列に実行することで並列化を達成している。しかし、複数の乱数生成器から生成される乱数が混ぜ合わされた乱数列が得られることとなるため、１つの乱数生成器を用いた場合とは本質的に異なる乱数列しか得られない。つまり、逐次処理であるｍｔ１９９３７ａｒと、並列処理であるＭＴＧＰとでは、得られる乱数が異なり、その結果としてモンテカルロ法の実行結果の再現性が失われてしまう。この再現性は、金融工学の分野などでバックテストを行うときや計算条件を変更した際の影響を検証する際に欠かせないものである。ＭＴＧＰにおいて再現性を保持するためには、実行する乱数生成器を１つとすればよいが、そうするとＧＰＵ内のある１ブロックしか用いないこととなるため、他のブロックの計算資源が無駄となる。これは図７に示した米ＮＶＩＤＩＡ社の態様でも同様である。また、１つ１つの乱数生成器の実行速度は速くはないという問題もある。
【００１７】
かかる状況に鑑みて、本発明は、１つの乱数生成器をＧＰＵ上で並列化し、逐次的に乱数生成処理を行った場合と同一の乱数列をより効率よく生成することを目的とする。
【００１８】
本発明は、更新処理の並列度と、調律処理及び変換処理の並列度とが異なることに着目してなされている。まず、更新処理については、ｎ≦Ｌとして、状態ビットＲＮ〜ＲＮ＋ｎ−１というｎ個の状態ビットの計算は、互いにデータ依存関係がないため同時に行うことができる。つまり、ＧＰＵ内のある１ブロックを用いて複数のスレッドを生成させて、更新処理を並列に行うことができる。このようにして、同一時刻に同一の演算しか実行できないというＧＰＵの制約をうまく利用して更新処理を並列化する。
【００１９】
他方、調律処理における整数値Ｓｉについては、異なる添え字ｉに対してデータ依存関係がない。そのため、任意の個数のＳｉを並列に計算することができる。すなわち、調律処理は、ＧＰＵ内の複数のブロックを用いて並列に行うことができる。
【００２０】
変換処理における乱数Ｔｉ＝Ｈ（Ｓｉ）についても、異なる添え字ｉに対してデータ依存関係がない。そのため、任意の個数のＴｉを並列に計算することができる。すなわち、変換処理も、ＧＰＵ内の複数のブロックを用いて並列に行うことができる。
【００２１】
このように、並列度の異なる更新処理と調律処理及び変換処理とを分離し、ＧＰＵにおいてそれぞれ異なる並列数で実行することで、計算の効率性を高めることを行う。つまり、更新処理は、ＧＰＵ内のある１ブロックを用いて整数Ｌに近い並列度ｎで実行する。その一方で、調律処理及び変換処理は、ＧＰＵ内の複数のブロックを用いることにより、ＧＰＵの計算資源をなるべく活用できる並列度で計算を行う。このようにすることで、単一の乱数生成器を用いて乱数を生成する際に計算資源を有効に利用することができる。
【課題を解決するための手段】
【００２２】
本発明は、まず、１個以上のコアからなるブロックを複数備えたＧＰＵを用いて、状態ベクトルを更新する更新処理と、更新された状態ベクトルを別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う方法であって、ある１つのブロックが前記更新処理を実行するステップと、前記更新処理の結果を受けて、複数のブロックの各々が前記変換処理を並列に実行するステップとを含む方法を提供する。
【００２３】
また、本発明は、１個以上のコアからなるブロックを複数備えたＧＰＵを用いて、状態ベクトルを更新する更新処理と、更新された状態ベクトルを整数乱数へと変換する調律処理と、変換の結果生成された整数乱数を別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う方法であって、ある１つのブロックが前記更新処理を実行するステップと、前記更新処理の結果を受けて、複数のブロックの各々が前記調律処理及び前記変換処理を並列に実行するステップとを含む方法を提供する。
【００２４】
ここで、前記変換処理あるいは前記調律処理及び前記変換処理を実行する複数のブロックが当該処理を実行する際に生成するスレッドの数が、前記更新処理を実行する１つのブロックが当該処理を実行する際に生成するスレッドの数よりも多くてもよい。
【００２５】
さらに、本発明は、１個以上のコアからなるブロックを複数備えたＧＰＵに対して、状態ベクトルを更新する更新処理と、更新された状態ベクトルを別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を実行させるコンピュータプログラムであって、ある１つのブロックが前記更新処理を実行するステップと、前記更新処理の結果を受けて、複数のブロックの各々が前記変換処理を並列に実行するステップとを含むコンピュータプログラムを提供する。
【００２６】
また、本発明は、１個以上のコアからなるブロックを複数備えたＧＰＵに対して、状態ベクトルを更新する更新処理と、更新された状態ベクトルを整数乱数へと変換する調律処理と、変換の結果生成された整数乱数を別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を実行させるコンピュータプログラムであって、ある１つのブロックが前記更新処理を実行するステップと、前記更新処理の結果を受けて、複数のブロックの各々が前記調律処理及び前記変換処理を並列に実行するステップとを含むコンピュータプログラムを提供する。
【００２７】
ここで、前記変換処理あるいは前記調律処理及び前記変換処理を実行する複数のブロックが当該処理を実行する際に生成するスレッドの数が、前記更新処理を実行する１つのブロックが当該処理を実行する際に生成するスレッドの数よりも多くてもよい。
【００２８】
加えて、本発明は、１個以上のコアからなるブロックを複数備えたＧＰＵを有しており、状態ベクトルを更新する更新処理と、更新された状態ベクトルを別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う計算装置であって、ある１つのブロックが前記更新処理を実行し、前記更新処理の結果を受けて、複数のブロックの各々が前記変換処理を並列に実行する、計算装置を提供する。
【００２９】
また、本発明は、１個以上のコアからなるブロックを複数備えたＧＰＵを有しており、状態ベクトルを更新する更新処理と、更新された状態ベクトルを整数乱数へと変換する調律処理と、変換の結果生成された整数乱数を別の分布を持った乱数へと変換する変換処理とを含む乱数生成処理を行う計算装置であって、ある１つのブロックが前記更新処理を実行し、前記更新処理の結果を受けて、複数のブロックの各々が前記調律処理及び前記変換処理を並列に実行する、計算装置を提供する。
【００３０】
ここで、前記変換処理あるいは前記調律処理及び前記変換処理を実行する複数のブロックが当該処理を実行する際に生成するスレッドの数が、前記更新処理を実行する１つのブロックが当該処理を実行する際に生成するスレッドの数よりも多くてもよい。
【００３１】
本明細書及び本願特許請求の範囲において、「乱数」とは、擬似乱数、準乱数を含む概念である。数学的に純粋な意味での乱数を生成することは、デジタルコンピュータでは不可能であることはよく知られており、ここでは、モンテカルロ法を実行するのに十分な擬似乱数や、数値積分に用いられる規則性の強い準乱数などを含む広い意味で「乱数」という用語を用いる。
【発明の効果】
【００３２】
本発明により、１つの乱数生成器をＧＰＵ上に実現することにより効率的に乱数列を発生させることができる。
【図面の簡単な説明】
【００３３】
【図１】本発明による乱数生成のフローチャートの一例である
【図２】更新処理を行う際のＧＰＵの構成を模式的に示したダイアグラムである。
【図３】調律処理及び変換処理を行う際のＧＰＵの構成を模式的に示したダイアグラムである。
【図４】ＧＰＵを備えた計算装置のブロック図である。
【図５】ＭＴＧＰにより乱数を生成する際のＧＰＵの構成を模式的に示したダイアグラムである。
【図６】ＭＴＧＰにおける乱数生成処理のフローチャートである。
【図７】米ＮＶＩＤＩＡ社の態様により乱数を生成する際のＧＰＵの構成を模式的に示したダイアグラムである。
【発明を実施するための形態】
【００３４】
図１は、本発明による乱数生成のフローチャートである。この図１に示した処理を行う際のＧＰＵの構成を図２及び図３のダイアグラムに模式的に示している。まず、ステップＳ１にて乱数生成処理を開始すると同時に、ＧＰＵ３は更新処理プログラムをＧＰＵメモリ４から読み込むなどの準備を行う。ステップＳ２において、変数Ｍの値を変数Ｎに代入する。この変数Ｍは、後述する変数Ｌとともに、メルセンヌ・ツイスター法の実行の態様によって定まるものである。ステップＳ３では、状態ベクトルの更新処理を行う。この更新処理は、単一のブロック内で生成したｎ個のスレッドを用いて、ｎ個の状態ビットＲＮ〜ＲＮ＋ｎ−１の次式に基づいた更新を並列に行う。
【数４】
2012094040.tif 000006
ただし、状態ビットＲｉはＫビットの変数であり、Ｆはビット演算及び四則演算などからなる関数であり、Ｌ及びＭはＬ＜Ｍを満たす正の整数である。これらＬ及びＭはメルセンヌ・ツイスター法の実行の態様によって定まる。
【００３５】
ステップＳ４では、ステップＳ３にて生成したｎ個の状態ビットＲＮ〜ＲＮ＋ｎ−１を次式に従ってグローバルメモリに保存する処理を並列に行う。
【数５】
2012094040.tif 000007
ただし、Ｍｅｍはグローバルメモリへの保存を示している。
【００３６】
ステップＳ５では、必要とされる個数の状態ビットＲｉが生成されたか否かを判断する。必要とされる個数の状態ビットが生成されていない場合には、ステップＳ６に進み、変数Ｎに並列実行数ｎを加えた値を変数Ｎに代入する。そして、ステップＳ３及びＳ４を実行する。このようにして、必要とされる個数の乱数が生成されるまでステップＳ３及びＳ４を繰り返す。
【００３７】
図１及び図２に示したように、ステップＳ３〜Ｓ６の各ステップは、ＧＰＵ３内のある単一のブロックで生成させたｎ個のスレッドにより並列に処理される。このとき、その他のブロックは使用しない。つまり、更新処理においてはスレッド間にデータの授受が発生するが、単一のブロック内で並列処理をすることにより、時間のかかるグローバルメモリへのアクセスを回避することができる。なお、更新処理を担当するのは、第１ブロックに限らず、いずれか１つのブロックとすることができる。
【００３８】
必要な個数の状態ビットＲｉが生成されると、ステップＳ７に進み、変数Ｍの値を変数Ｎに代入すると同時に、調律処理及び変換処理のプログラムをＧＰＵメモリ４から読み込む。
【００３９】
ステップＳ８では調律処理を行う。具体的には、グローバルメモリから状態ビットＲｉの値を取得し、その状態ビットを整数乱数Ｓｉへ変換する。これを次式に示す。
【数６】
2012094040.tif 000008
ただし、Ｇはビット演算からなる関数である。
【００４０】
ステップＳ９では変換処理を行う。具体的には、整数乱数Ｓｉを任意の分布を持つ乱数Ｈ（Ｓｉ）へと変換し、それらをグローバルメモリに保存する。これを次式に示す。
【数７】
2012094040.tif 000009
ただし、Ｈはそれぞれの分布を持つ乱数への変換関数である。例えば、連続一様分布や正規分布への変換関数が挙げられる。また、入力は複数必要な場合もある。加えて、ＭｅｍはＧＰＵメモリ４への保存を示している。
【００４１】
ステップＳ１０では、必要とされる個数の変換処理が行われたか否かを判断する。必要とされる個数の変換処理が行われていない場合には、ステップＳ１１に進み、変数Ｎに並列実行数ｋ×ｍを加えた値を変数Ｎに代入する。ただし、ｋはＧＰＵ内のブロックの個数であり、ｍは各ブロック内のスレッドの個数である。その後、ステップＳ８及びＳ９を実行する。このようにして、必要とされる個数の変換処理が行われるまでステップＳ８及びＳ９を繰り返す。そして、ステップＳ１２にて処理を終了する。
【００４２】
図１及び図３に示したように、ステップＳ８〜Ｓ１１の各ステップは、合計ｋ×ｍ個のスレッドにより並列に処理される。つまり、ＧＰＵ内の計算資源を効率的に使用することができる。
【００４３】
以上のように、更新処理においてはスレッド間にデータの授受が発生するため、図２に示したように、１つの乱数生成器をある１つのブロック内のスレッドで共有して１ブロックのみで計算を行い、その結果をグローバルメモリに保存する。その一方で、調律処理及び変換処理はデータの依存関係が存在しないため、図３に示したようにＧＰＵ内の全ての計算資源を用いることができる。これにより、モンテカルロ法の実行結果の再現性を保持することができるとともに、計算資源の無駄を減らして効率良く使用することができる。
【００４４】
［実施例］
実施例１として、Ｃ言語で書かれたメルセンヌ・ツイスター・２００２年標準コード「ｍｔ１９９３７ａｒ」を上記の実施形態に基づいてＧＰＵにより実行した。整数乱数、連続一様乱数、正規乱数をそれぞれ１億個生成して、その実行時間を測定した。このとき、更新処理と調律処理及び変換処理とを分離して、更新処理の並列実行数は１ブロック×２２４スレッドとし、調律処理及び変換処理の並列実行数は２８ブロック×６４スレッドとした。
【００４５】
一方、比較例１として、ｍｔ１９９３７ａｒをＣＰＵで実行した。さらに、比較例２として、更新処理、調律処理、変換処理の各処理を分離せずに、全ての処理をＧＰＵ内のある１つのブロックで実行した。比較例１及び２でも整数乱数、連続一様乱数、正規乱数をそれぞれ１億個生成した。
【００４６】
実施例１、比較例１、比較例２のいずれについても、正規乱数への変換処理においてはＭｏｒｏの逆関数法を用いた。詳細は非特許文献「Moro, B., “The Full Monte,” Risk, Vol. 8, No. 2, 1995」を参照されたい。また、検証環境として、米ＮＶＩＤＩＡ社のＧＰＵであるＴｅｓｌａＣ２０５０を取り付けた米Ｄｅｌｌ社のＰｒｅｃｉｓｉｏｎＴ５５００を用いた。その詳細を表１に示している。
【００４７】
【表１】
2012094040.tif 000010
【００４８】
以上のような条件下での実施例１、比較例１、比較例２の実行時間を表２に示している。整数乱数を生成した場合には、実施例１は比較例１と比べて実行時間を約８１．６％削減できている。また、正規乱数を生成した場合には、実施例１は比較例１と比べて実行時間を約９４．６％削減できている。同じく正規乱数を生成した場合の実施例１は、ＧＰＵを用いた比較例２と比べても、実行時間を約８４．１％削減できている。
【００４９】
【表２】
2012094040.tif 000011
【００５０】
なお、上記の実施形態及び実施例においては乱数生成アルゴリズムとしてメルセンヌ・ツイスター法を用いた。しかし、これに限られず、線形合同法といった他のアルゴリズムを用いることもできる。かかる他の乱数生成アルゴリズムの多くは、状態ベクトルを更新する更新処理と、更新された状態ベクトルそのものを別の分布を持った乱数へと変換する変換処理とを含んでおり、調律処理を含まない。これは、更新処理にて更新された状態ベクトルそれ自体が整数乱数であるために、調律処理が不要となるからである。このような、調律処理を含まないアルゴリズムを用いることもできる。その場合には、ＧＰＵ内のいずれか１つのブロックが状態ベクトルの更新処理を行い、続いてＧＰＵ内の全てのブロックが、更新された状態ベクトルの変換処理を行う。さらに別の形態では、少なくとも更新処理と変換処理とを含む準乱数生成アルゴリズムを用いることもできる。また、１つのスレッドを１個のコアが担当するように構成することもできるし、複数のスレッドを１個のコアが担当することもできる。
【００５１】
また、本発明の各実施の形態を実施するためのコンピュータのハードウエア構成は特段図示した構成に限定されるものではなく、数値演算の形式に併せて専用の数値演算ユニットを一つ以上有していたり、複数の筐体に分かれて互いにネットワークにより接続されているクラスタ構成にされていたりすることができる。なお、本実施例のみならず、本発明全般に、区別して記載された機能手段は、実質的にそのような区別された機能を果たす任意の構成要素によって実現される。このとき、その構成要素が、物理的にいくつあるのか、あるいは、複数あるのであれば互いにどのような位置関係にあるかといった、機能を果たす上で制限とならない属性によって、本発明が制限されることはない。例えば、複数の区別された機能が単一の構成要素によって経時的に異なるタイミングで実行されることも本発明の範囲に含まれる。
【００５２】
本発明の各実施の形態の機能処理を実現するコンピュータにおいて数値計算をするためのソフトウエア構成は、本発明の各実施の形態の数値情報処理を実現する限り任意の構成とすることができる。そのコンピュータは、基本入出力システム（ＢＩＯＳ）などのハードウエア制御のためのソフトウエアを搭載しており、これと連携して動作し、ファイル入出力やハードウエアリソースの割り振りを担当するオペレーティングシステム（ＯＳ）によって管理されている。当該ＯＳは、ＯＳやハードウエアと連携して動作するアプリケーションプログラムを、例えばユーザーからの明示の命令や、ユーザーからの間接的な命令や他のプログラムからの命令に基づいて実行することができる。アプリケーションプログラムは、このような動作を可能とし、ＯＳと関連して動作するように、ＯＳの規定する手続に依存して、あるいはＯＳに依存しないように適切にプログラムされている。本発明の各実施の形態を実現する場合には、一般に、専用のアプリケーションプログラムの形式で数値計算やファイル入出力等の処理を実行するが、本発明はそれのみに限定されるものではない。本発明の方法またはプログラムは、複数の専用または汎用アプリケーションプログラムを用いたり、既成の数値計算ライブラリを部分的に用いたり、他のコンピュータのハードウエアによって処理されるようにネットワークプログラミング手法によって実現されていたり、その他の任意の実施形態によって実現されうる。したがって、本発明の各実施の形態の計算手法をコンピュータ上において実行するための一連の命令を表現するソフトウエアを、単にプログラムと呼ぶ。プログラムは、コンピュータにより実行可能な任意の形式あるいはそのような形式に最終的に変換可能な任意の形式によって表現される。
【００５３】
本発明の各実施の形態のプログラムは、ハードウエア資源であるＧＰＵなどの演算手段が、ＯＳを介してあるいはＯＳを介することなく計算プログラムからの指令を受け、ハードウエア資源であるメインメモリや補助記憶装置などの記憶手段と協働して、ハードウエア資源である適当なバスなどを通じて演算処理を行うように構成される。つまり、本発明の各実施の形態の計算手法を実現するソフトウエアによる情報処理は、これらのハードウエア資源によって実現される。記憶手段あるいは記憶部は、任意の単位によって論理的に区分されている、コンピュータによって読み出すことができる情報記憶媒体の一部または全部またはそれらの組み合わせをいう。この記憶手段は、例えば、ＧＰＵ内のキャッシュメモリや、ＧＰＵと接続されたメモリや、ＧＰＵと適当なバスによって接続されたハードディスクドライブなどの不揮発性記憶媒体など、任意のハードウエア資源によって実現される。ここで、記憶手段は、ＧＰＵのアーキテクチャによって規定されるメモリ内の領域や、ＯＳが管理するファイルシステム上のファイルやフォルダ、同じコンピュータ内やネットワーク上のいずれかのコンピュータにあってアクセス可能なデータベースマネージメントシステム内のリストや、レコード、リレーショナルデータベースによって、相互に関連性がある複数のリストで管理されたレコードなど任意の形式によって実現され、論理的に他と区分され、少なくとも一時的に情報を識別可能に記憶または記録できる任意のものを含む。
【符号の説明】
【００５４】
１  ＣＰＵ
２  メインメモリ
３  ＧＰＵ
４  ＧＰＵメモリ
５  入力装置
６  表示装置
７  バス
【図１】
2012094040.tif 000012
【図２】
2012094040.tif 000013
【図３】
2012094040.tif 000014
【図４】
2012094040.tif 000015
【図５】
2012094040.tif 000016
【図６】
2012094040.tif 000017
【図７】
2012094040.tif 000018