C# Symmetric Matrix Example

using System;

using CenterSpace.NMath.Core;


namespace CenterSpace.NMath.Examples.CSharp
{
  /// <summary>
  /// A .NET example in C# demonstrating the features of the symmetric matrix classes.
  /// </summary>
  class SymmetricMatrixExample
  {

    static void Main( string[] args )
    {
      int order = 6;

      // Set up a symmetric matrix S as the transpose product of a general 
      // matrix with itself (which is symmetric).
      var rng = new RandGenUniform( -2, 2 );
      rng.Reset( 0x124 );
      var A = new FloatMatrix( order, order, rng );

      var S = new FloatSymmetricMatrix( NMathFunctions.TransposeProduct( A, A ) );

      Console.WriteLine();

      Console.WriteLine( "S =" );
      Console.WriteLine( S.ToTabDelimited( "G3" ) );

      // S =
      // 3.3     0.529   -1.24   0.262   2.76    1.44
      // 0.529   2.39    -0.864  0.315   1.59    -0.543
      // -1.24   -0.864  5.55    4.65    -2.57   -0.0683
      // 0.262   0.315   4.65    8.96    -4.57   -1.7
      // 2.76    1.59    -2.57   -4.57   11.3    5.73
      // 1.44    -0.543  -0.0683 -1.7    5.73    4.4

      // Indexer accessor works just like it does for general matrices. 
      Console.WriteLine( "S[2,2] = {0}", S[2, 2] );
      Console.WriteLine( "S[3,0] = {0}", S[3, 0] );

      // You can set the values of elements in a symmetric matrix using the 
      // indexer. Note that setting the element in row i and column j implicitly
      // sets the element in column j and row i to the same value
      S[2, 1] = 100;
      Console.WriteLine( "S[2,1] = {0}", S[2, 1] );
      Console.WriteLine( "S[1,2] = {0}", S[1, 2] );

      // Scalar addition/subtractions/multiplication and matrix 
      // addition/subtraction are supported.
      float a = -.123F;
      FloatSymmetricMatrix C2 = a * S;
      FloatSymmetricMatrix C = 3 - S;
      FloatSymmetricMatrix D = C2 + S;
      Console.WriteLine();
      Console.WriteLine( "D =" );
      Console.WriteLine( D.ToTabDelimited( "G3" ) );

      // Matrix/vector products too.
      rng = new RandGenUniform( -1, 1 );
      rng.Reset( 0x124 );
      var x = new FloatVector( S.Cols, rng ); // vector of random deviates
      FloatVector y = MatrixFunctions.Product( S, x );
      Console.WriteLine( "Sx = {0}", y.ToString() );

      // You can transform the elements of a symmetric matrix object by using
      // the Transform() method.
      C2.Transform( NMathFunctions.FloatPowFunc, 2 ); // Square every element of C2
      Console.WriteLine();
      Console.WriteLine( "C2^2 =" );
      Console.WriteLine( C2.ToTabDelimited( "G9" ) );

      // You can also solve linear systems.
      FloatVector x2 = MatrixFunctions.Solve( S, y );

      // x and x2 should be about the same. Lets look at the l2 norm of 
      // their difference.
      FloatVector residual = x - x2;
      float residualL2Norm = (float) Math.Sqrt( NMathFunctions.Dot( residual, residual ) );
      Console.WriteLine( "||x - x2|| = {0}", residualL2Norm );

      // Compute condition number
      float rcond = MatrixFunctions.ConditionNumber( S );
      Console.WriteLine();
      Console.WriteLine( "Reciprocal condition number = {0}", rcond );

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