0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 1 year has passed since last update.

無限に偶数(奇数)が連続する疑似乱数生成器

Last updated at Posted at 2022-05-06

はじめに

疑似乱数生成器が生成する数列のうち, その生成器の状態空間全体に比べてはるかに狭い区間で, 偶数と奇数の分布が期待したとおりになっていると思わない方がよいです. 偶奇を希望する分布で得たいなら, 疑似乱数生成器で値を生成した後で加工するとよいと思います.
もちろん, 例えば2値を均等な確率で得たい場合, 多くの標準ライブラリの実装では1ビットに依存するため難しい話ですが, 開発・運営コストとリスクを考えて選択した方がよいと思います.
追記, 疑似乱数生成器の質を議論することが目的ではありません, 特定のビット(最下位ビット)だけで目的の分布を得ようとしない方がよいという意見です.
追記, よい分布とは何か?には答えられないです. そのうえで, 特定のビットに依存した設計をしないほうがよいと言っています.

無限に奇数(偶数)が連続する疑似乱数生成器

簡単な"IBM RANDU"です. 初期値に応じて偶数か奇数しかでない疑似乱数生成器です.
今だ教育現場で話題にされるのでしょうか, 数年前に環境同梱の疑似乱数生成器に対して, 「使ってもいいですか?」と聞かれたことがあります. 学校か研修か以前の現場か, トラウマレベルで怒られたのかな?
下手に自作しようとすると余計に痛い目をみそうですが, うっかり簡単に自作できてしまうC++を使用します.

#include <cstdint>
#include <random>
#include <iostream>

using randu_ibm = std::linear_congruential_engine<std::uint_fast32_t, 65539, 0, 21474836478>;

int main(void)
{
    std::random_device device;
    randu_ibm engine(device()|0x01U);
    for(;;){
        uint32_t x = engine();
        if(0 == (x&0x01U)){
            break;
        }
        std::cout << x << std::endl;
    }
    return 0;
}

線形合同法の性質を使っているだけです.
JavaのRandomでも少し無理矢理ですが合法に改変できそうです, 暇ができたら.

624回奇数が連続するメルセンヌ・ツイスタ

生成された乱数列から内部状態を逆算します. 逆算できないときもあるため, 成功していなければカウントを打ち切っています. 流石に偏り過ぎなのでメルセンヌ・ツイスタが実際にこの状態にはならないと思われます. 逆算がどこか間違っていると思います.
これほど長く連続することはないにしても, 検査プログラムが正しければ40回程度連続する区間があることは確認しています.

Revert Mersenne Twister
#include <random>
#include <iostream>

class mt_random
{
public:
    static constexpr uint32_t N = 624;
    static constexpr uint32_t M = 397;
    static constexpr uint32_t MATRIX_A = 0x9908B0DFUL;
    static constexpr uint32_t UPPER_MASK = 0x80000000UL;
    static constexpr uint32_t LOWER_MASK = 0x7FFFFFFFUL;

    mt_random(uint32_t s=5489UL)
    {
        seed(s);
    }

    void seed(uint32_t s)
    {
        mt[0] = s & 0xFFFFFFFFUL;
        for(uint32_t i=1; i<N; ++i){
            mt[i] = (1812433253UL * (mt[i-1] ^ (mt[i-1] >> 30)) + i);
            mt[i] &= 0xFFFFFFFFUL;
        }
        mti = 0;
    }

    void set(const uint32_t history[N])
    {
        for(uint32_t i=0; i<N; ++i){
            mt[i] = inv_tempering(history[i]);
        }
    }

    void revert(const uint32_t history[N])
    {
        static const uint32_t mag0[2] = {0UL, MATRIX_A};
        set(history);
        for(int32_t i=static_cast<int32_t>(N-1); 0<=i; --i){
            uint32_t t = mt[(i+N-1)%N] ^ mt[(i+M-1)%N];
            uint32_t b = t ^ mag0[t>>31];
            uint32_t h = mt[(i+N)%N] ^ mt[(i+M)%N];
            mt[i] = ((h<<1) & UPPER_MASK) ^ ((b<<1) & LOWER_MASK) ^ (t>>31);
        }
        mti = 0;
    }

    uint32_t state(uint32_t index) const
    {
        return mt[index];
    }

    uint32_t operator()()
    {
        static const uint32_t mag0[2] = {0UL, MATRIX_A};
        uint32_t y = (mt[mti] & UPPER_MASK) | (mt[(mti+1)%N] & LOWER_MASK);
        mt[mti] = mt[(mti+M)%N] ^ (y>>1) ^ mag0[y&0x01U];
        uint32_t x = tempering(mt[mti]);
        mti = (mti+1)%N;
        return x;
    }

private:
    static uint32_t tempering(uint32_t x)
    {
        x ^= (x >> 11);
        x ^= (x << 7) & 0x9D2C5680UL;
        x ^= (x << 15) & 0xEFC60000UL;
        x ^= (x >> 18);
        return x;
    }

    static uint32_t inv_tempering(uint32_t x)
    {
        x ^= (x>>18);
        x ^= (x<<15) & 0xEFC60000UL;
        x ^= ((x <<  7) & 0x9D2C5680UL)
            ^ ((x << 14) & 0x94284000UL)
            ^ ((x << 21) & 0x14200000UL)
            ^ ((x << 28) & 0x10000000UL);
        x ^= (x>>11) ^ (x>>22);
        return x;
    }

    uint32_t mt[N];
    uint32_t mti = 0;
};

int main(void)
{
    std::random_device device;
    uint32_t seed = device();
    mt_random mt(seed);

    uint32_t history[mt_random::N];
    for(uint32_t i=0; i<mt_random::N; ++i){
        history[i] = mt() | 0x01UL;
    }
    mt.revert(history);
    uint32_t count = 0;
    for(uint32_t i=0; i<mt_random::N; ++i){
        uint32_t x = mt();
        if(x != history[i]){
            break;
        }
        ++count;
    }
    std::cout << "seed: " << seed << ", odd: " << count << std::endl;
    mt.revert(history);
    count = 0;
    for(uint32_t i=0; i<mt_random::N; ++i){
        std::cout << mt.state(i) << ",";
        ++count;
        if(0 == (count%10)){
            std::cout << std::endl;
        }
    }
    return 0;
}

次のように内部状態を設定すると624回奇数が連続するのではないでしょうか.

内部状態
1415792764,1591143757,2618347592,2576834658,3251543905,157795228,732642334,1141350247,1873486419,1400093550,
2134251629,1874313302,588213766,3253470429,820555899,3914459424,4250032066,456799647,1346811743,206928669,
4090247551,1014468791,175270048,1210337436,2627398968,1757056325,3167118687,2819605196,3113421382,2618228550,
566700899,3175908248,648217569,2190143718,2176026500,4019012103,4192646414,1894700821,2109381554,1716981773,
4266905327,160265715,2713796716,2826962181,2555584721,2816192425,2751467293,2173598437,1439310795,3260243917,
3514466655,43763434,2370745478,2118376521,3341378490,1121184195,1657770356,224884163,1688388109,4242433356,
1364156714,485684305,829186971,599670829,2397546754,4191595841,696001758,2783154569,1297652180,3810481105,
3989622969,3659460749,2218056749,1389882366,3092726522,4056876224,3791689744,2820877004,996811156,3081239073,
1530517431,3529145249,2160387619,662990209,3764403513,726807095,2274014875,2984366474,4260624687,3399560356,
3772953472,1874705995,2468192404,4051123005,2513636459,2424545118,3512876454,986418270,1127425119,1283592308,
1019678160,1319061290,4138054856,3837400799,512954812,3691786809,322010265,3837598032,2788451729,1745665644,
1331914059,2972105612,3683880905,2923452462,3235806540,1777490273,2084977325,3884694467,2393781570,61570023,
1876141612,2896575795,4233946495,3310882450,2915356963,1404781454,1776239343,1546407927,2767213378,1476085570,
123579476,839241349,1661404032,4076422685,3163239373,3229852292,2763311196,1269856308,4112925768,1024677928,
3677006562,363037389,1658778288,2695723364,4110410104,3472592287,237470382,2103898196,1747919446,2668079352,
3141869267,2620131071,515520208,1722903071,746192625,1379276061,2546591037,63346793,3544148610,441891017,
3501465788,1428011082,2188849336,1531282774,2900819663,2233131780,1402073319,321353465,2329711209,2755572455,
3267638969,2355362882,3130032261,670230484,2832667985,4138384226,4230804282,1139408139,2877620327,4166634455,
161252333,2503061007,302523335,1242271861,589285348,1572079715,2938950041,3691124999,2839462235,1881209598,
1051590372,1414635382,291276093,2444624669,2007982720,2692192887,2261129051,457299001,3139887755,2021166105,
3694557401,1852100953,3352109692,1074792784,3729941422,4014585027,2513180645,219460441,3495557467,2211175566,
1462122782,2757817225,2855142435,1162518608,4180568072,3520690562,1770131856,333659351,2459829129,931936787,
516283329,1961950466,1234136589,2571617471,2770560984,2819141831,1983189748,2730956013,692595247,3272623586,
755184965,720141234,3172640648,2884826203,369953031,3228485180,1058014509,3021383996,3969315018,2145882428,
398521089,2989648214,3116960388,1434872791,3902280290,2990056360,626534556,2241495279,979582593,2552850040,
1971501340,2641268332,3013450224,1363304271,4038230474,3596666444,511149355,3822847082,3041710491,1878075848,
1811281899,2572780485,2534086685,1579410088,4149696295,2084883781,2597039992,1469744685,1555887610,4219225978,
2597485913,3293162702,466130571,276875496,3975144712,15645296,1497495254,1735306716,4080507713,3985644979,
2102799632,1438550618,1042878567,1724618942,1145731959,3588904249,1257144052,3958256208,1193106757,39530741,
1302652761,2485082439,3034412279,1839809230,544582211,2252667804,1831916670,3106949078,421941594,4148042113,
1789390574,4005573612,1469413759,640406647,1349687689,2753059801,2915135897,2455422554,3730659820,4185809792,
787994487,3733681922,4044723621,2462546991,1489775765,4101295007,1946835034,4191966858,2277048281,1235367070,
3935276129,2339427819,2338133807,1801206532,1880880957,1690250851,2616998756,3159306123,3827931847,2382307749,
3452187297,4131818299,1490937444,3097162790,3196132880,2798440285,1292250275,2856670419,3804147283,3125183197,
3790862927,2655546705,371001223,1018545458,97539562,79412515,1761414589,3731573779,404324262,2900964315,
4222588667,81348919,1718178579,2193940713,3278528287,1872855919,1515466840,3566465404,957628289,2747568606,
552198076,3401817999,476307366,3416331241,3990855468,2022055678,3762109083,773455457,3662597927,4238066861,
827727756,1904036169,333052750,2476346945,1430158565,2365440761,2316335991,3784768526,1142360487,3053021695,
3876206413,696654033,911551469,332636818,1306685714,3580106086,4055007067,1556012713,570979686,297230673,
3726651765,3743490695,558188342,3991601111,258696238,3007129009,3698664889,2849133909,3108811379,264606860,
1779631604,2038334772,2762402939,3992160194,2825689793,3713857930,2511555633,2528094862,2657217134,4232954791,
2547724008,364786815,2303675535,1979407832,777267515,1440864486,561984451,3631431938,1746487091,3333599057,
3910798180,2339711450,3009534168,144521641,1237760967,1037850975,2795706251,1458246186,464351913,745123674,
3623836208,3664306894,3692607491,3655643723,3135138184,32230073,256148651,1441813950,996799607,132260778,
4071654488,36242182,1236459658,405068184,1378103988,2786683316,2965125728,1979346799,91346470,4102935743,
2746548860,3047260639,3221621607,1863144147,3977681202,83379332,1267939556,139637948,3963183716,3236151270,
3479121343,1344973049,3020514214,3620183057,1140618434,1351074388,118187605,2849069006,951257834,329515957,
4241681323,2302625067,2609175199,2863875786,2568222402,3542948763,1438137556,3495196532,1126471333,1209677869,
2476455740,3157359447,4124649323,148880505,2462089993,3451361488,4188803097,94788751,504630611,3658164026,
1352035937,3507154379,856152532,2587094787,1765911935,41506457,2086038997,365104669,2393476559,225612398,
127051098,3565871589,3598797754,4036485876,2636364907,2446100582,3308248316,1716487022,3026400607,895567528,
3431744468,2736600490,1975013090,2333840,616945362,4186352975,1495388680,138418476,4004705776,3799383372,
3192207667,1497484126,274994565,3114511492,886325624,4014868941,2690245210,2377903957,526480219,3304328652,
2332360643,3495607210,137162697,540057454,1056829314,602542689,4220691151,268625341,2230769573,3898373630,
1036876415,577901470,4103587046,2119682176,921754015,3096477006,3177570858,1789429015,3916404280,4031534078,
3173523299,4244640404,2087539123,173010779,1816439269,2564914987,1302424819,920647221,880686979,3378190308,
1275858771,1484154393,1930386852,3383771468,585270477,1514414000,3152655769,2969307990,4063334486,3260741968,
27314015,385681444,298851480,2128847037,1970951926,3747265418,2989731957,3360188058,1038341133,3883600126,
642052885,2694426396,1861331472,441030244,2499953750,1854235391,4181075572,2510550180,2690831448,663716775,
2134962085,1793670045,1349225294,3139397762,1569994048,2791892224,1838147364,712769774,2381168188,3037312863,
1778063505,1282982833,1982331478,2492734484,3461633176,3392145386,3704183769,2842984647,3577953561,3803483697,
304587906,2558521919,2587229437,2023074662,4023601307,2624816473,2507035215,3131879370,1051419730,3905483095,
561487498,3851149645,2545063080,3057942919,

まとめ

疑似乱数生成器が生成する値の偶奇に, 期待する分布があると思わない方がいいと思います.

参照

[1] 正逆双方向に計算できるメルセンヌ・ツイスタ

0
0
4

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
0
0

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?