特許ウォッチ

公開番号2025168621
公報種別公開特許公報(A)
公開日2025-11-11
出願番号2024082202
出願日2024-04-29
発明の名称2個以上の互いに素な正で奇数の単精度整数の積を法とする乗算合同法乱数の,孫子の定理の構造を利用する倍精度整数演算だけによる『再現可能性と移植可能性』を保つ生成と、実倍精度演算だけによる乱数出力計算と、が可能な高速MC乱数生成計算方式。
出願人個人
代理人
主分類G06F 7/58 20060101AFI20251104BHJP(計算;計数)
要約【課題】倍精度整数と倍精度実数の算法のみで、倍精度実数乱数を高速に生成する方法を提供する。
【解決手段】方法は、2個以上m個の共通素因数を持たない、即ち互いに素な、正の整数e1、e2、…、emの積である正の整数d=e1e2…emを法とし、『e1、e2、…、emは単精度整数』であると仮定し、dとは素な正の乗数z<d、dとは素な正の種n<dを用い、構成される乗算合同(MC)法乱数生成機構(d、z、n)において、与えられた初期値nから順次得られる漸化合同式による正でd未満の整数列{x^(t)|t=0、1、2、…}を得る
但し、
x⁽⁰⁾=mod(n、d)、x^(t+1)=mod(zx^(t)、d)、t=0、1、…、による一様独立乱数列{v^(t):=x^(t)/d|t=0、1、2、…}
【選択図】なし
特許請求の範囲【請求項１】
正の法ｄでのｍｏｄ関数ｍｏｄ（Ａ、ｄ）を０または正の整数Ａに対して
ｍｏｄ（Ａ、ｄ）＝『Ａをｄで割った整数の余りｒ、
ｒは０以上ｄ未満』
と定義し、２個以上ｍ個の共通素因数を持たない、即ち互いに素な、正の整数ｅ
１
、ｅ
２
、…、ｅ
ｍ
の積である正の整数ｄ＝ｅ
１
ｅ
２
…ｅ
ｍ
を法とし、『ｅ
１
、ｅ
２
、…、ｅ
ｍ
は単精度整数』であると仮定し、ｄとは素な正の乗数ｚ＜ｄ、ｄとは素な正の種ｎ＜ｄを用い、構成される乗算合同（ＭＣ）法乱数生成機構（ｄ、ｚ、ｎ）において、与えられた初期値ｎから順次得られる漸化合同式による正でｄ未満の整数列｛ｘ
（ｔ）
｜ｔ＝０、１、２、…｝、但し
ｘ
（０）
＝ｍｏｄ（ｎ、ｄ）、
ｘ
（ｔ＋１）
＝ｍｏｄ（ｚｘ
（ｔ）
、ｄ）、ｔ＝０、１、…、
による一様独立乱数列｛ｖ
（ｔ）
：＝ｘ
（ｔ）
／ｄ｜ｔ＝０、１、２、…｝を得る計算方法であって、孫子の定理の一つの形，部分法ｅ
１
、ｅ
２
、…、ｅ
ｍ
ぞれぞれだけに対する部分ＭＣ（ｅ
ｋ
、ｚ
ｋ
、ｎ
ｋ
）生成機構、即ち部分漸化合同式
ｘ
ｋ
（０）
＝ｍｏｄ（ｎ
ｋ
、ｅ
ｋ
）、ｘ
ｋ
（ｔ＋１）
＝ｍｏｄ（ｚ
ｋ
ｘ
ｋ
（ｔ）
、ｅ
ｋ
）、
ｋ＝１、２、…、ｍ、ｔ＝０、１、…、
ｚ
ｋ
：＝ｍｏｄ（ｚ，ｅ
ｋ
）、ｋ＝１、２、…、ｍ、
ｎ
ｋ
：＝ｍｏｄ（ｎ，ｅ
ｋ
）、ｋ＝１、２、…、ｍ、
が与える部分ＭＣ乱数列
｛ｘ
ｋ
（ｔ）
｜ｋ＝１、２、…ｍ、ｔ＝０、１、２、…｝
を用い、互いに素な部分法ｅ
１
、ｅ
２

発明の詳細な説明【技術分野】
【０００１】
互いに素な（共通素因数を持たない、最大公約数１の）正の奇数
ｍ個、ｍは２以上、の組ｅ
１
、ｅ
２
、…、ｅ
ｍ
の積ｄ、
ｄ＝ｅ
１
ｅ
２
…ｅ
ｍ
、（＄）
を法とする乗算合同ｍｕｌｔｉｐｌｉｃａｔｉｖｅｃｏｎｇｒｕｅｎｔｉａｌ（ＭＣ）法乱数生成機構（ｄ、ｚ、ｎ）、即ち法ｄ＞０、ｄとは素で『乗数』と呼ばれる正でｄ未満の整数ｚ、ｄとは素で『初期値或いは種或いはｓｅｅｄ』と呼ばれる正でｄ未満の整数ｎ、が与える乗算合同（ＭＣ）法乱数生成機構は、開区間（０、ｄ）内の整数の列｛ｚ
０
，ｚ
１
，ｚ
２
，…｝をｄを法とする再帰的合同関係
ｚ
０
＝ｍｏｄ（ｎ、ｄ）、
ｚ
ｋ＋１
＝ｍｏｄ（ｚ×ｚ
ｋ
、ｄ）、ｋ＝０、１、２、… （￥）
で生成する。よく知られたｍｏｄ関数の定義は、すぐ下に紛れなく定義し、乗法の記号×は今後省略する。ここに述べる発明は、ＭＣ法に必ず用いられる合同算法（￥）の形の法ｄの計算が孫子の定理の下で持つ重要な構造の指摘である。この利用方法の必要は乱数問題での素数の性質に関係する。シミュレーション問題が必要とする乱数列は見本過程だから、統計的一様独立性は『その様に見える』という事、仮説検定の問題である。本格的なＭＣ乱数生成方法はＦｉｓｈｍａｎとＭｏｏｒｅの有名な素数、整数単精度限界に近いメルセンヌ素数２
３１
－１の利用に始まる。彼らの行った検定が示す事は、『よい性能と見える』ＭＣ乱数生成機構の法を与え得る素数は整数単精度限界２
３２
に近い大きさが必要だという経験である。この詳細は（文献１）を見て頂きたい。大きさが小さいと検定で『適格』と看做さる乱数の発見は難しい。さらにシミュレーションの経験は、『倍精度実数計算』でなければ十分な精度は得られない、乱数周期もそれに相応しい巨大さでなければならない。これらに適合する乱数を構成する道は、２つ以上の単精度限界に近い素数の積を法とする事である。この、（￥）の形のＭＣ乱数にはシミュレーションの要請から単精度限界に近い２つ以上の部分法の積である法ｄが必要である事情、の説明から始める。
（文献１）ＮａｏｙａＮａｋａｚａｗａ／ＨｉｒｏｓｈｉＮａｋａｚａｗａ： “ＲａｎｄｏｍＮｕｍｂｅｒＧｅｎｅｒａｔｏｒｏｎＣｏｍｐｕｔｅｒｓ”、ＪＥＮＮＹＳＴＡＮＦＯＲＤ、校正刷段階のものを添付する。
まもなく出版予定。
続きを表示（約 8,300 文字）【０００２】
実倍精度ＭＣ乱数生成機構の技術設計の実際的可能性は
（１）１つの巨大な奇素数の法ｄ＝ｐの場合、
（２）２以上のｍ個の『互いに異なる大きい奇素数の積』ｄ＝ｐ
１
ｐ
２
…ｐ
ｍ
の場合、
の２つに大別される。計算実現の容易な場合として、上の（２）でｍ＝２の奇素数の積の法ｄの場合、単精度整数である奇素数ｐ
１
とｐ
２
が条件『積の２倍、２ｄ＝２ｐ
１
ｐ
２
、が倍精度用整数２
６４
限界未満を満たす場合』は特願２０２２－０２４２７７として出願され、審査中である。ここでは『積の個数ｍが２以上任意である場合』への一般化、
『互いに素な部分法ｅ
１
、ｅ
２
、…、ｅ
ｍ
には異なる単精度の
素数（従って互いに素）である以外の条件を付けない』
を可能とする新しい技術を述べる。まず（１）の場合が技術的に適さない理由を説明する。現在コンピュータシミュレーションで多く必要になるのは倍精度実数計算である。それに適合する倍精度実数乱数をＭＣ法で１つの素数の法ｄ＝ｐで与え、しかも乱数周期を通常のシミュレーションに不足しないほど巨大にするには、大きさが２
６４
程度、倍精度整数、の巨大な法ｄ＝ｐや乗数ｚが必要になる。しかしＭＣ乱数数列｛ｚ
０
、ｚ
１
、ｚ
２
、…｝の算出の各段階では法ｄの合同漸化式
ｚ
ｋ＋１
＝ｍｏｄ（ｚ×ｚ
ｋ
，ｄ）（＊）
の計算が必要で（乗法記号×は以後省略する）、ｚ
ｋ
が倍精度整数であるとこの掛け算部分はどうしても４倍精度実数による算術を経なければならず、それは大幅な（１０倍程度以上もの）計算時間の増大を招き、しかも現在（２０２４年）までに発見されたＭＣ乱数生成機構の優れた素数の法の最大周期は２
３４
程度で、実用シミュレーションには短かすぎ、効果に対して算出の負担は過大である。
【０００３】
前段落に指摘された困難を避ける道は、法ｄがｍ個の『単精度整数限界内』の大きさの奇素数ｍ個
（注１）
ｅ
１
、ｅ
２
、…、ｅ
ｍ
の積での形
ｄ＝ｅ
１
ｅ
２
…ｅ
ｍ
を採用する事であり、さらに『各部分法が与える単精度整数部分ＭＣ乱数生成機構に倍精度整数高速演算を可能にし、再現可能性も移植可能性も確保する』事である。前出願ではシミュレーション用の倍精度実数乱数出力が『ｍ＝２で、かつｄ＝ｅ
１
ｅ
２
の２倍が倍精度整数限界内に留まる』条件を満たす場合について、『倍精度整数演算だけを用いる方法』でこの『再現可能性と移植可能性とを持つ高速ＭＣ実倍精度乱数生成機構』を得た。ここでは倍精度『実数』演算も用いて再現可能性と移植可能性を確保する方法を『部分単精度整数の法の個数ｍが制限を持たない形』に提示したい。
【０００４】
まず正整数の法ｄと正の整数Ａを念頭において、以後多用されるｍｏｄ関数ｍｏｄ（Ａ、ｄ）の定義を紛れなく準備する。
（定義１）０以上の正の整数Ａと正整数の法ｄに対するｍｏｄ（Ａ、ｄ）は『Ａから法である正の整数ｄの整数商ｑ倍を除いて得られる０以上ｄ未満の整数の余りである』と定義する。割り算の等式で言えば、
Ａ＝ｑｄ＋ｒ、ｒは０以上ｄ未満の整数、
に基づいて、ｍｏｄ（Ａ、ｄ）＝ｒと定義する。（定義１終り）
実用コンピュータ計算では用いる事の可能な数全体に単精度、倍精度、４倍精度などの区別が厳しいが、議論の上ではいずれでも（定義１）を用いる事ができる。また計算諸言語では、変数Ａや法ｄが実数であっても（定義１）と同じ形に定められる定義でのｍｏｄ関数も既に実装されている。ここでは計算速度の考察が４倍精度計算を必要としない形の技術を必須とするが、それは計算プログラムの構成で述べるので、ｍｏｄ関数はすべて上の（定義１）による、とする。大切なのは現在の問題状況に見事に適合している孫子の定理の理解である。
（定理２．孫子の定理）整数の法ｄ＞０が正で互いに素な（最大公約数１の）整数の積ｄ＝ｅ
１
ｅ
２
…ｅ
ｍ
の形、ｍは２以上、だとする。『任意の正の整数Ａ』は
Ａ＝Ａ
１
Ｄ
１
Ｄ
１
－１
＋Ａ
２
Ｄ
２
Ｄ
２
－１
＋…＋Ａ
ｍ
Ｄ
ｍ
Ｄ
ｍ
－１
、
Ａ
ｋ
＝ｍｏｄ（Ａ、ｅ
ｋ
）、Ａ
ｋ
は０以上ｅ
ｋ
未満、
ｋ＝１、２、…、ｍ（＊）
の形に分解され、『法ｄで一意に』定まる。但しＤ
ｋ
＝ｄ／ｅ
ｋ
、Ｄ
ｋ
－１
は法ｅ
ｋ
でのＤ
ｋ
の正の逆数であって
ｍｏｄ（Ｄ
ｋ
Ｄ
ｋ
－１
、ｅ
ｊ
）＝δ
ｊｋ
（＃）
が成り立つ。ここでδ
ｊｋ
はクロネッカーのデルタでｊ＝ｋなら１、そうでなければ０である。
（証明）ｅ
１
、ｅ
２
、…、ｅ
ｍ
が互いに素、すなわち最大公約数が１だから、ｅ
ｋ
を考えるとｅ
ｋ
と残りの部分法の積を表すＤ
ｋ
＝ｅ／ｅ
ｋ
とは互いに素である。故にユークリッドの互除法の結論、整数ＭとＮがあって最大公約数１を次の『等式』で与えられ、
ＭＤ
ｋ
＋Ｎｅ
【０００５】
互いに素な正の部分法ｍ個の積ｄ＝ｅ
１
ｅ
２
…ｅ
ｍ
を法とするＭＣ乱数の構造を考える。部分法の表記のために（ｄ、ｚ、ｎ）が生成するＭＣ乱数列を一般に
｛ｚ
（ｔ）
＝ｍｏｄ（ｎｚ
ｔ
、ｄ）｜ｔ＝０，１、２、…、｝
と表す；ｊは見難いので、一般的に時刻と考えてｔを用いる。これは或いは漸化式形式で
ｚ
（０）
＝ｍｏｄ（ｎ、ｄ）、
ｚ
（ｔ＋１）
＝ｍｏｄ（ｚｚ
（ｔ）
、ｄ）、ｔ＝０、１、２、…（＃１
＊
）
とも記す。法ｄが互いに素な正の部分法による形ｄ＝ｅ
１
ｅ
２
…ｅ
ｍ
を持つとき、
部分乗数、部分出発値を
部分乗数：｛ｚ
ｋ
＝ｍｏｄ（ｚ、ｅ
ｋ
）、ｋ＝１、２、…、ｍ｝、
部分初期値：｛ｎ
ｋ
＝ｍｏｄ（ｎ、ｅ
ｋ
）、ｋ＝１、２、…、ｍ｝。
と定義する。これらが定める『部分ｅ
ｋ
のＭＣ乱数生成機構』を
｛（ｅ
ｋ
、ｚ
ｋ
、ｎ
ｋ
）｜ｋ＝１、２、…、ｍ｝
と記し、それらが生成する『ＭＣ部分整数乱数列』を
｛ｍｏｄ（ｎ
ｋ
（ｚ
ｋ
）
ｔ
、ｅ
ｋ
）、ｔ＝０，１、２、…、
ｋ＝１、２、…、ｍ｝（＃２）
或いは漸化式形式で
ｚ
ｋ
（０）
＝ｍｏｄ（ｎ
ｋ
、ｅ
ｋ
）、
ｚ
ｋ
（ｔ＋１）
＝ｍｏｄ（ｚ
ｋ
ｚ
ｋ
（ｔ）
、ｅ
ｋ
）、
ｔ＝０、１、２、…、ｋ＝１、２、…、ｍ（＃２
【０００７】
（ｄ、ｚ、ｎ）ＭＣ乱数の生成の技術として第一の鍵は、法ｄ＝ｅ
１
ｅ
２
…ｅ
ｍ
の奇数の部分法ｅ
１
、ｅ
２
、…、ｅ
ｍ
を整数単精度で互いに素に選ぶ事である。これによって『部分ＭＣ乱数生成機構
｛（ｅ
ｋ
、ｚ
ｋ
、ｎ
ｋ
）｜ｋ＝１，…、ｍ｝
の個数ｍには制限がなくなる』。この場合部分乱数の生成は『倍精度整数演算だけ』を用いて『再現可能かつ移植可能』に行う事ができる。それらの孫子の定理による（ｄ、ｚ、ｎ）ＭＣ乱数出力への合成は、再現可能性を考える必要がないので倍精度実数として行う事によって、幾つかの困難を減らし、しかも高速に実現される。これが現在の出願の骨子である。但し法ｄの大きさにはやはり制限が残り、個数ｍも現実的に考えて２ないし３が可能限界であろう。例えば４倍精度実数として計算するにしても、ｄは２
１２８
を越える使い方はできないし、限界に近い使い方では計算中に大きさが限界に達する不安は残る。実数でも規格を目いっぱい使う事は困難である。ただ漸化式のｚ×ｚ
（ｔ）
を小さい部分法で倍精度整数内で行う事は確かに可能性を広げている。次の段落でその実際プログラムを、ｍ＝３の場合で例示し、どの計算段階が新しい工夫であるかを明確に示す。
【０００８】
異なる３個の素数ｅ
１
、ｅ
２
、ｅ
３
を用いる法ｄ＝ｅ
１
ｅ
２
ｅ
３
による乱数ＭＣ計算プログラムをＦＯＲＴＲＡＮ形式で次ページ図１として掲げる。但しこのＭＣ乱数の統計性能は優れたものではないので御了解頂きたい。現在の視点から最も分りやすい表記を目指し、プログラムの中に変数の単、倍、４倍などの大きさの種別を詳しく記す。まず生成速度を実測する。
図１．ｍ＝３個の素数の積の法ｄ＝ｅ
１
ｅ
２
ｅ
３
で１０００万個の
乱数を発生するメインプログラム
ｐｒｏｇｒａｍｍａｉｎ
ｉｍｐｌｉｃｉｔｉｎｔｅｇｅｒ
＊
８（ｉ－ｐ），ｒｅａｌ
＊
８（ａ－ｈ，ｒ－ｚ）
ｃｏｍｍｏｎｐ１，ｐ２，ｐ３，ｄ，ｉｚ１，ｉｚ２，ｉｚ３，ｄ１，ｄ２，ｄ３，ｉｎｖｄ１，ｉｎｖｄ２，ｉｎｖｄ３
ｐ１＝１３４２６５０２３！ｐ１はおよそ２
∧
（２７）の素数単精度整数
ｐ２＝１０００９１９！ｐ２はおよそ２
∧
（２０）の素数単精度整数
ｐ３＝２５７！ｐ３はおよそ２
∧
（８）の素数単精度整数
ｎｕｍ＝１０００００００！発生させる乱数の個数を設定
！ｄ：＝ｐ１
＊
ｐ２
＊
ｐ３！法ｄはおよそ２
∧
（５４．９４）の整数
ｄ１＝ｐ２
＊
ｐ３！Ｄ１：＝ｄ／ｐ１（倍精度実数）
ｄ２＝ｐ３
＊
ｐ１！Ｄ２：＝ｄ／ｐ２（倍精度実数）
ｄ３＝ｐ１
＊
ｐ２！Ｄ３：＝ｄ／ｐ３（倍精度実数）
ｄ＝ｐ１
＊
ｐ２
＊
ｐ３！ｄ＝ｐ１
＊
ｐ２
＊
ｐ３（倍精度実数）
ｉｚ１＝１９０６１２５２！ｐ１の原始根ｉｚ１
ｉｚ２＝６４５４５！ｐ２の原始根ｉｚ２
ｉｚ３＝５３！ｐ３の原始根ｉｚ３
ｉｎｖｄ１＝３８７８９３６７！ｉｎｖｄ１：＝ｍｏｄ（Ｄ１
∧
（－１），ｐ１）
ｉｎｖｄ２＝３８８４９８！ｉｎｖｄ２：＝ｍｏｄ（Ｄ２
∧
（－１），ｐ２）
ｉｎｖｄ３＝８３！ｉｎｖｄ３：＝ｍｏｄ（Ｄ３
∧
（－１），ｐ３）
ｎ１＝１！法ｐ１での部分ＭＣ初期値
ｎ２＝１！法ｐ２での部分ＭＣ初期値
【０００９】
前段落［０００８］と同じ法ｄ＝ｅ
１
ｅ
２
ｅ
３
で、部分法に分解することなく４倍精度実数演算で直接乗算を行う原始的な方法のプログラムを次の図３に掲げる。
図３．法ｄでの４倍精度実数乗算による方法のプログラム
ｐｒｏｇｒａｍｍａｉｎ
ｉｍｐｌｉｃｉｔｉｎｔｅｒｇｅｒ
＊
８（ｉ－ｐ），ｒｅａｌ
＊
８（ａ－ｈ，ｒ－ｚ），ｒｅａｌ
＊
１６（ｑ）
ｎｕｍ＝１０００００００！発生させる乱数の個数を設定
ｐ１＝１３４２６５０２３！ｐ１はおよそ２
∧
（２７）の素数単精度整数
ｐ２＝１０００９１９！ｐ２はおよそ２
∧
（２０）の素数単精度整数
ｐ３＝２５７！ｐ３はおよそ２
∧
（８）の素数単精度整数
！ｄ：＝ｐ１
＊
ｐ２
＊
ｐ３！法ｄはおよそ２
∧
（５４．９４）の整数
ｑｄ＝ｐ１
＊
ｐ２
＊
ｐ３！ｑｄ＝ｐ１
＊
ｐ２
＊
ｐ３（４倍精度実数）
ｑｚ＝１９５０９０９６０５８４２８９９ｄ０！ｑｚ：乗数（４倍精度実数）
ｑｍｚ＝１！法ｄでのＭＣ乱数列（４倍精度実数）、初期値は１
！乱数をｎｕｍ＝１０００万回呼び発生出力させる
ｄｏｉ＝１，ｎｕｍ
ｑｍｚ＝ｍｏｄ（ｑｍｚ
＊
ｑｚ，ｑｄ）！法ｄでのＭＣ乱数列（４倍精度実数）
ｒａｎｄ＝ｑｍｚ／ｑｄ！出力
ｅｎｄｄｏ
ｅｎｄ
（図３．終り）
このプログラムで上の図１および図２と同じ出力結果が得られることが確認できる。図３は法ｄでの乗算を単純に繰り返すだけなので、プログラム構造としては一見簡単である。しかし、ＣＰＵ時間はおよそ２．７秒／１０００万個で、図２のサブルーチンと比べ１０倍近い時間を要する。このことは、本発明で提示される、孫子の定理の構造を利用して倍精度整数と倍精度実数のみによる生成方式を実用にすべき一つの重大な理由である。さらに、本発明の生成方式は、２
６４
を超える合成数の法ｄでのＭＣ乱数列出力を『倍精度実数演算のみ』で実現しうるものでもあり、単純な４倍精度乗算を用いるだけの方式と比べて、可能性と発展性が大きい。
【００１０】
乱数をシミュレーションやゲームなどのソフトウェアで実際に利用することを念頭に置いて、一つのプログラム例を挙げる。下の図４は、原始的で簡単ながら指定した個数だけ乱数を発生させるメインプログラムである。
図４．乱数をｎｕｍ個だけ発生させるメインプログラム
ｐｒｏｇｒａｍｍａｉｎ
ｉｍｐｌｉｃｉｔｉｎｔｅｇｅｒ
＊
８（ｉ－ｐ），ｒｅａｌ
＊
８（ａ－ｈ，ｒ－ｚ）
ｃｏｍｍｏｎｐ１，ｐ２，ｐ３，ｄ，ｉｚ１，ｉｚ２，ｉｚ３，ｄ１，ｄ２，ｄ３，ｉｎｖｄ１，ｉｎｖｄ２，ｉｎｖｄ３
ｎｕｍ＝１０！発生させる乱数の個数
ｐ１＝１３４２６５０２３！ｐ１はおよそ２
∧
（２７）の素数単精度整数
ｐ２＝１０００９１９！ｐ２はおよそ２
∧
（２０）の素数単精度整数
ｐ３＝２５７！ｐ３はおよそ２
∧
（８）の素数単精度整数
ｄ＝ｐ１
＊
ｐ２
＊
ｐ３！法ｄはおよそ２
∧
（５４．９４）の整数
ｄ１＝ｐ２
＊
ｐ３！Ｄ１：＝ｄ／ｐ１（倍精度実数）
ｄ２＝ｐ３
＊
ｐ１！Ｄ２：＝ｄ／ｐ２（倍精度実数）
ｄ３＝ｐ１
＊
ｐ２！Ｄ３：＝ｄ／ｐ３（倍精度実数）
ｉｚ１＝１９０６１２５２！ｐ１の原始根
ｉｚ２＝６４５４５！ｐ２の原始根
ｉｚ３＝５３！ｐ３の原始根
ｉｎｖｄ１＝３８７８９３６７！ｉｎｖｄ１：＝ｍｏｄ（Ｄ１
∧
（－１），ｉｐ１）
ｉｎｖｄ２＝３８８４９８！ｉｎｖｄ２：＝ｍｏｄ（Ｄ２
∧
（－１），ｉｐ２）
ｉｎｖｄ３＝８３！ｉｎｖｄ３：＝ｍｏｄ（Ｄ３
∧
（－１），ｉｐ３）
ｎ１＝１０！法ｐ１での部分ＭＣ初期値
ｎ２＝２２！法ｐ２での部分ＭＣ初期値
ｎ３＝４１！法ｐ３での部分ＭＣ初期値
ｍｚ１＝ｍｏｄ（ｎ１，ｐ１）！法ｐ１での部分ＭＣ乱数列の初期設定
ｍｚ２＝ｍｏｄ（ｎ２，ｐ２）！法ｐ２での部分ＭＣ乱数列の初期設定
ｍｚ３＝ｍｏｄ（ｎ３，ｐ３）！法ｐ３での部分ＭＣ乱数列の初期設定
！倍精度整数ｍｚ１，ｍｚ２，ｍｚ３は乱数発生の度に値が更新される
！ｒａｎｄｇｅｎｅｒａｔｅサブルーチンを呼び、乱数をｎｕｍ個発生させる
！引数ｍｚ１、ｍｚ２、ｍｚ３にはｎｕｍ番目の各部分法ＭＣ乱数（倍精度整数）が保存され
！次の乱数生成に引き継ぎが可能
ｃａｌｌｒａｎｄｇｅｎｅｒａｔｅ（ｍｚ１，ｍｚ２，ｍｚ３，ｎｕｍ）！ｒａｎｄｇｅｎｅｒａｔｅサブルーチン（図５）
ｅｎｄ
（図４．終り）
基本構造としては、前段落の図２のメインプログラムと同様である。ただし、シミュレーションやゲームなどのプログラム開発を想定して、『必要な個数だけ』の乱数生成をサブルーチンとして随時呼び出せるように手を加えている。乱数を任意のｎｕｍ個だけ発生させるサブルーチンを図５に掲げる。
図５．乱数をｎｕｍ個だけ発生させるサブルーチン
ｓｕｂｒｏｕｔｉｎｅｒａｎｄｇｅｎｅｒａｔｅ（ｍｚ１，ｍｚ２，ｍｚ３，ｎｕｍ）
ｉｍｐｌｉｃｉｔｉｎｔｅｇｅｒ
＊

関連特許