Starting to add documentation for all methods and public variables. Trying to make sure advanced concepts are explained in an easy to understand method so let me know if you like the for beginner sections

Author: ooples

src/DecompositionMethods/MatrixDecomposition/BidiagonalDecomposition.cs

@@ -1,26 +1,61 @@
-using AiDotNet.Enums.AlgorithmTypes;
-
 namespace AiDotNet.DecompositionMethods.MatrixDecomposition;
 
 public class BidiagonalDecomposition<T> : IMatrixDecomposition<T>
 {
-    private readonly INumericOperations<T> NumOps;
+    private readonly INumericOperations<T> _numOps;
 
+    /// <summary>
+    /// Gets the original matrix that was decomposed.
+    /// </summary>
     public Matrix<T> A { get; }
+
+    /// <summary>
+    /// Gets the left orthogonal matrix in the decomposition.
+    /// In simpler terms, this matrix helps transform the original matrix's columns.
+    /// </summary>
     public Matrix<T> U { get; private set; }
+
+    /// <summary>
+    /// Gets the bidiagonal matrix in the decomposition.
+    /// A bidiagonal matrix is a special matrix where non-zero values appear only on the main diagonal
+    /// and the diagonal immediately above it (called the superdiagonal).
+    /// </summary>
     public Matrix<T> B { get; private set; }
+
+    /// <summary>
+    /// Gets the right orthogonal matrix in the decomposition.
+    /// In simpler terms, this matrix helps transform the original matrix's rows.
+    /// </summary>
     public Matrix<T> V { get; private set; }
 
+    /// <summary>
+    /// Creates a new bidiagonal decomposition of the specified matrix.
+    /// </summary>
+    /// <param name="matrix">The matrix to decompose.</param>
+    /// <param name="algorithm">The algorithm to use for decomposition (default is Householder).</param>
+    /// <remarks>
+    /// Bidiagonal decomposition breaks down a matrix A into three simpler matrices: U, B, and V,
+    /// where A = U*B*V^T. This makes many matrix operations easier to perform.
+    /// </remarks>
     public BidiagonalDecomposition(Matrix<T> matrix, BidiagonalAlgorithmType algorithm = BidiagonalAlgorithmType.Householder)
     {
         A = matrix;
-        NumOps = MathHelper.GetNumericOperations<T>();
+        _numOps = MathHelper.GetNumericOperations<T>();
         U = new Matrix<T>(matrix.Rows, matrix.Rows);
         B = new Matrix<T>(matrix.Columns, matrix.Columns);
         V = new Matrix<T>(matrix.Columns, matrix.Columns);
         Decompose(algorithm);
     }
 
+    /// <summary>
+    /// Performs the bidiagonal decomposition using the specified algorithm.
+    /// </summary>
+    /// <param name="algorithm">The algorithm to use for decomposition.</param>
+    /// <exception cref="ArgumentException">Thrown when an unsupported algorithm is specified.</exception>
+    /// <remarks>
+    /// Different algorithms have different performance characteristics and numerical stability.
+    /// Householder is generally the most stable and commonly used method.
+    /// </remarks>
     public void Decompose(BidiagonalAlgorithmType algorithm = BidiagonalAlgorithmType.Householder)
     {
         switch (algorithm)
@@ -39,6 +74,13 @@ public void Decompose(BidiagonalAlgorithmType algorithm = BidiagonalAlgorithmTyp
         }
     }
 
+    /// <summary>
+    /// Performs bidiagonal decomposition using Householder reflections.
+    /// </summary>
+    /// <remarks>
+    /// Householder reflections are a way to transform vectors by reflecting them across a plane.
+    /// This method is numerically stable and commonly used for matrix decompositions.
+    /// </remarks>
     private void DecomposeHouseholder()
     {
         int m = A.Rows;
@@ -54,7 +96,7 @@ private void DecomposeHouseholder()
             Vector<T> v = HouseholderVector(x);
 
             // Apply Householder reflection to B
-            Matrix<T> P = Matrix<T>.CreateIdentity(m - k).Subtract(v.OuterProduct(v).Multiply(NumOps.FromDouble(2)));
+            Matrix<T> P = Matrix<T>.CreateIdentity(m - k).Subtract(v.OuterProduct(v).Multiply(_numOps.FromDouble(2)));
             Matrix<T> subB = B.GetSubMatrix(k, k, m - k, n - k);
             B.SetSubMatrix(k, k, P.Multiply(subB));
 
@@ -69,7 +111,7 @@ private void DecomposeHouseholder()
                 v = HouseholderVector(x);
 
                 // Apply Householder reflection to B
-                P = Matrix<T>.CreateIdentity(n - k - 1).Subtract(v.OuterProduct(v).Multiply(NumOps.FromDouble(2)));
+                P = Matrix<T>.CreateIdentity(n - k - 1).Subtract(v.OuterProduct(v).Multiply(_numOps.FromDouble(2)));
                 subB = B.GetSubMatrix(k, k + 1, m - k, n - k - 1);
                 B.SetSubMatrix(k, k + 1, subB.Multiply(P));
 
@@ -80,6 +122,13 @@ private void DecomposeHouseholder()
         }
     }
 
+    /// <summary>
+    /// Performs bidiagonal decomposition using Givens rotations.
+    /// </summary>
+    /// <remarks>
+    /// Givens rotations are a way to zero out specific elements in a matrix by rotating
+    /// two rows or columns. This method is useful for creating bidiagonal matrices.
+    /// </remarks>
     private void DecomposeGivens()
     {
         int m = A.Rows;
@@ -105,6 +154,14 @@ private void DecomposeGivens()
         }
     }
 
+    /// <summary>
+    /// Performs bidiagonal decomposition using the Lanczos algorithm.
+    /// </summary>
+    /// <remarks>
+    /// The Lanczos algorithm is an iterative method that can efficiently compute
+    /// partial decompositions of large matrices. It's particularly useful for
+    /// sparse matrices (matrices with mostly zero values).
+    /// </remarks>
     private void DecomposeLanczos()
     {
         int m = A.Rows;
@@ -120,7 +177,7 @@ private void DecomposeLanczos()
         Random rand = new();
         for (int i = 0; i < n; i++)
         {
-            v[i] = NumOps.FromDouble(rand.NextDouble());
+            v[i] = _numOps.FromDouble(rand.NextDouble());
         }
         v = v.Divide(v.Norm());
 
@@ -142,6 +199,16 @@ private void DecomposeLanczos()
         }
     }
 
+    /// <summary>
+    /// Solves the linear system Ax = b using the bidiagonal decomposition.
+    /// </summary>
+    /// <param name="b">The right-hand side vector of the equation Ax = b.</param>
+    /// <returns>The solution vector x.</returns>
+    /// <exception cref="ArgumentException">Thrown when the length of vector b doesn't match the number of rows in matrix A.</exception>
+    /// <remarks>
+    /// This method uses the decomposition to efficiently solve the system without
+    /// directly inverting the matrix, which is more numerically stable.
+    /// </remarks>
     public Vector<T> Solve(Vector<T> b)
     {
         if (b.Length != A.Rows)
@@ -153,6 +220,14 @@ public Vector<T> Solve(Vector<T> b)
         return V.Multiply(z);
     }
 
+    /// <summary>
+    /// Computes the inverse of the original matrix A using the bidiagonal decomposition.
+    /// </summary>
+    /// <returns>The inverse matrix of A.</returns>
+    /// <remarks>
+    /// Matrix inversion is computationally expensive and can be numerically unstable.
+    /// When possible, use the Solve method instead of explicitly computing the inverse.
+    /// </remarks>
     public Matrix<T> Invert()
     {
         int n = A.Columns;
@@ -161,46 +236,69 @@ public Matrix<T> Invert()
         for (int i = 0; i < n; i++)
         {
             Vector<T> ei = new Vector<T>(n);
-            ei[i] = NumOps.One;
+            ei[i] = _numOps.One;
             inverse.SetColumn(i, Solve(ei));
         }
 
         return inverse;
     }
 
+    /// <summary>
+    /// Computes a Householder reflection vector for the given input vector.
+    /// </summary>
+    /// <param name="x">The input vector.</param>
+    /// <returns>A Householder vector that can be used to create a reflection matrix.</returns>
+    /// <remarks>
+    /// A Householder reflection is a transformation that reflects a vector across a plane.
+    /// It's used to zero out specific elements in a matrix during decomposition.
+    /// </remarks>
     private Vector<T> HouseholderVector(Vector<T> x)
     {
         T norm = x.Norm();
         Vector<T> v = x.Copy();
-        v[0] = NumOps.Add(v[0], NumOps.Multiply(NumOps.SignOrZero(x[0]), norm));
+        v[0] = _numOps.Add(v[0], _numOps.Multiply(_numOps.SignOrZero(x[0]), norm));
 
         return v.Divide(v.Norm());
     }
 
+        /// <summary>
+    /// Applies a Givens rotation to the specified matrices.
+    /// </summary>
+    /// <param name="M">The matrix to which the rotation is applied.</param>
+    /// <param name="Q">The orthogonal matrix that accumulates the rotations.</param>
+    /// <param name="i">The first row/column index for the rotation.</param>
+    /// <param name="k">The second row/column index for the rotation.</param>
+    /// <param name="j">The first column/row index for the rotation.</param>
+    /// <param name="l">The second column/row index for the rotation.</param>
+    /// <param name="isLeft">If true, applies a left rotation (row operation); otherwise, applies a right rotation (column operation).</param>
+    /// <remarks>
+    /// A Givens rotation is a simple rotation in a plane spanned by two coordinate axes.
+    /// It's used to selectively zero out specific elements in a matrix during decomposition.
+    /// </remarks>
     private void GivensRotation(Matrix<T> M, Matrix<T> Q, int i, int k, int j, int l, bool isLeft)
     {
         T a = M[i, j];
         T b = M[k, l];
-        T r = NumOps.Sqrt(NumOps.Add(NumOps.Multiply(a, a), NumOps.Multiply(b, b)));
-        T c = NumOps.Divide(a, r);
-        T s = NumOps.Divide(b, r);
+        T r = _numOps.Sqrt(_numOps.Add(_numOps.Multiply(a, a), _numOps.Multiply(b, b)));
+        T c = _numOps.Divide(a, r);
+        T s = _numOps.Divide(b, r);
 
         if (isLeft)
         {
             for (int j2 = 0; j2 < M.Columns; j2++)
             {
                 T temp1 = M[i, j2];
                 T temp2 = M[k, j2];
-                M[i, j2] = NumOps.Add(NumOps.Multiply(c, temp1), NumOps.Multiply(s, temp2));
-                M[k, j2] = NumOps.Subtract(NumOps.Multiply(NumOps.Negate(s), temp1), NumOps.Multiply(c, temp2));
+                M[i, j2] = _numOps.Add(_numOps.Multiply(c, temp1), _numOps.Multiply(s, temp2));
+                M[k, j2] = _numOps.Subtract(_numOps.Multiply(_numOps.Negate(s), temp1), _numOps.Multiply(c, temp2));
             }
 
             for (int i2 = 0; i2 < Q.Rows; i2++)
             {
                 T temp1 = Q[i2, i];
                 T temp2 = Q[i2, k];
-                Q[i2, i] = NumOps.Add(NumOps.Multiply(c, temp1), NumOps.Multiply(s, temp2));
-                Q[i2, k] = NumOps.Subtract(NumOps.Multiply(NumOps.Negate(s), temp1), NumOps.Multiply(c, temp2));
+                Q[i2, i] = _numOps.Add(_numOps.Multiply(c, temp1), _numOps.Multiply(s, temp2));
+                Q[i2, k] = _numOps.Subtract(_numOps.Multiply(_numOps.Negate(s), temp1), _numOps.Multiply(c, temp2));
             }
         }
         else
@@ -209,20 +307,30 @@ private void GivensRotation(Matrix<T> M, Matrix<T> Q, int i, int k, int j, int l
             {
                 T temp1 = M[i2, j];
                 T temp2 = M[i2, l];
-                M[i2, j] = NumOps.Add(NumOps.Multiply(c, temp1), NumOps.Multiply(s, temp2));
-                M[i2, l] = NumOps.Subtract(NumOps.Multiply(NumOps.Negate(s), temp1), NumOps.Multiply(c, temp2));
+                M[i2, j] = _numOps.Add(_numOps.Multiply(c, temp1), _numOps.Multiply(s, temp2));
+                M[i2, l] = _numOps.Subtract(_numOps.Multiply(_numOps.Negate(s), temp1), _numOps.Multiply(c, temp2));
             }
 
             for (int j2 = 0; j2 < Q.Columns; j2++)
             {
                 T temp1 = Q[j, j2];
                 T temp2 = Q[l, j2];
-                Q[j, j2] = NumOps.Add(NumOps.Multiply(c, temp1), NumOps.Multiply(s, temp2));
-                Q[l, j2] = NumOps.Subtract(NumOps.Multiply(NumOps.Negate(s), temp1), NumOps.Multiply(c, temp2));
+                Q[j, j2] = _numOps.Add(_numOps.Multiply(c, temp1), _numOps.Multiply(s, temp2));
+                Q[l, j2] = _numOps.Subtract(_numOps.Multiply(_numOps.Negate(s), temp1), _numOps.Multiply(c, temp2));
             }
         }
     }
 
+    /// <summary>
+    /// Solves a bidiagonal system of linear equations.
+    /// </summary>
+    /// <param name="y">The right-hand side vector of the equation Bx = y.</param>
+    /// <returns>The solution vector x.</returns>
+    /// <remarks>
+    /// This method efficiently solves a system where the coefficient matrix is bidiagonal
+    /// (has non-zero elements only on the main diagonal and the superdiagonal).
+    /// It uses back-substitution, which is much faster than general matrix inversion.
+    /// </remarks>
     private Vector<T> SolveBidiagonal(Vector<T> y)
     {
         int n = B.Columns;
@@ -232,13 +340,26 @@ private Vector<T> SolveBidiagonal(Vector<T> y)
         {
             T sum = y[i];
             if (i < n - 1)
-                sum = NumOps.Subtract(sum, NumOps.Multiply(B[i, i + 1], x[i + 1]));
-            x[i] = NumOps.Divide(sum, B[i, i]);
+                sum = _numOps.Subtract(sum, _numOps.Multiply(B[i, i + 1], x[i + 1]));
+            x[i] = _numOps.Divide(sum, B[i, i]);
         }
 
         return x;
     }
 
+    /// <summary>
+    /// Returns the three factor matrices of the bidiagonal decomposition.
+    /// </summary>
+    /// <returns>
+    /// A tuple containing:
+    /// - U: The left orthogonal matrix
+    /// - B: The bidiagonal matrix
+    /// - V: The right orthogonal matrix
+    /// </returns>
+    /// <remarks>
+    /// The original matrix A can be reconstructed as A = U * B * V^T.
+    /// This decomposition is useful for solving linear systems and computing singular values.
+    /// </remarks>
     public (Matrix<T> U, Matrix<T> B, Matrix<T> V) GetFactors()
     {
         return (U, B, V);