// testBoostJan2014.cpp : Defines the entry point for the console application.
//

#include "targetver.h"

#include <stdio.h>
#include <tchar.h>

#include <fstream>
#include <iostream>
#include <iomanip>
#include <cmath>

#include <stdio.h>
#include <string.h>

#include <boost/format.hpp>

#define BIND_FORTRAN_LOWERCASE_UNDERSCORE
#define BOOST_UBLAS_ENABLE_PROXY_SHORTCUTS

#include <boost/numeric/bindings/lapack/workspace.hpp>
#include <boost/numeric/bindings/ublas/matrix.hpp> 
#include <boost/numeric/bindings/ublas/vector.hpp> 
#include <boost/numeric/ublas/matrix.hpp> 
#include <boost/numeric/ublas/matrix_proxy.hpp>
#include <boost/numeric/ublas/vector.hpp> 
#include <boost/numeric/ublas/io.hpp>
//#include <boost/numeric/bindings/lapack/geev.hpp>
#include <boost/numeric/bindings/lapack/driver/heev.hpp>
#include <boost/numeric/bindings/lapack/driver/gesdd.hpp>
// New options to test:
#include <boost/numeric/bindings/lapack/driver/gesvd.hpp>
// #include <boost/numeric/bindings/lapack/driver/heevr.hpp> // We would like this, but it isn't in the version of the numeric bindings in Boost 1.36.


using namespace boost::numeric::bindings::lapack; 
using namespace boost::numeric::ublas; 

typedef matrix<std::complex<double>, column_major> cmat;
typedef vector< std::complex<double> > cvec; // Must be > > within a nested template so as not to be confused with >> operator.
typedef vector<double> vec;
typedef std::complex<double> cdouble;


int _tmain(int argc, _TCHAR* argv[])
{
	using namespace std;

	cout << "TESTING..." << endl;
	
	cmat eigVec (4,4);
	vec eigVal (4); 

	int idx,jdx;

	for (idx=0;idx<4;++idx)
		for (jdx=0;jdx<4;++jdx)
			eigVec(idx,jdx) = 10-idx-jdx; // Make this a Hermitian matrix.

	cout << "eigVec input = " << eigVec << endl;

	// Try using heev or heev because this is a Hermitian matrix:
	heev('V',eigVec,eigVal,boost::numeric::bindings::lapack::optimal_workspace());

	cout << "eigVec output = " << eigVec << endl;

	// Prompt before return.
	cout << "Type a letter and press return to finish..." << endl;
	string tmp;
	cin >> tmp;

	return 0;
}