#include <vector>
#include "moab/Matrix3.hpp"
#include "MBCore.hpp"
#include "MBCartVect.hpp"

#include <armadillo>

int main (int argc, char** argv)
{

  //Create a matrix
  moab::Matrix3 mat;

  mat(0) = 0.093662468394669018;
  mat(1) = -0.10086327887079793;
  mat(2) = -0.13677885066200801;
  mat(3) = -0.10086327887079793;
  mat(4) = 0.13531200872876126;
  mat(5) = 0.16300669153639191;
  mat(6) = -0.13677885066200801;
  mat(7) = 0.16300669153639191;
  mat(8) = 555.5555555555552;


  // mat(0) = 0.0043698935017632001;
  // mat(1) = -0.014255577039788442;
  // mat(2) = -0.060767582383949303;
  // mat(3) = -0.014255577039788442;
  // mat(4) = 0.057226424277920775;
  // mat(5) = 0.22678772541701031;
  // mat(6) = -0.060767582383949303;
  // mat(7) = 0.22678772541701031;
  // mat(8) = 555.55555555555566;


  std::cout << mat << std::endl;

  //now do the Eigen Decomposition of this Matrix

  double lamda[3]; 
  MBCartVect vectors[3];
  moab::Matrix::EigenDecomp( mat, lamda, vectors);

  std::cout << "----------MOAB----------" << std::endl;

  std::cout << "Eigenvalues: " 
	    << lamda[0] << " "
	    << lamda[1] << " "
	    << lamda[2] << " "
	    << std::endl;

  std::cout << "Eigenvectors: " << std::endl
	    << vectors[0] << std::endl
	    << vectors[1] << std::endl
	    << vectors[2] << std::endl;


  std::cout << "----------ARMADILLO----------" << std::endl;

  //duplicate this matrix in Armadillo
  arma::Mat<double> M;

  M << mat(0) << mat(1) << mat(2) << arma::endr
    << mat(3) << mat(4) << mat(5) << arma::endr
    << mat(5) << mat(7) << mat(8) << arma::endr;


  std::cout << M << std::endl;

  arma::cx_vec eigenval;
  arma::cx_mat eigenvec;

  arma::eig_gen(eigenval, eigenvec, M);


  std::cout << eigenval << std::endl;

  std::cout << eigenvec << std::endl;


  return 0; 

}