---

Random.h

Go to the documentation of this file.

/* SimData: Data Infrastructure for Simulations
 * Copyright (C) 2002, 2003 Mark Rose <tm2@stm.lbl.gov>
 * 
 * This file is part of SimData.
 * 
 * This program 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.
 * 
 * This program 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 program; if not, write to the Free Software
 * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.
 */


#ifndef __SIMDATA_RANDOM_H__
#define __SIMDATA_RANDOM_H__


#include <cmath>
#include <SimData/Ref.h>
#include <SimData/Export.h>


NAMESPACE_SIMDATA


// Random Number Generators

namespace rng { // random number generators

class SIMDATA_EXPORT MT19937 {
        static const int N = 624;       /* Period parameters */
        static const int M = 397;

        // most significant w-r bits
        static const unsigned long UPPER_MASK = 0x80000000UL;   

        // least significant r bits 
        static const unsigned long LOWER_MASK = 0x7fffffffUL;   

        // state
        unsigned long _mt[N];
        int _mti;

        static inline unsigned long MAGIC(unsigned long y) {
                return (((y)&0x1) ? 0x9908b0dfUL : 0);
        }

        void update();
        
        inline unsigned long generate() {
                if (_mti >= N) update(); // generate N words at one time
                unsigned long k = _mt[_mti];
                k ^= (k >> 11);
                k ^= (k << 7) & 0x9d2c5680UL;
                k ^= (k << 15) & 0xefc60000UL;
                k ^= (k >> 18);
                _mti++;
                return k;
        }

public:
        struct State {
                unsigned long _mt[N];
                int _mti;
        };

        inline std::string getName() const { return "Mersenne Twister (MT19937)"; }

        MT19937() { setSeed(0); }

        template <typename T>
        inline T uniformInt(T lower, T upper) {
                return lower + static_cast<T>(unit() * (upper-lower));
        }

        template <typename T>
        inline T uniformInt(T upper) {
                return uniformInt(static_cast<T>(0), upper);
        }

        inline double unit() {
                return generate() / 4294967296.0 ;
        };

        inline double uniform(double lower, double upper) {
                return lower + unit() * (upper-lower);
        }

        void setSeed(unsigned long int s);

        void getState(State &state) const;

        void setState(State const &state);
};




class SIMDATA_EXPORT Taus2 {
        // state
        unsigned long int _s1, _s2, _s3;

        static const unsigned long int MASK = 0xffffffffUL;

        inline static unsigned long int LCG(unsigned long int n) {
                return ((69069 * n) & 0xffffffffUL);
        }

        inline static unsigned long int TAUSWORTHE(unsigned long int s,
                                                   unsigned long int a,
                                                   unsigned long int b,
                                                   unsigned long int c,
                                                   unsigned long int d) {
                return (((s &c) <<d) &MASK) ^ ((((s <<a) &MASK)^s) >>b);
        }

        inline unsigned long generate() {
                _s1 = TAUSWORTHE(_s1, 13, 19, 4294967294UL, 12);
                _s2 = TAUSWORTHE(_s2, 2, 25, 4294967288UL, 4);
                _s3 = TAUSWORTHE(_s3, 3, 11, 4294967280UL, 17);
                return (_s1 ^ _s2 ^ _s3);
        }
public:
        struct State {
                unsigned long int _s1, _s2, _s3;
        };

        inline std::string getName() const { return "Tausworthe (Taus2)"; }

        Taus2() { setSeed(0); }

        template <typename T>
        inline T uniformInt(T lower, T upper) {
                return lower + static_cast<T>(unit() * (upper-lower));
        }

        template <typename T>
        inline T uniformInt(T upper) {
                return uniformInt(static_cast<T>(0), upper);
        }

        inline double unit() {
                return generate() / 4294967296.0 ;
        };

        inline double uniform(double lower, double upper) {
                return lower + unit() * (upper-lower);
        }

        void setSeed(unsigned long int s);

        void getState(State &state) const;

        void setState(State const &state);
};


} // namespace rng



// Random Number Generator Wrappers

class SIMDATA_EXPORT RandomInterface {
protected:
        struct _State: Referenced { };
public:
        typedef Ref<_State> State;
        typedef unsigned long int SeedType;

        virtual ~RandomInterface();

        virtual void setSeed(SeedType seed)=0;
        
        virtual State getState() const=0;

        virtual void setState(State const &state)=0;

        virtual std::string getName() const=0;
};


class SIMDATA_EXPORT RandomNumberGeneratorInterface: public RandomInterface {
public:
        virtual double unit()=0;

        virtual double uniform(double lower, double upper)=0;
        
        virtual long uniformInt(long lower, long upper)=0;

        virtual long uniformInt(long upper)=0;

        virtual unsigned long uniformUInt(unsigned long lower, unsigned long upper)=0;

        virtual unsigned long uniformUInt(unsigned long upper)=0;
};

template <class RNG>
class RandomNumberGenerator: public RandomNumberGeneratorInterface {

        struct RNGState: _State {
                typename RNG::State _state;
        };
        RNG _gen;
public:
        virtual void setSeed(SeedType seed);

        virtual State getState() const;

        virtual void setState(State const &);

        virtual double unit();

        virtual double uniform(double, double);

        virtual long uniformInt(long lower, long upper);

        virtual long uniformInt(long upper);

        virtual unsigned long uniformUInt(unsigned long lower, unsigned long upper);

        virtual unsigned long uniformUInt(unsigned long upper);

        virtual std::string getName() const;

        inline RNG *operator->() { return &_gen; }
};


template <class RNG>
RandomInterface::State RandomNumberGenerator<RNG>::getState() const {
        RNGState *rng_state = new RNGState;
        _gen.getState(rng_state->_state);
        return rng_state;
}

template <class RNG>
void RandomNumberGenerator<RNG>::setState(RandomInterface::State const &state) {
        Ref<RNGState> rng_state = state;
        _gen.setState(rng_state->_state);
}

template <class RNG>
void RandomNumberGenerator<RNG>::setSeed(RandomInterface::SeedType seed) {
        _gen.setSeed(seed);
}

template <class RNG>
double RandomNumberGenerator<RNG>::unit() {
        return _gen.unit();
}

template <class RNG>
double RandomNumberGenerator<RNG>::uniform(double low, double high) {
        return _gen.uniform(low, high);
}

template <class RNG>
long RandomNumberGenerator<RNG>::uniformInt(long low, long high) {
        return _gen.uniformInt(low, high);
}

template <class RNG>
long RandomNumberGenerator<RNG>::uniformInt(long high) {
        return _gen.uniformInt(high);
}

template <class RNG>
unsigned long RandomNumberGenerator<RNG>::uniformUInt(unsigned long low, unsigned long high) {
        return _gen.uniformInt(low, high);
}

template <class RNG>
unsigned long RandomNumberGenerator<RNG>::uniformUInt(unsigned long high) {
        return _gen.uniformInt(high);
}

template <class RNG>
std::string RandomNumberGenerator<RNG>::getName() const { 
        return _gen.getName(); 
}


// Random Distributions

namespace rd { // random distributions

class SIMDATA_EXPORT Gauss {
        rng::Taus2 _gen;
        double _mean, _sigma;
        double _x;
        bool _odd;
public:

        struct State {
                rng::Taus2::State _state;
                double _mean, _sigma, _x;
                bool _odd;
        };

        std::string getName() const { return "Gaussian (Gauss)"; }

        Gauss(double mean=0.0, double sigma=0.0): _mean(mean), _sigma(sigma), _odd(true) {}

        inline void setDistribution(double mean, double sigma) {
                _mean = mean;
                _sigma = sigma;
        }

        inline double getMean() const { return _mean; }

        inline double getSigma() const { return _sigma; }

        double sample();

        inline void setSeed(unsigned long int seed) { _gen.setSeed(seed); _odd = true; }

        void getState(State &) const;

        void setState(State const &);
};


extern SIMDATA_EXPORT double BoxMueller(RandomNumberGeneratorInterface &_gen, double _mean, double _sigma);

} // namespace rd



// Random Distribution Wrappers


class SIMDATA_EXPORT RandomDistributionInterface: public RandomInterface {
public:
        virtual double sample()=0;
};



template <class RD>
class RandomDistribution: public RandomDistributionInterface {
        struct RDState: _State {
                typename RD::State _state;
        };
        RD _dist;
public:
        double sample();

        virtual void setSeed(SeedType seed);

        virtual State getState() const;

        virtual void setState(State const &);

        virtual std::string getName() const;

        inline RD *operator->() { return &_dist; }
};


template <class RD>
double RandomDistribution<RD>::sample() {
        return _dist.sample();
}

template <class RD>
RandomInterface::State RandomDistribution<RD>::getState() const { 
        RDState *rd_state = new RDState;
        _dist.getState(rd_state->_state);
        return rd_state;
}

template <class RD>
void RandomDistribution<RD>::setState(RandomInterface::State const &state) { 
        Ref<RDState> rd_state = state;
        _dist.setState(rd_state->_state);
}

template <class RD>
void RandomDistribution<RD>::setSeed(RandomInterface::SeedType s) {
        _dist.setSeed(s);
}

template <class RD>
std::string RandomDistribution<RD>::getName() const { 
        return _dist.getName(); 
}




namespace random { // generators and distributions
        typedef RandomDistribution<rd::Gauss> Gauss;
        typedef RandomNumberGenerator<rng::MT19937> MersenneTwister;
        typedef RandomNumberGenerator<rng::Taus2> Taus2;
        typedef Taus2 Standard;
}


extern SIMDATA_EXPORT random::Taus2 g_Random;


NAMESPACE_SIMDATA_END // simdata

#endif // __SIMDATA_RANDOM_H__

---

SimData version pre-0.4.0. For more information on SimData, visit the SimData Homepage. Generated on Tue Oct 14 12:06:39 2003, using Doxygen 1.2.18.