Subversion Repositories rng_lib

----------------------------------------------------------------------
----                                                              ----
---- Rand number generator library.                               ----
----                                                              ----
---- This file is part of the Random Number Generator project     ----
---- http://www.opencores.org/cores/rng_lib/                      ----
----                                                              ----
---- Description                                                  ----
---- This library has function for generation random numbers with ----
---- the following distributions:                                 ----
---- - Uniform (continous)                                        ----
---- - Exponential (continous)                                    ----
---- - Gaussian (continous)                                       ----
----                                                              ----
---- Random numbers are produced with a combination of 3          ----
---- Tausworthe generators which gives very good statistical      ----
---- properties.                                                  ----
----                                                              ----
---- NOTE! These functions will NOT synthesize. They are for test ----
----       bench use only!                                        ----
----                                                              ----
---- To Do:                                                       ----
---- -                                                            ----
----                                                              ----
---- Author(s):                                                   ----
---- - Geir Drange, gedra@opencores.org                           ----
----                                                              ----
----------------------------------------------------------------------
----                                                              ----
---- Copyright (C) 2004 Authors and OPENCORES.ORG                 ----
----                                                              ----
---- This source file may be used and distributed without         ----
---- restriction provided that this copyright statement is not    ----
---- removed from the file and that any derivative work contains  ----
---- the original copyright notice and the associated disclaimer. ----
----                                                              ----
---- This source file 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.0 of the License, or (at your option) any   ----
---- later version.                                               ----
----                                                              ----
---- This source 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 this source; if not, download it   ----
---- from http://www.gnu.org/licenses/gpl.txt                     ----
----                                                              ----
----------------------------------------------------------------------
--
-- CVS Revision History
--
-- $Log: not supported by cvs2svn $
-- Revision 1.1  2004/09/28 15:12:28  gedra
-- Random number library functions.
--
--
--
 
library ieee;
use ieee.std_logic_1164.all;
use ieee.numeric_std.all;
use std.textio.all;
use work.math_lib.all;
 
package rng_lib is
 
   type distribution is (UNIFORM, GAUSSIAN, EXPONENTIAL);
   type rand_var is record                          -- random variable record
      rnd                 : real;                   -- random number
      rnd_v               : unsigned(31 downto 0);  -- random number vector
      dist                : distribution;           -- distribution type
      y, z                : real;                   -- distribution parameters
      s1, s2, s3          : unsigned(31 downto 0);  -- seeds
      mask1, mask2, mask3 : unsigned(31 downto 0);
      shft1, shft2, shft3 : natural;
   end record;
 
   function rand (rnd : rand_var) return rand_var;
   function init_uniform(constant a, b, c : natural;
                         constant lo, hi  : real) return rand_var;
   function init_gaussian(constant a, b, c     : natural;
                          constant mean, stdev : real) return rand_var;
   function init_exponential(constant a, b, c : natural;
                             constant mean    : real) return rand_var;
 
   constant q1 : natural := 13;
   constant q2 : natural := 2;
   constant q3 : natural := 3;
   constant p1 : natural := 12;
   constant p2 : natural := 4;
   constant p3 : natural := 17;
 
end rng_lib;
 
package body rng_lib is
 
-- Function to convert 32bit unsigned vector to real
-- Integers only go to 2**31 (VHDL'87), so do it clever
   function unsigned_2_real (constant a : unsigned(31 downto 0)) return real is
      variable r : real;
   begin
      r := 2.0*real(to_integer(a(31 downto 1)));
      if a(0) = '1' then
         r := r + 1.0;
      end if;
      return(r);
   end unsigned_2_real;
 
-- Generate random number using a combination of 3 tausworthe generators 
-- Source: Pierre L'Ecuyer, "Maximally Equidistributed Combined Tausworthe
-- Generators". Mathematics of Computation, vol.65, no.213(1996), pp203--213.
   function rng (rnd : rand_var) return rand_var is
      variable new_rnd : rand_var;
      variable b       : unsigned(31 downto 0);
   begin
      new_rnd       := rnd;
      b             := ((new_rnd.s1 sll q1) xor new_rnd.s1) srl new_rnd.shft1;
      new_rnd.s1    := ((new_rnd.s1 and new_rnd.mask1) sll p1) xor b;
      b             := ((new_rnd.s2 sll q2) xor new_rnd.s2) srl new_rnd.shft2;
      new_rnd.s2    := ((new_rnd.s2 and new_rnd.mask2) sll p2) xor b;
      b             := ((new_rnd.s3 sll q3) xor new_rnd.s3) srl new_rnd.shft3;
      new_rnd.s3    := ((new_rnd.s3 and new_rnd.mask3) sll p3) xor b;
      new_rnd.rnd_v := new_rnd.s1 xor new_rnd.s2 xor new_rnd.s3;
      -- normalize to range [0,1)
      new_rnd.rnd   := unsigned_2_real(new_rnd.rnd_v) / 65536.0;
      new_rnd.rnd   := new_rnd.rnd / 65536.0;
      return (new_rnd);
   end rng;
 
-- rand function generates a random variable with different distributions
   function rand (rnd : rand_var) return rand_var is
      variable rnd_out : rand_var;
      variable x, y, z : real;
      variable t       : real := 0.0;
   begin
      case rnd.dist is
         -- Uniform distribution
         when UNIFORM =>
            rnd_out     := rng(rnd);
            rnd_out.rnd := rnd.y + (rnd_out.rnd * (rnd.z - rnd.y));
            -- Gaussian distribution
         when GAUSSIAN =>               -- Box-Mueller method
            z       := 2.0;
            rnd_out := rnd;
            while z > 1.0 or z = 0.0 loop
               -- choose x,y in uniform square (-1,-1) to (+1,+1)
               rnd_out := rng(rnd_out);
               x       := -1.0 + 2.0 * rnd_out.rnd;
               rnd_out := rng(rnd_out);
               y       := -1.0 + 2.0 * rnd_out.rnd;
               z       := (x * x) + (y * y);
            end loop;
            -- Box-Mueller transform
            rnd_out.rnd := rnd_out.y + rnd_out.z * y * sqrt(-2.0 * ln(z)/z);
            -- Exponential distribution
         when EXPONENTIAL =>
            rnd_out     := rng(rnd);
            rnd_out.rnd := -rnd_out.y * log(1.0 - rnd_out.rnd);
         when others =>
            report "rand() function encountered an error!"
               severity failure;
      end case;
      return (rnd_out);
   end rand;
 
-- Initialize seeds, used by all init_ functions
   function gen_seed (constant a, b, c : natural) return rand_var is
      variable seeded : rand_var;
      variable x      : unsigned(31 downto 0) := "11111111111111111111111111111111";
      constant k1     : natural               := 31;
      constant k2     : natural               := 29;
      constant k3     : natural               := 28;
   begin
      seeded.shft1 := k1-p1;
      seeded.shft2 := k2-p2;
      seeded.shft3 := k3-p3;
      seeded.mask1 := x sll (32-k1);
      seeded.mask2 := x sll (32-k2);
      seeded.mask3 := x sll (32-k3);
      seeded.s1    := to_unsigned(390451501, 32);
      seeded.s2    := to_unsigned(613566701, 32);
      seeded.s3    := to_unsigned(858993401, 32);
      if to_unsigned(a, 32) > (to_unsigned(1, 32) sll (32-k1)) then
         seeded.s1 := to_unsigned(a, 32);
      end if;
      if to_unsigned(b, 32) > (to_unsigned(1, 32) sll (32-k2)) then
         seeded.s2 := to_unsigned(b, 32);
      end if;
      if to_unsigned(c, 32) > (to_unsigned(1, 32) sll (32-k3)) then
         seeded.s3 := to_unsigned(c, 32);
      end if;
      return(seeded);
   end gen_seed;
 
-- Uniform distribution random variable initialization
-- a,b,c are seeds
-- lo,hi is the range for the uniform distribution
   function init_uniform(constant a, b, c : natural;
                         constant lo, hi  : real) return rand_var is
      variable rnd, rout : rand_var;
   begin
      if lo >= hi then
         report "Uniform parameter error: 'hi' must be > 'lo'!"
            severity failure;
      end if;
      rnd      := gen_seed(a, b, c);
      rnd.dist := UNIFORM;
      rnd.y    := lo;
      rnd.z    := hi;
      rout     := rand(rnd);
      return(rout);
   end init_uniform;
 
-- Gaussian distribution random variable initialization
-- a,b,c are seeds
-- mean,stdev is mean and standard deviation
   function init_gaussian(constant a, b, c     : natural;
                          constant mean, stdev : real) return rand_var is
      variable rnd, rout : rand_var;
   begin
      if stdev = 0.0 then
         report "Gaussian parameter error: 'stdev' must be non-zero!"
            severity failure;
      end if;
      rnd      := gen_seed(a, b, c);
      rnd.dist := GAUSSIAN;
      rnd.y    := mean;
      rnd.z    := stdev;
      rout     := rand(rnd);
      return(rout);
   end init_gaussian;
 
-- Exponential distribution random variable initialization
-- a,b,c are seeds
-- mean: mean value
   function init_exponential(constant a, b, c : natural;
                             constant mean    : real) return rand_var is
      variable rnd, rout : rand_var;
   begin
      if mean <= 0.0 then
         report "Exponential parameter error: 'mean' must be > 0!"
            severity failure;
      end if;
      rnd      := gen_seed(a, b, c);
      rnd.dist := EXPONENTIAL;
      rnd.y    := mean;
      rout     := rand(rnd);
      return(rout);
   end init_exponential;
 
end rng_lib;