// Copyright 2004 The Trustees of Indiana University.

// Use, modification and distribution is subject to the Boost Software
// License, Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
// http://www.boost.org/LICENSE_1_0.txt)

//  Authors: Douglas Gregor
//           Andrew Lumsdaine
#ifndef BOOST_GRAPH_PLOD_GENERATOR_HPP
#define BOOST_GRAPH_PLOD_GENERATOR_HPP

#include <iterator>
#include <utility>
#include <boost/random/uniform_int.hpp>
#include <boost/shared_ptr.hpp>
#include <boost/graph/graph_traits.hpp>
#include <vector>
#include <cmath>

namespace boost {

  template<typename RandomGenerator, typename Graph>
  class plod_iterator
  {
    typedef std::vector<std::pair<std::size_t, std::size_t> > out_degrees_t;
    typedef typename graph_traits<Graph>::directed_category directed_category;

  public:
    typedef std::input_iterator_tag iterator_category;
    typedef std::pair<std::size_t, std::size_t> value_type;
    typedef const value_type& reference;
    typedef const value_type* pointer;
    typedef void difference_type;

    plod_iterator() 
      : gen(0), out_degrees(), degrees_left(0), allow_self_loops(false) { }

    plod_iterator(RandomGenerator& gen, std::size_t n,  
                  double alpha, double beta, bool allow_self_loops = false)
      : gen(&gen), n(n), out_degrees(new out_degrees_t(n)),
        degrees_left(0), allow_self_loops(allow_self_loops)
    {
      using std::pow;

      uniform_int<std::size_t> x(0, n-1);
      for (out_degrees_t::iterator i = out_degrees->begin();
           i != out_degrees->end(); ++i) {
        std::size_t xv = x(gen);
        i->first = i - out_degrees->begin();
        i->second = (xv == 0? 0 : std::size_t(beta * pow(xv, -alpha)));
        degrees_left += i->second;
      }

      next(directed_category());
    }

    reference operator*() const { return current; }
    pointer operator->() const { return &current; }
    
    plod_iterator& operator++()
    { 
      next(directed_category());
      return *this;
    }

    plod_iterator operator++(int)
    {
      plod_iterator temp(*this);
      ++(*this);
      return temp;
    }

    bool operator==(const plod_iterator& other) const
    { 
      return degrees_left == other.degrees_left; 
    }

    bool operator!=(const plod_iterator& other) const
    { return !(*this == other); }

  private:
    void next(directed_tag)
    {
      uniform_int<std::size_t> x(0, out_degrees->size()-1);
      std::size_t source;
      do {
        source = x(*gen);
      } while ((*out_degrees)[source].second == 0);
      current.first = (*out_degrees)[source].first;
      current.second = x(*gen);
      --degrees_left;
      if (--(*out_degrees)[source].second == 0) {
        (*out_degrees)[source] = out_degrees->back();
        out_degrees->pop_back();
      }
    }

    void next(undirected_tag)
    {
      std::size_t source, target;
      while (true) {
        /* We may get to the point where we can't actually find any
           new edges, so we just add some random edge and set the
           degrees left = 0 to signal termination. */
        if (out_degrees->size() < 2) {
          uniform_int<std::size_t> x(0, n);
          current.first  = x(*gen);
          do {
            current.second = x(*gen);
          } while (current.first == current.second && !allow_self_loops);
          degrees_left = 0;
          out_degrees->clear();
          return;
        }

        uniform_int<std::size_t> x(0, out_degrees->size()-1);

        // Select source vertex
        source = x(*gen);
        if ((*out_degrees)[source].second == 0) {
          (*out_degrees)[source] = out_degrees->back();
          out_degrees->pop_back();
          continue;
        } 

        // Select target vertex
        target = x(*gen);
        if ((*out_degrees)[target].second == 0) {
          (*out_degrees)[target] = out_degrees->back();
          out_degrees->pop_back();
          continue;
        } else if (source != target 
                   || (allow_self_loops && (*out_degrees)[source].second > 2)) {
          break;
        }
      }

      // Update degree counts
      --(*out_degrees)[source].second;
      --degrees_left;
      --(*out_degrees)[target].second;
      --degrees_left;
      current.first  = (*out_degrees)[source].first;
      current.second = (*out_degrees)[target].first;
    }

    RandomGenerator* gen;
    std::size_t n;
    shared_ptr<out_degrees_t> out_degrees;
    std::size_t degrees_left;
    bool allow_self_loops;
    value_type current;
  };

} // end namespace boost

#endif // BOOST_GRAPH_PLOD_GENERATOR_HPP