home *** CD-ROM | disk | FTP | other *** search
- /**************************************************************************
- * *
- * Author : Dr. Thomas Brandes, GMD, I1.HR *
- * Copyright : GMD St. Augustin, Germany *
- * Date : Feb 92 *
- * Last Update : Apr 93 *
- * *
- * This Module is part of the DALIB *
- * *
- * Module : random.c *
- * *
- * Function : Getting random numbers (integer, float, double) *
- * *
- * Export : Internal Use *
- * *
- * void random_block_init () *
- * *
- * Export : FORTRAN Interface *
- * *
- * dalib_random_init (n) : initialization of generator by seed *
- * int *n; *
- * *
- * dalib_get_real_randoms (a, size, limit) *
- * float a[]; *
- * int *size; *
- * float *limit; *
- * *
- * dalib_get_double_randoms (a, size, limit) *
- * double a[]; *
- * int *size; *
- * double *limit; *
- * *
- * dalib_get_int_randoms (a, size, limit) *
- * int a[]; *
- * int *size; *
- * int *limit; ! integers between 0 .. limit-1 *
- * *
- **************************************************************************/
-
- # include "system.h" /* needs pcb.i */
-
- /*******************************************************************
- * *
- * Local Block for random numbers on each process *
- * *
- *******************************************************************/
-
- #define RANDOM_SIZE 20 /* 250 */
- #define RANDOM_K1 17 /* 103 */
- #define RANDOM_K2 5
-
- struct
- { int state[RANDOM_SIZE];
- int pos;
- } random_block;
-
- void random_block_init ()
- /* This routine is called concurrently */
- { int j;
- #if defined(SUN4) || defined(ALLIANT) || defined(KSR1) || defined(OS2)
- srand (pcb.i*(2*3*5*7 - 1));
- #endif
- #if defined(IBM) || defined (GC) || defined(IPSC) || defined(SRM) || defined(MEIKO) || defined(SGI)
- srand48 (pcb.i*(2*3*5*7 - 1));
- #endif
- for (j=0;j<RANDOM_SIZE;j++)
- #if defined(SUN4) || defined(ALLIANT) || defined(KSR1) || defined(OS2)
- random_block.state[j] = rand();
- #endif
- #if defined(IBM) || defined (GC) || defined(IPSC) || defined(SRM) || defined(MEIKO) || defined(SGI)
- random_block.state[j] = mrand48 ();
- #endif
- random_block.pos = 0;
- }
-
- int random_block_get ()
- /* get a new 32-bit random value from random_block
- this routine can be called concurrently */
- { int k1, k2, pos, z;
- pos = random_block.pos;
- k1 = pos - RANDOM_K1;
- k2 = pos - RANDOM_K2;
- if (k1 < 0) k1+=RANDOM_SIZE;
- if (k2 < 0) k2+=RANDOM_SIZE;
- z = random_block.state[k1] - random_block.state[k2]; /* xor */
- random_block.state[pos] = z;
- pos += 1;
- if (pos >= RANDOM_SIZE) pos = 0;
- random_block.pos = pos;
- return(z);
- }
-
- /*******************************************************************
- * *
- * FORTRAN - Interface *
- * *
- * dalib_random_init (n) : initialization of generator by node *
- * *
- * dalib_get_int_randoms (a, size, limit) *
- * *
- *******************************************************************/
-
- void dalib_random_init__ (n)
- /* new initialization of the number generator with n */
- int *n;
- { int j, p;
- #if defined(SUN4) || defined(ALLIANT) || defined(KSR1) || defined(OS2)
- srand(*n * pcb.i); /* set value */
- #endif
- #if defined(IBM) || defined (GC) || defined(IPSC) || defined(SRM) || defined(MEIKO) || defined(SGI)
- srand48(*n * pcb.i); /* set value */
- #endif
- random_block_init ();
- }
-
- void dalib_get_int_randoms__ (a, size, limit)
- int a[], *size, *limit;
- { int j, hz;
- unsigned char *ptr;
-
- for (j=0;j<*size;j++)
- { hz = random_block_get ();
- ptr = (unsigned char *) &hz;
- if (*limit != 0)
- { /* make value positive , big endian , otherwise replace 3 with 0 */
- #if defined(SUN4) || defined(IBM) || defined(KSR1) || defined(SGI)
- ptr[0] = ptr[0] & 127;
- #endif
- #if defined(ALLIANT) || defined(GC) || defined(IPSC) || defined(SRM) || defined(MEIKO) || defined(OS2)
- ptr[3] = ptr[3] & 127;
- #endif
- a[j] = hz % *limit;
- }
- else
- a[j] = hz;
- }
- }
-
- void dalib_get_real_randoms__ (a, size, limit)
- float a[];
- int *size;
- float *limit;
-
- { int j, hz;
- unsigned char *ptr;
-
- /****************************************************************
- * *
- * float byte1 byte2 byte3 byte4 *
- * ============================================= *
- * *
- * 0.5 3f 00 00 00 *
- * 1.0 3f 80 00 00 *
- * 1.+ 3f 80 00 01 *
- * 2.- 3f ff ff ff *
- * 2.0 40 00 00 00 *
- * *
- ****************************************************************/
-
- for (j=0;j<*size;j++)
- { hz = random_block_get ();
- ptr = (unsigned char *) &hz;
- /* little endian, otherwise replace 3 with 0, 2 with 1 */
- #if defined(SUN4) || defined(IBM) || defined(KSR1) || defined(SGI)
- ptr[0] = 63; /* 3f */
- ptr[1] = ptr[1] | 128; /* 80 */
- #endif
- #if defined(ALLIANT) || defined(GC) || defined(IPSC) || defined(SRM) || defined(MEIKO) || defined(OS2)
- ptr[3] = 63; /* 3f */
- ptr[2] = ptr[1] | 128; /* 80 */
- #endif
- /* hz is now between 1.0 and 2.0 */
- a[j] = *((float *) &hz) - 1.0;
- }
- }
-
- void dalib_get_double_randoms__ (a, size, limit)
- double a[], limit;
- int *size;
- { int i, j;
- char *ptr;
- int hz[2];
-
- /****************************************************************
- * *
- * double byte1 byte2 byte3 byte4 ... byte8 *
- * ========================================================== *
- * *
- * 0.5 3f E0 00 00 .... 00 *
- * 1.0 3f F0 00 00 .... 00 *
- * 1.+ 3f F0 00 00 .... 01 *
- * 2.- 3f ff ff ff .... ff *
- * 2.0 40 00 00 00 .... 00 *
- * *
- ****************************************************************/
-
- for (j=0;j<*size;j++)
- { hz[0] = random_block_get ();
- hz[1] = random_block_get ();
- ptr = (char *) hz;
- /* big endian, otherwise replace 7 with 0, 6 with 1 */
- #if defined(SUN4) || defined(IBM) || defined(KSR1) || defined(SGI)
- ptr[0] = 63; /* 3f */
- ptr[1] = ptr[1] | 240; /* F0 */
- #endif
- #if defined(ALLIANT) || defined(GC) || defined(IPSC) || defined(SRM) || defined(MEIKO) || defined(OS2)
- ptr[7] = 63; /* 3f */
- ptr[6] = ptr[1] | 240; /* F0 */
- #endif
- /* hz is now between 1.0 and 2.0 */
- a[j] = *((double *) ptr) - (double)1.0;
- }
- }