C# Sparse Factorization Example

using System;
using System.Collections.Generic;

using CenterSpace.NMath.Core;
using CenterSpace.NMath.Matrix;

namespace CenterSpace.NMath.Matrix.Examples.CSharp
{
  /// <summary>
  /// A .NET example in C# showing how to solve a symmetric sparse linear system using the 
  /// DoubleSparseSymFact class.
  /// </summary>
  class SparseFactorizationExample
  {

    static void Main(string[] args)
    {
      Dictionary<IntPair, double> values = new Dictionary<IntPair, double>();
      values.Add(new IntPair(0, 0), 7);
      values.Add(new IntPair(0, 2), 1);
      values.Add(new IntPair(0, 5), 2);
      values.Add(new IntPair(0, 6), 7);

      values.Add(new IntPair(1, 1), -4);
      values.Add(new IntPair(1, 2), 8);
      values.Add(new IntPair(1, 4), 2);

      values.Add(new IntPair(2, 2), 1);
      values.Add(new IntPair(2, 7), 5);

      values.Add(new IntPair(3, 3), 7);
      values.Add(new IntPair(3, 6), 9);

      values.Add(new IntPair(4, 4), 5);
      values.Add(new IntPair(4, 5), 1);
      values.Add(new IntPair(4, 6), 5);

      values.Add(new IntPair(5, 5), -1);
      values.Add(new IntPair(5, 7), 5);

      values.Add(new IntPair(6, 6), 11);

      values.Add(new IntPair(7, 7), 5);

      // Create a symmetric sparse matrix with 8 columns from the above data.
      DoubleSymCsrSparseMatrix S = new DoubleSymCsrSparseMatrix(values, 8);

      // Factor the matrix. The factorization can then be used to solve for various
      // right hand sides.
      DoubleSparseSymFact fact = new DoubleSparseSymFact(S);
      // Check that factorization succeeded.
      if (fact.ErrorStatus != SparseMatrixFact<double>.Error.NoError)
      {
        Console.WriteLine("Error {0} in factoring sparse matrix S.", fact.ErrorStatus);
        return;
      }

      // Right hand side a vector containing all ones.
      DoubleVector b = new DoubleVector(8, 1);
      DoubleVector x = fact.Solve(b);
      // Check that solve succeeded.
      if (fact.ErrorStatus != SparseMatrixFact<double>.Error.NoError)
      {
        Console.WriteLine("Error {0} in solving sparse system Sx = b.", fact.ErrorStatus);
        return;
      }
      Console.WriteLine("Solution for one right hand side = " + x);

      // Solve for several right hand sides.
      int nrhs = 3;
      DoubleMatrix B = new DoubleMatrix(8, nrhs, new RandGenBeta());
      DoubleMatrix X = fact.Solve(B);
      // Check that solve succeeded.
      if (fact.ErrorStatus != SparseMatrixFact<double>.Error.NoError)
      {
        Console.WriteLine("Error {0} in solving sparse system Sx = b.", fact.ErrorStatus);
        return;
      }
      Console.WriteLine("\nSolving for {0} right hand sides.", nrhs);
      for (int i = 0; i < nrhs; i++)
      {
        Console.WriteLine("Solution for rhs {0} \nis\n {1}\n", B.Col(i), X.Col(i));
      }

      Console.WriteLine();
      Console.WriteLine("Press Enter Key");
      Console.Read();
    }
  }
}