// Author:  Douglas Wilhelm Harder
// Copyright (c) 2009 by Douglas Wilhelm Harder.  All rights reserved.

#include <iostream>
#include <cmath>
#include <string>
#include <cstdio>

/***************************************************
 * Print the matrix in a Matlab style format:
 *   Row by row
 *   Rows separated by a ';'
 *   Surrounded by '[' and ']'
 *
 * Run time:   O( M + na )
 * Memory:     O( 1 )
 ***************************************************/

template <int M, int N>
std::ostream &operator<< ( std::ostream &out, Matrix<M, N> const &A ) {
	char str[100];

	int min10 = 0;
	int max10 = 0;

	for ( int i = 0; i < Matrix<M,N>::minMN; ++i ) {
		if ( A.diagonal[i] != 0.0 ) {
			int log10x = int( std::floor( std::log10( std::fabs( A.diagonal[i] ) ) ) );

			min10 = std::min( min10, log10x );
			max10 = std::max( max10, log10x );
		}
	}

	for ( int j = 0; j < A.row_index[M]; ++j ) {
		int log10x = int( std::floor( std::log10( std::fabs( A.off_diagonal[j] ) ) ) );

		min10 = std::min( min10, log10x );
		max10 = std::max( max10, log10x );
	}

	++max10;

	cout << min10 << ", " << max10 << endl;

	if ( M != 0 && N != 0 ) {
		for ( int i = 0; i < M; ++i ) {
			int k = A.row_index[i];

			for ( int j = 0; j < N; ++j ) {
				if ( i == j ) {
					if ( A.diagonal[i] == 0.0 ) {
						std::sprintf( str, " % 15.*f", 10 - max10, 0.0 );
						out << str;
					} else {
						std::sprintf( str, " % 15.*f", 10 - max10, A.diagonal[i] );
						out << str;
					}
				} else if ( k < A.row_index[i + 1] && A.column_index[k] == j ) {
					std::sprintf( str, " % 15.*f", 10 - max10, A.off_diagonal[k] );
					out << str;
					++k;
				} else {
					std::sprintf( str, " % 15.*f", 10 - max10, 0.0 );
					out << str;
				}
			}

			out << std::endl;
		}
	}

	return out;
}

/***************************************************
 * Print information about the matrx.
 *
 * Run time:   O( M + na )
 * Memory:     O( 1 )
 ***************************************************/

template <int M, int N>
void Matrix<M, N>::details() const {
	std::cout << "Rows:                 " << M << std::endl;
	std::cout << "Columns:              " << N << std::endl;
	std::cout << "Minimum:              " << minMN << std::endl;
	std::cout << "Capacity:             " << capacity << std::endl;
	std::cout << "Size:                 " << row_index[M] << std::endl;

	std::cout << "Diagonal entries:     ";
	for ( int i = 0; i < minMN; ++i ) {
		std::cout << diagonal[i] << " ";
	}
	std::cout << std::endl;

	std::cout << "Row index:            ";
	for ( int i = 0; i <= M; ++i ) {
		std::cout << row_index[i] << " ";
	}
	std::cout << std::endl;

	std::cout << "Column index and off-diagonal entries:" << std::endl;

	for ( int j = 0; j < row_index[M]; ++j ) {
		std::cout << '\t' << j << '\t' << column_index[j] << '\t' << off_diagonal[j] << std::endl;
	}

	for ( int i = row_index[M]; i < capacity; ++i ) {
		std::cout << "\t-  \t-" << std::endl;
	}

}

/***************************************************
 * Print the matrix in such a format as to allow
 * it to be entered into Matlab.
 *
 * Run time:   O( M + na )
 * Memory:     O( M + na )
 ***************************************************/

template <int M, int N>
void Matrix<M, N>::matlab( std::string const &s ) const {
	std::cout << s << " = spalloc(" << M << "," << N << "," << capacity << ");" << std::endl;
	std::cout << s << " = " << s << " + spdiag([" << diagonal[0];

	for ( int i = 1; i < minMN; ++i ) {
		std::cout << " " << diagonal[i];
	}

	std::cout << "]);" << std::endl;

	for ( int i = 0; i < M; ++i ) {
		for ( int j = row_index[i]; j < row_index[i + 1]; ++j ) {
			std::cout << s << "(" << i + 1 << "," << column_index[j] + 1 << ") = " << off_diagonal[j] << ";" << std::endl;
		}
	}
}

/***************************************************
 * Print the matrix in such a format as to allow
 * it to be entered into Matlab.
 *
 * Run time:   O( M*N )
 * Memory:     O( M*N )
 ***************************************************/

template <int M, int N>
void Matrix<M, N>::matlab_dense( std::string const &s ) const {
	std::cout << s << " = [";

	if ( M != 0 && N != 0 ) {
		for ( int i = 0; i < M; ++i ) {
			if ( i != 0 ) {
				std::cout << "; ";
			}

			int k = row_index[i];

			for ( int j = 0; j < N; ++j ) {
				if ( j != 0 ) {
					std::cout << " ";
				}

				if ( i == j ) {
					std::cout << diagonal[i];
				} else if ( k < row_index[i + 1] && column_index[k] == j ) {
					std::cout << off_diagonal[k];
					++k;
				} else {
					std::cout << 0;
				}
			}
		}
	}

	std::cout << ']' << endl;
}