# Applying the result of Cuthill-McKee in SciPy

I have applied SciPy's implementation of the Cuthill-McKee algorithm to a $48 \times 48$ sparse non-symmetric matrix in Compressed Sparse Row (CSR) format and the output is an array of length $48$ which the documentation calls a "permutation array".

How do I use this "permutation" array?

Obviously, the array is in some way representing how my matrix is to be reordered, but a permutation array of length $48$ is only enough to reorder one axis. What am I missing?

• Is your matrix symmetric and positive definite? In that case you would normally apply the permutation symmetrically to both rows and columns of the matrix. – Brian Borchers Aug 24 '16 at 20:50
• Positive definite, but not symmetric, but apparently scipy adds the matrix and its transpose to make it symmetric before applying RCM. – Johan Falkenjack Aug 25 '16 at 11:10
• @BrianBorchers, I'm not aware that there is any requirement on positive definiteness to use RCM. Is there a reason why you say this? I thought it just needs to be symmetric in structure for the algorithm to "make sense". – Nick C. Aug 25 '16 at 14:55
• If you're going to use RCM with a sparse Cholesky factorization then the matrix needs to be positive definite. That's the most common use-case. – Brian Borchers Aug 25 '16 at 16:02
• Thanks, I see. My only experience with RCM is in the context of reordering the nodes of a finite element mesh, and then you are really just thinking about how those nodes are "connected" on a graph and the values (in the context of the final assembled FE system) don't matter, from what I recall. – Nick C. Aug 26 '16 at 3:17

perm = reverse_cuthill_mckee(arr_orig)

• For larger sparse matrices it is much faster to apply the permutation directly to the sparse matrix index array: mat.indices = perm.take(mat.indices) than convert the matrix to csc: mat = mat.tocsc() and apply the permutation again: mat.indices = perm.take(mat.indices). If you need the matrix in csr you can convert it back again. Using the advanced indexing is much slower than this - even when the matrix is converted twice between the formats. Tested with 25k x 25k matrices. – w1th0utnam3 Mar 27 '17 at 15:09