\documentclass[12pt]{article}


\marginparwidth 0pt
\marginparsep 0pt
\topskip 0pt
\headsep 0pt
\headheight 0pt
\oddsidemargin 0pt
\evensidemargin 0pt
\textwidth 6.5in
\topmargin 0.0in
\textheight 9.0in


% Macros pour inserer des bouts de code (programmes).
% Faire  \code ...  \endcode
{\obeyspaces\gdef {\ }}
\def\setverbatim{\def\par{\leavevmode\endgraf}
            \parskip=0pt\parindent=0pt\obeylines\obeyspaces }
\chardef\other=12
\def\ttverbatim{\setverbatim\tt
       \catcode`\{=\other \catcode`\}=\other \catcode`\_=\other
       \catcode`\^=\other \catcode`\$=\other \catcode`\%=\other
       \catcode`\#=\other \catcode`\&=\other \baselineskip=11pt
       }  %$
    % Reproduit tel quel ce qui est ecrit, en caracteres \tt.
    % On doit faire  \begingroup\ttverbatim   ....  \endgroup
\def\code {\vfil\vfilneg\vbox\bgroup\ttverbatim}
\def\longcode {\vfil\vfilneg\bgroup\ttverbatim}
\def\smallc {\small\tt\baselineskip=9.5pt}
\let\endcode=\egroup
\def\parup{\nobreak\vskip -2pt\nobreak}
\def\tab{\small\parindent=0pt\advance\leftskip by 1.5em\parup}
\def\tabb{\small\parindent=0pt\advance\leftskip by 3.0em\parup}
\def\endtab{\vskip 0.01pt\advance\leftskip by -1.5em\normalsize}
\def\endtabb{\vskip 0.01pt\advance\leftskip by -3.0em\normalsize}
%
\def\hide{\iffalse}
\let\endhide=\fi
\catcode`\@=11  %%%%%%%%%%%

% For the bibliography.

\newenvironment {refer}{
  \begin {list}{}{
  \def\newblock {\hskip .11em plus .33em minus .07em}
  \setlength {\itemsep}{2pt} 
  \setlength {\parsep}{0pt}
  \setlength {\parindent}{20pt}
  \setlength {\leftmargin}{+\parindent}
  \setlength {\itemindent}{-\parindent}}}{\end{list}}
\def\@biblabel#1{}

%\def\cite{\@ifnextchar [{\@tempswatrue\@citex}{\@tempswafalse\@citex[]}}

\def\citep{\@ifnextchar [{\@tempswatrue\@citexp}{\@tempswafalse\@citexp[]}}

\def\@cite#1#2{\null{#1\if@tempswa , #2\fi}\relax }

\def\@citex[#1]#2{\if@filesw\immediate\write\@auxout{\string\citation{#2}}\fi
  \def\@citea{}\@cite{\@for\@citeb:=#2\do
    {\@citea\def\@citea{, }\@ifundefined
       {b@\@citeb}{{\bf ?}\@warning
       {Citation `\@citeb' on page \thepage \space undefined}}%
  {\csname b@\@citeb\endcsname}}}{#1}}

\def\@citexp[#1]#2{\if@filesw\immediate\write\@auxout{\string\citation{#2}}\fi
  \def\@citea{}\@cite{\@for\@citeb:=#2\do
    {\@citea\def\@citea{, }\@ifundefined
       {p@\@citeb}{{\bf ?}\@warning
       {Citation `\@citeb' on page \thepage \space undefined}}%
  {\csname p@\@citeb\endcsname}}}{#1}}

\def\bibcite#1#2#3{\global\@namedef{b@#1}{#2\ #3}}

\def\bibcitep#1#2#3{\global\@namedef{p@#1}{#2\ (#3)}}

\def\@lbibitem[#1][#2]#3{\item[\@biblabel{#1}]\if@filesw
        { \def\protect##1{\string##1\space}
          \immediate\write\@auxout{\string\bibcite{#3}{#1}{#2}}
          \immediate\write\@auxout{\string\bibcitep{#3}{#1}{#2}}
          \fi\ignorespaces}}


\begin{document}

\begin{center}
{\large\bf AN OBJECT-ORIENTED RANDOM-NUMBER PACKAGE IN JAVA
 WITH MANY LONG STREAMS AND SUBSTREAMS}
\end{center}


\bigskip


This class implements the {\tt RngStream} class, with a few
additional tools.
The backbone generator is the combined multiple recursive 
generator (CMRG) {\tt Mrg32k3a} proposed in (\cite{rLEC99b}),
implemented in 64-bit floating-point arithmetic.
This backbone generator has period length $\rho\approx 2^{191}$.
The values of $V$, $W$, and $Z$ are $2^{51}$, $2^{76}$, and $2^{127}$,
respectively. % (See {\tt RandomStream} for their definition.)
The seed of the RNG, and the state of a stream at any given step,
are 6-dimensional vectors of 32-bit integers.
The default initial seed of the package is 
$(12345, 12345, 12345, 12345, 12345,$  $12345)$.


\iffalse  %%%%%%%%
The procedure {\tt getState} returns the state of a stream.
One can change the state of a given stream, without modifying
the state of the other streams, by calling {\tt setSeed} 
or {\tt advanceState}.
However, after calling {\tt setSeed} for a given stream,
the initial states of the different streams are no longer spaced
$Z$ values apart.  Therefore, this procedure should be called
very sparingly, if at all.  The procedure {\tt resetStream}
suffices for almost all applications.
\fi  %%%%%%%%%%


%%%%%%%%%%%%%%%%%%%%%%%
\bigskip\hrule

\code\hide
/**
 * Title:          RngStream.java
 * Description:    Multiple Streams and Substreams of Random Numbers
 * Copyright:      Pierre L'Ecuyer, University of Montreal
 * Notice:         This code can be used freely for personal, academic,
 *                 or non-commercial purposes. For commercial purposes, 
 *                 please contact P. L'Ecuyer at: lecuyer@iro.UMontreal.ca
 * Version         1.0
 * Date:           14 August 2001
 */
\endhide


public class RngStream {\hide

// %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Private constants.


   private static final double norm   = 2.328306549295727688e-10;
   private static final double m1     = 4294967087.0;
   private static final double m2     = 4294944443.0;
   private static final double a12    =  1403580.0;
   private static final double a13n   =   810728.0;
   private static final double a21    =   527612.0;
   private static final double a23n   =   1370589.0;
   private static final double two17    =  131072.0;
   private static final double two53    =  9007199254740992.0;
   private static final double invtwo24 = 5.9604644775390625e-8;

   private static final double InvA1[][] = {   // Inverse of A1p0
      { 184888585.0, 0.0, 1945170933.0 },
      {         1.0, 0.0,          0.0 },
      {         0.0, 1.0,          0.0 }
      };
   private static final double InvA2[][] = {   // Inverse of A2p0
      { 0.0, 360363334.0, 4225571728.0 },
      { 1.0,         0.0,          0.0 },
      { 0.0,         1.0,          0.0 }
      };
   private static final double A1p0[][]  =  {
      {       0.0,       1.0,      0.0 },
      {       0.0,       0.0,      1.0 },
      { -810728.0, 1403580.0,      0.0 }
      };
   private static final double A2p0[][]  =  {
      {        0.0,   1.0,         0.0 },
      {        0.0,   0.0,         1.0 },
      { -1370589.0,   0.0,    527612.0 }
      };
   private static final double A1p76[][] = {
      {      82758667.0, 1871391091.0, 4127413238.0 },
      {    3672831523.0,   69195019.0, 1871391091.0 },
      {    3672091415.0, 3528743235.0,   69195019.0 }
      };
   private static final double A2p76[][] = {
      {    1511326704.0, 3759209742.0, 1610795712.0 },
      {    4292754251.0, 1511326704.0, 3889917532.0 },
      {    3859662829.0, 4292754251.0, 3708466080.0 }
      };
   private static final double A1p127[][] = {
      {    2427906178.0, 3580155704.0,  949770784.0 },
      {     226153695.0, 1230515664.0, 3580155704.0 },
      {    1988835001.0,  986791581.0, 1230515664.0 }
      };
   private static final double A2p127[][] = {
      {    1464411153.0,  277697599.0, 1610723613.0 },
      {      32183930.0, 1464411153.0, 1022607788.0 },
      {    2824425944.0,   32183930.0, 2093834863.0 }
      };


// %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Private variables (fields) for each stream.

   private static double nextSeed[] = {12345, 12345, 12345, 12345, 12345, 12345};
   // Default seed of the package and seed for the next stream to be created.

   private double Cg[] = new double[6];
   private double Bg[] = new double[6];
   private double Ig[] = new double[6];
   // The arrays {\tt Cg}, {\tt Bg}, and {\tt Ig} contain the current state, 
   // the starting point of the current substream,
   // and the starting point of the stream, respectively.

   private boolean anti;
   // This stream generates antithetic variates if 
   // and only if {\tt anti = true}.

   private boolean prec53;
   // The precision of the output numbers is ``increased'' (see
   // {\tt increasedPrecis}) if and only if {\tt prec53 = true}.

   private String descriptor;
   // Describes the stream (for writing the state, error messages, etc.).

// %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
// Private methods

    //--------------------------------------------------------------
    /* Compute (a*s + c) MOD m ; m must be < 2^35 */
    /* Works also for s, c < 0.                   */
    private static double multModM 
           (double a, double s, double c, double m) {
        double v;
        int a1;
        v = a * s + c;
        if (v >= two53 || v <= -two53 ) {
                a1 = (int)(a / two17);  a -= a1 * two17;
                v  = a1 * s;
                a1 = (int)(v / m);    v -= a1 * m;
                v  = v * two17 + a * s + c;
        }
        a1 = (int)(v / m);
        if ((v -= a1 * m) < 0.0) return v += m; else return v;
    }


    //-----------------------------------------------------------
    /* Returns v = A*s MOD m.  Assumes that -m < s[i] < m. */
    /* Works even if v = s.                                */
    private static void matVecModM (double A[][], double s[], 
                                    double v[], double m) {
            int i;
            double x[] = new double[3];
            for (i = 0; i < 3;  ++i) {
                x[i] = multModM (A[i][0], s[0], 0.0, m);
                x[i] = multModM (A[i][1], s[1], x[i], m);
                x[i] = multModM (A[i][2], s[2], x[i], m);
            }
            for (i = 0; i < 3;  ++i)  v[i] = x[i];
    }

    //------------------------------------------------------------
    /* Returns C = A*B MOD m */
    /* Note: work even if A = C or B = C or A = B = C.         */
    private static void matMatModM (double A[][], double B[][], 
                                    double C[][], double m){
             int i, j;
             double V[] = new double[3], W[][] = new double[3][3];
             for (i = 0; i < 3;  ++i) {
                    for (j = 0; j < 3;  ++j)  V[j] = B[j][i];
                            matVecModM (A, V, V, m);
                    for (j = 0; j < 3;  ++j)  W[j][i] = V[j];
             }
             for (i = 0; i < 3;  ++i) {
                    for (j = 0; j < 3;  ++j)
                        C[i][j] = W[i][j];
             }
    }

    //-------------------------------------------------------------
    /* Compute matrix B = (A^(2^e) Mod m);  works even if A = B */
    private static void matTwoPowModM (double A[][], double B[][], 
                                       double m, int e) {
            int i, j;
            /* initialize: B = A */
            if (A != B) {
                for (i = 0; i < 3; i++) {
                    for (j = 0; j < 3;  ++j)  B[i][j] = A[i][j];
                }
            }
            /* Compute B = A^{2^e} */
            for (i = 0; i < e; i++) matMatModM (B, B, B, m);
    }

    //-------------------------------------------------------------
    /* Compute matrix D = A^c Mod m ;  works even if A = B */
    private static void matPowModM (double A[][], double B[][], 
                                    double m, int c){
            int i, j;
            int n = c;
            double W[][] = new double[3][3];

            /* initialize: W = A; B = I */
            for (i = 0; i < 3; i++) {
                for (j = 0; j < 3;  ++j)  {
                    W[i][j] = A[i][j];
                    B[i][j] = 0.0;
                }
            }
            for (j = 0; j < 3;  ++j)   B[j][j] = 1.0;

            /* Compute B = A^c mod m using the binary decomp. of c */
            while (n > 0) {
                if ((n % 2)==1) matMatModM (W, B, B, m);
                matMatModM (W, W, W, m);
                n /= 2;
            }
    } 

   //-------------------------------------------------------------
   // Generate a uniform random number, with 32 bits of resolution.
   private double U01 () {
        int k;
        double p1, p2, u;
        /* Component 1 */
        p1 = a12 * Cg[1] - a13n * Cg[0];
        k = (int)(p1 / m1);
        p1 -= k * m1;
        if (p1 < 0.0) p1 += m1;
        Cg[0] = Cg[1];   Cg[1] = Cg[2];   Cg[2] = p1;
        /* Component 2 */
        p2 = a21 * Cg[5] - a23n * Cg[3];
        k  = (int)(p2 / m2);
        p2 -= k * m2;
        if (p2 < 0.0) p2 += m2;
        Cg[3] = Cg[4];   Cg[4] = Cg[5];   Cg[5] = p2;
        /* Combination */
        u = ((p1 > p2) ? (p1 - p2) * norm : (p1 - p2 + m1) * norm);
        return (anti) ? (1 - u) : u;
   }

   //-------------------------------------------------------------
   // Generate a uniform random number, with 52 bits of resolution.
   private double U01d () {
        double u = U01();
        if (anti) {
            // Antithetic case: note that U01 already returns 1-u.
            u += (U01() - 1.0) * invtwo24;
            return (u < 0.0) ? u + 1.0 : u;
        } else {
            u += U01() * invtwo24;
            return (u < 1.0) ? u : (u - 1.0);
        }
   }
\endhide
\endcode

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%\subsubsection* {Constructors}

\code

   public RngStream () \hide{
      anti = false;
      prec53 = false;
      for (int i = 0; i < 6; ++i)  
         Bg[i] = Cg[i] = Ig[i] = nextSeed[i];
      matVecModM (A1p127, nextSeed, nextSeed, m1);
      double temp[] = new double[3];
      for (int i = 0; i < 3; ++i)  
         temp[i] = nextSeed[i + 3];
      matVecModM (A2p127, temp, temp, m2);
      for (int i = 0; i < 3; ++i)  
         nextSeed[i + 3] = temp[i];
   }\endhide
\endcode
 \tabb Constructs a new stream, initializes its seed $I_g$,
   sets $B_g$ and $C_g$ equal to $I_g$.
   The seed $I_g$ is equal to the initial seed of the package given by 
   {\tt setPackageSeed} if this is the first created stream,
   otherwise it is $Z$ steps ahead of the seed of the most recently
   created stream. Also sets its antithetic and
   increased precision switches to {\tt false}.
 \endtabb
\code

   public RngStream (String name) \hide{
      this ();
      descriptor = name;
   }\endhide
\endcode
 \tabb  Constructs a new stream with identifier {\tt name}
   (used when printing the full state of the stream).
 \endtabb
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%\subsubsection* {Methods}
\code

   public static boolean setPackageSeed (long seed[]) \hide{
      // Must use long because there is no unsigned int type.
        if (CheckSeed (seed) != 0)
            return false;                   // FAILURE     
      for (int i = 0; i < 6;  ++i)  nextSeed[i] = seed[i];
          return true;                     // SUCCESS
   }\endhide
\endcode
  \tabb  Sets the initial seed for the class {\tt RngStream} to the 
   six integers in the vector {\tt seed[0..5]}.
   This will be the seed (initial state) of the first stream.
   If this procedure is not called, the default initial seed
   is $(12345, 12345, 12345, 12345, 12345, 12345)$.
   If it is called, the first 3 values of the seed must all be
   less than $m_1 = 4294967087$, and not all 0;
   and the last 3 values 
   must all be less than $m_2 = 4294944443$, and not all 0.
   Returns {\tt false} for invalid seeds, and {\tt true} otherwise.
 \endtabb
\code

   public void resetStartStream () \hide{
      for (int i = 0; i < 6;  ++i)  Cg[i] = Bg[i] = Ig[i];
   } \endhide
\endcode
 \tabb Reinitializes the stream to its initial state:
   $C_g$ and $B_g$ are set to $I_g$.
 \endtabb
\code

   public void resetStartSubstream () \hide{
      for (int i = 0; i < 6;  ++i)  Cg[i] = Bg[i];
   } \endhide
\endcode
 \tabb Reinitializes the stream to the beginning of its current
   substream: $C_g$ is set to $B_g$.
 \endtabb
\code

   public void resetNextSubstream () \hide{
      int i;
      matVecModM (A1p76, Bg, Bg, m1);
      double temp[] = new double[3];
      for (i = 0; i < 3; ++i) temp[i] = Bg[i + 3];
      matVecModM (A2p76, temp, temp, m2);
      for (i = 0; i < 3; ++i) Bg[i + 3] = temp[i];
      for (i = 0; i < 6;  ++i) Cg[i] = Bg[i];
   } \endhide
\endcode
 \tabb Reinitializes the stream to the beginning of its next
   substream: $N_g$ is computed, and
   $C_g$ and $B_g$ are set to $N_g$.
 \endtabb
\code

   public void setAntithetic (boolean a) \hide{
      anti = a;
   }\endhide
\endcode
  \tabb After this procedure is called with {\tt a = true}, the stream 
   starts generating antithetic variates, 
   i.e., $1-U$ instead of $U$,
   until the procedure is called again with {\tt a = false}.
 \endtabb
\code

   public void increasedPrecis (boolean incp) \hide{
      prec53 = incp;
   }\endhide
\endcode
 \tabb After calling this procedure with {\tt incp = true}, each call to 
  the generator (direct or indirect) for this stream 
  will return a uniform random number with  more bits of resolution
  (53 bits if machine follows IEEE-754 floating-point standard)
  instead of 32 bits,
  and will advance the state of the stream by 2 steps instead of 1.  
  More precisely, if {\tt s} is a stream of the class {\tt RngStream},
  in the non-antithetic case, the instruction
  ``{\tt u = s.randU01()}'' will be equivalent to
  ``{\tt u = (s.randU01() + s.randU01(g) * fact) \%\ 1.0}'' 
  where the constant {\tt fact} is equal to $2^{-24}$.
  This also applies when calling {\tt randU01} indirectly
  (e.g., via {\tt randInt}, etc.).
  By default, or if this procedure is called again with {\tt incp = false}, 
  each call to {\tt randU01} for this stream advances the state by 1 step
  and returns a number with 32 bits of resolution.
 \endtabb
\code

   public void advanceState (int e, int c) \hide{
      double B1[][]= new double[3][3], C1[][]=new double[3][3];
      double B2[][]= new double[3][3], C2[][]=new double[3][3];

      if (e > 0) {
          matTwoPowModM (A1p0, B1, m1, e);
          matTwoPowModM (A2p0, B2, m2, e);
      }
      else if (e < 0) {
          matTwoPowModM (InvA1, B1, m1, -e);
          matTwoPowModM (InvA2, B2, m2, -e);
      }

      if (c >= 0) {
          matPowModM (A1p0, C1, m1, c);
          matPowModM (A2p0, C2, m2, c);
      }
      else if (c < 0) {
          matPowModM (InvA1, C1, m1, -c);
          matPowModM (InvA2, C2, m2, -c);
      }

      if (e!=0) {
          matMatModM (B1, C1, C1, m1);
          matMatModM (B2, C2, C2, m2);
      }

      matVecModM (C1, Cg, Cg, m1);
      double[] cg3 = new double[3];
      for (int i = 0; i < 3; i++)  cg3[i] = Cg[i+3];
      matVecModM (C2, cg3, cg3, m2);
      for (int i = 0; i < 3; i++)  Cg[i+3] = cg3[i];
   }\endhide
\endcode
 \tabb Advances the state of this stream by $k$ values,
  without modifying the states of other streams (as in {\tt setSeed}),
  nor the values of $B_g$ and $I_g$ associated with this stream.
  If $e > 0$, then $k=2^e + c$; 
  if $e < 0$,  then $k=-2^{-e} + c$; and if $e = 0$,  then $k=c$.
  Note: $c$ is allowed to take negative values. 
  We discourage the use of this method.
 \endtabb
\code\hide

    private static int CheckSeed (long seed[]) {
    /* 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 || seed[i] < 0) {
              System.out.println("****************************************\n"
                  + "ERROR: Seed[" + i + "],   Seed is not set."
                  + "\n****************************************\n");
               return -1;
           }
       }
       for (i = 3; i < 6; ++i) {
          if (seed[i] >= m2 || seed[i] < 0) {
             System.out.println("*****************************************\n"
                + "ERROR: Seed[" + i + "],   Seed is not set."
                + "\n*****************************************\n");
               return -1;
           }
       }
       if (seed[0] == 0 && seed[1] == 0 && seed[2] == 0) {
           System.out.println ("****************************\n"
             + "ERROR: First 3 seeds = 0.\n"
             +  "****************************\n");
            return -1;
       }
       if (seed[3] == 0 && seed[4] == 0 && seed[5] == 0) {
          System.out.println ("****************************\n"
             + "ERROR: Last 3 seeds = 0.\n"
             + "****************************\n");
            return -1;
       }
       return 0;
   }
\endhide

   public boolean setSeed (long seed[]) \hide{
        int i;
        if (CheckSeed (seed) != 0)
            return false;                   // FAILURE     
        for (i = 0; i < 6;  ++i)
           Cg[i] = Bg[i] = Ig[i] = seed[i];
       return true;                        // SUCCESS
    }\endhide
\endcode
 \tabb  Sets the initial seed $I_g$ of this stream 
  to the vector {\tt seed[0..5]}.  This vector should contain valid
  seed values as described in {\tt SetPackageSeed}.
  The state of the stream is then reset to this initial seed.
  The states and seeds of the other streams are not modified.
  As a result, after calling this
  method, the initial seeds of the streams are no longer spaced
  $Z$ values apart. We discourage the use of this method; proper
  use of the {\tt Reset} methods is preferable.
   Returns {\tt false} for invalid seeds, and {\tt true} otherwise.
 \endtabb
\code

   public double[] getState() \hide{
      return Cg;
   }\endhide
\endcode
 \tabb Returns the current state $C_g$ of this stream.
  This is a vector of 6 integers represented in floating point.
  This method is convenient if we want to save the state for 
  subsequent use.  
 \endtabb
\code

   public void writeState () \hide {
       System.out.print ("The current state of the RngStream");
       if (descriptor.length() > 0)
          System.out.print (" " + descriptor);
       System.out.print (":\n   Cg = { ");
       for (int i = 0; i < 5; i++)
          System.out.print ((long) Cg[i] + ", ");
       System.out.println ((long) Cg[5] + " }\n");
   }\endhide 
\endcode
 \tabb Prints (to standard output) the name and current state $C_g$ 
   of this stream.
 \endtabb
\code

   public void writeStateFull () \hide{
       System.out.print ("The RngStream");
       if (descriptor.length() > 0)
          System.out.print (" " + descriptor);
       System.out.println (":\n   anti = " + (anti ? "true" : "false"));
       System.out.println ("   prec53 = " + (prec53 ? "true" : "false"));

       System.out.print ("   Ig = { ");
       for (int i = 0; i < 5; i++)
          System.out.print ((long) Ig[i] + ", ");
       System.out.println ((long) Ig[5] + " }");

       System.out.print ("   Bg = { ");
       for (int i = 0; i < 5; i++)
          System.out.print ((long) Bg[i] + ", ");
       System.out.println ((long) Bg[5] + " }");

       System.out.print ("   Cg = { ");
       for (int i = 0; i < 5; i++)
          System.out.print ((long) Cg[i] + ", ");
       System.out.println ((long) Cg[5] + " }\n");
   }\endhide
\endcode
 \tabb Prints (to standard output) the name of this stream and the 
   values of all its internal variables.
 \endtabb
\code

   public double randU01 () \hide{
      if (prec53) return this.U01d();
      else return this.U01();
   }\endhide
\endcode
 \tabb Returns a (pseudo)random number from the uniform distribution
   over the interval $(0,1)$, using this stream,
   after advancing its state by one step.  
   Normally, the returned number has 32 bits of resolution,
   in the sense that it is always a multiple of $1/(2^{32}-208)$.
   However, if the precision has been increased by calling 
   {\tt increasedPrecis} for this stream, the resolution is higher
   and the stream state advances by two steps.
 \endtabb
\code

   public int randInt (int i, int j) \hide{
      return (i + (int)(randU01() * (j - i + 1.0)));
   }\endhide
\endcode
 \tabb Returns a (pseudo)random number from the discrete uniform 
   distribution over the integers $\{i,i+1,\dots,j\}$,
   using this stream. Makes one call to {\tt randU01}.
 \endtabb
\code

}
\endcode

\bigskip
\hrule

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\bibliography {random,simul,math}
\bibliographystyle {plain}

\end{document}
