ClassicMC-MoRiBs/mc_randg.cc at master · Lechengwang/ClassicMC-MoRiBs · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
#include <math.h>
#include <stdio.h>
#include <ctime>
#include <cstdlib>

int SEED;

#include "sprng.h"

#include "mc_input.h"
#include "mc_randg.h"
#include "mc_utils.h"

const int MAXGENS=15;  // max number of rnd generators [gauss #8,9 N,M #>=10]
int *STREAM[MAXGENS];  // rnd numbers streams


void RandomInit(int mpi_id, int max_mpi_procs)  // sprng init
{
  int max_streams=MAXGENS*max_mpi_procs;        // number of streams

  srand(time(0));       //seed
  SEED = rand()%1000000000;        //seed for PIMC,seed is (1,1000000000)

  cout<<"The random seed for this work is: "<<SEED<<endl;

  for (int ip=0;ip<MAXGENS;ip++)
  {
     int streamnum = mpi_id*MAXGENS+ip;
     STREAM[ip] = init_sprng(streamnum,max_streams,SEED,SPRNG_DEFAULT); // initialize stream
//   print_sprng(stream);
  }
}

void RandomFree(void) // free memory used by sprng
{
  for (int ip=0;ip<MAXGENS;ip++)
  free_sprng(STREAM[ip]);
}

double rnd(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[0]);
}

double rnd1(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[1]);
}

double rnd2(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[2]);
}

double rnd3(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[3]);
}

double rnd4(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[4]);
}

double rnd5(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[5]);
}

double rnd6(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[6]);
}

double rnd7(void)
// double random number generator in the range [0,1]
{
  return sprng(STREAM[7]);
}

double gauss(double alpha)
// random numbers according Gaussian distribution exp{-alpha*x^2}
{
    double r1 = sprng(STREAM[8]);
    double r2 = sprng(STREAM[9]);

    double x1 = sqrt(-log(r1))*cos(2.0*M_PI*r2);

//  double x1 = sqrt(-log(r1))*cos(2.0*M_PI*r2);
//  double x2 = sqrt(-log(r1))*sin(2.0*M_PI*r2);

    return (x1/sqrt(alpha));
}

int nrnd1(int n)
// integer random number generator in the range [0,n-1]
{
  return (int)floor(n*sprng(STREAM[10]));
}

int nrnd2(int n)
// integer random number generator in the range [0,n-1]
{
  return (int)floor(n*sprng(STREAM[11]));
}

int nrnd3(int n)
// integer random number generator in the range [0,n-1]
{
  return (int)floor(n*sprng(STREAM[12]));
}

int nrnd4(int n)
// integer random number generator in the range [0,n-1]
{
  return (int)floor(n*sprng(STREAM[13]));
}