LCPQ/dft_tools
3
0
Fork 0
mirror of https://github.com/triqs/dft_tools synced 2024-12-25 05:43:40 +01:00
Code Issues Releases Wiki Activity

Change implementation of random_generator & python wrapper.

Browse Source 
- Use a new buffered_function to replace the complicated generator code from ALPS.
- Clean the implementation of the random_generator
- update the documentation
- update to the new python wrapper (could not be done with the previous
  version, because of lack of move constructor).
This commit is contained in:
Olivier Parcollet 2014-05-29 21:35:00 +02:00
parent 2b5dd322ce
commit e6529b608e
12 changed files with 270 additions and 248 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
doc/reference/python/random/random.py
View File

@ -1,7 +1,7 @@
from pytriqs.random import *
from pytriqs.random_generator import *
from pytriqs.plot.mpl_interface import *
r = RandomGenerator("lagged_fibonacci607", 237489)
l = []
for i in range(10000): l += [r.rand(),]
plt.hist(l, 30, normed=True);
for i in range(10000): l += [r(),]
plt.hist(l, 30, normed=True)

29
doc/reference/python/random/random.rst
View File

@ -1,33 +1,34 @@
.. index:: Random number generator
.. module:: pytriqs.random
.. module:: pytriqs.random_generator
.. _random_generator:
Random generator class
======================
TRIQS has a simple random number generator class called ``RandomGenerator``. It is based on boost
just like the C++ random generator provided by TRIQS.
TRIQS exposes to python the random number generators used in C++,
in the module ``RandomGenerator``.
The generators are the boost random generators.
Usage
-----
The generator is constructed from a name and a seed::
The generator is constructed from a name (the name of the boost generator) and a seed::
from pytriqs.random import *
from pytriqs.random_generator import *
r = RandomGenerator("mt19937", 237849)
A list of possible random generator names is obtained with::
A list of available random generators is obtained with::
print available_generator_names()
print random_generator_names_list()
Then you can either generate float number on the interval :math:`[0,1[` using
the ``rand()`` method or integer numbers in the inverval :math:`[0,N-1]` using
``int_rand(N)``::
Then you can either generate float number on the interval :math:`[0,1[`
simply by calling the generator, or integer numbers in the inverval :math:`[0,N-1]` by calling
it with `N` ::
print r.rand()
print r.int_rand(10)
print r()
print r(10)
Example
-------
@ -42,7 +43,7 @@ Here's a simple example showing how to use the generator.
Complete reference
------------------
.. autoclass:: pytriqs.random.RandomGenerator
.. autoclass:: pytriqs.random_generator.RandomGenerator
:members:
.. autofunction:: pytriqs.random.available_generator_names
.. autofunction:: pytriqs.random_generator.random_generator_names_list

5
pytriqs/random/__init__.py
View File

@ -1,5 +0,0 @@
from random_generator import RandomGenerator, available_generator_names
__all__ = ['RandomGenerator','available_generator_names']

49
pytriqs/random/random_generator.pyx
View File

@ -1,49 +0,0 @@
from dcomplex cimport *
from libcpp.string cimport string
# The c++ random generator class
# We wrap two members only
cdef extern from "triqs/mc_tools/random_generator.hpp":
cdef cppclass random_generator "triqs::mc_tools::random_generator":
random_generator(string, long) except +
double operator() () except +
int operator() (int) except +
# This is the wrapping of the static member random_generator_names
# It is a bit of a hack but there is no notion of ststic members in cython
cdef extern from "triqs/mc_tools/random_generator.hpp" namespace "triqs::mc_tools::random_generator":
string random_generator_names(string) except +
# The python RandomGenerator class
cdef class RandomGenerator:
cdef random_generator * _c
def __init__(self, name, seed):
"""This is a random number generator class based on boost.
:param name: (string) Name of the random number generator
:param seed: (int) Random number seed
"""
self._c = new random_generator(name, seed)
def __dealloc__(self):
del self._c
def rand(self):
"""Generate a float random number in [0,1["""
return self._c[0]()
def int_rand(self, N):
"""Generate an integer random number in [0,N-1]"""
return self._c[0](N)
# The list of generator names accessed through the static member
def available_generator_names():
"""Get a list of available random generator names"""
return random_generator_names(",").split(',')

7
pytriqs/random/CMakeLists.txt → pytriqs/random_generator/CMakeLists.txt
View File

@ -2,7 +2,8 @@ SET(PYTHON_SOURCES
${CMAKE_CURRENT_SOURCE_DIR}/__init__.py
)
install (FILES ${PYTHON_SOURCES} DESTINATION ${TRIQS_PYTHON_LIB_DEST}/random)
install (FILES ${PYTHON_SOURCES} DESTINATION ${TRIQS_PYTHON_LIB_DEST}/random_generator)
# Build C extension module
triqs_python_extension(random_generator random_generator)
# the c++ code
cython_module(random random_generator random)

5
pytriqs/random_generator/__init__.py Normal file
View File

@ -0,0 +1,5 @@
from random_generator import RandomGenerator, random_generator_names_list
__all__ = ['RandomGenerator','random_generator_names_list']

46
pytriqs/random_generator/random_generator_desc.py Normal file
View File

@ -0,0 +1,46 @@
from wrap_generator import *
module = module_(full_name = "pytriqs.random_generator", doc = "")
module.add_include("<triqs/mc_tools/random_generator.hpp>")
module.add_using("namespace triqs::mc_tools")
# Not needed. Reorganize the hpp wrapper tool
module.add_include("<triqs/h5.hpp>")
module.add_include("<triqs/arrays.hpp>")
# --------- RandomGenerator ----------------------------------
r = class_(py_type = "RandomGenerator",
c_type = "random_generator",
c_type_absolute = "triqs::mc_tools::random_generator",
)
r.add_constructor(signature = "(std::string name, int seed)",
doc =
"""
This is a random number generator class based on boost.
name Name of the random number generator
seed Random number seed
""")
r.add_call(signature = "int(int N)", doc = """Generate an integer random number in [0,N-1]""")
r.add_call(signature = "double()", doc = """Generate a float random number in [0,1[""")
module.add_class(r)
# --------- Module functions ----------------------------------
module.add_function(name = "random_generator_names_list",
signature = "std::vector<std::string>()",
doc = """Get a list of available random generator names"""
)
########################
## Code generation
########################
if __name__ == '__main__' :
module.generate_code(mako_template = sys.argv[1], wrap_file = sys.argv[2])
module.generate_py_converter_header(mako_template = sys.argv[3], wrap_file = sys.argv[4])

37
test/triqs/utility/buffered_function.cpp Normal file
View File

@ -0,0 +1,37 @@
/*******************************************************************************
*
* TRIQS: a Toolbox for Research in Interacting Quantum Systems
*
* Copyright (C) 2014 by O. Parcollet
*
* TRIQS 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.
*
* TRIQS 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
* TRIQS. If not, see <http://www.gnu.org/licenses/>.
*
******************************************************************************/
#include <triqs/utility/buffered_function.hpp>
#include <iostream>
int main(int argc, char **argv) {
// a function that generates all the square ....
int x = 0;
auto f = [x]() mutable { return x * (x++); };
// C++14 : init-capture
// auto f = [x=0]() mutable { return x*(x++);};
auto gen = triqs::utility::buffered_function<double>(f, 5);
for (int u = 0; u < 22; ++u)
if (gen() != u * u) throw "Error";
}

96
triqs/mc_tools/generator.hpp
View File

@ -1,96 +0,0 @@
// (C) Copyright 2008-10 Lukas Gamper <gamperl -at- gmail.com>
// Brigitte Surer <surerb -at- phys.ethz.ch>
// Bela Bauer <bauerb -at- phys.ethz.ch>
//
// downgraded by O. Parcollet for compilation with icc11 & boost 1.45
//
// Use, modification, and distribution are subject to the Boost Software
// License, Version 1.0. (See at <http://www.boost.org/LICENSE_1_0.txt>.)
#include <vector>
#include <algorithm>
#include <boost/function.hpp>
#include <boost/ref.hpp>
#include <boost/shared_ptr.hpp>
#include <boost/make_shared.hpp>
#include <boost/lambda/lambda.hpp>
#include <boost/lambda/bind.hpp>
#include <boost/lambda/construct.hpp>
#include <boost/lambda/exceptions.hpp>
namespace boost{
template<typename R> class generator {
public:
static std::size_t const buffer_size = 1000;
template<typename Generator> generator(Generator const & gen, std::size_t size = buffer_size)
:
buf(make_shared<buffer>(gen, size))
, cur(buf->cur)
, end(buf->end())
{}
inline R const operator()() {
if (cur == end)
buf->fill();
return *cur++;
}
inline R const preview() {
if (cur == end)
buf->fill();
return *cur;
}
generator(generator<R> const & rhs)
:
buf(new buffer (*rhs.buf))
, cur(buf->cur)
, end(buf->end())
{}
generator<R> & operator=(generator<R> rhs) {
swap(*this, rhs);
return *this;
}
private:
class buffer: public std::vector<R> {
public:
friend class generator<R>;
template<typename Generator> buffer(Generator const & gen, std::size_t size)
: std::vector<R>(size)
, engine(make_shared<Generator>(gen))
, cur(this->end())
, fill_helper(lambda::bind(std::generate<typename buffer::iterator, Generator &>, lambda::_1, lambda::_2, ref(*reinterpret_cast<Generator *>(engine.get()))))
, clone(lambda::bind(buffer::template type_keeper<Generator>, lambda::_1, lambda::_2))
{}
template<typename Generator> buffer(reference_wrapper<Generator> gen, std::size_t size)
: std::vector<R>(size)
, cur(this->end())
, fill_helper(lambda::bind(std::generate<typename buffer::iterator, Generator &>, lambda::_1, lambda::_2, gen))
{}
buffer(buffer const & rhs)
: std::vector<R>(rhs)
, cur(this->begin() + (rhs.cur - rhs.begin()))
, fill_helper(rhs.fill_helper)
, clone(rhs.clone)
{
if (rhs.engine)
clone(*this, rhs);
}
inline void fill() {
fill_helper(this->begin(), this->end());
cur = this->begin();
}
private:
template<typename Generator> static void type_keeper(buffer & lhs, buffer const & rhs) {
lhs.engine = make_shared<Generator>(*reinterpret_cast<Generator *>(rhs.engine.get()));
lhs.fill_helper = lambda::bind(std::generate<typename buffer::iterator, Generator &>, lambda::_1, lambda::_2, ref(*reinterpret_cast<Generator *>(lhs.engine.get())));
}
shared_ptr<void> engine;
typename std::vector<R>::const_iterator cur;
function<void(typename std::vector<R>::iterator, typename std::vector<R>::iterator)> fill_helper;
function<void(buffer &, buffer const &)> clone;
};
shared_ptr<buffer> buf;
typename buffer::const_iterator & cur;
typename buffer::const_iterator end;
};
}

74
triqs/mc_tools/random_generator.cpp
View File

@ -1,9 +1,8 @@
/*******************************************************************************
*
* TRIQS: a Toolbox for Research in Interacting Quantum Systems
*
* Copyright (C) 2011 by M. Ferrero, O. Parcollet
* Copyright (C) 2011-2014 by M. Ferrero, O. Parcollet
*
* TRIQS 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
@ -19,7 +18,6 @@
* TRIQS. If not, see <http://www.gnu.org/licenses/>.
*
******************************************************************************/
#include "random_generator.hpp"
#include "./MersenneRNG.hpp"
#include <boost/random.hpp>
@ -32,50 +30,52 @@
#include <sstream>
#include <boost/preprocessor/seq.hpp>
#include <boost/preprocessor/control/if.hpp>
// List of All available Boost random number generator
#define RNG_LIST \
(mt19937)(mt11213b)(lagged_fibonacci607)(lagged_fibonacci1279)(lagged_fibonacci2281)(lagged_fibonacci3217)( \
lagged_fibonacci4423)(lagged_fibonacci9689)(lagged_fibonacci19937)(lagged_fibonacci23209)(lagged_fibonacci44497)(ranlux3)
namespace triqs {
namespace mc_tools {
random_generator::random_generator(std::string const& RandomGeneratorName, uint32_t seed_) {
_name = RandomGeneratorName;
if (RandomGeneratorName == "") {
gen = utility::buffered_function<double>(mc_tools::RandomGenerators::RandMT(seed_));
return;
}
boost::uniform_real<> dis;
#define AS_STRING(X) AS_STRING2(X)
#define AS_STRING2(X) #X
// List of All available Boost random number generator
#define RNG_LIST (mt19937)(mt11213b)\
(lagged_fibonacci607) (lagged_fibonacci1279) (lagged_fibonacci2281) (lagged_fibonacci3217) (lagged_fibonacci4423)\
(lagged_fibonacci9689) (lagged_fibonacci19937) (lagged_fibonacci23209) (lagged_fibonacci44497) (ranlux3)
// now boost random number generators
#define DRNG(r, data, XX) \
if (RandomGeneratorName == AS_STRING(XX)) { \
gen = utility::buffered_function<double>(boost::variate_generator<boost::XX, boost::uniform_real<>>(boost::XX(seed_), dis)); \
return; \
}
BOOST_PP_SEQ_FOR_EACH(DRNG, ~, RNG_LIST)
namespace triqs {
namespace mc_tools {
typedef boost::generator <double> gen_type;
inline gen_type * choose_gen(std::string const & RandomGeneratorName, std::size_t seed_, boost::shared_ptr<void> &ptr) {
if (RandomGeneratorName=="") return new gen_type( mc_tools::RandomGenerators::RandMT (seed_));
// now boost random number generators
#define DRNG(r,data,XX) if (RandomGeneratorName==AS_STRING(XX)) { \
boost::shared_ptr<boost::XX> localptr = boost::make_shared<boost::XX>(seed_);\
ptr = localptr; boost::uniform_real<> dis;\
return new gen_type( boost::variate_generator<boost::XX&, boost::uniform_real<> >(*localptr,dis));}
BOOST_PP_SEQ_FOR_EACH(DRNG,~,RNG_LIST)
TRIQS_RUNTIME_ERROR<<"The random generator "<<RandomGeneratorName<<" is not recognized";
TRIQS_RUNTIME_ERROR << "The random generator " << RandomGeneratorName << " is not recognized";
}
//---------------------------------------------
random_generator::random_generator(std::string const & RandomGeneratorName, std::size_t seed_ ) :
rng_ptr(), gen (choose_gen(RandomGeneratorName,seed_,rng_ptr)), name(RandomGeneratorName),seed(seed_) {}
std::string random_generator_names(std::string const &sep) {
#define PR(r, sep, p, XX) BOOST_PP_IF(p, +sep +, ) std::string(AS_STRING(XX))
return BOOST_PP_SEQ_FOR_EACH_I(PR, sep, RNG_LIST);
}
//---------------------------------------------
random_generator::random_generator( random_generator const & p) :
rng_ptr(), gen (choose_gen(p.name,p.seed,rng_ptr)), name (p.name),seed(p.seed) {}
//---------------------------------------------
std::string random_generator::random_generator_names(std::string const & sep) {
#define PR(r,sep,p,XX) BOOST_PP_IF(p,+ sep +,) std::string(AS_STRING(XX))
return BOOST_PP_SEQ_FOR_EACH_I (PR,sep,RNG_LIST);
std::vector<std::string> random_generator_names_list() {
std::vector<std::string> res;
#define PR2(r, sep, p, XX) res.push_back(AS_STRING(XX));
BOOST_PP_SEQ_FOR_EACH_I(PR2, sep, RNG_LIST);
return res;
}
}
}

74
triqs/mc_tools/random_generator.hpp
View File

@ -2,7 +2,7 @@
*
* TRIQS: a Toolbox for Research in Interacting Quantum Systems
*
* Copyright (C) 2011 by M. Ferrero, O. Parcollet
* Copyright (C) 2011-2014 by M. Ferrero, O. Parcollet
*
* TRIQS 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
@ -18,62 +18,70 @@
* TRIQS. If not, see <http://www.gnu.org/licenses/>.
*
******************************************************************************/
#ifndef POLYMORPH_RANDOM_GENERATOR_H
#define POLYMORPH_RANDOM_GENERATOR_H
#pragma once
#include <triqs/utility/first_include.hpp>
#include <boost/scoped_ptr.hpp>
#include "../utility/exceptions.hpp"
#include "./generator.hpp"
#include "../utility/buffered_function.hpp"
#include "math.h"
#include <string>
#include <assert.h>
#include <type_traits>
namespace triqs {
namespace mc_tools {
/// Return a list of the names of available generators, with separator sep
std::string random_generator_names(std::string const & sep=" ");
std::vector<std::string> random_generator_names_list();
/**
* This is a random generator that take the name of the random generator at construct
* and serve random numbers, buffered with a version boost::generator ,
* to hide the cost of dynamical polymorphism.
* Random generator, adapting the boost random generator.
*
* The name of the generator is given at construction, and its type is erased in this class.
* For performance, the call to the generator is bufferized, with chunks of 1000 numbers.
*/
class random_generator {
boost::shared_ptr<void> rng_ptr;
boost::scoped_ptr< boost::generator <double> > gen;
std::string name;
void operator = ( random_generator const & p); //forbid
std::size_t seed;
public:
utility::buffered_function<double> gen;
std::string _name;
/// Takes the boost name of the generator e.g. mt19937,...
random_generator(std::string const & RandomGeneratorName, std::size_t seed_ );
public:
/** Constructor
* @param RandomGeneratorName : Name of a boost generator e.g. mt19937, or "" (another Mersenne Twister).
* @param seed : The seed of the random generator
*/
random_generator(std::string const& RandomGeneratorName, uint32_t seed_);
random_generator() : random_generator("mt19937", 198) {}
///
random_generator( random_generator const & p);
random_generator(random_generator const& p);
/// Return a list of the names of available generators, with separator sep
static std::string random_generator_names(std::string const & sep=" ");
random_generator(random_generator&&) = default;
//
random_generator & operator=(random_generator&&) = default;
/// Name of the random generator
std::string name() const { return _name; }
/// Returns a integer in [0,i-1] with flat distribution
#define INTEGER_OVERLOAD(T) T operator()(T i) { return (i==1 ? 0 : T(floor(i*((*gen)()))));}
INTEGER_OVERLOAD(int)
INTEGER_OVERLOAD(size_t)
#undef INTEGER_OVERLOAD
template <typename T> typename std::enable_if<std::is_integral<T>::value, T>::type operator()(T i) {
return (i == 1 ? 0 : T(floor(i * (gen()))));
}
/// Returns a double in [0,1[ with flat distribution
double preview() { return gen->preview();}
double preview() { return gen.preview();}
/// Returns a double in [0,1[ with flat distribution
double operator()() { return ((*gen)());}
double operator()() { return gen(); }
/// Returns a double in [0,x[ with flat distribution
double operator()(double x) { return x*((*gen)());}
double operator()(double x) { return x * (gen()); }
/// Returns a double in [a,b[ with flat distribution
double operator()(double a, double b) { assert (b>a); return a + (b-a)*((*gen)());}
double operator()(double a, double b) {
assert(b > a);
return a + (b - a) * (gen());
}
};
}
}
#endif

74
triqs/utility/buffered_function.hpp Normal file
View File

@ -0,0 +1,74 @@
/*******************************************************************************
*
* TRIQS: a Toolbox for Research in Interacting Quantum Systems
*
* Copyright (C) 2014 by O. Parcollet
*
* TRIQS 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.
*
* TRIQS 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
* TRIQS. If not, see <http://www.gnu.org/licenses/>.
*
******************************************************************************/
#pragma once
#include <vector>
#include <functional>
namespace triqs {
namespace utility {
/**
* A simple buffer for a generator.
* Given a function, it provides a buffer of this function
* Advantage :
* - do not pay the indirection cost at each call, but once every size call.
* - erase the function type
* It is a semi-regular type.
*/
template <typename R> struct buffered_function {
/// Default constructor : no function bufferized. () will throw in this state
buffered_function() = default;
/** Constructor
*
* @tparam Function : type of the function to bufferize
* @param f : function to bufferize
* @param size : size of the buffer [optional]
*/
template <typename Function> buffered_function(Function f, size_t size = 1000) : buffer(size) {
refill = [f](buffered_function *bf) mutable { // without the mutable, the () of the lambda object is const, hence f
for (auto &x : bf->buffer) x = f();
bf->index = 0;
};
refill(this); // first filling of the buffer
}
/// Returns the next element. Refills the buffer if necessary.
R operator()() {
if (index > buffer.size() - 1) refill(this);
return buffer[index++];
}
/// Returns the future next element, without increasing the index. Refills the buffer if necessary.
R preview() {
if (index > buffer.size() - 1) refill(this);
return buffer[index];
}
private:
size_t index;
std::vector<R> buffer;
std::function<void(buffered_function *)> refill; // this refills the buffer and reset index of a buffered_function.
// NB : can not capture this in refill because we want the object to be copyable and movable
};
}
}