//***********************************************************
//
//  NXV4 "Math-on-Paper" software
//
//  package : MTX class _____________________________________
//            
//            supports: move semantix
//                      operators + - * / L,R-val context
//                                += -= *= /=
//                      fabs()
//                      ~  -      unary operators L,R-val
//                      mtx<arithmetic> i.e.     int
//                                               float 
//                                               double
//                                               real
//                                           cpx<int> 
//                                           cpx<float> 
//                                           cpx<double> 
//                                           cpx<real> 
//                      type_casts between the types 
//                      cout << overload 
//
//                      TOTAL = 5084 instantiated op/f's
//                      _____________________________________
//
//
//  authors : Mihai-Tiberiu Dima (1)
//            Maria Dima         (2)
//
//             Meshcheryakov Laboratory of Information Techn.
//             Joint Institute for Nuclear Research
//             Dubna-Russia
//
//            on leave from:
//
//             Hyperion Univ. Bucharest-Romania
//
//
//  date    : MTX4/v1.0.alfa   / Mon Mar 07 12:45:54 CET 2022
//
//***********************************************************
#include <cmath>
#include <iostream>

#include "cpx.hh"
#include "vec.hh"
#include "mtx.hh"

using namespace std                                         ;
using namespace cpx4                                        ;
using namespace vec4                                        ;
using namespace mtx4                                        ;


int main()
{

// vec4::utf8x = true                                         ;
// mtx4::utf8x = true                                         ;

// this is a matrix diagonalisation
//
//   M_d = (U^-1) M (U)
//
// since U+ = U^-1, we can use
//
//   M_d = (U+) M (U)
//
// matrix eigen-vectors = columns of U
//
//  U = eigen(M)
//
//
// to use this example do the following:
//   rm test.cc
//   ln -s e3.cc test.cc
//   make test
//   make run

  using fltx = float                                        ;
  using dblx = double                                       ;
  using real = long double                                  ;

  double e   = exp(1)                                       ;
  double pi  = 3.14159265853                                ;

  cpx<int >  i(0,1)                                         ;

  mtx<dblx>  M( 0.807, 0.807, 0,
               -0.807, 0.807, 0,
                    0,     0, 1 )                           ;

  auto E = eigen(M)                                         ;

  vec<cpx<dblx>> v1(E.x11, E.x21, E.x31)                    ;
  vec<cpx<dblx>> v2(E.x12, E.x22, E.x32)                    ;
  vec<cpx<dblx>> v3(E.x13, E.x23, E.x33)                    ;

  cout                                             << endl  ;
  cout << noboolalpha                              << endl  ;

  cout                                             << endl  ;
  cout << " _____________________________________" << endl  ;
  cout                                             << endl  ;
  cout                                             << endl  ;
  cout                                             << endl  ;
  cout << "   " << M << "_diag = "
                << (~eigen(M)) * M * eigen(M) <<endl<<endl  ;
  cout                                             << endl  ;
  cout << " v1 = " << v1                           << endl  ;
  cout << " v2 = " << v2                           << endl  ;
  cout << " v3 = " << v3                           << endl  ;
  cout << " _____________________________________" << endl  ;
  cout                                             << endl  ;
  cout                                             << endl  ;
  cout                                             << endl  ;

}