+/* Copyright (c) 2012, 2014. The SimGrid Team.
+ * All rights reserved. */
+
+/* This program is free software; you can redistribute it and/or modify it
+ * under the terms of the license (GNU LGPL) which comes with this package. */
+
/***********************************************************************\
*
* File: RngStream.c for multiple streams of Random Numbers
* Language: ANSI C
* Copyright: Pierre L'Ecuyer, University of Montreal
* Date: 14 August 2001
- * License: GPL version 2 or later
+ * License: GPL version 2 or later
*
* Notice: Please contact P. L'Ecuyer at <lecuyer@iro.UMontreal.ca>
* for commercial purposes.
*
\***********************************************************************/
-
#include "xbt/RngStream.h"
+#include "xbt/sysdep.h"
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
-
/*---------------------------------------------------------------------*/
/* Private part. */
/*---------------------------------------------------------------------*/
-
#define norm 2.328306549295727688e-10
#define m1 4294967087.0
#define m2 4294944443.0
#define two53 9007199254740992.0
#define fact 5.9604644775390625e-8 /* 1 / 2^24 */
-
-
-
-/* Default initial seed of the package. Will be updated to become
- the seed of the next created stream. */
+/* Default initial seed of the package. Will be updated to become the seed of the next created stream. */
static double nextSeed[6] = { 12345, 12345, 12345, 12345, 12345, 12345 };
-
/* The following are the transition matrices of the two MRG components */
/* (in matrix form), raised to the powers -1, 1, 2^76, and 2^127, resp.*/
static double InvA1[3][3] = { /* Inverse of A1p0 */
};
static double A1p76[3][3] = {
- { 82758667.0, 1871391091.0, 4127413238.0 },
- { 3672831523.0, 69195019.0, 1871391091.0 },
+ { 82758667.0, 1871391091.0, 4127413238.0 },
+ { 3672831523.0, 69195019.0, 1871391091.0 },
{ 3672091415.0, 3528743235.0, 69195019.0 }
};
static double A2p76[3][3] = {
- { 1511326704.0, 3759209742.0, 1610795712.0 },
- { 4292754251.0, 1511326704.0, 3889917532.0 },
+ { 1511326704.0, 3759209742.0, 1610795712.0 },
+ { 4292754251.0, 1511326704.0, 3889917532.0 },
{ 3859662829.0, 4292754251.0, 3708466080.0 }
};
static double A1p127[3][3] = {
- { 2427906178.0, 3580155704.0, 949770784.0 },
+ { 2427906178.0, 3580155704.0, 949770784.0 },
{ 226153695.0, 1230515664.0, 3580155704.0 },
{ 1988835001.0, 986791581.0, 1230515664.0 }
};
{ 2824425944.0, 32183930.0, 2093834863.0 }
};
-
-
-
-
-/*-------------------------------------------------------------------------*/
-
static double MultModM (double a, double s, double c, double m)
/* Compute (a*s + c) % m. m must be < 2^35. Works also for s, c < 0 */
{
return v;
}
-
-/*-------------------------------------------------------------------------*/
-
static void MatVecModM (double A[3][3], double s[3], double v[3], double m)
/* Returns v = A*s % m. Assumes that -m < s[i] < m. */
/* Works even if v = s. */
v[i] = x[i];
}
-
-/*-------------------------------------------------------------------------*/
-
-static void MatMatModM (double A[3][3], double B[3][3], double C[3][3],
- double m)
+static void MatMatModM (double A[3][3], double B[3][3], double C[3][3], double m)
/* Returns C = A*B % m. Work even if A = C or B = C or A = B = C. */
{
int i, j;
}
}
-
-/*-------------------------------------------------------------------------*/
-
static void MatTwoPowModM (double A[3][3], double B[3][3], double m, long e)
/* Compute matrix B = (A^(2^e) % m); works even if A = B */
{
MatMatModM (B, B, B, m);
}
-
-/*-------------------------------------------------------------------------*/
-
static void MatPowModM (double A[3][3], double B[3][3], double m, long n)
/* Compute matrix B = A^n % m ; works even if A = B */
{
}
}
-
-/*-------------------------------------------------------------------------*/
-
static double U01 (RngStream g)
{
long k;
return (g->Anti) ? (1 - u) : u;
}
-
-/*-------------------------------------------------------------------------*/
-
static double U01d (RngStream g)
{
double u;
}
}
-
-/*-------------------------------------------------------------------------*/
-
static int CheckSeed (unsigned long seed[6])
{
- /* Check that the seeds are legitimate values. Returns 0 if legal seeds,
- -1 otherwise */
+ /* Check that the seeds are legitimate values. Returns 0 if legal seeds, -1 otherwise */
int i;
for (i = 0; i < 3; ++i) {
if (seed[i] >= m1) {
- fprintf (stderr, "****************************************\n"
- "ERROR: Seed[%1d] >= m1, Seed is not set.\n"
- "****************************************\n\n", i);
- return (-1);
+ fprintf (stderr, "****************************************\n"
+ "ERROR: Seed[%1d] >= m1, Seed is not set.\n"
+ "****************************************\n\n", i);
+ return (-1);
}
}
for (i = 3; i < 6; ++i) {
if (seed[i] >= m2) {
- fprintf (stderr, "****************************************\n"
- "ERROR: Seed[%1d] >= m1, Seed is not set.\n"
- "****************************************\n\n", i);
- return (-1);
+ fprintf (stderr, "****************************************\n"
+ "ERROR: Seed[%1d] >= m1, Seed is not set.\n"
+ "****************************************\n\n", i);
+ return (-1);
}
}
if (seed[0] == 0 && seed[1] == 0 && seed[2] == 0) {
fprintf (stderr, "****************************\n"
- "ERROR: First 3 seeds = 0.\n"
- "****************************\n\n");
+ "ERROR: First 3 seeds = 0.\n"
+ "****************************\n\n");
return (-1);
}
if (seed[3] == 0 && seed[4] == 0 && seed[5] == 0) {
fprintf (stderr, "****************************\n"
- "ERROR: Last 3 seeds = 0.\n"
- "****************************\n\n");
+ "ERROR: Last 3 seeds = 0.\n"
+ "****************************\n\n");
return (-1);
}
return 0;
}
-
/*---------------------------------------------------------------------*/
/* Public part. */
/*---------------------------------------------------------------------*/
-
RngStream RngStream_CreateStream (const char name[])
{
int i;
RngStream g;
size_t len;
- g = (RngStream) malloc (sizeof (struct RngStream_InfoState));
+ g = (RngStream) xbt_malloc (sizeof (struct RngStream_InfoState));
if (g == NULL) {
printf ("RngStream_CreateStream: No more memory\n\n");
exit (EXIT_FAILURE);
}
if (name) {
len = strlen (name);
- g->name = (char *) malloc ((len + 1) * sizeof (char));
+ g->name = (char *) xbt_malloc ((len + 1) * sizeof (char));
strncpy (g->name, name, len + 1);
} else
g->name = 0;
return g;
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_DeleteStream (RngStream * p)
{
if (*p == NULL)
*p = NULL;
}
-/*-------------------------------------------------------------------------*/
+RngStream RngStream_CopyStream (const RngStream src)
+{
+ RngStream g;
+
+ if(src == NULL) {
+ printf ("RngStream_CopyStream: 'src' not initialized\n\n");
+ exit (EXIT_FAILURE);
+ }
+
+ g = (RngStream) xbt_malloc (sizeof (struct RngStream_InfoState));
+ if (g == NULL) {
+ printf ("RngStream_CopyStream: No more memory\n\n");
+ exit (EXIT_FAILURE);
+ }
+ memcpy((void*) g, (void*) src, sizeof (struct RngStream_InfoState));
+
+ return g;
+}
void RngStream_ResetStartStream (RngStream g)
{
g->Cg[i] = g->Bg[i] = g->Ig[i];
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_ResetNextSubstream (RngStream g)
{
int i;
g->Cg[i] = g->Bg[i];
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_ResetStartSubstream (RngStream g)
{
int i;
g->Cg[i] = g->Bg[i];
}
-/*-------------------------------------------------------------------------*/
-
int RngStream_SetPackageSeed (unsigned long seed[6])
{
int i;
return 0; /* SUCCESS */
}
-/*-------------------------------------------------------------------------*/
-
int RngStream_SetSeed (RngStream g, unsigned long seed[6])
{
int i;
return -1; /* FAILURE */
for (i = 0; i < 6; ++i)
g->Cg[i] = g->Bg[i] = g->Ig[i] = seed[i];
- return 0; /* SUCCESS */
+ return 0; /* SUCCESS */
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_AdvanceState (RngStream g, long e, long c)
{
double B1[3][3], C1[3][3], B2[3][3], C2[3][3];
MatVecModM (C2, &g->Cg[3], &g->Cg[3], m2);
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_GetState (RngStream g, unsigned long seed[6])
{
int i;
seed[i] = g->Cg[i];
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_WriteState (RngStream g)
{
int i;
printf ("%lu }\n\n", (unsigned long) g->Cg[5]);
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_WriteStateFull (RngStream g)
{
int i;
printf ("%lu }\n\n", (unsigned long) g->Cg[5]);
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_IncreasedPrecis (RngStream g, int incp)
{
g->IncPrec = incp;
}
-/*-------------------------------------------------------------------------*/
-
void RngStream_SetAntithetic (RngStream g, int a)
{
g->Anti = a;
}
-/*-------------------------------------------------------------------------*/
-
double RngStream_RandU01 (RngStream g)
{
if (g->IncPrec)
return U01 (g);
}
-/*-------------------------------------------------------------------------*/
-
int RngStream_RandInt (RngStream g, int i, int j)
{
return i + (int) ((j - i + 1.0) * RngStream_RandU01 (g));
}
-
-/* Undefine this terms, or supernovae build will fail. */
-#undef norm
-#undef m1
-#undef m2
-#undef a12
-#undef a13n
-#undef a21
-#undef a23n
-#undef two17
-#undef two53
-#undef fact
