C# Hermitian Matrix Example

using System;

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

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

    static void Main(string[] args)
    {
      int order = 5;
      string numberFormatString = "F4"; // Format number strings as fixed, 4 digits.

      // Set up a Hermitian matrix S as the conjugate transpose product of a general 
      // matrix with itself (which is Hermitian).
      RandGenUniform rng = new RandGenUniform( -1, 1 );
      rng.Reset( 0x124 );
      DoubleComplexMatrix A = new DoubleComplexMatrix( order, order, rng );

      DoubleHermitianMatrix S = 
        new DoubleHermitianMatrix( NMathFunctions.ConjTransposeProduct(A,A) );

      Console.WriteLine();
      Console.WriteLine( "S = {0}", S.ToString(numberFormatString) );

      // S = 5x5 [ (3.1219,0.0000)  (0.1293,0.7632)   (-0.5926,-0.5191) (1.0169,-0.4854) (-0.6211,-0.7439)  
      //           (0.1293,-0.7632) (1.0186,0.0000)   (-0.6158,0.5822)  (0.3471,-1.1798) (-0.3765,0.3526)  
      //           (-0.5926,0.5191) (-0.6158,-0.5822) (2.6691,0.0000)   (-0.7861,2.2311) (0.0942,0.1853)  
      //           (1.0169,0.4854)  (0.3471,1.1798)   (-0.7861,-2.2311) (4.1041,0.0000)  (0.5963,0.6960)  
      //           (-0.6211,0.7439) (-0.3765,-0.3526) (0.0942,-0.1853)  (0.5963,-0.6960) (3.9763,0.0000) ]

      // Indexer accessor works just like it does for general matrices. 
      Console.WriteLine();
      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 Hermitian matrix using the 
      // indexer. Note that setting the element in row i and column j to
      // a value implicitly sets the element in column j and row i to the 
      // complex conjugate of that value.
      S[2,1] = new DoubleComplex( 100, -99 );
      Console.WriteLine( "S[2,1] = {0}", S[2,1] ); // (100, -99)
      Console.WriteLine( "S[1,2] = {0}", S[1,2] ); // (100, 99)
      
      // Scalar multiplication and matrix addition/subtraction are supported.
      Double a = -.123;
      DoubleHermitianMatrix C2 = a*S;
      DoubleHermitianMatrix D = C2 + S;
      Console.WriteLine();
      Console.WriteLine( "D = {0}", D.ToString(numberFormatString) );

      // Matrix/vector products too.
      DoubleComplexVector x = new DoubleComplexVector( S.Cols, rng ); // vector of random deviates
      DoubleComplexVector y = MatrixFunctions.Product( S, x );
      Console.WriteLine();
      Console.WriteLine( "Sx = {0}", y.ToString(numberFormatString) );

      // You can also solve linear systems.
      DoubleComplexVector x2 = MatrixFunctions.Solve( S, y );
      
      // x and x2 should be about the same. Let's look at the l2 norm of 
      // their difference.
      DoubleComplexVector residual = x - x2;
      double residualL2Norm = Math.Sqrt( NMathFunctions.ConjDot(residual, residual).Real );
      Console.WriteLine();
      Console.WriteLine( "||x - x2|| = {0}", residualL2Norm );

      // You can transform the elements of a Hermitian matrix object by using
      // the Transform() method.
      C2.DataVector.Transform( NMathFunctions.DoubleComplexCoshFunction );
      Console.WriteLine();
      Console.WriteLine( "cosh(C2) = {0}", C2.ToString(numberFormatString) );

      // For a matrix to satisfy the strict definition of a Hermitian matrix,
      // its diagonal elements must be real. The Hermitian matrix classes provide
      // a MakeDigaonalReal() method to ensure that your matrix satisfies the
      // the strict definition of Hermitian.
      C2.MakeDiagonalReal();
      Console.WriteLine();
      Console.WriteLine( "Diagonal element is real: {0}", C2[3,3].Imag == 0.0 ); // True

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

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