/* Implementation of the RANDOM intrinsics
- Copyright 2002, 2004, 2005, 2006, 2007 Free Software Foundation, Inc.
+ Copyright 2002, 2004, 2005, 2006, 2007, 2009 Free Software Foundation, Inc.
Contributed by Lars Segerlund <seger@linuxmail.org>
and Steve Kargl.
Libgfortran is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public
License as published by the Free Software Foundation; either
-version 2 of the License, or (at your option) any later version.
-
-In addition to the permissions in the GNU General Public License, the
-Free Software Foundation gives you unlimited permission to link the
-compiled version of this file into combinations with other programs,
-and to distribute those combinations without any restriction coming
-from the use of this file. (The General Public License restrictions
-do apply in other respects; for example, they cover modification of
-the file, and distribution when not linked into a combine
-executable.)
+version 3 of the License, or (at your option) any later version.
Ligbfortran is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
-You should have received a copy of the GNU General Public
-License along with libgfortran; see the file COPYING. If not,
-write to the Free Software Foundation, Inc., 51 Franklin Street, Fifth Floor,
-Boston, MA 02110-1301, USA. */
+Under Section 7 of GPL version 3, you are granted additional
+permissions described in the GCC Runtime Library Exception, version
+3.1, as published by the Free Software Foundation.
+
+You should have received a copy of the GNU General Public License and
+a copy of the GCC Runtime Library Exception along with this program;
+see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
+<http://www.gnu.org/licenses/>. */
#include "libgfortran.h"
#include <gthr.h>
GFC_REAL_* types in the range of [0,1). If GFC_REAL_*_RADIX are 2
or 16, respectively, we mask off the bits that don't fit into the
correct GFC_REAL_*, convert to the real type, then multiply by the
- correct offset.
-*/
+ correct offset. */
static inline void
We do this by using three generators with different seeds, the
first one always for the most significant bits, the second one
for bits 33..64 (if present in the REAL kind), and the third one
- (called twice) for REAL(16).
-*/
+ (called twice) for REAL(16). */
#define GFC_SL(k, n) ((k)^((k)<<(n)))
#define GFC_SR(k, n) ((k)^((k)>>(n)))
with 0<=x<2^32, 0<y<2^32, 0<=z<2^32, 0<=c<698769069
except that the two pairs
z=0,c=0 and z=2^32-1,c=698769068
- should be avoided.
-*/
+ should be avoided. */
+
+/* Any modifications to the seeds that change kiss_size below need to be
+ reflected in check.c (gfc_check_random_seed) to enable correct
+ compile-time checking of PUT size for the RANDOM_SEED intrinsic. */
#define KISS_DEFAULT_SEED_1 123456789, 362436069, 521288629, 316191069
#define KISS_DEFAULT_SEED_2 987654321, 458629013, 582859209, 438195021
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;
+ stride[n] = GFC_DESCRIPTOR_STRIDE(x,n);
+ extent[n] = GFC_DESCRIPTOR_EXTENT(x,n);
if (extent[n] <= 0)
return;
}
while (dest)
{
- /* random_r4 (dest); */
+ /* random_r4 (dest); */
kiss = kiss_random_kernel (kiss_seed_1);
rnumber_4 (dest, kiss);
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;
+ stride[n] = GFC_DESCRIPTOR_STRIDE(x,n);
+ extent[n] = GFC_DESCRIPTOR_EXTENT(x,n);
if (extent[n] <= 0)
return;
}
while (dest)
{
- /* random_r8 (dest); */
+ /* random_r8 (dest); */
kiss = ((GFC_UINTEGER_8) kiss_random_kernel (kiss_seed_1)) << 32;
kiss += kiss_random_kernel (kiss_seed_2);
rnumber_8 (dest, kiss);
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;
+ stride[n] = GFC_DESCRIPTOR_STRIDE(x,n);
+ extent[n] = GFC_DESCRIPTOR_EXTENT(x,n);
if (extent[n] <= 0)
return;
}
while (dest)
{
- /* random_r10 (dest); */
+ /* random_r10 (dest); */
kiss = ((GFC_UINTEGER_8) kiss_random_kernel (kiss_seed_1)) << 32;
kiss += kiss_random_kernel (kiss_seed_2);
rnumber_10 (dest, kiss);
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;
+ stride[n] = GFC_DESCRIPTOR_STRIDE(x,n);
+ extent[n] = GFC_DESCRIPTOR_EXTENT(x,n);
if (extent[n] <= 0)
return;
}
while (dest)
{
- /* random_r16 (dest); */
+ /* random_r16 (dest); */
kiss1 = ((GFC_UINTEGER_8) kiss_random_kernel (kiss_seed_1)) << 32;
kiss1 += kiss_random_kernel (kiss_seed_2);
#endif
+
+
+static void
+scramble_seed (unsigned char *dest, unsigned char *src, int size)
+{
+ int i;
+
+ for (i = 0; i < size; i++)
+ dest[(i % 2) * (size / 2) + i / 2] = src[i];
+}
+
+
+static void
+unscramble_seed (unsigned char *dest, unsigned char *src, int size)
+{
+ int i;
+
+ for (i = 0; i < size; i++)
+ dest[i] = src[(i % 2) * (size / 2) + i / 2];
+}
+
+
+
/* random_seed is used to seed the PRNG with either a default
set of seeds or user specified set of seeds. random_seed
must be called with no argument or exactly one argument. */
random_seed_i4 (GFC_INTEGER_4 *size, gfc_array_i4 *put, gfc_array_i4 *get)
{
int i;
+ unsigned char seed[4*kiss_size];
__gthread_mutex_lock (&random_lock);
runtime_error ("Array rank of PUT is not 1.");
/* If the array is too small, abort. */
- if (((put->dim[0].ubound + 1 - put->dim[0].lbound)) < kiss_size)
+ if (GFC_DESCRIPTOR_EXTENT(put,0) < kiss_size)
runtime_error ("Array size of PUT is too small.");
- /* This code now should do correct strides. */
+ /* We copy the seed given by the user. */
for (i = 0; i < kiss_size; i++)
- kiss_seed[i] = (GFC_UINTEGER_4) put->data[i * put->dim[0].stride];
+ memcpy (seed + i * sizeof(GFC_UINTEGER_4),
+ &(put->data[(kiss_size - 1 - i) * GFC_DESCRIPTOR_STRIDE(put,0)]),
+ sizeof(GFC_UINTEGER_4));
+
+ /* We put it after scrambling the bytes, to paper around users who
+ provide seeds with quality only in the lower or upper part. */
+ scramble_seed ((unsigned char *) kiss_seed, seed, 4*kiss_size);
}
/* Return the seed to GET data. */
runtime_error ("Array rank of GET is not 1.");
/* If the array is too small, abort. */
- if (((get->dim[0].ubound + 1 - get->dim[0].lbound)) < kiss_size)
+ if (GFC_DESCRIPTOR_EXTENT(get,0) < kiss_size)
runtime_error ("Array size of GET is too small.");
- /* This code now should do correct strides. */
+ /* Unscramble the seed. */
+ unscramble_seed (seed, (unsigned char *) kiss_seed, 4*kiss_size);
+
+ /* Then copy it back to the user variable. */
for (i = 0; i < kiss_size; i++)
- get->data[i * get->dim[0].stride] = (GFC_INTEGER_4) kiss_seed[i];
+ memcpy (&(get->data[(kiss_size - 1 - i) * GFC_DESCRIPTOR_STRIDE(get,0)]),
+ seed + i * sizeof(GFC_UINTEGER_4),
+ sizeof(GFC_UINTEGER_4));
}
__gthread_mutex_unlock (&random_lock);
runtime_error ("Array rank of PUT is not 1.");
/* If the array is too small, abort. */
- if (((put->dim[0].ubound + 1 - put->dim[0].lbound)) < kiss_size / 2)
+ if (GFC_DESCRIPTOR_EXTENT(put,0) < kiss_size / 2)
runtime_error ("Array size of PUT is too small.");
/* This code now should do correct strides. */
for (i = 0; i < kiss_size / 2; i++)
- memcpy (&kiss_seed[2*i], &(put->data[i * put->dim[0].stride]),
+ memcpy (&kiss_seed[2*i], &(put->data[i * GFC_DESCRIPTOR_STRIDE(put,0)]),
sizeof (GFC_UINTEGER_8));
}
runtime_error ("Array rank of GET is not 1.");
/* If the array is too small, abort. */
- if (((get->dim[0].ubound + 1 - get->dim[0].lbound)) < kiss_size / 2)
+ if (GFC_DESCRIPTOR_EXTENT(get,0) < kiss_size / 2)
runtime_error ("Array size of GET is too small.");
/* This code now should do correct strides. */
for (i = 0; i < kiss_size / 2; i++)
- memcpy (&(get->data[i * get->dim[0].stride]), &kiss_seed[2*i],
+ memcpy (&(get->data[i * GFC_DESCRIPTOR_STRIDE(get,0)]), &kiss_seed[2*i],
sizeof (GFC_UINTEGER_8));
}
OSDN Git Service
