Re: [PATCH] Fix biased shuffle by avoiding already "struck" elements

[Sergei Trofimovich]

On Wed, 19 Jun 2024 10:55:53 -0400
Dmitry Goncharov <dgoncha...@users.sf.net> wrote:

> On Tue, Jun 18, 2024 at 5:38 PM Sergei Trofimovich <sly...@gmail.com> wrote:
> > After the change probabilities are not as biased:
> >
> >             0     1     2     3
> >         _____ _____ _____ _____
> >     0 | 24.99 24.99 25.01 25.01
> >     1 | 24.99 25.04 24.99 24.99
> >     2 | 25.01 25.00 25.00 24.99
> >     3 | 25.01 24.98 25.00 25.01
> >  
> 
> Looks like you wrote a program to figure out these probabilities.
> Should such a program be included as a part of make test suite?

I wrote a one-off `a.c` printer (attached).

It might be tricky to write safe ranges for probabilities: 4% is
"obviously" biased, but 0.02% is trickier (PRNG / modulo bias,
truncation, or expected deviation).

Glancing at tests/ all the tests exercise user-facing `make` API. What
would be the best way to validate probabilities?

-- 

  Sergei

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>

#define LEN 4
static int a[LEN];

static void
random_shuffle_array (void)
{
  size_t i;
  for (i = 0; i < LEN; i++)
    {
      int t;

      /* Pick random element and swap. */
      unsigned int j = rand () % LEN;
      if (i == j)
        continue;

      /* Swap. */
      t = a[i];
      a[i] = a[j];
      a[j] = t;
    }
}

static void
random_shuffle_array_fixed (void)
{
  size_t i;
  for (i = LEN - 1; i >= 1; i--)
    {
      int t;

      /* Pick random element and swap. */
      unsigned int j = rand () % (i + 1);
      if (i == j)
        continue;

      /* Swap. */
      t = a[i];
      a[i] = a[j];
      a[j] = t;
    }
}

int main() {
    size_t hist[LEN][LEN];
    memset(hist, 0, sizeof(hist));

    srand(time(NULL));

    size_t niters = 10000000;
    for (size_t iter = 0; iter < niters; ++iter) {
        for (size_t i = 0; i < LEN; ++i) a[i] = i;
//        random_shuffle_array ();
        random_shuffle_array_fixed ();
        for (size_t i = 0; i < LEN; ++i) hist[i][a[i]] += 1;
    }

    printf("       ");
    for (size_t j = 0; j < LEN; ++j) {
        printf(" %5zu", j);
    }
    printf("\n");
    printf("       ");
    for (size_t j = 0; j < LEN; ++j) {
        printf(" _____");
    }
    printf("\n");
    for (size_t i = 0; i < LEN; ++i) {
        printf("%5zu |", i);
        for (size_t j = 0; j < LEN; ++j) {
            printf(" %5.2f", (double)(hist[i][j]) / (double)(niters) * 100.0);
        }
        printf("\n");
    }
}