Na pierwszy ogień – Randomizery

Przy okazji odgruzowywania repozytoriów z Flaksatorem – 2-3 wersje ;-), ta na Githubie może się okazać, że nie jest najlepszą… (Nic to, po drodze się naprawi.) – odgrzebałem rozmaite wersje tzw. “wspólnych bibliotek”, do których w(y)rzucałem przydatne we współdzieleniu klasy i komponenty.

W pierwszej kolejności zająłem się ważnym komponentem generującym losowe wartości. Zwłaszcza, że w pewnych częściach przydatny jest / będzie w kilku moich innych pet-projektach. Ktoś powie – a na grzyba, przecież jest klasa Random?!
I w zasadzie słusznie. Problemy są trzy:

  • Podobno nie jest zbyt szybka
  • Podobno jej rozkład pozostawia wiele do życzenia
  • Nie da się spersystować / serializować jej stanu

Co do ostatniego – to BinaryFormatter powinien sobie z tym łatwo poradzić (Random implementuje interfejs ISerializable, mimo, że publicznych pól nie posiada…) – o tym w kolejnym poście.

W międzyczasie – myślę – warto spróbować zweryfikować dwie pierwsze tezy. Najlepiej móc porównać z czym innym. Idealnym kandydatem wydaje się szeroko rozpowszechniony algorytm Mersenne Twister, którego mam kilka wersji:

W zasadzie dwie pierwsze implementacje są niemal identyczne – i w zasadzie w moim 64bitowym systemie dają identyczne rezultaty. Trzecia – wyniki z tego samego ziarna daje już inne.
Zanim zacznę rozkminiać jak to jest z tą bitowością tych implementacji postanowiłem sobie porównać…

(Uwaga – biedny kod testujący…)

[Test]
public void CompareSequencesOfMersennes()
{
    Random rnd = new Random();
    for (int i = 0; i < 10; i++)
    {
        int seed = rnd.Next();
        var m1 = new MT19937();
        m1.init_genrand((ulong)seed);
        var m2 = new RandomMT((ulong)seed);
        var m3 = new SharpDevs.Randomization.MersenneTwister(seed);
        var m4 = new Random(seed);

        for (int j = 0; j < 10; j++)
        {
            var r1 = m1.genrand_int32();
            var r2 = m2.RandomInt();
            var r3 = m3.GetNextULong();
            var r4 = m4.Next();

            Console.WriteLine("{0}\t{1}\t{2}\t{3}", r1, r2, r3, r4);

            //Assert.AreEqual(r1, r2);
            //Assert.AreEqual(r2, r3);
        }

        Console.WriteLine();
    }
}

I wyniki (dwa pierwsze sety):

2658104388	1306749176	2658104388	1140326036
4271057295	3682971234	4271057295	1810704752
2445844545	1710876972	2445844545	1713401550
2995924299	563814871	2995924299	201792233
2670607760	3380877015	2670607760	186337359
3177152157	3093296799	3177152157	1465407615
3047290646	2006867027	3047290646	429790269
3928896115	4237850699	3928896115	720743159
167834357	3811397742	167834357	651523369
1853511702	428212234	1853511702	1657147808

1993119808	2224621582	1993119808	1362312584
1841636498	4085140072	1841636498	1189751224
1387250596	1132592719	1387250596	1453857008
1163256458	1057422393	1163256458	1502410392
493414620	499985414	493414620	2075395322
570257101	2608426189	570257101	131272236
2017910186	2871847751	2017910186	2100987751
1355161007	2389705801	1355161007	1056760113
3399534222	1066188073	3399534222	2080094558
3204281763	3048615178	3204281763	1990654215

Nuda, Panie.
Na podstawie tego na pewno nie wywnioskuję jaki jest rozkład. Muszę rzutować wyniki na mniejszy zakres. Teraz się zacznie. Klasy nie posiadają jednolitego interfejsu, by jakoś ładnie ten test napisać. Trochę szkoda mi czasu na dekorowanie implementacji, które potem odrzucę. Pozostaje napisać brzydki kod testujący 😉 I od razu zbadam czas wykonania.
Oraz… Uch… statystyka… Badanie rozkładu zmiennej losowej. Nic nie pamiętam. Chyba czas spytać Googla. Choć z drugiej strony… Na razie roboczo przyjmę sobie oczekiwane prawdopodobieństwo jako równe i policzę coś co wydaje mi się sensownym średnim odchyleniem standardowym…

[Test]
public void CompareDistributionOfMersennes()
{
    Random rnd = new Random();
    int baseSeed = rnd.Next();

    var m1 = new MT19937();
    m1.init_genrand((ulong)baseSeed);
    var m2 = new RandomMT((ulong)baseSeed);
    var m3 = new SharpDevs.Randomization.MersenneTwister(baseSeed);
    var m4 = new Random(baseSeed);
    var maxRange = 100;
    var iterations = 100000;

    var dist1 = CalculateDistribution(() => m1.RandomRange(0, maxRange - 1), maxRange, iterations);
    var dist2 = CalculateDistribution(() => m2.RandomRange(0, maxRange - 1), maxRange, iterations);
    var dist3 = CalculateDistribution(() => m3.GetNext(0, maxRange), maxRange, iterations);
    var dist4 = CalculateDistribution(() => m4.Next(0, maxRange), maxRange, iterations);

    var distribution1 = AnalyzeDistribution(dist1, iterations);
    var distribution2 = AnalyzeDistribution(dist2, iterations);
    var distribution3 = AnalyzeDistribution(dist3, iterations);
    var distribution4 = AnalyzeDistribution(dist4, iterations);

    Console.WriteLine("{0}\t{1}\t{2}\t{3}",
        distribution1, distribution2, distribution3, distribution4);
}

private decimal AnalyzeDistribution(int[] distribution, int iterations)
{
    int expectedScores = distribution.Length;
    decimal expectedScore = (decimal)iterations/ expectedScores;
    var deviation = distribution.Average(i => Math.Abs(i - expectedScore));
    return deviation;
}

private int[] CalculateDistribution(Func<int> generator, int maxRange, int iterations)
{
    Stopwatch sw = Stopwatch.StartNew();

    int[] resultScore = new int[maxRange];

    for (int i = 0; i < iterations; i++)
    {
        int score = generator.Invoke();
        resultScore[score]++;
    }

    sw.Stop();
    Console.WriteLine(@"Czas wykonania: " + sw.ElapsedMilliseconds + @"ms");

    return resultScore;
}

Wyniki są bezlitosne, ale dość zaskakujące:

Czas wykonania: 11ms
Czas wykonania: 9ms
Czas wykonania: 14ms
Czas wykonania: 2ms
25,44	248,88	 25,98 	25,46
Przy oczekiwanym napełnieniu: 1000
Co daje procentowo:
2,83400  	34,4500  	2,41200	   2,44800

Po zwiększeniu maxRange do 1000 a iterations do 10 000 000 jeszcze bardziej druzgocące:

Czas wykonania: 941ms
Czas wykonania: 867ms
Czas wykonania: 1297ms
Czas wykonania: 220ms
81,448    10819,288	79,22	81,946
Przy oczekiwanym napełnieniu: 10000
Co daje procentowo:
0,7937600	107,1931200	0,80400	   0,7739400

Wydaje się, że mogę swobodnie zignorować istnienie implementacji MT i pozostać przy klasie Random – rozkład jest pomijalnie gorszy niż dwóch implementacji MT, za to znacznie szybszy.
O ile oczywiście nie rąbnąłem się gdzieś w kodzie (ale jak pobieżnie przeglądałem zawartość tabel dist* to by się pokrywały z wynikami…). W takim przypadku, dalszą analizę klas MT póki co sobie odpuszczę.

Do następnego.

PS. Swoją drogą – “moja” implementacja daje najlepsze wyniki ale jest najwolniejsza… Grrr.
EDIT: Podczas poszukiwania zagubionego artykuły z MT, w komentarzach wpadłem też na tę bibliotekę. Czasy wywołania prostego Next() wyglądają wielce obiecująco. Pytanie jak się zachowa Next(min, max)…