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

/*Headers:*/
/*{{{*/
#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"

/*}}}*/

/*SeqVec constructors and destructor*/
/*FUNCTION SeqVec::SeqVec(){{{*/
SeqVec::SeqVec(){

	this->M=0;
	this->vector=NULL;
}
/*}}}*/
/*FUNCTION SeqVec::SeqVec(int M,bool fromlocalsize){{{*/
SeqVec::SeqVec(int pM,bool fromlocalsize){

	this->M=pM;
	this->vector=NULL;
	if(this->M) this->vector=xNewZeroInit<IssmDouble>(pM);
}
/*}}}*/
/*FUNCTION SeqVec::SeqVec(IssmDouble* serial_vec,int M){{{*/
SeqVec::SeqVec(IssmDouble* buffer,int pM){

	int i,j;

	this->M=pM;
	this->vector=NULL;
	if(this->M){
		this->vector=xNew<IssmDouble>(pM);
		xMemCpy<IssmDouble>(this->vector,buffer,pM);
	}
}
/*}}}*/
		/*FUNCTION SeqVec::~SeqVec(){{{*/
SeqVec::~SeqVec(){
	xDelete<IssmDouble>(this->vector);
	M=0;
}
/*}}}*/

/*SeqVec specific routines: */
/*FUNCTION SeqVec::Echo{{{*/
void SeqVec::Echo(void){

	int i;
	_printLine_("SeqVec size " << this->M);
	for(i=0;i<M;i++){
		_printString_(vector[i] << "\n ");
	}
}
/*}}}*/

/*FUNCTION SeqVec::Assemble{{{*/
void SeqVec::Assemble(void){
		
	/*do nothing*/

}
/*}}}*/
/*FUNCTION SeqVec::SetValues{{{*/
void SeqVec::SetValues(int ssize, int* list, IssmDouble* values, InsMode mode){
	
	int i;
	switch(mode){
		case ADD_VAL:
			for(i=0;i<ssize;i++) this->vector[list[i]]+=values[i];
			break;
		case INS_VAL:
			for(i=0;i<ssize;i++) this->vector[list[i]]=values[i];
			break;
		default:
			_error2_("unknown insert mode!");
			break;
	}

}
/*}}}*/
/*FUNCTION SeqVec::SetValue{{{*/
void SeqVec::SetValue(int dof, IssmDouble value, InsMode mode){

	switch(mode){
		case ADD_VAL:
			this->vector[dof]+=value;
			break;
		case INS_VAL:
			this->vector[dof]=value;
			break;
		default:
			_error2_("unknown insert mode!");
			break;
	}
}
/*}}}*/
/*FUNCTION SeqVec::GetValue{{{*/
void SeqVec::GetValue(IssmDouble* pvalue,int dof){

	*pvalue=this->vector[dof];

}
/*}}}*/
		
/*FUNCTION SeqVec::GetSize{{{*/
void SeqVec::GetSize(int* pM){

	*pM=this->M;

}
/*}}}*/
/*FUNCTION SeqVec::GetLocalSize{{{*/
void SeqVec::GetLocalSize(int* pM){
	
	*pM=this->M;

}
/*}}}*/
/*FUNCTION SeqVec::Duplicate{{{*/
SeqVec* SeqVec::Duplicate(void){
	
	return new SeqVec(this->vector,this->M);

}
/*}}}*/
/*FUNCTION SeqVec::Set{{{*/
void SeqVec::Set(IssmDouble value){

	int i;
	for(i=0;i<this->M;i++)this->vector[i]=value;

}
/*}}}*/
/*FUNCTION SeqVec::AXPY{{{*/
void SeqVec::AXPY(SeqVec* X, IssmDouble a){

	int i;

	/*y=a*x+y where this->vector is y*/
	for(i=0;i<this->M;i++)this->vector[i]=a*X->vector[i]+this->vector[i];
	
}
/*}}}*/
/*FUNCTION SeqVec::AYPX{{{*/
void SeqVec::AYPX(SeqVec* X, IssmDouble a){
	
	int i;

	/*y=x+a*y where this->vector is y*/
	for(i=0;i<this->M;i++)this->vector[i]=X->vector[i]+a*this->vector[i];

}
/*}}}*/
/*FUNCTION SeqVec::ToMPISerial{{{*/
IssmDouble* SeqVec::ToMPISerial(void){

	IssmDouble* buffer=NULL;

	if(this->M){
		buffer=xNew<IssmDouble>(this->M);
		xMemCpy<IssmDouble>(buffer,this->vector,this->M);
	}
	return buffer;

}
/*}}}*/
/*FUNCTION SeqVec::Copy{{{*/
void SeqVec::Copy(SeqVec* to){

	int i;

	to->M=this->M;
	for(i=0;i<this->M;i++)to->vector[i]=this->vector[i];

}
/*}}}*/
/*FUNCTION SeqVec::Norm{{{*/
IssmDouble SeqVec::Norm(NormMode mode){

	IssmDouble norm;
	int i;

	switch(mode){
		case NORM_INF:
			norm=0; for(i=0;i<this->M;i++)norm=max(norm,fabs(this->vector[i]));
			return norm;
			break;
		case NORM_TWO:
			norm=0; 
			for(i=0;i<this->M;i++)norm+=pow(this->vector[i],2);
			return sqrt(norm);
			break;
		default:
			_error2_("unknown norm !");
			break;
	}
}
/*}}}*/
/*FUNCTION SeqVec::Scale{{{*/
void SeqVec::Scale(IssmDouble scale_factor){

	int i;
	for(i=0;i<this->M;i++)this->vector[i]=scale_factor*this->vector[i];
	
}
/*}}}*/
/*FUNCTION SeqVec::Dot{{{*/
IssmDouble SeqVec::Dot(SeqVec* input){

	int i;

	IssmDouble dot=0;
	for(i=0;i<this->M;i++)dot+=this->vector[i]*input->vector[i];
	return dot;

}
/*}}}*/
/*FUNCTION SeqVec::PointwiseDivide{{{*/
void SeqVec::PointwiseDivide(SeqVec* x,SeqVec* y){

	int i;
	/*pointwise w=x/y where this->vector is w: */
	for(i=0;i<this->M;i++)this->vector[i]=x->vector[i]/y->vector[i];
}
/*}}}*/
