/*!\file SeqMat.cpp
 * \brief: implementation of the SeqMat object
 */

/*Headers:*/
/*{{{1*/
#ifdef HAVE_CONFIG_H
	#include <config.h>
#else
#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
#endif

#include <stdio.h>
#include <string.h>
#include "../toolkits.h"
#include "../../shared/shared.h"

/*}}}*/

/*SeqMat constructors and destructor*/
/*FUNCTION SeqMat::SeqMat(){{{1*/
SeqMat::SeqMat(){

	this->M=0;
	this->N=0;
	this->matrix=NULL;
}
/*}}}*/
/*FUNCTION SeqMat::SeqMat(int M,int N){{{1*/
SeqMat::SeqMat(int pM,int pN){

	this->M=pM;
	this->N=pN;
	this->matrix=NULL;
	if(M*N) this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
}
/*}}}*/
/*FUNCTION SeqMat::SeqMat(int M,int N, double sparsity){{{1*/
SeqMat::SeqMat(int pM,int pN, double sparsity){

	this->M=pM;
	this->N=pN;
	this->matrix=NULL;
	if(M*N) this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
}
/*}}}*/
/*FUNCTION SeqMat(double* serial_mat,int M,int N,double sparsity){{{1*/
SeqMat::SeqMat(double* serial_mat,int pM,int pN,double sparsity){

	int i,j;

	this->M=pM;
	this->N=pN;
	this->matrix=NULL;
	if(M*N){
		this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
		memcpy(this->matrix,serial_mat,pM*pN*sizeof(double));
	}

}
/*}}}*/
/*FUNCTION SeqMat::SeqMat(int M,int N, int connectivity, int numberofdofspernode){{{1*/
SeqMat::SeqMat(int pM,int pN, int connectivity,int numberofdofspernode){

	this->M=pM;
	this->N=pN;
	this->matrix=NULL;
	if(M*N)this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
}
/*}}}*/
/*FUNCTION SeqMat::~SeqMat(){{{1*/
SeqMat::~SeqMat(){

	xfree((void**)&this->matrix);
	M=0;
	N=0;
}
/*}}}*/

/*SeqMat specific routines: */
/*FUNCTION SeqMat::Echo{{{1*/
void SeqMat::Echo(void){

	int i,j;
	printf("SeqMat size %i-%i\n",this->M,this->N);
	for(i=0;i<M;i++){
		for(j=0;j<N;j++){
			printf("%g ",this->matrix[N*i+j]);
		}
		printf("\n");
	}
}
/*}}}*/
/*FUNCTION SeqMat::Assemble{{{1*/
void SeqMat::Assemble(void){
		
	/*do nothing*/

}
/*}}}*/
/*FUNCTION SeqMat::Norm{{{1*/
double SeqMat::Norm(NormMode mode){

	double norm;
	double absolute;
	int i,j;

	switch(mode){
		case NORM_INF:
			norm=0;
			for(i=0;i<this->M;i++){
				absolute=0;
				for(j=0;j<this->N;j++){
					absolute+=fabs(this->matrix[N*i+j]);
				}
				norm=max(norm,absolute);
			}
			return norm;
			break;
		default:
			_error_("unknown norm !");
			break;
	}
}
/*}}}*/
/*FUNCTION SeqMat::GetSize{{{1*/
void SeqMat::GetSize(int* pM,int* pN){

	*pM=this->M;
	*pN=this->N;

}
/*}}}*/
/*FUNCTION SeqMat::GetLocalSize{{{1*/
void SeqMat::GetLocalSize(int* pM,int* pN){
	
	*pM=this->M;
	*pN=this->N;

}
/*}}}*/
/*FUNCTION SeqMat::MatMult{{{1*/
void SeqMat::MatMult(SeqVec* X,SeqVec* AX){

	int i,j;
	int XM,AXM;
	double dummy;

	X->GetSize(&XM);
	AX->GetSize(&AXM);

	if(M!=AXM)_error_("A and AX should have the same number of rows!");
	if(N!=XM)_error_("A and X should have the same number of columns!");

	for(i=0;i<M;i++){
		dummy=0;
		for(j=0;j<N;j++){
			dummy+= this->matrix[N*i+j]*X->vector[j];
		}
		AX->vector[i]=dummy;
	}

}
/*}}}*/
/*FUNCTION SeqMat::Duplicate{{{1*/
SeqMat* SeqMat::Duplicate(void){

	double dummy=0;

	return new SeqMat(this->matrix,this->M,this->N,dummy);

}
/*}}}*/
/*FUNCTION SeqMat::ToSerial{{{1*/
double* SeqMat::ToSerial(void){

	double* buffer=NULL;

	if(this->M*this->N){
		buffer=(double*)xmalloc(this->M*this->N*sizeof(double));
		memcpy(buffer,this->matrix,this->M*this->N*sizeof(double));
	}
	return buffer;

}
/*}}}*/
/*FUNCTION SeqMat::SetValues{{{1*/
void SeqMat::SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode){
	
	int i,j;
	switch(mode){
		case ADD_VAL:
			for(i=0;i<m;i++) for(j=0;j<n;j++) this->matrix[N*idxm[i]+idxn[j]]+=values[n*i+j];
			break;
		case INS_VAL:
			for(i=0;i<m;i++) for(j=0;j<n;j++) this->matrix[N*idxm[i]+idxn[j]]=values[n*i+j];
			break;
		default:
			_error_("unknown insert mode!");
			break;
	}

}
/*}}}*/
/*FUNCTION SeqMat::Convert{{{1*/
void SeqMat::Convert(MatrixType type){

	/*do nothing*/

}
/*}}}*/
