Tell me more ×
Code Review Stack Exchange is a question and answer site for peer programmer code reviews. It's 100% free, no registration required.

I am very much interested in the Reverse Cuthil McKee Algorithm. I have seen Fortran and C or C++ implementations of it, and I decided that it would be a nice exercise to implement it in Python. I know this algorithm is quite domain specific, but I would still be happy to see what kind of comments I get regarding:

  • Correctness - I am not sure my only test case works for others, although I did some comparison to the Octave and Matlab version.
  • Speed - Of course a C version would be faster. However, is there some Python improvements which can be done?
  • Readability - Is this code clear enough to other peer programmers?

The code:

import numpy as np

def getDegree(Graph):
    """
    find the degree of each node. That is the number
    of neighbours or connections.
    (number of non-zero elements) in each row minus 1.
    Graph is a Cubic Matrix.
    """
    degree = [0]*Graph.shape[0]
    for row in range(Graph.shape[0]):
        degree[row] = len(np.flatnonzero(Graph[row]))-1
    return degree

def getAdjcncy(Mat):
    """
    return the adjacncy matrix for each node
    """
    adj = [0]*Mat.shape[0]
    for i in xrange(Mat.shape[0]):
        q=np.flatnonzero(Mat[i])
        q=list(q)
        q.pop(q.index(i))
        adj[i] = q
    return adj

def RCM_loop(deg,start, adj,pivots,R):
    """
    Reverse Cuthil McKee ordering of an adjacency Matrix
    """
    digar=np.array(deg)
    # use np.where here to get indecies of minimums
    if start not in R:
        R.append(start)
    Q=adj[start]
    for idx, item in enumerate(Q):
        if item not in R:
            R.append(item)
    Q=adj[R[-1]]
    if set(Q).issubset(set(R)) and len(R) < len(deg) :
         p = pivots[0]
         pivots.pop(0)
         return RCM_loop(deg,p,adj,pivots,R)
    elif len(R) < len(deg):
         return RCM_loop(deg,R[-1],adj,pivots,R)
    else:
         R.reverse() 
         return R      

def test():
    """
    test the RCM loop
    """
    A = np.diag(np.ones(8))
    print A
    nzc=[[4],[2,5,7],[1,4],[6],[0,2],[1,7],[3],[1,5]]

    for i in range(len(nzc)):
        for j in nzc[i]:
            A[i,j]=1
    # define the Result queue
    R = ["C"]*A.shape[0]
    adj = getAdjcncy(A)
    degree = getDegree(A)
    digar=np.array(degree)
    pivots = list(np.where(digar == digar.min())[0])
    inl=[]
    Res = RCM_loop(degree,0, adj,pivots,inl)
    print degree
    print adj
    print "solution:", Res
    print "correct:", [6,3,7,5,1,2,4,0] 

if __name__ == '__main__':
    test()
share|improve this question

Know someone who can answer? Share a link to this question via email, Google+, Twitter, or Facebook.

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.