ROOT  6.06/09
Reference Guide
mixmax.cxx
Go to the documentation of this file.
1 // $Id:$
2 // -*- C++ -*-
3 //
4 // -----------------------------------------------------------------------
5 // MixMax Matrix PseudoRandom Number Generator
6 // --- MixMax ---
7 // class header file
8 // -----------------------------------------------------------------------
9 //
10 //
11 // Created by Konstantin Savvidy on Sun Feb 22 2004.
12 // As of version 0.99 and later, the code is being released under
13 // GNU Lesser General Public License v3
14 //
15 // Generator described in
16 // N.Z.Akopov, G.K.Savvidy and N.G.Ter-Arutyunian, Matrix Generator of Pseudorandom Numbers,
17 // J.Comput.Phys. 97, 573 (1991);
18 // Preprint EPI-867(18)-86, Yerevan Jun.1986;
19 //
20 // and
21 //
22 // K.Savvidy
23 // The MIXMAX random number generator
24 // Comp. Phys. Commun. (2015)
25 // http://dx.doi.org/10.1016/j.cpc.2015.06.003
26 //
27 // -----------------------------------------------------------------------
28 
29 #include <stdio.h>
30 #include <stdlib.h>
31 #include <stdint.h>
32 
33 #define __MIXMAX_C // do NOT define it in your own program, just include mixmax.h
34 
35 //#define USE_INLINE_ASM //LM: uncomment if want to use inline asm
36 
37 #include "mixmax.h"
38 
39 int nskip = 2; // number of iterations we want to skip (default is 2)
40 
41 int firstElement = 3; // first element we want to return of the vector and skipping the previous ones
42 
43 void set_skip_number(int n){
44  nskip = n;
45 }
47  firstElement = n;
48 }
49 
51  return nskip;
52 }
54  return firstElement;
55 }
56 
58  int niter = nskip+1;
59  for (int i = 0; i < niter; ++i)
60  X->sumtot = iterate_raw_vec(X->V, X->sumtot);
61  return 0;
62 }
63 
64 #if (SPECIALMUL!=0)
65 inline uint64_t MULWU (uint64_t k){ return (( (k)<<(SPECIALMUL) & M61) | ( (k) >> (BITS-SPECIALMUL)) ) ;}
66 #elif (SPECIALMUL==0)
67 inline uint64_t MULWU (uint64_t k){ (void)k; return 0;}
68 #else
69 #error SPECIALMUL not undefined
70 #endif
71 
73  // operates with a raw vector, uses known sum of elements of Y
74  int i;
75 #ifdef SPECIAL
76  myuint temp2 = Y[1];
77 #endif
78  myuint tempP, tempV;
79  Y[0] = ( tempV = sumtotOld);
80  myuint sumtot = Y[0], ovflow = 0; // will keep a running sum of all new elements (except Y[0])
81  tempP = 0; // will keep a partial sum of all old elements (except Y[0])
82  for (i=1; i<N; i++){
83 #if (SPECIALMUL!=0)
84  myuint tempPO = MULWU(tempP);
85  tempP = modadd(tempP,Y[i]);
86  tempV = MOD_MERSENNE(tempV + tempP + tempPO); // edge cases ?
87 #else
88  tempP = modadd(tempP , Y[i]);
89  tempV = modadd(tempV , tempP);
90 #endif
91  Y[i] = tempV;
92  sumtot += tempV; if (sumtot < tempV) {ovflow++;}
93  }
94 #ifdef SPECIAL
95  temp2 = MOD_MULSPEC(temp2);
96  Y[2] = modadd( Y[2] , temp2 );
97  sumtot += temp2; if (sumtot < temp2) {ovflow++;}
98 #endif
99  return MOD_MERSENNE(MOD_MERSENNE(sumtot) + (ovflow <<3 ));
100 }
101 
103  int i;
104  i=X->counter;
105 
106  if (i<(N) ){
107  X->counter++;
108  return X->V[i];
109  }else{
110  int niter = nskip + 1;
111  for (int iter = 0; iter < niter; ++iter) {
112  X->sumtot = iterate_raw_vec(X->V, X->sumtot);
113  }
114 
115  X->counter=firstElement+1;
116  return X->V[firstElement];
117  }
118 }
119 
121  /* cast to signed int trick suggested by Andrzej Gòˆrlich */
122  // times for sknuth_Gap test: N = 1, n = 500000000, r = 0, Alpha = 0, Beta = 0.0625
123  //return ((int64_t)get_next(X)) * INV_MERSBASE; // 00:01:17.23
124  //return ( (double)get_next(X)) * INV_MERSBASE; // 00:01:57.89
125  int64_t Z=(int64_t)get_next(X);
126 #if defined(__SSE__) && defined(USE_INLINE_ASM)
127  double F;
128  __asm__ ("pxor %0, %0; "
129  "cvtsi2sdq %1, %0; "
130  :"=x"(F)
131  :"r"(Z)
132  );
133  return F*INV_MERSBASE;
134 #else
135  return Z*INV_MERSBASE;
136 #endif
137 
138 }
139 
140 // LM: this should be fixed for skipping the first number (firstElement(
141 void fill_array(rng_state_t* X, unsigned int n, double *array)
142 {
143  // Return an array of n random numbers uniformly distributed in (0,1]
144  unsigned int i,j;
145  const int M=N-1;
146  for (i=0; i<(n/M); i++){
147  int niter = nskip+1;
148  for (int iter = 0; iter < niter; ++iter)
149  iterate_and_fill_array(X, array+i*M);
150  }
151  unsigned int rem=(n % M);
152  if (rem) {
153  iterate(X);
154  for (j=0; j< (rem); j++){
155  array[M*i+j] = (int64_t)X->V[j] * (double)(INV_MERSBASE);
156  }
157  X->counter = j; // needed to continue with single fetches from the exact spot, but if you only use fill_array to get numbers then it is not necessary
158  }else{
159  X->counter = N;
160  }
161 }
162 
163 void iterate_and_fill_array(rng_state_t* X, double *array){
164  myuint* Y=X->V;
165  int i;
166  myuint tempP, tempV;
167 #if (SPECIAL != 0)
168  myuint temp2 = Y[1];
169 #endif
170  Y[0] = (tempV = modadd(Y[0] , X->sumtot));
171  //array[0] = (double)tempV * (double)(INV_MERSBASE);
172  myuint sumtot = 0, ovflow = 0; // will keep a running sum of all new elements (except Y[0])
173  tempP = 0; // will keep a partial sum of all old elements (except Y[0])
174  for (i=1; i<N; i++){
175  tempP = modadd(tempP,Y[i]);
176  Y[i] = ( tempV = modadd(tempV,tempP) );
177  sumtot += tempV; if (sumtot < tempV) {ovflow++;}
178  array[i-1] = (int64_t)tempV * (double)(INV_MERSBASE);
179  }
180 #if (SPECIAL != 0)
181  temp2 = MOD_MULSPEC(temp2);
182  Y[2] = modadd( Y[2] , temp2 );
183  sumtot += temp2; if (sumtot < temp2) {ovflow++;}
184 #endif
185  X->sumtot = MOD_MERSENNE(MOD_MERSENNE(sumtot) + (ovflow <<3 ));
186 }
187 
189 #if defined(__x86_64__) && defined(USE_INLINE_ASM)
190  myuint out;
191  /* Assembler trick suggested by Andrzej Gòˆrlich */
192  __asm__ ("addq %2, %0; "
193  "btrq $61, %0; "
194  "adcq $0, %0; "
195  :"=r"(out)
196  :"0"(foo), "r"(bar)
197  );
198  return out;
199 #else
200  return MOD_MERSENNE(foo+bar);
201 #endif
202 }
203 
205 {
206 /* allocate the state */
207  rng_state_t *p = (rng_state_t*)malloc(sizeof(rng_state_t));
208  p->fh=NULL; // by default, set the output file handle to stdout
209  return p;
210 }
211 
212 int rng_free(rng_state_t* X) /* free the memory occupied by the state */
213 {
214  free(X);
215  return 0;
216 }
217 
219 {
220  /* copy the vector stored at Y, and return pointer to the newly allocated and initialized state.
221  It is the user's responsibility to make sure that Y is properly allocated with rng_alloc,
222  then pass Y->V or it can also be an array -- such as myuint Y[N+1] and Y[1]...Y[N] have been set to legal values [0 .. MERSBASE-1]
223  Partial sums on this new state are recalculated, and counter set to zero, so that when get_next is called,
224  it will output the initial vector before any new numbers are produced, call iterate(X) if you want to advance right away */
225  rng_state_t* X = rng_alloc();
226  myuint sumtot=0,ovflow=0;
227  X->counter = 2;
228  int i;
229  for ( i=0; i < N; i++){
230  X->V[i] = Y[i];
231  sumtot += X->V[(i)]; if (sumtot < X->V[(i)]) {ovflow++;}
232 
233  }
234  X->sumtot = MOD_MERSENNE(MOD_MERSENNE(sumtot) + (ovflow <<3 ));
235  return X;
236 }
237 
238 void seed_vielbein(rng_state_t* X, unsigned int index)
239 {
240 int i;
241  if (index<N){
242  for (i=0; i < N; i++){
243  X->V[i] = 0;
244  }
245  X->V[index] = 1;
246  }else{
247  fprintf(stderr, "Out of bounds index, is not ( 0 <= index < N )\n"); exit(ARRAY_INDEX_OUT_OF_BOUNDS);
248  }
249  X->counter = N; // set the counter to N if iteration should happen right away
250  //precalc(X);
251  X->sumtot = 1; //(index ? 1:0);
252  if (X->fh==NULL){X->fh=stdout;}
253 }
254 
256 { // a 64-bit LCG from Knuth line 26, in combination with a bit swap is used to seed
257  const myuint MULT64=6364136223846793005ULL;
258  int i;
259  myuint sumtot=0,ovflow=0;
260  if (seed == 0){
261  fprintf(stderr, " try seeding with nonzero seed next time!\n");
262  exit(SEED_WAS_ZERO);
263  }
264 
265  myuint l = seed;
266 
267  //X->V[0] = l & MERSBASE;
268  if (X->fh==NULL){X->fh=stdout;} // if the filehandle is not yet set, make it stdout
269  for (i=0; i < N; i++){
270  l*=MULT64; l = (l << 32) ^ (l>>32);
271  X->V[i] = l & MERSBASE;
272  sumtot += X->V[(i)]; if (sumtot < X->V[(i)]) {ovflow++;}
273  }
274  X->counter = N; // set the counter to N if iteration should happen right away
275  X->sumtot = MOD_MERSENNE(MOD_MERSENNE(sumtot) + (ovflow <<3 ));
276 }
277 
279  int i;
280  myuint temp;
281  temp = 0;
282  for (i=0; i < N; i++){
283  temp = MOD_MERSENNE(temp + X->V[i]);
284  }
285  X->sumtot = temp;
286  return temp;
287 }
288 
289 
290 int rng_get_N(void){return N;}
291 
292 //#define MASK32 0xFFFFFFFFULL
293 
294 
295 //inline myuint modmulM61(myuint a, myuint b){
296 // // my best modmul so far
297 // __uint128_t temp;
298 // temp = (__uint128_t)a*(__uint128_t)b;
299 // return mod128(temp);
300 //}
301 
302 #if defined(__x86_64__)
303 inline myuint mod128(__uint128_t s){
304  myuint s1;
305  s1 = ( ( ((myuint)s)&MERSBASE ) + ( ((myuint)(s>>64)) * 8 ) + ( ((myuint)s) >>BITS) );
306  return MOD_MERSENNE(s1);
307 }
308 
309 inline myuint fmodmulM61(myuint cum, myuint a, myuint b){
310  __uint128_t temp;
311  temp = (__uint128_t)a*(__uint128_t)b + cum;
312  return mod128(temp);
313 }
314 
315 #else // on all other platforms, including 32-bit linux, PPC and PPC64 and all Windows
316 #define MASK32 0xFFFFFFFFULL
317 
319 {
320  register myuint o,ph,pl,ah,al;
321  o=(s)*a;
322  ph = ((s)>>32);
323  pl = (s) & MASK32;
324  ah = a>>32;
325  al = a & MASK32;
326  o = (o & M61) + ((ph*ah)<<3) + ((ah*pl+al*ph + ((al*pl)>>32))>>29) ;
327  o += cum;
328  o = (o & M61) + ((o>>61));
329  return o;
330 }
331 #endif
332 
334  int j;
335  fprintf(X->fh, "mixmax state, file version 1.0\n" );
336  fprintf(X->fh, "N=%u; V[N]={", rng_get_N() );
337  for (j=0; (j< (rng_get_N()-1) ); j++) {
338  fprintf(X->fh, "%llu, ", X->V[j] );
339  }
340  fprintf(X->fh, "%llu", X->V[rng_get_N()-1] );
341  fprintf(X->fh, "}; " );
342  fprintf(X->fh, "counter=%u; ", X->counter );
343  fprintf(X->fh, "sumtot=%llu;\n", X->sumtot );
344 }
345 
346 void read_state(rng_state_t* X, const char filename[] ){
347  // a function for reading the state from a file, after J. Apostolakis
348  FILE* fin;
349  if( ( fin = fopen(filename, "r") ) ){
350  char l=0;
351  while ( l != '{' ) { // 0x7B = "{"
352  l=fgetc(fin); // proceed until hitting opening bracket
353  }
354  ungetc(' ', fin);
355  }else{
356  fprintf(stderr, "mixmax -> read_state: error reading file %s\n", filename);
358  }
359 
360  myuint vecVal;
361  //printf("mixmax -> read_state: starting to read state from file\n");
362  if (!fscanf(fin, "%llu", &X->V[0]) ) {fprintf(stderr, "mixmax -> read_state: error reading file %s\n", filename); exit(ERROR_READING_STATE_FILE);}
363  //printf("V[%d] = %llu\n",0, X->V[0]);
364  int i;
365  for( i = 1; i < rng_get_N(); i++){
366  if (!fscanf(fin, ", %llu", &vecVal) ) {fprintf(stderr, "mixmax -> read_state: error reading vector component i=%d from file %s\n", i, filename); exit(ERROR_READING_STATE_FILE);}
367  //printf("V[%d] = %llu\n",i, vecVal);
368  if( vecVal <= MERSBASE ){
369  X->V[i] = vecVal;
370  }else{
371  fprintf(stderr, "mixmax -> read_state: Invalid state vector value= %llu"
372  " ( must be less than %llu ) "
373  " obtained from reading file %s\n"
374  , vecVal, MERSBASE, filename);
375 
376  }
377  }
378 
379  unsigned int counter;
380  if (!fscanf( fin, "}; counter=%u; ", &counter)){fprintf(stderr, "mixmax -> read_state: error reading counter from file %s\n", filename); exit(ERROR_READING_STATE_FILE);}
381  if( counter <= N ) {
382  X->counter= counter;
383  }else{
384  fprintf(stderr, "mixmax -> read_state: Invalid counter = %d"
385  " Must be 0 <= counter < %u\n" , counter, N);
386  print_state(X);
388  }
389  precalc(X);
390  myuint sumtot;
391  if (!fscanf( fin, "sumtot=%llu\n", &sumtot)){fprintf(stderr, "mixmax -> read_state: error reading checksum from file %s\n", filename); exit(ERROR_READING_STATE_FILE);}
392 
393  if (X->sumtot != sumtot) {
394  fprintf(stderr, "mixmax -> checksum error while reading state from file %s - corrupted?\n", filename);
396  }
397 // else{fprintf(stderr, "mixmax -> read_state: checksum ok: %llu == %llu\n",X->sumtot, sumtot);}
398  fclose(fin);
399 }
400 
401 
402 #define FUSEDMODMULVEC \
403 { for (i =0; i<N; i++){ \
404 cum[i] = fmodmulM61( cum[i], coeff , Y[i] ) ; \
405 } }
406 
407 
408 #define SKIPISON 1
409 
410 #if ( ( (N==88)||(N==256) ||(N==1000) || (N==3150) ) && BITS==61 && SKIPISON!=0)
411 void seed_uniquestream( rng_state_t* Xin, myID_t clusterID, myID_t machineID, myID_t runID, myID_t streamID ){
412  seed_vielbein(Xin,0);
413  Xin->sumtot = apply_bigskip(Xin->V, Xin->V, clusterID, machineID, runID, streamID );
414  if (Xin->fh==NULL){Xin->fh=stdout;} // if the filehandle is not yet set, make it stdout
415 }
416 
417 void branch_inplace( rng_state_t* Xin, myID_t* IDvec ){
418  Xin->sumtot = apply_bigskip(Xin->V, Xin->V, IDvec[3], IDvec[2], IDvec[1], IDvec[0] );
419 }
420 
421 myuint apply_bigskip(myuint* Vout, myuint* Vin, myID_t clusterID, myID_t machineID, myID_t runID, myID_t streamID ){
422  /*
423  makes a derived state vector, Vout, from the mother state vector Vin
424  by skipping a large number of steps, determined by the given seeding ID's
425 
426  it is mathematically guaranteed that the substreams derived in this way from the SAME (!!!) Vin will not collide provided
427  1) at least one bit of ID is different
428  2) less than 10^100 numbers are drawn from the stream
429  (this is good enough : a single CPU will not exceed this in the lifetime of the universe, 10^19 sec,
430  even if it had a clock cycle of Planch time, 10^44 Hz )
431 
432  Caution: never apply this to a derived vector, just choose some mother vector Vin, for example the unit vector by seed_vielbein(X,0),
433  and use it in all your runs, just change runID to get completely nonoverlapping streams of random numbers on a different day.
434 
435  clusterID and machineID are provided for the benefit of large organizations who wish to ensure that a simulation
436  which is running in parallel on a large number of clusters and machines will have non-colliding source of random numbers.
437 
438  did i repeat it enough times? the non-collision guarantee is absolute, not probabilistic
439 
440  */
441 
442 
443  const myuint skipMat[128][N] =
444 
445 #if (N==88)
446 #include "mixmax_skip_N88.icc" // to make this file, delete all except some chosen 128 rows of the coefficients table
447 #elif (N==256)
448 #include "mixmax_skip_N256.icc"
449 #elif (N==1000)
450 #include "mixmax_skip_N1000.icc"
451 #elif (N==3150)
452 #include "mixmax_skip_N3150.icc"
453 #endif
454  ;
455 
456  myID_t IDvec[4] = {streamID, runID, machineID, clusterID};
457  int r,i,j, IDindex;
458  myID_t id;
459  myuint Y[N], cum[N];
460  myuint coeff;
461  myuint* rowPtr;
462  myuint sumtot=0;
463 
464 
465  for (i=0; i<N; i++) { Y[i] = Vin[i]; sumtot = modadd( sumtot, Vin[i]); } ;
466  for (IDindex=0; IDindex<4; IDindex++) { // go from lower order to higher order ID
467  id=IDvec[IDindex];
468  //printf("now doing ID at level %d, with ID = %d\n", IDindex, id);
469  r = 0;
470  while (id){
471  if (id & 1) {
472  rowPtr = (myuint*)skipMat[r + IDindex*8*sizeof(myID_t)];
473  //printf("free coeff for row %d is %llu\n", r, rowPtr[0]);
474  for (i=0; i<N; i++){ cum[i] = 0; }
475  for (j=0; j<N; j++){ // j is lag, enumerates terms of the poly
476  // for zero lag Y is already given
477  coeff = rowPtr[j]; // same coeff for all i
479  sumtot = iterate_raw_vec(Y, sumtot);
480  }
481  sumtot=0;
482  for (i=0; i<N; i++){ Y[i] = cum[i]; sumtot = modadd( sumtot, cum[i]); } ;
483  }
484  id = (id >> 1); r++; // bring up the r-th bit in the ID
485  }
486  }
487  sumtot=0;
488  for (i=0; i<N; i++){ Vout[i] = Y[i]; sumtot = modadd( sumtot, Y[i]); } ; // returns sumtot, and copy the vector over to Vout
489  return (sumtot) ;
490 }
491 #else
492 #warning For this N, we dont have the skipping coefficients yet, using alternative method to seed
493 
494 void seed_uniquestream( rng_state_t* Xin, myID_t clusterID, myID_t machineID, myID_t runID, myID_t streamID ){
495  Xin->V[0] = (myuint)clusterID;
496  Xin->V[1] = (myuint)machineID;
497  Xin->V[2] = (myuint)runID;
498  Xin->V[3] = (myuint)streamID;
499  Xin->V[4] = (myuint)clusterID << 5;
500  Xin->V[5] = (myuint)machineID << 7;
501  Xin->V[6] = (myuint)runID << 11;
502  Xin->V[7] = (myuint)streamID << 13;
503  precalc(Xin);
504  Xin->sumtot = iterate_raw_vec(Xin->V, Xin->sumtot);
505  Xin->sumtot = iterate_raw_vec(Xin->V, Xin->sumtot);
506 }
507 #endif // SKIPISON
508 
509 
510 
511 
int firstElement
Definition: mixmax.cxx:41
rng_state_t * rng_copy(myuint *Y)
Definition: mixmax.cxx:218
myuint V[N]
Definition: mixmax.h:61
#define SPECIALMUL
Definition: mixmax.h:143
#define ERROR_READING_STATE_FILE
Definition: mixmax.h:221
void set_skip_number(int n)
Definition: mixmax.cxx:43
void branch_inplace(rng_state_t *Xin, myID_t *IDvec)
Definition: mixmax.cxx:417
static const char * filename()
#define MOD_MULSPEC(k)
Definition: mixmax.h:145
#define N
int rng_free(rng_state_t *X)
Definition: mixmax.cxx:212
#define SEED_WAS_ZERO
Definition: mixmax.h:220
rng_state_t * rng_alloc()
Definition: mixmax.cxx:204
#define FUSEDMODMULVEC
Definition: mixmax.cxx:402
int get_skip_number()
Definition: mixmax.cxx:50
double get_next_float(rng_state_t *X)
Definition: mixmax.cxx:120
#define MOD_MERSENNE(k)
Definition: mixmax.h:135
Vc_ALWAYS_INLINE void free(T *p)
Frees memory that was allocated with Vc::malloc.
Definition: memory.h:94
#define ARRAY_INDEX_OUT_OF_BOUNDS
Definition: mixmax.h:219
std::map< std::string, std::string >::const_iterator iter
Definition: TAlienJob.cxx:54
uint64_t myuint
Definition: mixmax.h:52
XFontStruct * id
Definition: TGX11.cxx:108
char * out
Definition: TBase64.cxx:29
void set_first_return_element(int n)
Definition: mixmax.cxx:46
#define F(x, y, z)
ROOT::R::TRInterface & r
Definition: Object.C:4
myuint get_next(rng_state_t *X)
Definition: mixmax.cxx:102
void seed_uniquestream(rng_state_t *Xin, myID_t clusterID, myID_t machineID, myID_t runID, myID_t streamID)
Definition: mixmax.cxx:411
FILE * fh
Definition: mixmax.h:64
uint32_t myID_t
Definition: mixmax.h:88
int counter
Definition: mixmax.h:63
void seed_vielbein(rng_state_t *X, unsigned int index)
Definition: mixmax.cxx:238
myuint fmodmulM61(myuint cum, myuint s, myuint a)
Definition: mixmax.cxx:318
TLine * l
Definition: textangle.C:4
myuint sumtot
Definition: mixmax.h:62
#define M61
Definition: mixmax.h:126
#define BITS
Definition: gifdecode.c:7
#define MERSBASE
Definition: mixmax.h:132
void read_state(rng_state_t *X, const char filename[])
Definition: mixmax.cxx:346
#define ERROR_READING_STATE_CHECKSUM
Definition: mixmax.h:223
void seed_spbox(rng_state_t *X, myuint seed)
Definition: mixmax.cxx:255
void print_state(rng_state_t *X)
Definition: mixmax.cxx:333
myuint apply_bigskip(myuint *Vout, myuint *Vin, myID_t clusterID, myID_t machineID, myID_t runID, myID_t streamID)
Definition: mixmax.cxx:421
uint64_t MULWU(uint64_t k)
Definition: mixmax.cxx:67
void iterate_and_fill_array(rng_state_t *X, double *array)
Definition: mixmax.cxx:163
#define INV_MERSBASE
#define INV_MERSBASE (0x1p-61)
Definition: mixmax.h:138
int nskip
Definition: mixmax.cxx:39
myuint precalc(rng_state_t *X)
Definition: mixmax.cxx:278
int get_first_return_element()
Definition: mixmax.cxx:53
myuint iterate_raw_vec(myuint *Y, myuint sumtotOld)
Definition: mixmax.cxx:72
typedef void((*Func_t)())
int rng_get_N(void)
Definition: mixmax.cxx:290
#define NULL
Definition: Rtypes.h:82
int iterate(rng_state_t *X)
Definition: mixmax.cxx:57
std::complex< float_v > Z
Definition: main.cpp:120
#define ERROR_READING_STATE_COUNTER
Definition: mixmax.h:222
Vc_ALWAYS_INLINE_L T *Vc_ALWAYS_INLINE_R malloc(size_t n)
Allocates memory on the Heap with alignment and padding suitable for vectorized access.
Definition: memory.h:67
const Int_t n
Definition: legend1.C:16
#define MASK32
Definition: mixmax.cxx:316
void fill_array(rng_state_t *X, unsigned int n, double *array)
Definition: mixmax.cxx:141
myuint modadd(myuint foo, myuint bar)
Definition: mixmax.cxx:188