import numpy as np
from fractions import Fraction
from scipy.sparse import isspmatrix, dok_matrix, find
from .mcl import sparse_allclose

def is_undirected(matrix):
    """
    Determine if the matrix reprensents a directed graph

    :param matrix: The matrix to tested
    :returns: boolean
    """
    if isspmatrix(matrix):
        return sparse_allclose(matrix, matrix.transpose())
    
    return np.allclose(matrix, matrix.T)


def convert_to_adjacency_matrix(matrix):
    """
    Converts transition matrix into adjacency matrix

    :param matrix: The matrix to be converted
    :returns: adjacency matrix
    """
    for i in range(matrix.shape[0]):
        
        if isspmatrix(matrix):
            col = find(matrix[:,i])[2]
        else:
            col = matrix[:,i]

        coeff = max( Fraction(c).limit_denominator().denominator for c in col )
        matrix[:,i] *= coeff

    return matrix


def modularity(matrix, result):
    """
    Compute the modularity

    :param matrix: The adjacency matrix
    :param result: The matrix result of mcl
    :returns: modularity value
    """
    matrix = convert_to_adjacency_matrix(matrix)

    # makes sure result[i,j]>0 if i and j belong to the same cluster
    result = result + result.T
    
    m = matrix.sum()

    if isspmatrix(matrix):
        matrix_2 = matrix.tocsr(copy=True)
    else :
        matrix_2 = matrix

    if is_undirected(matrix):
        expected = lambda i,j : (( matrix_2[i,:].sum() + matrix[:,i].sum() )*
                                 ( matrix[:,j].sum() + matrix_2[j,:].sum() ))
    else:
        expected = lambda i,j : ( matrix_2[i,:].sum()*matrix[:,j].sum() )
    
    indices = np.array(result.nonzero())
    Q = sum( matrix[i, j] - expected(i, j)/m for i, j in indices.T if i != j )/m
    
    return Q
    

@@ -43,16 +43,52 @@
 raction%0A
+from itertools import combinations%0A%0A
 from sci
@@ -943,16 +943,277 @@
 atrix%0A%0A%0A
+def delta_matrix(matrix, clusters):%0A    if isspmatrix(matrix):%0A        delta = dok_matrix(matrix.shape)%0A    else :%0A        delta = np.zeros(matrix.shape)%0A%0A    for i in clusters :%0A        for j in combinations(l, 2):%0A            delta%5Bj%5D = 1%0A%0A    return delta%0A%0A%0A
 def modu
@@ -1227,22 +1227,24 @@
 matrix, 
-result
+clusters
 ):%0A    %22
@@ -1329,40 +1329,55 @@
 ram 
-result: The matrix result of mcl
+clusters: The clusters returned by get_clusters
 %0A   
@@ -1465,108 +1465,43 @@
 ix)%0A
-%0A
     
-# makes sure result%5Bi,j%5D%3E0 if i and j belong to the same cluster%0A    result = result + result.T
+delta  = delta_matrix(clusters)
 %0A   
@@ -1934,14 +1934,13 @@
 ray(
-resu
+de
 lt
+a
 .non
@@ -2018,18 +2018,8 @@
 s.T 
-if i != j 
 )/m%0A
@@ -2036,13 +2036,8 @@
 eturn Q%0A
-    %0A

Add delta_matrix where delta[i,j] = 1 if i,j belong to the same cluster