Documentation updates and minor changes.

Author: wo80

CSparse.Tests/Ordering/TestDulmageMendelsohn.cs

@@ -1,6 +1,8 @@
 namespace CSparse.Tests.Ordering
 {
+    using CSparse.Double;
     using CSparse.Ordering;
+    using CSparse.Storage;
     using NUnit.Framework;
     using System;
 
@@ -33,5 +35,46 @@ public void TestGenerate2()
 
             Assert.That(dm.StructuralRank == n, Is.True);
         }
+
+        [Test]
+        public void TestGenerate3()
+        {
+            var A = SparseMatrix.OfRowMajor(8, 8,
+            [
+                11, 12,  0,  0,  0,  0,  0,  0,
+                 0, 22, 23,  0, 25, 26,  0,  0,
+                 0,  0, 33, 34,  0,  0, 37,  0,
+                 0,  0, 43, 44,  0,  0,  0, 48,
+                51,  0,  0,  0, 55, 56,  0,  0,
+                 0,  0,  0,  0,  0, 66, 67,  0,
+                 0,  0,  0,  0,  0, 76, 77,  0,
+                 0,  0,  0, 84,  0,  0, 87, 88
+            ]);
+
+            var S = SymbolicColumnStorage.Create(A);
+
+            var dm = DulmageMendelsohn.Generate(S, 1);
+
+            Assert.That(dm.StructuralRank, Is.EqualTo(8));
+            Assert.That(dm.Blocks, Is.EqualTo(3));
+
+            int[] expected = [0, 1, 4, 2, 3, 7, 5, 6];
+
+            Assert.That(dm.RowPermutation, Is.EqualTo(expected).AsCollection);
+            Assert.That(dm.ColumnPermutation, Is.EqualTo(expected).AsCollection);
+
+            expected = [0, 3, 6, 8];
+
+            Assert.That(dm.BlockRowPointers, Is.EqualTo(expected).AsCollection);
+            Assert.That(dm.BlockColumnPointers, Is.EqualTo(expected).AsCollection);
+
+            expected = [0, 0, 8, 8, 8];
+
+            Assert.That(dm.CoarseRowDecomposition, Is.EqualTo(expected).AsCollection);
+
+            expected = [0, 0, 0, 8, 8];
+
+            Assert.That(dm.CoarseColumnDecomposition, Is.EqualTo(expected).AsCollection);
+        }
     }
 }

CSparse/Ordering/AMD.cs

@@ -25,6 +25,22 @@ public static class AMD
         /// </remarks>
         public static int[] Generate<T>(CompressedColumnStorage<T> A, ColumnOrdering order)
             where T : struct, IEquatable<T>, IFormattable
+        {
+            return Generate(SymbolicColumnStorage.Create(A), order);
+        }
+
+        /// <summary>
+        /// Generate minimum degree ordering of A+A' (if A is symmetric) or A'A.
+        /// </summary>
+        /// <param name="A">Column-compressed matrix</param>
+        /// <param name="order">Column ordering method</param>
+        /// <returns>amd(A+A') if A is symmetric, or amd(A'A) otherwise, null on 
+        /// error or for natural ordering</returns>
+        /// <remarks>
+        /// See Chapter 7.1 (Fill-reducing orderings: Minimum degree ordering) in 
+        /// "Direct Methods for Sparse Linear Systems" by Tim Davis.
+        /// </remarks>
+        public static int[] Generate(SymbolicColumnStorage A, ColumnOrdering order)
         {
             int[] Cp, Ci, P, W, nv, next, head, elen, degree, w, hhead;
 
@@ -42,7 +58,7 @@ public static int[] Generate<T>(CompressedColumnStorage<T> A, ColumnOrdering ord
                 return Permutation.Create(n);
             }
 
-            var C = ConstructMatrix(SymbolicColumnStorage.Create(A), order);
+            var C = ConstructMatrix(A, order);
 
             Cp = C.ColumnPointers;
             cnz = Cp[n];
@@ -139,7 +155,7 @@ public static int[] Generate<T>(CompressedColumnStorage<T> A, ColumnOrdering ord
                             Ci[p] = -(j + 2); // first entry is now CS_FLIP(j)
                         }
                     }
-                    for (q = 0, p = 0; p < cnz; ) // scan all of memory
+                    for (q = 0, p = 0; p < cnz;) // scan all of memory
                     {
                         if ((j = FLIP(Ci[p++])) >= 0) // found object j
                         {
@@ -303,7 +319,7 @@ public static int[] Generate<T>(CompressedColumnStorage<T> A, ColumnOrdering ord
                         eln = elen[i];
                         for (p = Cp[i] + 1; p <= Cp[i] + ln - 1; p++) w[Ci[p]] = mark;
                         jlast = i;
-                        for (j = next[i]; j != -1; ) // compare i with all j
+                        for (j = next[i]; j != -1;) // compare i with all j
                         {
                             ok = (W[j] == ln) && (elen[j] == eln);
                             for (p = Cp[j] + 1; ok && p <= Cp[j] + ln - 1; p++)
@@ -411,13 +427,13 @@ private static bool KeepOffDiag(int i, int j)
         private static SymbolicColumnStorage ConstructMatrix(SymbolicColumnStorage A, ColumnOrdering order)
         {
             SymbolicColumnStorage result = null;
-            
+
             // Compute A'
             var AT = A.Transpose();
 
             int m = A.RowCount;
             int n = A.ColumnCount;
-            
+
             if (order == ColumnOrdering.MinimumDegreeAtPlusA)
             {
                 if (n != m)
@@ -444,7 +460,7 @@ private static SymbolicColumnStorage ConstructMatrix(SymbolicColumnStorage A, Co
                 {
                     // Column j of AT starts here.
                     p = colptr[j];
-                    
+
                     // New column j starts here.
                     colptr[j] = p2;
 

CSparse/Ordering/DulmageMendelsohn.cs

@@ -22,7 +22,7 @@ public sealed class DulmageMendelsohn
         private int nb;    // number of blocks in fine dmperm decomposition
 
         /// <summary>
-        /// Create a new Decomposition instance. 
+        /// Create a new decomposition instance. 
         /// </summary>
         private DulmageMendelsohn(int m, int n)
         {
@@ -97,11 +97,10 @@ public int Singletons
         public int[] CoarseColumnDecomposition => cc;
 
         /// <summary>
-        /// Compute coarse and then fine Dulmage-Mendelsohn decomposition. seed
-        /// optionally selects a randomized algorithm.
+        /// Compute coarse and fine Dulmage-Mendelsohn decomposition.
         /// </summary>
         /// <param name="matrix">column-compressed matrix</param>
-        /// <param name="seed">0: natural, -1: reverse, random order otherwise</param>
+        /// <param name="seed"> The seed optionally selects a randomized algorithm (0 = default (natural), -1 = reverse, random order otherwise).</param>
         /// <returns>Dulmage-Mendelsohn analysis</returns>
         public static DulmageMendelsohn Generate<T>(CompressedColumnStorage<T> matrix, int seed = 0)
              where T : struct, IEquatable<T>, IFormattable
@@ -110,11 +109,10 @@ public static DulmageMendelsohn Generate<T>(CompressedColumnStorage<T> matrix, i
         }
 
         /// <summary>
-        /// Compute coarse and then fine Dulmage-Mendelsohn decomposition. seed
-        /// optionally selects a randomized algorithm.
+        /// Compute coarse and fine Dulmage-Mendelsohn decomposition.
         /// </summary>
         /// <param name="A">The matrix represented by <see cref="SymbolicColumnStorage"/>.</param>
-        /// <param name="seed">0: natural, -1: reverse, random order otherwise</param>
+        /// <param name="seed"> The seed optionally selects a randomized algorithm (0 = default (natural), -1 = reverse, random order otherwise).</param>
         /// <returns>Dulmage-Mendelsohn analysis</returns>
         public static DulmageMendelsohn Generate(SymbolicColumnStorage A, int seed = 0)
         {

CSparse/Storage/SymbolicColumnStorage.cs

@@ -94,12 +94,12 @@ public SymbolicColumnStorage(int rowCount, int columnCount, int[] columnPointers
 
             if (columnPointers.Length < columnCount + 1)
             {
-                throw new ArgumentOutOfRangeException("Column pointers array size don't match given column count argument.");
+                throw new ArgumentOutOfRangeException(nameof(columnPointers), "Column pointers array size doesn't match given column count argument.");
             }
 
             if (rowIndices.Length < columnPointers[columnCount])
             {
-                throw new ArgumentOutOfRangeException("Row indices array size don't match non-zeros count.");
+                throw new ArgumentOutOfRangeException(nameof(rowIndices), "Row indices array size doesn't match non-zeros count.");
             }
 
             this.rowCount = rowCount;
@@ -187,13 +187,16 @@ public void Sort()
             }
         }
 
+        /// <summary>
+        /// Returns a copy of the <see cref="SymbolicColumnStorage"/>.
+        /// </summary>
+        /// <returns></returns>
         public SymbolicColumnStorage Clone()
         {
-            int m = RowCount;
             int n = ColumnCount;
             int nnz = NonZerosCount;
 
-            var result = new SymbolicColumnStorage(m, n, nnz, true);
+            var result = new SymbolicColumnStorage(RowCount, n, nnz, true);
 
             Buffer.BlockCopy(ColumnPointers, 0, result.ColumnPointers, 0, (n + 1) * Constants.SizeOfInt);
             Buffer.BlockCopy(RowIndices, 0, result.RowIndices, 0, nnz * Constants.SizeOfInt);