"""
Sequential Graph Generation
"""
__author__ = 'Brian Cloteaux (brian.cloteaux@nist.gov)'

# Mathematical and Computational Sciences Division
# National Institute of Standards and Technology,
# Gaithersburg, MD USA
#
# This software was developed at the National Institute of Standards and
# Technology (NIST) by employees of the Federal Government in the course
# of their official duties. Pursuant to title 17 Section 105 of the
# United States Code, this software is not subject to copyright protection
# and is in the public domain. NIST assumes no responsibility whatsoever for
# its use by other parties, and makes no guarantees, expressed or implied,
# about its quality, reliability, or any other characteristic.

import copy
import random

import networkx as nx
import numpy as np

import ordered_sequence

def blitzstein_diaconis_random_graph(deg_sequence, create_using=None,
                                     name=None):
    """ 
    Returns a random graph with the given degree sequence.

    Notes
    -----
    This algorithm is based on Blitzstein and Diaconis [1]_.

    References
    ----------
    ..[1] J. Blitzstein and P. Diaconis, "A Sequential Importance Sampling
          Algorithm for Generating Random Graphs with Prescribed Degrees,"
          Internet Mathematics, vol. 6, no. 4, p. 489, 2010.
          http://dx.doi.org/10.1080/15427951.2010.557277
    """ 

    sequence = ordered_sequence.ordered_sequence(deg_sequence)
    degs = list(deg_sequence)
    degs.sort(reverse=True)

    graph = nx.Graph()
    if name is not None:
        graph.name = name
    graph.add_nodes_from(xrange(sequence.get_length()))

    while sequence.get_largest_val()>0:
        # Choose the least i with  d_i a minimal positive entry
        minVal = len(degs)
        indx = -1
        for i in reversed(xrange(len(degs))):
            if degs[i]>0 and degs[i] < minVal:
                minVal = degs[i]
                indx = i

        while degs[indx] > 0:
            candList = []

            # Compute candidate list
            for j in xrange(len(degs)):
                if degs[j]>0 and j!=indx and not graph.has_edge(indx,j):
                    val1 = degs[indx]
                    val2 = degs[j]
                    sequence.decrement( val1 )
                    sequence.decrement( val2 )
                    if sequence.is_graphical():
                        candList.append(j)
                    sequence.increment( val1-1 )
                    sequence.increment( val2-1 )

            # Pick candidate edge with probability proportional to its degree
            weights = []
            indices = [] 
            counter = 0
            for val in candList:
                indices.append(val)
                weights.append(degs[val])
                counter += 1
            s = sum(weights)
            weights = np.cumsum(weights)
            rv =  random.randint(1,s)
            x = np.searchsorted(weights,rv) 
            n = indices[x]

            # Add edge {i,j} to E and update d
            graph.add_edge(indx,n)
            sequence.decrement(degs[indx])
            sequence.decrement(degs[n])
            degs[indx] -= 1
            degs[n] -= 1

    return graph