Save and Load Random Number Generator State in C++11

Question

This question has been asked before (stackoverflow) but the (accepted) answer is not satisfactory.

The following example saves and loads the state but depending on the number of generated values it works or it doesn't:

#include <fstream>
#include <iostream>
#include <random>
#include <cassert>

int main()
{
  const int preN = 4;
  const int middleN = 0;

  // initialize randGen
  std::mt19937 randGen1;
  std::normal_distribution<double> distribution1;


  // print some initial random numbers
  for (int i=0;i<preN;++i)
    std::cout << distribution1(randGen1)<<" ";

  // save state
  std::cout << std::endl << "Saving...\n";
  {
    std::ofstream fout("seed.dat");
    fout << randGen1;
  }

  // maybe advance randGen1
  for (int i=0;i<middleN;++i)
    std::cout << distribution1(randGen1)<<" ";

  // load saved state into randGen2 
  std::cout << std::endl << "Loading...\n";
  std::ifstream fin("seed.dat");
  std::mt19937 randGen2;
  fin >> randGen2;
  std::normal_distribution<double> distribution2;

  // are both randGen equal?
  assert(randGen1 == randGen2);

  // print numbers from both generators
  std::cout << "Generator1\tGenerator2\n";
  std::cout << distribution1(randGen1) << "\t"
            << distribution2(randGen2) << "\n";

  return 0;

}

With these parameters it works like expected. However, if I set preN=3 the output looks like:

0.13453 -0.146382 0.46065 
Saving...

Loading...
Generator1  Generator2
-1.87138    0.163712

Why did the assert not apply? Now I set preN=3 and middleN=1 and the output is

0.13453 -0.146382 0.46065 
Saving...
-1.87138 
Loading...
Generator1  Generator2
0.163712    0.163712

If I set middleN to anything larger than 1 the assert applies. Can anyone explain what is going on? What am I doing wrong or not understanding?

Tested with GCC5.4.0 and CLANG3.8.0 on Linux

Nicol Bolas · Accepted Answer · 2017-02-22T02:10:35.087

8

The problem is not your random number generator's state. The problem is your distribution's state. Yes, distributions can have state too.

You can get the same values by resetting the normal distribution's state with reset. Alternatively, you can preserve and reconstitute the distribution's state too, using << and >>.

edited Feb 22 '17 at 02:10

answered Feb 21 '17 at 21:23

Nicol Bolas

449,505
63
781
982

Tried it too (because that behaviour looked fascinating and i did not try C++ for a long time) and can say for sure, that a uniform-dist behaves differently! – sascha Feb 21 '17 at 21:25
@sascha: I'm not sure what you mean that it behaves differently. – Nicol Bolas Feb 21 '17 at 21:26
That it probably does not have an internal-state which matters (much simpler than normal-dist), but that's just empirical observation. It was just my first trial. Maybe it does too. Ignore my comment. Nice answer! – sascha Feb 21 '17 at 21:26
@sascha: Whether a distribution has state that affects its output is implementation dependent. – Nicol Bolas Feb 21 '17 at 21:27
You got any idea on the why of this kind of free specification? I can't think of any good reason not to assume these are state-less except for caching (build a random-pool). I also don't think that scipy does that (at first look even GSL). – sascha Feb 21 '17 at 21:35
`normal_distribution` may or may not have state. It all depends on the implementation's choice of algorithm. – T.C. Feb 22 '17 at 02:02
Thank you very much. I considered this to be the case and something made me think to circumvent the problem by introducing `distribution2`. In any case, I'll mark the answer as the solution but add the corrected code below. – Pieter Bruegel the Elder Feb 22 '17 at 08:58

score 4 · Answer 2 · edited May 23 '17 at 12:17

Thanks to the answer from Nicol Bolas above, I can add the corrected code below:

#include <fstream>
#include <iostream>
#include <random>
#include <cassert>

int main()
{
  const int preN = 7;
  const int middleN = 0;

  // initialize another randGen
  std::mt19937 randGen1;
  std::normal_distribution<double> distribution1;

  // print some initial random numbers
  for (int i=0;i<preN;++i)
    std::cout << distribution1(randGen1)<<" ";

  // save state
  std::cout << std::endl << "Saving...\n";
  {
    std::ofstream fout("seed.dat");
    fout << randGen1;
    fout.close();
    std::ofstream fout2("distribution.dat");
    fout2 << distribution1;
    fout2.close();
  }

  // maybe advance randGen
  for (int i=0;i<middleN;++i)
    std::cout << distribution1(randGen1)<<" ";

  // load saved state into randGen2
  std::cout << std::endl << "Loading...\n";
  std::mt19937 randGen2;
  std::normal_distribution<double> distribution2;
  {
    std::ifstream fin("seed.dat");
    fin >> randGen2;
    fin.close();
    std::ifstream fin2("distribution.dat");
    fin2 >> distribution2;
    fin2.close();
  }

  // are both randGen equal?
  assert(randGen1 == randGen2);
  assert(distribution1 == distribution2);

  // print numbers from both generators
  std::cout << "Generator1\tGenerator2\n";
  std::cout << distribution1(randGen1) << "\t"
            << distribution2(randGen2) << "\n";

  return 0;
}

score 1 · Answer 3 · edited Feb 15 '20 at 17:06

Here is how one can save and restore seed forrandom numbers of double float precision. It should be similar for integers - use jrand48 instead of erand48.

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

int main(int argc, char **argv)
{

     unsigned short savedseed[3];
     unsigned short currentseed[3];
     double x;

     /*-- initialize ramdom seed to whatever --*/
     currentseed[0]= 23;
     currentseed[1]= 45;
     currentseed[2]= 67;

     printf("\n");

     /*-- generate three random numbers  --*/
     x =  erand48(currentseed);     printf("%g\n", x);
     x =  erand48(currentseed);     printf("%g\n", x);
     x =  erand48(currentseed);     printf("%g\n", x);
     printf("\n");

     /*-- save seed  --*/
     memcpy(savedseed, currentseed, 3*sizeof(unsigned short));

     /*-- generate next three random numbers  --*/     
     x =  erand48(currentseed);     printf("%g\n", x);
     x =  erand48(currentseed);     printf("%g\n", x);
     x =  erand48(currentseed);     printf("%g\n", x);
     printf("\n", x);

     /*-- restore seed  --*/
     memcpy(currentseed, savedseed,  3*sizeof(unsigned short));

     /*-- generate the same three random numbers again --*/
     x =  erand48(currentseed);     printf("%g\n", x);
     x =  erand48(currentseed);     printf("%g\n", x);
     x =  erand48(currentseed);     printf("%g\n", x);
     printf("\n");  
}

Save and Load Random Number Generator State in C++11

3 Answers3