williamr@2: // Copyright 2004, 2005 The Trustees of Indiana University. williamr@2: williamr@2: // Distributed under the Boost Software License, Version 1.0. williamr@2: // (See accompanying file LICENSE_1_0.txt or copy at williamr@2: // http://www.boost.org/LICENSE_1_0.txt) williamr@2: williamr@2: // Authors: Jeremiah Willcock williamr@2: // Douglas Gregor williamr@2: // Andrew Lumsdaine williamr@2: #ifndef BOOST_GRAPH_ERDOS_RENYI_GENERATOR_HPP williamr@2: #define BOOST_GRAPH_ERDOS_RENYI_GENERATOR_HPP williamr@2: williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: #include williamr@2: williamr@2: namespace boost { williamr@2: williamr@2: template williamr@2: class erdos_renyi_iterator williamr@2: { williamr@2: typedef typename graph_traits::directed_category directed_category; williamr@2: typedef typename graph_traits::vertices_size_type vertices_size_type; williamr@2: typedef typename graph_traits::edges_size_type edges_size_type; williamr@2: williamr@2: BOOST_STATIC_CONSTANT williamr@2: (bool, williamr@2: is_undirected = (is_base_and_derived::value williamr@2: || is_same::value)); williamr@2: williamr@2: public: williamr@2: typedef std::input_iterator_tag iterator_category; williamr@2: typedef std::pair value_type; williamr@2: typedef const value_type& reference; williamr@2: typedef const value_type* pointer; williamr@2: typedef void difference_type; williamr@2: williamr@2: erdos_renyi_iterator() : gen(), n(0), edges(0), allow_self_loops(false) {} williamr@2: erdos_renyi_iterator(RandomGenerator& gen, vertices_size_type n, williamr@2: double fraction = 0.0, bool allow_self_loops = false) williamr@2: : gen(&gen), n(n), edges(edges_size_type(fraction * n * n)), williamr@2: allow_self_loops(allow_self_loops) williamr@2: { williamr@2: if (is_undirected) edges = edges / 2; williamr@2: next(); williamr@2: } williamr@2: williamr@2: erdos_renyi_iterator(RandomGenerator& gen, vertices_size_type n, williamr@2: edges_size_type m, bool allow_self_loops = false) williamr@2: : gen(&gen), n(n), edges(m), williamr@2: allow_self_loops(allow_self_loops) williamr@2: { williamr@2: next(); williamr@2: } williamr@2: williamr@2: reference operator*() const { return current; } williamr@2: pointer operator->() const { return ¤t; } williamr@2: williamr@2: erdos_renyi_iterator& operator++() williamr@2: { williamr@2: --edges; williamr@2: next(); williamr@2: return *this; williamr@2: } williamr@2: williamr@2: erdos_renyi_iterator operator++(int) williamr@2: { williamr@2: erdos_renyi_iterator temp(*this); williamr@2: ++(*this); williamr@2: return temp; williamr@2: } williamr@2: williamr@2: bool operator==(const erdos_renyi_iterator& other) const williamr@2: { return edges == other.edges; } williamr@2: williamr@2: bool operator!=(const erdos_renyi_iterator& other) const williamr@2: { return !(*this == other); } williamr@2: williamr@2: private: williamr@2: void next() williamr@2: { williamr@2: uniform_int rand_vertex(0, n-1); williamr@2: current.first = rand_vertex(*gen); williamr@2: do { williamr@2: current.second = rand_vertex(*gen); williamr@2: } while (current.first == current.second && !allow_self_loops); williamr@2: } williamr@2: williamr@2: RandomGenerator* gen; williamr@2: vertices_size_type n; williamr@2: edges_size_type edges; williamr@2: bool allow_self_loops; williamr@2: value_type current; williamr@2: }; williamr@2: williamr@2: template williamr@2: class sorted_erdos_renyi_iterator williamr@2: { williamr@2: typedef typename graph_traits::directed_category directed_category; williamr@2: typedef typename graph_traits::vertices_size_type vertices_size_type; williamr@2: typedef typename graph_traits::edges_size_type edges_size_type; williamr@2: williamr@2: BOOST_STATIC_CONSTANT williamr@2: (bool, williamr@2: is_undirected = (is_base_and_derived::value williamr@2: || is_same::value)); williamr@2: williamr@2: public: williamr@2: typedef std::input_iterator_tag iterator_category; williamr@2: typedef std::pair value_type; williamr@2: typedef const value_type& reference; williamr@2: typedef const value_type* pointer; williamr@2: typedef void difference_type; williamr@2: williamr@2: sorted_erdos_renyi_iterator() williamr@2: : gen(), rand_vertex(0.5), n(0), allow_self_loops(false), williamr@2: src((std::numeric_limits::max)()), tgt(0), prob(0) {} williamr@2: sorted_erdos_renyi_iterator(RandomGenerator& gen, vertices_size_type n, williamr@2: double prob = 0.0, williamr@2: bool allow_self_loops = false) williamr@2: : gen(), williamr@2: // The "1.0 - prob" in the next line is to work around a Boost.Random williamr@2: // (and TR1) bug in the specification of geometric_distribution. It williamr@2: // should be replaced by "prob" when the issue is fixed. williamr@2: rand_vertex(1.0 - prob), williamr@2: n(n), allow_self_loops(allow_self_loops), src(0), tgt(0), prob(prob) williamr@2: { williamr@2: this->gen.reset(new uniform_01(gen)); williamr@2: williamr@2: if (prob == 0.0) {src = (std::numeric_limits::max)(); return;} williamr@2: next(); williamr@2: } williamr@2: williamr@2: reference operator*() const { return current; } williamr@2: pointer operator->() const { return ¤t; } williamr@2: williamr@2: sorted_erdos_renyi_iterator& operator++() williamr@2: { williamr@2: next(); williamr@2: return *this; williamr@2: } williamr@2: williamr@2: sorted_erdos_renyi_iterator operator++(int) williamr@2: { williamr@2: sorted_erdos_renyi_iterator temp(*this); williamr@2: ++(*this); williamr@2: return temp; williamr@2: } williamr@2: williamr@2: bool operator==(const sorted_erdos_renyi_iterator& other) const williamr@2: { return src == other.src && tgt == other.tgt; } williamr@2: williamr@2: bool operator!=(const sorted_erdos_renyi_iterator& other) const williamr@2: { return !(*this == other); } williamr@2: williamr@2: private: williamr@2: void next() williamr@2: { williamr@2: using std::sqrt; williamr@2: using std::floor; williamr@2: williamr@2: // In order to get the edges from the generator in sorted order, one williamr@2: // effective (but slow) procedure would be to use a williamr@2: // bernoulli_distribution for each legal (src, tgt) pair. Because of the williamr@2: // O(n^2) cost of that, a geometric distribution is used. The geometric williamr@2: // distribution tells how many times the bernoulli_distribution would williamr@2: // need to be run until it returns true. Thus, this distribution can be williamr@2: // used to step through the edges which are actually present. Everything williamr@2: // beyond "tgt += increment" is done to effectively convert linear williamr@2: // indexing (the partial sums of the geometric distribution output) into williamr@2: // graph edges. williamr@2: assert (src != (std::numeric_limits::max)()); williamr@2: vertices_size_type increment = rand_vertex(*gen); williamr@2: tgt += increment; williamr@2: if (is_undirected) { williamr@2: // Update src and tgt based on position of tgt williamr@2: // Basically, we want the greatest src_increment such that (in \bbQ): williamr@2: // src_increment * (src + allow_self_loops + src_increment - 1/2) <= tgt williamr@2: // The result of the LHS of this, evaluated with the computed williamr@2: // src_increment, is then subtracted from tgt williamr@2: double src_minus_half = (src + allow_self_loops) - 0.5; williamr@2: double disc = src_minus_half * src_minus_half + 2 * tgt; williamr@2: double src_increment_fp = floor(sqrt(disc) - src_minus_half); williamr@2: vertices_size_type src_increment = vertices_size_type(src_increment_fp); williamr@2: if (src + src_increment >= n) { williamr@2: src = n; williamr@2: } else { williamr@2: tgt -= (src + allow_self_loops) * src_increment + williamr@2: src_increment * (src_increment - 1) / 2; williamr@2: src += src_increment; williamr@2: } williamr@2: } else { williamr@2: // Number of out edge positions possible from each vertex in this graph williamr@2: vertices_size_type possible_out_edges = n - (allow_self_loops ? 0 : 1); williamr@2: src += (std::min)(n - src, tgt / possible_out_edges); williamr@2: tgt %= possible_out_edges; williamr@2: } williamr@2: // Set end of graph code so (src, tgt) will be the same as for the end williamr@2: // sorted_erdos_renyi_iterator williamr@2: if (src >= n) {src = (std::numeric_limits::max)(); tgt = 0;} williamr@2: // Copy (src, tgt) into current williamr@2: current.first = src; williamr@2: current.second = tgt; williamr@2: // Adjust for (src, src) edge being forbidden williamr@2: if (!allow_self_loops && tgt >= src) ++current.second; williamr@2: } williamr@2: williamr@2: shared_ptr > gen; williamr@2: geometric_distribution rand_vertex; williamr@2: vertices_size_type n; williamr@2: bool allow_self_loops; williamr@2: vertices_size_type src, tgt; williamr@2: value_type current; williamr@2: double prob; williamr@2: }; williamr@2: williamr@2: } // end namespace boost williamr@2: williamr@2: #endif // BOOST_GRAPH_ERDOS_RENYI_GENERATOR_HPP