/*------------------------------------------------------*/
/* matrix.cpp											*/
/* Author: Aaron Pulley									*/
/* (801) 687-0496 - awpulley@cavtel.net					*/
/* Created: 2 Mar 2010									*/
/* With matrix.h, supplies various functions.			*/
/* History of modifications:							*/
/* -	23 Feb 2011 Modified while re-taking class,		*/
/*		fixed bug.										*/
/* -	2 Mar 2011 Modified for lab 7 while re-taking	*/
/*		class.											*/
/*------------------------------------------------------*/
#include "matrix.h"

//public:
matrix::matrix()
{
	resize(1,1);
}

matrix::matrix(ifstream& filein)
{
	readfile(filein);
}

matrix::~matrix()
{

}

void matrix::readfile(ifstream& filein)
{
	int rows, cols;
	filein>>rows>>cols;
	resize(rows, cols);

	for(int i=0; i<nRows; i++)
	{
		for(int j=0; j<nCols; j++)
		{
			filein>>mat[i][j];
		}
	}
}

void matrix::print()
{
	cout<<endl;
	for(int i=0; i<nRows; i++)
	{
		for(int j=0; j<nCols; j++)
		{
			cout<<mat[i][j]<<"\t";
		}
		cout<<endl;
	}
}

int matrix::getrows()
{
	return nRows;
}

int matrix::getcols()
{
	return nCols;
}

vector<vector<double> > matrix::getmat()
{
	return mat;
}

int matrix::addtest(matrix B)
{
	if((nRows==B.nRows) && (nCols==B.nCols))
	{
			return 1;
	}
	else
	{
		cout<<"Error: cannot add matrices that are not same dimensions."<<endl<<endl;
		return 2;
	}
}

void matrix::add(matrix A, matrix B)
{
	resize(A.nRows, A.nCols);
	for(int i=0; i<nRows; i++)
	{
		for(int j=0; j<nCols; j++)
		{
			mat[i][j]=A.mat[i][j]+B.mat[i][j];
		}
	}
}

int matrix::multtest(matrix B)
{
	if(nCols==B.nRows)
	{
		return 1;
	}
	else
	{
		cout<<"Error: matrix dimensions are not compatible for purposes of multiplying."<<endl<<endl;
		return 2;
	}
}

void matrix::multiply(matrix A, matrix B)
{
	double x;
	resize(A.nRows, B.nCols);
	for(int i=0; i<nRows; i++)
	{
		for(int j=0; j<nCols; j++)
		{
			x=0;
			for(int k=0; k<A.nCols; k++)
			{
				x=x+(A.mat[i][k]*B.mat[k][j]);
			}
			mat[i][j]=x;
		}
	}
}

int matrix::gausstest(vector<double> b) //tests for solvability of AX=B
{
	if(nRows==nCols)
	{
		if(nCols==b.size())
		{
			return 1;
		}
		else
		{
			cout<<"Error: cannot solve. Number of rows in b and columns in A must be equal."<<endl<<endl;
			return 2;
		}
	}
	else
	{
		cout<<"Error: cannot solve. A is not a square matrix."<<endl<<endl;
		return 2;
	}
}

vector<double> matrix::gauss(vector<double> b) //solves AX=B for X
{
	//declarations
	vector<vector<double> > A = mat;
	vector<double> x;
	x.resize(b.size());

	//forward elimination
	double fac;
	int i, j, k;

	for(k=0; k<(nCols-1); k++)
	{
		for(i=(k+1); i<nRows; i++)
		{
			fac=A[i][k]/A[k][k];

			for(j=k; j<nCols; j++)
			{
				A[i][j]=A[i][j]-fac*A[k][j];
			}

			b[i]=b[i]-fac*b[k];
		}
	}

	//back substitution
	double sum;

	x[x.size()-1]=b[b.size()-1]/A[nRows-1][nCols-1];

	for(i=(nCols-2); i>-1; i--)
	{
		sum=b[i];

		for(j=(i+1); j<nCols; j++)
		{
			sum=sum-(A[i][j]*x[j]);
		}
		x[i]=sum/(A[i][i]);
	}

	return x;
}

//private:
void matrix::resize(int rows, int cols)
{
	nRows=rows;
	nCols=cols;
	mat.resize(nRows);

	for(int i=0; i<nRows; i++)
	{
		mat[i].resize(nCols);
	}
}