/home/docs/checkouts/readthedocs.org/user_builds/advanced-micro-devices-rocrand/checkouts/develop/projects/rocrand/library/include/rocrand/rocrand_mrg32k3a.h Source File

• library
• include
• rocrand

rocrand_mrg32k3a.h

1 // Copyright (c) 2017-2025 Advanced Micro Devices, Inc. All rights reserved.

2 //

3 // Permission is hereby granted, free of charge, to any person obtaining a copy

4 // of this software and associated documentation files (the "Software"), to deal

5 // in the Software without restriction, including without limitation the rights

6 // to use, copy, modify, merge, publish, distribute, sublicense, and/or sell

7 // copies of the Software, and to permit persons to whom the Software is

8 // furnished to do so, subject to the following conditions:

9 //

10 // The above copyright notice and this permission notice shall be included in

11 // all copies or substantial portions of the Software.

12 //

13 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR

14 // IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,

15 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE

16 // AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER

17 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,

18 // OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN

19 // THE SOFTWARE.

20 

21 #ifndef ROCRAND_MRG32K3A_H_

22 #define ROCRAND_MRG32K3A_H_

23 

24 #include "rocrand/rocrand_common.h"

25 #include "rocrand/rocrand_mrg32k3a_precomputed.h"

26 

27 #include <hip/hip_runtime.h>

28 

29 #define ROCRAND_MRG32K3A_POW32 4294967296U

30 #define ROCRAND_MRG32K3A_M1 4294967087U

31 #define ROCRAND_MRG32K3A_M1C 209U

32 #define ROCRAND_MRG32K3A_M2 4294944443U

33 #define ROCRAND_MRG32K3A_M2C 22853U

34 #define ROCRAND_MRG32K3A_A12 1403580U

35 #define ROCRAND_MRG32K3A_A13 (4294967087U - 810728U)

36 #define ROCRAND_MRG32K3A_A13N 810728U

37 #define ROCRAND_MRG32K3A_A21 527612U

38 #define ROCRAND_MRG32K3A_A23 (4294944443U - 1370589U)

39 #define ROCRAND_MRG32K3A_A23N 1370589U

40 #define ROCRAND_MRG32K3A_NORM_DOUBLE (2.3283065498378288e-10) // 1/ROCRAND_MRG32K3A_M1

41 #define ROCRAND_MRG32K3A_UINT_NORM \

42  (1.000000048661607) // (ROCRAND_MRG32K3A_POW32 - 1)/(ROCRAND_MRG32K3A_M1 - 1)

43 

52  #define ROCRAND_MRG32K3A_DEFAULT_SEED 12345ULL // end of group rocranddevice

54 

55 namespace rocrand_device {

56 

57 class mrg32k3a_engine

58 {

59 public:

60  struct mrg32k3a_state

61  {

62  #ifndef ROCRAND_DETAIL_BM_NOT_IN_STATE

63  // The Box–Muller transform requires two inputs to convert uniformly

64  // distributed real values [0; 1] to normally distributed real values

65  // (with mean = 0, and stddev = 1). Often user wants only one

66  // normally distributed number, to save performance and random

67  // numbers the 2nd value is saved for future requests.

68  double boxmuller_double; // normally distributed double

69  float boxmuller_float; // normally distributed float

70  #endif

71  unsigned int g1[3];

72  unsigned int g2[3];

73  };

74 

75  __forceinline__ __device__ __host__ mrg32k3a_engine()

76  {

77  this->seed(ROCRAND_MRG32K3A_DEFAULT_SEED, 0, 0);

78  }

79 

88  __forceinline__ __device__ __host__ mrg32k3a_engine(const unsigned long long seed,

89  const unsigned long long subsequence,

90  const unsigned long long offset)

91  {

92  this->seed(seed, subsequence, offset);

93  }

94 

103  __forceinline__ __device__ __host__ void seed(unsigned long long seed_value,

104  const unsigned long long subsequence,

105  const unsigned long long offset)

106  {

107  if(seed_value == 0)

108  {

109  seed_value = ROCRAND_MRG32K3A_DEFAULT_SEED;

110  }

111  unsigned int x = (unsigned int) seed_value ^ 0x55555555U;

112  unsigned int y = (unsigned int) ((seed_value >> 32) ^ 0xAAAAAAAAU);

113  m_state.g1[0] = mod_mul_m1(x, seed_value);

114  m_state.g1[1] = mod_mul_m1(y, seed_value);

115  m_state.g1[2] = mod_mul_m1(x, seed_value);

116  m_state.g2[0] = mod_mul_m2(y, seed_value);

117  m_state.g2[1] = mod_mul_m2(x, seed_value);

118  m_state.g2[2] = mod_mul_m2(y, seed_value);

119  this->restart(subsequence, offset);

120  }

121 

123  __forceinline__ __device__ __host__ void discard(unsigned long long offset)

124  {

125  this->discard_impl(offset);

126  }

127 

130  __forceinline__ __device__ __host__ void discard_subsequence(unsigned long long subsequence)

131  {

132  this->discard_subsequence_impl(subsequence);

133  }

134 

137  __forceinline__ __device__ __host__ void discard_sequence(unsigned long long sequence)

138  {

139  this->discard_sequence_impl(sequence);

140  }

141 

142  __forceinline__ __device__ __host__ void restart(const unsigned long long subsequence,

143  const unsigned long long offset)

144  {

145  #ifndef ROCRAND_DETAIL_BM_NOT_IN_STATE

146  m_state.boxmuller_float = ROCRAND_NAN_FLOAT;

147  m_state.boxmuller_double = ROCRAND_NAN_DOUBLE;

148  #endif

149  this->discard_subsequence_impl(subsequence);

150  this->discard_impl(offset);

151  }

152 

153  __forceinline__ __device__ __host__ unsigned int operator()()

154  {

155  return this->next();

156  }

157 

158  // Returned value is in range [1, ROCRAND_MRG32K3A_M1],

159  // where ROCRAND_MRG32K3A_M1 < UINT_MAX

160  __forceinline__ __device__ __host__

161  unsigned int next()

162  {

163  const unsigned int p1 = mod_m1(detail::mad_u64_u32(

164  ROCRAND_MRG32K3A_A12,

165  m_state.g1[1],

166  detail::mul_u64_u32(ROCRAND_MRG32K3A_A13N, (ROCRAND_MRG32K3A_M1 - m_state.g1[0]))));

167 

168  m_state.g1[0] = m_state.g1[1];

169  m_state.g1[1] = m_state.g1[2];

170  m_state.g1[2] = p1;

171 

172  const unsigned int p2 = mod_m2(detail::mad_u64_u32(

173  ROCRAND_MRG32K3A_A21,

174  m_state.g2[2],

175  detail::mul_u64_u32(ROCRAND_MRG32K3A_A23N, (ROCRAND_MRG32K3A_M2 - m_state.g2[0]))));

176 

177  m_state.g2[0] = m_state.g2[1];

178  m_state.g2[1] = m_state.g2[2];

179  m_state.g2[2] = p2;

180 

181  return (p1 - p2) + (p1 <= p2 ? ROCRAND_MRG32K3A_M1 : 0);

182  }

183 

184 protected:

185  // Advances the internal state to skip \p offset numbers.

186  // DOES NOT CALCULATE NEW ULONGLONG

187  __forceinline__ __device__ __host__ void discard_impl(unsigned long long offset)

188  {

189  discard_state(offset);

190  }

191 

192  // DOES NOT CALCULATE NEW ULONGLONG

193  __forceinline__ __device__ __host__ void

194  discard_subsequence_impl(unsigned long long subsequence)

195  {

196  int i = 0;

197 

198  while(subsequence > 0) {

199  if (subsequence & 1) {

200  #if defined(__HIP_DEVICE_COMPILE__)

201  mod_mat_vec_m1(d_A1P76 + i, m_state.g1);

202  mod_mat_vec_m2(d_A2P76 + i, m_state.g2);

203  #else

204  mod_mat_vec_m1(h_A1P76 + i, m_state.g1);

205  mod_mat_vec_m2(h_A2P76 + i, m_state.g2);

206  #endif

207  }

208  subsequence >>= 1;

209  i += 9;

210  }

211  }

212 

213  // DOES NOT CALCULATE NEW ULONGLONG

214  __forceinline__ __device__ __host__ void discard_sequence_impl(unsigned long long sequence)

215  {

216  int i = 0;

217 

218  while(sequence > 0) {

219  if (sequence & 1) {

220  #if defined(__HIP_DEVICE_COMPILE__)

221  mod_mat_vec_m1(d_A1P127 + i, m_state.g1);

222  mod_mat_vec_m2(d_A2P127 + i, m_state.g2);

223  #else

224  mod_mat_vec_m1(h_A1P127 + i, m_state.g1);

225  mod_mat_vec_m2(h_A2P127 + i, m_state.g2);

226  #endif

227  }

228  sequence >>= 1;

229  i += 9;

230  }

231  }

232 

233  // Advances the internal state by offset times.

234  // DOES NOT CALCULATE NEW ULONGLONG

235  __forceinline__ __device__ __host__ void discard_state(unsigned long long offset)

236  {

237  int i = 0;

238 

239  while(offset > 0) {

240  if (offset & 1) {

241  #if defined(__HIP_DEVICE_COMPILE__)

242  mod_mat_vec_m1(d_A1 + i, m_state.g1);

243  mod_mat_vec_m2(d_A2 + i, m_state.g2);

244  #else

245  mod_mat_vec_m1(h_A1 + i, m_state.g1);

246  mod_mat_vec_m2(h_A2 + i, m_state.g2);

247  #endif

248  }

249  offset >>= 1;

250  i += 9;

251  }

252  }

253 

254  // Advances the internal state to the next state

255  // DOES NOT CALCULATE NEW ULONGLONG

256  __forceinline__ __device__ __host__ void discard_state()

257  {

258  discard_state(1);

259  }

260 

261 private:

262  __forceinline__ __device__ __host__

263  static void mod_mat_vec_m1(const unsigned int* A, unsigned int* s)

264  {

265  unsigned long long x[3] = {s[0], s[1], s[2]};

266 

267  s[0] = mod_m1(mod_m1(A[0] * x[0]) + mod_m1(A[1] * x[1]) + mod_m1(A[2] * x[2]));

268 

269  s[1] = mod_m1(mod_m1(A[3] * x[0]) + mod_m1(A[4] * x[1]) + mod_m1(A[5] * x[2]));

270 

271  s[2] = mod_m1(mod_m1(A[6] * x[0]) + mod_m1(A[7] * x[1]) + mod_m1(A[8] * x[2]));

272  }

273 

274  __forceinline__ __device__ __host__

275  static void mod_mat_vec_m2(const unsigned int* A, unsigned int* s)

276  {

277  unsigned long long x[3] = {s[0], s[1], s[2]};

278 

279  s[0] = mod_m2(mod_m2(A[0] * x[0]) + mod_m2(A[1] * x[1]) + mod_m2(A[2] * x[2]));

280 

281  s[1] = mod_m2(mod_m2(A[3] * x[0]) + mod_m2(A[4] * x[1]) + mod_m2(A[5] * x[2]));

282 

283  s[2] = mod_m2(mod_m2(A[6] * x[0]) + mod_m2(A[7] * x[1]) + mod_m2(A[8] * x[2]));

284  }

285 

286  __forceinline__ __device__ __host__ static unsigned long long mod_mul_m1(unsigned int i,

287  unsigned long long j)

288  {

289  long long hi, lo, temp1, temp2;

290 

291  hi = i / 131072;

292  lo = i - (hi * 131072);

293  temp1 = mod_m1(hi * j) * 131072;

294  temp2 = mod_m1(lo * j);

295  lo = mod_m1(temp1 + temp2);

296 

297  if (lo < 0)

298  lo += ROCRAND_MRG32K3A_M1;

299  return lo;

300  }

301 

302  __forceinline__ __device__ __host__

303  static unsigned long long mod_m1(unsigned long long p)

304  {

305  p = detail::mad_u64_u32(ROCRAND_MRG32K3A_M1C,

306  static_cast<unsigned int>(p >> 32),

307  static_cast<unsigned int>(p));

308  if(p >= ROCRAND_MRG32K3A_M1)

309  p -= ROCRAND_MRG32K3A_M1;

310 

311  return p;

312  }

313 

314  __forceinline__ __device__ __host__

315  static unsigned long long mod_mul_m2(unsigned int i, unsigned long long j)

316  {

317  long long hi, lo, temp1, temp2;

318 

319  hi = i / 131072;

320  lo = i - (hi * 131072);

321  temp1 = mod_m2(hi * j) * 131072;

322  temp2 = mod_m2(lo * j);

323  lo = mod_m2(temp1 + temp2);

324 

325  if (lo < 0)

326  lo += ROCRAND_MRG32K3A_M2;

327  return lo;

328  }

329 

330  __forceinline__ __device__ __host__

331  static unsigned long long mod_m2(unsigned long long p)

332  {

333  p = detail::mad_u64_u32(ROCRAND_MRG32K3A_M2C,

334  static_cast<unsigned int>(p >> 32),

335  static_cast<unsigned int>(p));

336  p = detail::mad_u64_u32(ROCRAND_MRG32K3A_M2C,

337  static_cast<unsigned int>(p >> 32),

338  static_cast<unsigned int>(p));

339  if(p >= ROCRAND_MRG32K3A_M2)

340  p -= ROCRAND_MRG32K3A_M2;

341 

342  return p;

343  }

344 

345 protected:

346  // State

347  mrg32k3a_state m_state;

348 

349  #ifndef ROCRAND_DETAIL_BM_NOT_IN_STATE

350  friend struct detail::engine_boxmuller_helper<mrg32k3a_engine>;

351  #endif

352 

353 }; // mrg32k3a_engine class

354 

355 } // end namespace rocrand_device

356 

363 typedef rocrand_device::mrg32k3a_engine rocrand_state_mrg32k3a;

365 

377 __forceinline__ __device__ __host__

378 void rocrand_init(const unsigned long long seed,

379  const unsigned long long subsequence,

380  const unsigned long long offset,

381  rocrand_state_mrg32k3a* state)

382 {

383  *state = rocrand_state_mrg32k3a(seed, subsequence, offset);

384 }

385 

398 __forceinline__ __device__ __host__

399 unsigned int rocrand(rocrand_state_mrg32k3a* state)

400 {

401  // next() in [1, ROCRAND_MRG32K3A_M1]

402  return static_cast<unsigned int>((state->next() - 1) * ROCRAND_MRG32K3A_UINT_NORM);

403 }

404 

413 __forceinline__ __device__ __host__

414 void skipahead(unsigned long long offset, rocrand_state_mrg32k3a* state)

415 {

416  return state->discard(offset);

417 }

418 

428 __forceinline__ __device__ __host__

429 void skipahead_subsequence(unsigned long long subsequence, rocrand_state_mrg32k3a* state)

430 {

431  return state->discard_subsequence(subsequence);

432 }

433 

443 __forceinline__ __device__ __host__

444 void skipahead_sequence(unsigned long long sequence, rocrand_state_mrg32k3a* state)

445 {

446  return state->discard_sequence(sequence);

447 }

448  // end of group rocranddevice

450 

451 #endif // ROCRAND_MRG32K3A_H_

ROCRAND_MRG32K3A_DEFAULT_SEED

#define ROCRAND_MRG32K3A_DEFAULT_SEED

Default seed for MRG32K3A PRNG.

Definition: rocrand_mrg32k3a.h:52

rocrand_init

__forceinline__ __device__ __host__ void rocrand_init(const unsigned long long seed, const unsigned long long subsequence, const unsigned long long offset, rocrand_state_mrg32k3a *state)

Initializes MRG32K3A state.

Definition: rocrand_mrg32k3a.h:378

skipahead_subsequence

__forceinline__ __device__ __host__ void skipahead_subsequence(unsigned long long subsequence, rocrand_state_mrg32k3a *state)

Updates MRG32K3A state to skip ahead by subsequence subsequences.

Definition: rocrand_mrg32k3a.h:429

skipahead

__forceinline__ __device__ __host__ void skipahead(unsigned long long offset, rocrand_state_mrg32k3a *state)

Updates MRG32K3A state to skip ahead by offset elements.

Definition: rocrand_mrg32k3a.h:414

rocrand

__forceinline__ __device__ __host__ unsigned int rocrand(rocrand_state_mrg32k3a *state)

Returns uniformly distributed random unsigned int value from [0; 2^32 - 1] range.

Definition: rocrand_mrg32k3a.h:399

skipahead_sequence

__forceinline__ __device__ __host__ void skipahead_sequence(unsigned long long sequence, rocrand_state_mrg32k3a *state)

Updates MRG32K3A state to skip ahead by sequence sequences.

Definition: rocrand_mrg32k3a.h:444