// $Id$

// Simple test application to demonstrate the free unary function count_structural_elements.
// The function is not expected to be particularly useful, but is used to experiment with
// efficient methods for iterating over containers.


#include <iostream>
#include <complex>

// Boost numeric uBlas library
#include <boost/numeric/ublas/traits.hpp>
#include <boost/numeric/ublas/vector.hpp>
#include <boost/numeric/ublas/vector_proxy.hpp>
#include <boost/numeric/ublas/matrix.hpp>
#include <boost/numeric/ublas/matrix_proxy.hpp>
#include <boost/numeric/ublas/matrix_sparse.hpp>
#include <boost/numeric/ublas/triangular.hpp>
#include <boost/numeric/ublas/symmetric.hpp>
#include <boost/numeric/ublas/io.hpp>
#include <boost/numeric/ublas/storage.hpp>
#include <boost/numeric/ublas/functional.hpp>

// Example free function header
#include "count_structural_elements.hpp"

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

namespace ublas = boost::numeric::ublas;

int main( int argc, char **argv )
{
	std::cout << "Testing count_structural_elements..." << std::endl << std::endl;
	
	// ------ Create some vectors and matrices ------

	// Define size of vectors and matrices
	std::size_t n(5),m(5);

	int c = 0;

	// Create dense vector
  ublas::vector<double> V(n*m);
	c = 0;
	for (ublas::vector<double>::size_type i = 0; i < V.size(); i++)
      V(i) = static_cast<double>((c++) % 4);
  std::cout << "V = " << V << std::endl;

	// Create sparse vector and initialise elements using copy constructor
  ublas::compressed_vector<double> CV(V);

	// Create dense matrix
  ublas::matrix<double> A(n, m);
	// Fill matrix A
	c = 0;
	for (ublas::matrix<double>::size_type i = 0; i < A.size1(); i++)
    for (ublas::matrix<double>::size_type j = 0; j < A.size2(); j++)
      A(i,j) = static_cast<double>((c++) % 4);
  std::cout << "A = " << A << std::endl;

	// Create sparse matrix and initialise elements using copy constructor
  ublas::compressed_matrix<double> CM(A); 
  std::cout << "CM = " << CM << std::endl;

	// Create lower triangular matrix (uses packed storage)
	ublas::triangular_matrix<double, ublas::lower> LT(n, m);
  for (unsigned i = 0; i < LT.size1(); ++ i)
    for (unsigned j = 0; j <= i; ++ j)
      LT(i, j) = 3 * i + j;
  std::cout << "LT = " << LT << std::endl;


	// Create symmetrical matrix (uses packed storage)
	ublas::symmetric_matrix<double, ublas::lower> SYM(n, m);
  for (unsigned i = 0; i < SYM.size1 (); ++ i)
    for (unsigned j = 0; j <= i; ++ j)
      SYM(i, j) = 3 * i + j;
  std::cout << "SYM = " << SYM << std::endl;

	// ----- Count structural elements -----

	std::cout << "count_structural_elements(V) = " << count_structural_elements(V) << std::endl;
	std::cout << "count_structural_elements(V+V) = " << count_structural_elements(V+V) << std::endl;
	std::cout << "count_structural_elements(CV) = " << count_structural_elements(CV) << std::endl;
	std::cout << "count_structural_elements(CV+CV) = " << count_structural_elements(CV+CV) << std::endl;
	std::cout << "count_structural_elements(A) = " << count_structural_elements(A) << std::endl;			// Matrix container
	std::cout << "count_structural_elements(A+A) = " << count_structural_elements(A+A) << std::endl;	// matrix expression
	std::cout << "count_structural_elements(CM) = " << count_structural_elements(CM) << std::endl;
	std::cout << "count_structural_elements(CM+CM) = " << count_structural_elements(CM+CM) << std::endl;
	std::cout << "count_structural_elements(LT) = " << count_structural_elements(LT) << std::endl;
	std::cout << "count_structural_elements(SYM) = " << count_structural_elements(SYM) << std::endl;

	return 0;
}