Subversion Repositories public iLand

/********************************************************************************************

**    iLand - an individual based forest landscape and disturbance model

**    http://iland.boku.ac.at

**    Copyright (C) 2009-  Werner Rammer, Rupert Seidl

**

**    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 3 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, see <http://www.gnu.org/licenses/>.

********************************************************************************************/

#include "randomgenerator.h"

#include <QtGlobal>

#include "../3rdparty/MersenneTwister.h"

#include <QMutex>

class RGenerators

{

public:

    RGenerators() {}

    void seed();

    void seed(unsigned int oneSeed);

    inline unsigned int random_function(const int type) {   if (type==0) return WELLRNG512();

                                                            if (type==1) return xorshf96();

                                                            if (type==2) return fastrand();

                                                            return 0;

                                                        }

private:

    // see  http://www.lomont.org/Math/Papers/2008/Lomont_PRNG_2008.pdf

    // for details on the WellRNG512 algorithm

    inline unsigned long WELLRNG512(void);

    /* initialize state to random bits  */

    unsigned long state[16];

    /* init should also reset this to 0 */

    unsigned int index;

    /* return 32 bit random number      */

    int g_seed;

    inline unsigned int fastrand()

    {

        g_seed = (214013*g_seed+2531011);

        return g_seed;

    }

    unsigned long x,y,z;

    unsigned long xorshf96(void);

};

inline unsigned long RGenerators::WELLRNG512(void)

{

unsigned long a, b, c, d;

a  = state[index];

c  = state[(index+13)&15];

b  = a^c^(a<<16)^(c<<15);

c  = state[(index+9)&15];

c ^= (c>>11);

a  = state[index] = b^c;

d  = a^((a<<5)&0xDA442D24UL);

index = (index + 15)&15;

a  = state[index];

state[index] = a^b^d^(a<<2)^(b<<18)^(c<<28);

return state[index];

}

// The Marsaglia's xorshf generator:

// see: http://stackoverflow.com/questions/1640258/need-a-fast-random-generator-for-c and

// http://www.cse.yorku.ca/~oz/marsaglia-rng.html

//unsigned long x=123456789, y=362436069, z=521288629;

unsigned long RGenerators::xorshf96(void) {          //period 2^96-1

unsigned long t;

    x ^= x << 16;

    x ^= x >> 5;

    x ^= x << 1;

   t = x;

   x = y;

   y = z;

   z = t ^ x ^ y;

  return z;

}

inline void RGenerators::seed()

{

    srand ( time(NULL) );

    for (int i=0;i<16;i++)

        state[i] = std::rand();

    index = 0;

    // inits for the xorshift algorithm...

    x=123456789, y=362436069, z=521288629;

    // inits for the fast rand....

    g_seed = std::rand();

}

inline void RGenerators::seed(unsigned int oneSeed)

{

    srand ( oneSeed );

    for (int i=0;i<16;i++)

        state[i] = std::rand();

    index = 0;

    x=123456789, y=362436069, z=521288629;

    g_seed = oneSeed;

}

// static variables

unsigned int RandomGenerator::mBuffer[RANDOMGENERATORSIZE+5];

int RandomGenerator::mIndex = 0;

int RandomGenerator::mRotationCount = RANDOMGENERATORROTATIONS + 1;

int RandomGenerator::mRefillCounter = 0;

RandomGenerator::ERandomGenerators RandomGenerator::mGeneratorType = RandomGenerator::ergFastRandom;

/// fill the internal buffer with random numbers choosen from the defined random generator

QMutex random_generator_refill_mutex;

void RandomGenerator::refill() {

    QMutexLocker lock(&random_generator_refill_mutex); // serialize access

    if (mRotationCount<RANDOMGENERATORROTATIONS) // another thread might already succeeded in refilling....

        return;

    RGenerators gen;

    gen.seed(mBuffer[RANDOMGENERATORSIZE+4]); // use the last value as seed for the next round....

    switch (mGeneratorType) {

    case ergMersenneTwister: {

        MTRand mersenne;

        // qDebug() << "refill random numbers. seed" <<mBuffer[RANDOMGENERATORSIZE+4];

        mersenne.seed(mBuffer[RANDOMGENERATORSIZE+4]);

        //mersenne.seed(); // use random number seed from mersenne twister (hash of time and clock)

        for (int i=0;i<RANDOMGENERATORSIZE+5;++i)

            mBuffer[i] = mersenne.randInt();

        break;

    }

    case ergWellRNG512: {

        for (int i=0;i<RANDOMGENERATORSIZE+5;++i)

            mBuffer[i] = gen.random_function(0);

        break;

    }

    case ergXORShift96: {

        for (int i=0;i<RANDOMGENERATORSIZE+5;++i)

            mBuffer[i] = gen.random_function(1);

        break;

    }

    case ergFastRandom: {

        for (int i=0;i<RANDOMGENERATORSIZE+5;++i)

            mBuffer[i] = gen.random_function(2);

        break;

    }

    } // switch

    mIndex = 0; // reset the index

    mRotationCount=0;

    mRefillCounter++;

}

void RandomGenerator::seed(const unsigned oneSeed)

{

    if (oneSeed==0) {

        srand ( time(NULL) );

        mBuffer[RANDOMGENERATORSIZE+4] = std::rand();

    } else {

        mBuffer[RANDOMGENERATORSIZE+4] = oneSeed; // set a specific seed as seed for the next round

    }

}