#if 0
#include "config.h"
#include <stdio.h>
#include <stdlib.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <fcntl.h>
#ifdef HAVE_UNISTD_H
#include <unistd.h>
#endif
#include "libgfortran.h"
#define MT19937
#ifdef MT19937
enum constants
{
N = 624, M = 397, R = 19, TU = 11, TS = 7, TT = 15, TL = 17
};
#define M_A 0x9908B0DF
#define T_B 0x9D2C5680
#define T_C 0xEFC60000
#else
enum constants
{
N = 351, M = 175, R = 19, TU = 11, TS = 7, TT = 15, TL = 17
};
#define M_A 0xE4BD75F5
#define T_B 0x655E5280
#define T_C 0xFFD58000
#endif
static int i = N;
static unsigned int seed[N];
void
random_seed (GFC_INTEGER_4 * size, const gfc_array_i4 * put,
const gfc_array_i4 * get)
{
if (get == NULL && put == NULL && size == NULL)
{
int fd;
fd = open ("/dev/urandom", O_RDONLY);
if (fd == 0)
{
GFC_UINTEGER_4 s = (GFC_UINTEGER_4) seed;
for (i = 0; i < N; i++)
{
s = s * 29943829 - 1;
seed[i] = s;
}
}
else
{
read (fd, &seed[0], sizeof (GFC_UINTEGER_4) * N);
close (fd);
}
return;
}
if (size != NULL)
{
*size = N;
return;
}
if (put != NULL)
{
if (GFC_DESCRIPTOR_RANK (put) != 1)
abort ();
if (((put->dim[0].ubound + 1 - put->dim[0].lbound)) < N)
abort ();
if (put->dim[0].stride == 0)
return;
for (i = 0; i < N; i++)
seed[i] = put->data[i * put->dim[0].stride];
}
if (get != NULL)
{
if (GFC_DESCRIPTOR_RANK (get) != 1)
abort ();
if (((get->dim[0].ubound + 1 - get->dim[0].lbound)) < N)
abort ();
if (get->dim[0].stride == 0)
return;
for (i = 0; i < N; i++)
get->data[i * get->dim[0].stride] = seed[i];
}
}
static void
random_generate (void)
{
GFC_UINTEGER_4 y;
int k, m;
for (k = 0, m = M; k < N - 1; k++)
{
y = (seed[k] & (-1 << R)) | (seed[k + 1] & ((1u << R) - 1));
seed[k] = seed[m] ^ (y >> 1) ^ (-(GFC_INTEGER_4) (y & 1) & M_A);
if (++m >= N)
m = 0;
}
y = (seed[N - 1] & (-1 << R)) | (seed[0] & ((1u << R) - 1));
seed[N - 1] = seed[M - 1] ^ (y >> 1) ^ (-(GFC_INTEGER_4) (y & 1) & M_A);
i = 0;
}
#define random_r4 prefix(random_r4)
void
random_r4 (GFC_REAL_4 * harv)
{
if (i >= N)
random_generate ();
*harv = (GFC_REAL_4) ((GFC_REAL_4) (GFC_UINTEGER_4) seed[i++] /
(GFC_REAL_4) (~(GFC_UINTEGER_4) 0));
}
#define random_r8 prefix(random_r8)
void
random_r8 (GFC_REAL_8 * harv)
{
if ((i + 1) >= N)
random_generate ();
*harv = ((GFC_REAL_8) ((((GFC_UINTEGER_8) seed[i+1]) << 32) + seed[i])) /
(GFC_REAL_8) (~(GFC_UINTEGER_8) 0);
i += 2;
}
#define arandom_r4 prefix(arandom_r4)
void
arandom_r4 (gfc_array_r4 * harv)
{
index_type count[GFC_MAX_DIMENSIONS - 1];
index_type extent[GFC_MAX_DIMENSIONS - 1];
index_type stride[GFC_MAX_DIMENSIONS - 1];
index_type stride0;
index_type dim;
GFC_REAL_4 *dest;
int n;
dest = harv->data;
if (harv->dim[0].stride == 0)
harv->dim[0].stride = 1;
dim = GFC_DESCRIPTOR_RANK (harv);
for (n = 0; n < dim; n++)
{
count[n] = 0;
stride[n] = harv->dim[n].stride;
extent[n] = harv->dim[n].ubound + 1 - harv->dim[n].lbound;
if (extent[n] <= 0)
return;
}
stride0 = stride[0];
while (dest)
{
if (i >= N)
random_generate ();
*dest = (GFC_REAL_4) ((GFC_REAL_4) (GFC_UINTEGER_4) seed[i++] /
(GFC_REAL_4) (~(GFC_UINTEGER_4) 0));
dest += stride0;
count[0]++;
n = 0;
while (count[n] == extent[n])
{
count[n] = 0;
dest -= stride[n] * extent[n];
n++;
if (n == dim)
{
dest = NULL;
break;
}
else
{
count[n]++;
dest += stride[n];
}
}
}
}
#define arandom_r8 prefix(arandom_r8)
void
arandom_r8 (gfc_array_r8 * harv)
{
index_type count[GFC_MAX_DIMENSIONS - 1];
index_type extent[GFC_MAX_DIMENSIONS - 1];
index_type stride[GFC_MAX_DIMENSIONS - 1];
index_type stride0;
index_type dim;
GFC_REAL_8 *dest;
int n;
dest = harv->data;
if (harv->dim[0].stride == 0)
harv->dim[0].stride = 1;
dim = GFC_DESCRIPTOR_RANK (harv);
for (n = 0; n < dim; n++)
{
count[n] = 0;
stride[n] = harv->dim[n].stride;
extent[n] = harv->dim[n].ubound + 1 - harv->dim[n].lbound;
if (extent[n] <= 0)
return;
}
stride0 = stride[0];
while (dest)
{
if ((i + 1) >= N)
random_generate ();
*dest = ((GFC_REAL_8) ((((GFC_UINTEGER_8) seed[i+1]) << 32) + seed[i])) /
(GFC_REAL_8) (~(GFC_UINTEGER_8) 0);
i += 2;
dest += stride0;
count[0]++;
n = 0;
while (count[n] == extent[n])
{
count[n] = 0;
dest -= stride[n] * extent[n];
n++;
if (n == dim)
{
dest = NULL;
break;
}
else
{
count[n]++;
dest += stride[n];
}
}
}
}
#endif
#include "config.h"
#include "libgfortran.h"
#define GFC_SL(k, n) ((k)^((k)<<(n)))
#define GFC_SR(k, n) ((k)^((k)>>(n)))
static const GFC_INTEGER_4 kiss_size = 4;
#define KISS_DEFAULT_SEED {123456789, 362436069, 521288629, 916191069};
static const GFC_UINTEGER_4 kiss_default_seed[4] = KISS_DEFAULT_SEED;
static GFC_UINTEGER_4 kiss_seed[4] = KISS_DEFAULT_SEED;
static GFC_UINTEGER_4
kiss_random_kernel(void)
{
GFC_UINTEGER_4 kiss;
kiss_seed[0] = 69069 * kiss_seed[0] + 1327217885;
kiss_seed[1] = GFC_SL(GFC_SR(GFC_SL(kiss_seed[1],13),17),5);
kiss_seed[2] = 18000 * (kiss_seed[2] & 65535) + (kiss_seed[2] >> 16);
kiss_seed[3] = 30903 * (kiss_seed[3] & 65535) + (kiss_seed[3] >> 16);
kiss = kiss_seed[0] + kiss_seed[1] + (kiss_seed[2] << 16) + kiss_seed[3];
return kiss;
}
void
prefix(random_r4) (GFC_REAL_4 *x)
{
GFC_UINTEGER_4 kiss;
kiss = kiss_random_kernel ();
kiss_random_kernel ();
*x = normalize_r4_i4 (kiss, ~(GFC_UINTEGER_4) 0);
}
void
prefix(random_r8) (GFC_REAL_8 *x)
{
GFC_UINTEGER_8 kiss;
kiss = ((GFC_UINTEGER_8)kiss_random_kernel ()) << 32;
kiss += kiss_random_kernel ();
*x = normalize_r8_i8 (kiss, ~(GFC_UINTEGER_8) 0);
}
void
prefix(arandom_r4) (gfc_array_r4 *x)
{
index_type count[GFC_MAX_DIMENSIONS - 1];
index_type extent[GFC_MAX_DIMENSIONS - 1];
index_type stride[GFC_MAX_DIMENSIONS - 1];
index_type stride0;
index_type dim;
GFC_REAL_4 *dest;
int n;
dest = x->data;
if (x->dim[0].stride == 0)
x->dim[0].stride = 1;
dim = GFC_DESCRIPTOR_RANK (x);
for (n = 0; n < dim; n++)
{
count[n] = 0;
stride[n] = x->dim[n].stride;
extent[n] = x->dim[n].ubound + 1 - x->dim[n].lbound;
if (extent[n] <= 0)
return;
}
stride0 = stride[0];
while (dest)
{
prefix(random_r4) (dest);
dest += stride0;
count[0]++;
n = 0;
while (count[n] == extent[n])
{
count[n] = 0;
dest -= stride[n] * extent[n];
n++;
if (n == dim)
{
dest = NULL;
break;
}
else
{
count[n]++;
dest += stride[n];
}
}
}
}
void
prefix(arandom_r8) (gfc_array_r8 *x)
{
index_type count[GFC_MAX_DIMENSIONS - 1];
index_type extent[GFC_MAX_DIMENSIONS - 1];
index_type stride[GFC_MAX_DIMENSIONS - 1];
index_type stride0;
index_type dim;
GFC_REAL_8 *dest;
int n;
dest = x->data;
if (x->dim[0].stride == 0)
x->dim[0].stride = 1;
dim = GFC_DESCRIPTOR_RANK (x);
for (n = 0; n < dim; n++)
{
count[n] = 0;
stride[n] = x->dim[n].stride;
extent[n] = x->dim[n].ubound + 1 - x->dim[n].lbound;
if (extent[n] <= 0)
return;
}
stride0 = stride[0];
while (dest)
{
prefix(random_r8) (dest);
dest += stride0;
count[0]++;
n = 0;
while (count[n] == extent[n])
{
count[n] = 0;
dest -= stride[n] * extent[n];
n++;
if (n == dim)
{
dest = NULL;
break;
}
else
{
count[n]++;
dest += stride[n];
}
}
}
}
void
random_seed (GFC_INTEGER_4 *size, gfc_array_i4 * put,
gfc_array_i4 * get)
{
int i;
if (size == NULL && put == NULL && get == NULL)
{
kiss_seed[0] = kiss_default_seed[0];
kiss_seed[1] = kiss_default_seed[1];
kiss_seed[2] = kiss_default_seed[2];
kiss_seed[3] = kiss_default_seed[3];
}
if (size != NULL)
*size = kiss_size;
if (put != NULL)
{
if (GFC_DESCRIPTOR_RANK (put) != 1)
runtime_error ("Array rank of PUT is not 1.");
if (((put->dim[0].ubound + 1 - put->dim[0].lbound)) < kiss_size)
runtime_error ("Array size of PUT is too small.");
if (put->dim[0].stride == 0)
put->dim[0].stride = 1;
for (i = 0; i < kiss_size; i++)
kiss_seed[i] =(GFC_UINTEGER_4) put->data[i * put->dim[0].stride];
}
if (get != NULL)
{
if (GFC_DESCRIPTOR_RANK (get) != 1)
runtime_error ("Array rank of GET is not 1.");
if (((get->dim[0].ubound + 1 - get->dim[0].lbound)) < kiss_size)
runtime_error ("Array size of GET is too small.");
if (get->dim[0].stride == 0)
get->dim[0].stride = 1;
for (i = 0; i < kiss_size; i++)
get->data[i * get->dim[0].stride] = (GFC_INTEGER_4) kiss_seed[i];
}
}