[Lemon-commits] [lemon_svn] deba: r2969 - in hugo/trunk: lemon lemon/bits test
Lemon SVN
svn at lemon.cs.elte.hu
Mon Nov 6 21:51:32 CET 2006
Author: deba
Date: Mon Oct 2 18:11:00 2006
New Revision: 2969
Added:
hugo/trunk/lemon/random.cc
hugo/trunk/lemon/random.h
Removed:
hugo/trunk/lemon/bits/mingw32_rand.cc
hugo/trunk/lemon/bits/mingw32_rand.h
Modified:
hugo/trunk/lemon/Makefile.am
hugo/trunk/lemon/hypercube_graph.h
hugo/trunk/lemon/simann.h
hugo/trunk/test/graph_utils_test.h
hugo/trunk/test/heap_test.h
hugo/trunk/test/test_tools.h
Log:
Mersenne Twister random number generator
The code is based on the official MT19937 implementation
It is fully rewritten:
http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/emt.html
todo: fixing copyright information
Modified: hugo/trunk/lemon/Makefile.am
==============================================================================
--- hugo/trunk/lemon/Makefile.am (original)
+++ hugo/trunk/lemon/Makefile.am Mon Oct 2 18:11:00 2006
@@ -12,8 +12,8 @@
lemon/base.cc \
lemon/color.cc \
lemon/eps.cc \
- lemon/bits/mingw32_rand.cc \
- lemon/bits/mingw32_time.cc
+ lemon/bits/mingw32_time.cc \
+ lemon/random.cc
lemon_libemon_la_CXXFLAGS = $(GLPK_CFLAGS) $(CPLEX_CFLAGS)
lemon_libemon_la_LDFLAGS = $(GLPK_LIBS) $(CPLEX_LIBS)
@@ -86,6 +86,7 @@
lemon/prim.h \
lemon/radix_heap.h \
lemon/radix_sort.h \
+ lemon/random.h \
lemon/refptr.h \
lemon/simann.h \
lemon/smart_graph.h \
@@ -112,7 +113,6 @@
lemon/bits/item_reader.h \
lemon/bits/item_writer.h \
lemon/bits/map_extender.h \
- lemon/bits/mingw32_rand.h \
lemon/bits/mingw32_time.h \
lemon/bits/traits.h \
lemon/bits/utility.h \
Modified: hugo/trunk/lemon/hypercube_graph.h
==============================================================================
--- hugo/trunk/lemon/hypercube_graph.h (original)
+++ hugo/trunk/lemon/hypercube_graph.h Mon Oct 2 18:11:00 2006
@@ -260,8 +260,8 @@
/// HyperCubeGraph graph(DIM);
/// dim2::Point<double> base[DIM];
/// for (int k = 0; k < DIM; ++k) {
- /// base[k].x = rand() / (RAND_MAX + 1.0);
- /// base[k].y = rand() / (RAND_MAX + 1.0);
+ /// base[k].x = random.getReal();
+ /// base[k].y = random.getReal();
/// }
/// HyperCubeGraph::HyperMap<dim2::Point<double> >
/// pos(graph, base, base + DIM, dim2::Point<double>(0.0, 0.0));
Added: hugo/trunk/lemon/random.cc
==============================================================================
--- (empty file)
+++ hugo/trunk/lemon/random.cc Mon Oct 2 18:11:00 2006
@@ -0,0 +1,5 @@
+#include <lemon/random.h>
+
+namespace lemon {
+ Random rnd;
+}
Added: hugo/trunk/lemon/random.h
==============================================================================
--- (empty file)
+++ hugo/trunk/lemon/random.h Mon Oct 2 18:11:00 2006
@@ -0,0 +1,511 @@
+/* -*- C++ -*-
+ *
+ * This file is a part of LEMON, a generic C++ optimization library
+ *
+ * Copyright (C) 2003-2006
+ * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
+ * (Egervary Research Group on Combinatorial Optimization, EGRES).
+ *
+ * Permission to use, modify and distribute this software is granted
+ * provided that this copyright notice appears in all copies. For
+ * precise terms see the accompanying LICENSE file.
+ *
+ * This software is provided "AS IS" with no warranty of any kind,
+ * express or implied, and with no claim as to its suitability for any
+ * purpose.
+ *
+ */
+
+#ifndef LEMON_RANDOM_H
+#define LEMON_RANDOM_H
+
+#include <algorithm>
+
+#include <ctime>
+#include <cmath>
+#include <cmath>
+
+///\ingroup misc
+///\file
+///\brief Mersenne Twister random number generator
+///
+///\author Balazs Dezso
+
+namespace lemon {
+
+#if WORD_BIT == 32
+
+ /// \ingroup misc
+ ///
+ /// \brief Mersenne Twister random number generator
+ ///
+ /// The Mersenne Twister is a twisted generalized feedback
+ /// shift-register generator of Matsumoto and Nishimura. The period of
+ /// this generator is \f$ 2^{19937} - 1 \f$ and it is equi-distributed in
+ /// 623 dimensions. The time performance of this generator is comparable
+ /// to the commonly used generators.
+ ///
+ /// \author Balazs Dezso
+ class Random {
+
+ static const int length = 624;
+ static const int shift = 397;
+
+ public:
+
+ static const unsigned long min = 0ul;
+ static const unsigned long max = ~0ul;
+
+ /// \brief Constructor
+ ///
+ /// Constructor with time dependent seeding.
+ Random() { initState(std::time(0)); }
+
+ /// \brief Constructor
+ ///
+ /// Constructor
+ Random(unsigned long seed) { initState(seed); }
+
+ /// \brief Copy constructor
+ ///
+ /// Copy constructor. The generated sequence will be identical to
+ /// the other sequence.
+ Random(const Random& other) {
+ std::copy(other.state, other.state + length, state);
+ current = state + (other.current - other.state);
+ }
+
+ /// \brief Assign operator
+ ///
+ /// Assign operator. The generated sequence will be identical to
+ /// the other sequence.
+ Random& operator=(const Random& other) {
+ if (&other != this) {
+ std::copy(other.state, other.state + length, state);
+ current = state + (other.current - other.state);
+ }
+ return *this;
+ }
+
+ /// \brief Returns an unsigned random number
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, 2^{32} - 1\} \f$.
+ unsigned long getUnsigned() {
+ if (current == state) fillState();
+ --current;
+ unsigned long rnd = *current;
+ rnd ^= (rnd >> 11);
+ rnd ^= (rnd << 7) & 0x9D2C5680ul;
+ rnd ^= (rnd << 15) & 0xEFC60000ul;
+ rnd ^= (rnd >> 18);
+ return rnd;
+ }
+
+ /// \brief Returns a random number
+ ///
+ /// It returns an integer random number from the range
+ /// \f$ \{-2^{31}, \dots, 2^{31} - 1\} \f$.
+ long getInt() {
+ return (long)getUnsigned();
+ }
+
+ /// \brief Returns an unsigned integer random number
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, 2^{31} - 1\} \f$.
+ long getNatural() {
+ return (long)(getUnsigned() >> 1);
+ }
+
+ /// \brief Returns a random bool
+ ///
+ /// It returns a random bool.
+ bool getBool() {
+ return (bool)(getUnsigned() & 1);
+ }
+
+ /// \brief Returns a real random number
+ ///
+ /// It returns a real random number from the range
+ /// \f$ [0, 1) \f$. The double will have 32 significant bits.
+ double getReal() {
+ return std::ldexp((double)getUnsigned(), -32);
+ }
+
+ /// \brief Returns a real random number
+ ///
+ /// It returns a real random number from the range
+ /// \f$ [0, 1) \f$. The double will have 53 significant bits,
+ /// but it is slower than the \c getReal().
+ double getPrecReal() {
+ return std::ldexp((double)(getUnsigned() >> 5), -27) +
+ std::ldexp((double)(getUnsigned() >> 6), -53);
+ }
+
+ /// \brief Returns an unsigned random number from a given range
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, n - 1\} \f$.
+ unsigned long getUnsigned(unsigned long n) {
+ unsigned long mask = n - 1, rnd;
+ mask |= mask >> 1;
+ mask |= mask >> 2;
+ mask |= mask >> 4;
+ mask |= mask >> 8;
+ mask |= mask >> 16;
+ do rnd = getUnsigned() & mask; while (rnd >= n);
+ return rnd;
+ }
+
+ /// \brief Returns a random number from a given range
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, n - 1\} \f$.
+ long getInt(long n) {
+ long mask = n - 1, rnd;
+ mask |= mask >> 1;
+ mask |= mask >> 2;
+ mask |= mask >> 4;
+ mask |= mask >> 8;
+ mask |= mask >> 16;
+ do rnd = getUnsigned() & mask; while (rnd >= n);
+ return rnd;
+ }
+
+ private:
+
+ void initState(unsigned long seed) {
+ static const unsigned long mul = 0x6c078965ul;
+
+ current = state;
+
+ unsigned long *curr = state + length - 1;
+ curr[0] = seed; --curr;
+ for (int i = 1; i < length; ++i) {
+ curr[0] = (mul * ( curr[1] ^ (curr[1] >> 30) ) + i);
+ --curr;
+ }
+ }
+
+ void fillState() {
+ static const unsigned long mask[2] = { 0x0ul, 0x9908B0DFul };
+ static const unsigned long loMask = (1ul << 31) - 1;
+ static const unsigned long hiMask = ~loMask;
+
+ current = state + length;
+
+ register unsigned long *curr = state + length - 1;
+ register long num;
+
+ num = length - shift;
+ while (num--) {
+ curr[0] = (((curr[0] & hiMask) | (curr[-1] & loMask)) >> 1) ^
+ curr[- shift] ^ mask[curr[-1] & 1ul];
+ --curr;
+ }
+ num = shift - 1;
+ while (num--) {
+ curr[0] = (((curr[0] & hiMask) | (curr[-1] & loMask)) >> 1) ^
+ curr[length - shift] ^ mask[curr[-1] & 1ul];
+ --curr;
+ }
+ curr[0] = (((curr[0] & hiMask) | (curr[length - 1] & loMask)) >> 1) ^
+ curr[length - shift] ^ mask[curr[length - 1] & 1ul];
+
+ }
+
+ unsigned long *current;
+ unsigned long state[length];
+
+ };
+
+#elif WORD_BIT == 64
+
+ /// \brief Mersenne Twister random number generator
+ ///
+ /// The Mersenne Twister is a twisted generalized feedback
+ /// shift-register generator of Matsumoto and Nishimura. The period of
+ /// this generator is \f$ 2^{19937} - 1 \f$ and it is equi-distributed in
+ /// 311 dimensions. The time performance of this generator is comparable
+ /// to the commonly used generators.
+ class Random {
+
+ static const int length = 312;
+ static const int shift = 156;
+
+ public:
+
+ static const unsigned long min = 0ul;
+ static const unsigned long max = ~0ul;
+
+ /// \brief Constructor
+ ///
+ /// Constructor with time dependent seeding.
+ Random() { initState(std::time(0)); }
+
+ /// \brief Constructor
+ ///
+ /// Constructor
+ Random(unsigned long seed) { initState(seed); }
+
+ /// \brief Copy constructor
+ ///
+ /// Copy constructor. The generated sequence will be identical to
+ /// the other sequence.
+ Random(const Random& other) {
+ std::copy(other.state, other.state + length, state);
+ current = state + (other.current - other.state);
+ }
+
+ /// \brief Assign operator
+ ///
+ /// Assign operator. The generated sequence will be identical to
+ /// the other sequence.
+ Random& operator=(const Random& other) {
+ if (&other != this) {
+ std::copy(other.state, other.state + length, state);
+ current = state + (other.current - other.state);
+ }
+ return *this;
+ }
+
+ /// \brief Returns an unsigned random number
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, 2^{64} - 1\} \f$.
+ unsigned long getUnsigned() {
+ if (current == state) fillState();
+ --current;
+ unsigned long rnd = *current;
+ rnd ^= (rnd >> 29) & 0x5555555555555555ul;
+ rnd ^= (rnd << 17) & 0x71D67FFFEDA60000ul;
+ rnd ^= (rnd << 37) & 0xFFF7EEE000000000ul;
+ rnd ^= (rnd >> 43);
+ return rnd;
+ }
+
+ /// \brief Returns a random number
+ ///
+ /// It returns an integer random number from the range
+ /// \f$ \{-2^{63}, \dots, 2^{63} - 1\} \f$.
+ long getInt() {
+ return (long)getUnsigned();
+ }
+
+ /// \brief Returns an unsigned integer random number
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, 2^{63} - 1\} \f$.
+ long getNatural() {
+ return (long)(getUnsigned() >> 1);
+ }
+
+ /// \brief Returns a random bool
+ ///
+ /// It returns a random bool.
+ bool getBool() {
+ return (bool)(getUnsigned() & 1);
+ }
+
+ /// \brief Returns a real random number
+ ///
+ /// It returns a real random number from the range
+ /// \f$ [0, 1) \f$.
+ double getReal() {
+ return std::ldexp((double)(getUnsigned() >> 11), -53);
+ }
+
+ /// \brief Returns a real random number
+ ///
+ /// It returns a real random number from the range
+ /// \f$ [0, 1) \f$. This function is identical to the \c getReal().
+ double getPrecReal() {
+ return getReal();
+ }
+
+ /// \brief Returns an unsigned random number from a given range
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, n - 1\} \f$.
+ unsigned long getUnsigned(unsigned long n) {
+ unsigned long mask = n - 1, rnd;
+ mask |= mask >> 1;
+ mask |= mask >> 2;
+ mask |= mask >> 4;
+ mask |= mask >> 8;
+ mask |= mask >> 16;
+ mask |= mask >> 32;
+ do rnd = getUnsigned() & mask; while (rnd >= n);
+ return rnd;
+ }
+
+ /// \brief Returns a random number from a given range
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, n - 1\} \f$.
+ long getInt(long n) {
+ long mask = n - 1, rnd;
+ mask |= mask >> 1;
+ mask |= mask >> 2;
+ mask |= mask >> 4;
+ mask |= mask >> 8;
+ mask |= mask >> 16;
+ mask |= mask >> 32;
+ do rnd = getUnsigned() & mask; while (rnd >= n);
+ return rnd;
+ }
+
+ private:
+
+ void initState(unsigned long seed) {
+
+ static const unsigned long mul = 0x5851F42D4C957F2Dul;
+
+ current = state;
+
+ unsigned long *curr = state + length - 1;
+ curr[0] = seed; --curr;
+ for (int i = 1; i < length; ++i) {
+ curr[0] = (mul * ( curr[1] ^ (curr[1] >> 62) ) + i);
+ --curr;
+ }
+ }
+
+ void fillState() {
+ static const unsigned long mask[2] = { 0x0ul, 0xB5026F5AA96619E9ul };
+ static const unsigned long loMask = (1ul << 31) - 1;
+ static const unsigned long hiMask = ~loMask;
+
+ current = state + length;
+
+ register unsigned long *curr = state + length - 1;
+ register int num;
+
+ num = length - shift;
+ while (num--) {
+ curr[0] = (((curr[0] & hiMask) | (curr[-1] & loMask)) >> 1) ^
+ curr[- shift] ^ mask[curr[-1] & 1ul];
+ --curr;
+ }
+ num = shift - 1;
+ while (num--) {
+ curr[0] = (((curr[0] & hiMask) | (curr[-1] & loMask)) >> 1) ^
+ curr[length - shift] ^ mask[curr[-1] & 1ul];
+ --curr;
+ }
+ curr[0] = (((curr[0] & hiMask) | (curr[length - 1] & loMask)) >> 1) ^
+ curr[length - shift] ^ mask[curr[length - 1] & 1ul];
+
+ }
+
+ unsigned long *current;
+ unsigned long state[length];
+
+ };
+
+#else
+
+ /// \brief Mersenne Twister random number generator
+ ///
+ /// The Mersenne Twister is a twisted generalized feedback
+ /// shift-register generator of Matsumoto and Nishimura. There is
+ /// two different implementation in the Lemon library, one for
+ /// 32-bit processors and one for 64-bit processors. The period of
+ /// the generated sequence is \f$ 2^{19937} - 1 \f$, the generated
+ /// sequence of 32-bit random numbers is equi-distributed in 623
+ /// dimensions. The time performance of this generator is comparable
+ /// to the commonly used generators.
+ class Random {
+ public:
+
+ static const unsigned long min = 0ul;
+ static const unsigned long max = ~0ul;
+
+ /// \brief Constructor
+ ///
+ /// Constructor with time dependent seeding.
+ Random() {}
+
+ /// \brief Constructor
+ ///
+ /// Constructor
+ Random(unsigned long seed) {}
+
+ /// \brief Copy constructor
+ ///
+ /// Copy constructor. The generated sequence will be identical to
+ /// the other sequence.
+ Random(const Random& other) {}
+
+ /// \brief Assign operator
+ ///
+ /// Assign operator. The generated sequence will be identical to
+ /// the other sequence.
+ Random& operator=(const Random& other) { return *this; }
+
+ /// \brief Returns an unsigned random number
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, 2^{64} - 1\} \f$ for 64-bit processors and from
+ /// \f$ \{0, 1, \dots, 2^{32} - 1\} \f$ for 32-bit processors.
+ unsigned long getUnsigned() { return 0ul; }
+
+ /// \brief Returns a random number
+ ///
+ /// It returns an integer random number from the range
+ /// \f$ \{-2^{63}, \dots, 2^{63} - 1\} \f$ for 64-bit processors and from
+ /// \f$ \{-2^{31}, \dots, 2^{31} - 1\} \f$ for 32-bit processors.
+ long getInt() { return 0l; }
+
+ /// \brief Returns an unsigned integer random number
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, 2^{63} - 1\} \f$ for 64-bit processors and
+ /// from \f$ \{0, 1, \dots, 2^{31} - 1\} \f$ for 32-bit processors.
+ long getNatural() { return 0l; }
+
+ /// \brief Returns a random bool
+ ///
+ /// It returns a random bool.
+ bool getBool() { return false; }
+
+ /// \brief Returns a real random number
+ ///
+ /// It returns a real random number from the range
+ /// \f$ [0, 1) \f$. For 32-bit processors the generated random
+ /// number will just have 32 significant bits.
+ double getReal() { return 0.0; }
+
+ /// \brief Returns a real random number
+ ///
+ /// It returns a real random number from the range
+ /// \f$ [0, 1) \f$. This function returns random numbers with 53
+ /// significant bits for 32-bit processors. For 64-bit processors
+ /// it is identical to the \c getReal().
+ double getPrecReal() { return 0.0; }
+
+ /// \brief Returns an unsigned random number from a given range
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, n - 1\} \f$.
+ unsigned long getUnsigned(unsigned long n) { return 0; }
+
+ /// \brief Returns a random number from a given range
+ ///
+ /// It returns an unsigned integer random number from the range
+ /// \f$ \{0, 1, \dots, n - 1\} \f$.
+ long getInt(long n) { return 0; }
+
+ };
+
+#endif
+
+ /// \brief Global random number generator instance
+ ///
+ /// A global mersenne twister random number generator instance
+ extern Random rnd;
+
+}
+
+#endif
Modified: hugo/trunk/lemon/simann.h
==============================================================================
--- hugo/trunk/lemon/simann.h (original)
+++ hugo/trunk/lemon/simann.h Mon Oct 2 18:11:00 2006
@@ -31,10 +31,7 @@
#include <cmath>
#include <limits>
#include <lemon/time_measure.h>
-
-#ifdef WIN32
-#include <lemon/bits/mingw32_rand.h>
-#endif
+#include <lemon/random.h>
namespace lemon {
@@ -241,7 +238,6 @@
double _temp = 1000.0, double _ann_fact = 0.9999) : max_iter(_max_iter),
max_no_impr(_max_no_impr), temp(_temp), ann_fact(_ann_fact)
{
- srand48(time(0));
}
/// \brief This is called when a neighbouring state gets accepted.
void acceptEvent() {}
@@ -261,7 +257,7 @@
/// \brief Decides whether to accept the current solution or not.
bool accept() {
double cost_diff = simann->getCurrCost() - simann->getPrevCost();
- return (drand48() <= exp(-(cost_diff / temp)));
+ return (rnd.getReal() <= exp(-(cost_diff / temp)));
}
/// \brief Destructor.
virtual ~SimpleController() {}
@@ -321,7 +317,6 @@
alpha(_alpha), beta(_beta), gamma(_gamma), end_time(_end_time),
ann_fact(_ann_fact), init_ann_fact(_ann_fact), start(false)
{
- srand48(time(0));
}
/// \brief Does initializations before each run.
void init() {
@@ -370,7 +365,7 @@
}
else {
double cost_diff = simann->getCurrCost() - simann->getPrevCost();
- return (drand48() <= exp(-(cost_diff / temp)));
+ return (rnd.getReal() <= exp(-(cost_diff / temp)));
}
}
/// \brief Destructor.
Modified: hugo/trunk/test/graph_utils_test.h
==============================================================================
--- hugo/trunk/test/graph_utils_test.h (original)
+++ hugo/trunk/test/graph_utils_test.h Mon Oct 2 18:11:00 2006
@@ -50,15 +50,14 @@
typedef typename Graph::NodeIt NodeIt;
typedef typename Graph::EdgeIt EdgeIt;
Graph graph;
- srand(time(0));
for (int i = 0; i < 10; ++i) {
graph.addNode();
}
DescriptorMap<Graph, Node> nodes(graph);
typename DescriptorMap<Graph, Node>::InverseMap invNodes(nodes);
for (int i = 0; i < 100; ++i) {
- int src = (int)(rand() / (RAND_MAX + 1.0) * invNodes.size());
- int trg = (int)(rand() / (RAND_MAX + 1.0) * invNodes.size());
+ int src = rnd.getInt(invNodes.size());
+ int trg = rnd.getInt(invNodes.size());
graph.addEdge(invNodes[src], invNodes[trg]);
}
typename Graph::template EdgeMap<bool> found(graph, false);
Modified: hugo/trunk/test/heap_test.h
==============================================================================
--- hugo/trunk/test/heap_test.h (original)
+++ hugo/trunk/test/heap_test.h Mon Oct 2 18:11:00 2006
@@ -45,7 +45,7 @@
std::vector<int> v(n);
for (int i = 0; i < n; ++i) {
- v[i] = rand() % 1000;
+ v[i] = rnd.getInt(1000);
heap.push(i, v[i]);
}
std::sort(v.begin(), v.end());
@@ -65,11 +65,11 @@
std::vector<int> v(n);
for (int i = 0; i < n; ++i) {
- v[i] = rand() % 1000;
+ v[i] = rnd.getInt(1000);
heap.push(i, v[i]);
}
for (int i = 0; i < n; ++i) {
- v[i] += rand() % 1000;
+ v[i] += rnd.getInt(1000);
heap.increase(i, v[i]);
}
std::sort(v.begin(), v.end());
Modified: hugo/trunk/test/test_tools.h
==============================================================================
--- hugo/trunk/test/test_tools.h (original)
+++ hugo/trunk/test/test_tools.h Mon Oct 2 18:11:00 2006
@@ -28,6 +28,8 @@
#include <lemon/concept_check.h>
#include <lemon/concept/graph.h>
+#include <lemon/random.h>
+
using namespace lemon;
//! \ingroup misc
@@ -176,16 +178,8 @@
return n;
}
-int _urandom_init() {
- int seed = time(0);
- srand(seed);
- return seed;
-}
-
int urandom(int n) {
- static int seed = _urandom_init();
- ignore_unused_variable_warning(seed);
- return (int)(rand() / (1.0 + RAND_MAX) * n);
+ return rnd.getInt(n);
}
#endif
More information about the Lemon-commits
mailing list