/*!\file PetscVec.cpp
 * \brief: implementation of the PetscVec 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 "../petscincludes.h"
#include "../../../shared/shared.h"

/*}}}*/

/*PetscVec constructors and destructor*/
/*FUNCTION PetscVec::PetscVec(){{{*/
PetscVec::PetscVec(){
	this->vector=NULL;
	#ifdef _HAVE_ADOLC_
	this->avector=NULL;
	#endif
}
/*}}}*/
/*FUNCTION PetscVec::PetscVec(int M,bool fromlocalsize){{{*/
PetscVec::PetscVec(int M,bool fromlocalsize){
	
	this->vector=NewVec(M,fromlocalsize);

}
/*}}}*/
/*FUNCTION PetscVec::PetscVec(Vec petsc_vec){{{*/
PetscVec::PetscVec(Vec petsc_vec){

	if(petsc_vec==NULL){
		this->vector=NewVec(0);
	}
	else{
		/*copy vector*/
		VecDuplicate(petsc_vec,&this->vector);
		VecCopy(petsc_vec,this->vector);
	}

}
/*}}}*/
/*FUNCTION PetscVec::PetscVec(IssmDouble* serial_vec,int M){{{*/
PetscVec::PetscVec(IssmDouble* serial_vec,int M){

	int* idxm=NULL;
	if(M)idxm=xNew<int>(M);
	for(int i=0;i<M;i++) idxm[i]=i;

	this->vector=NewVec(M);
	VecSetValues(this->vector,M,idxm,serial_vec,INSERT_VALUES);
	VecAssemblyBegin(this->vector);
	VecAssemblyEnd(this->vector);

	xDelete<int>(idxm);
}
/*}}}*/
/*FUNCTION PetscVec::~PetscVec(){{{*/
PetscVec::~PetscVec(){
    VecFree(&this->vector);
}
/*}}}*/

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

	_assert_(this->vector);
	VecView(this->vector,PETSC_VIEWER_STDOUT_WORLD);
}
/*}}}*/
/*FUNCTION PetscVec::Assemble{{{*/
void PetscVec::Assemble(void){
		
	_assert_(this->vector);
	VecAssemblyBegin(this->vector); 
	VecAssemblyEnd(this->vector);

}
/*}}}*/
/*FUNCTION PetscVec::SetValues{{{*/
void PetscVec::SetValues(int ssize, int* list, IssmDouble* values, InsMode mode){
	
	_assert_(this->vector);
	VecSetValues(this->vector,ssize,list,values,ISSMToPetscInsertMode(mode));

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

	_assert_(this->vector);
	VecSetValues(this->vector,1,&dof,&value,ISSMToPetscInsertMode(mode));

}
/*}}}*/
/*FUNCTION PetscVec::GetValue{{{*/
void PetscVec::GetValue(IssmDouble* pvalue,int dof){

	_assert_(this->vector);
	VecGetValues(this->vector,1,&dof,pvalue);

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

	_assert_(this->vector);
	VecGetSize(this->vector,pM);

}
/*}}}*/
/*FUNCTION PetscVec::GetLocalSize{{{*/
void PetscVec::GetLocalSize(int* pM){

	_assert_(this->vector);
	VecGetLocalSize(this->vector,pM);

}
/*}}}*/
/*FUNCTION PetscVec::Duplicate{{{*/
PetscVec* PetscVec::Duplicate(void){
	
	PetscVec* output=NULL;
	_assert_(this->vector);
	Vec vec_output=NULL;
	VecDuplicate(this->vector,&vec_output);
	output=new PetscVec(vec_output);
	VecFree(&vec_output);

	return output;
}
/*}}}*/
/*FUNCTION PetscVec::Set{{{*/
void PetscVec::Set(IssmDouble value){

	_assert_(this->vector);
	VecSet(this->vector,value);

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

	_assert_(this->vector);
	VecAXPY(this->vector,a,X->vector);

}
/*}}}*/
/*FUNCTION PetscVec::AYPX{{{*/
void PetscVec::AYPX(PetscVec* X, IssmDouble a){

	_assert_(this->vector);
	VecAYPX(this->vector,a,X->vector);

}
/*}}}*/
/*FUNCTION PetscVec::ToMPISerial{{{*/
IssmDouble* PetscVec::ToMPISerial(void){
	
	IssmDouble* vec_serial=NULL;
	VecToMPISerial(&vec_serial, this->vector);
	return vec_serial;

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

	if(this->vector) VecCopy(this->vector,to->vector);

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

	IssmDouble norm=0;
	_assert_(this->vector);
	VecNorm(this->vector,ISSMToPetscNormMode(mode),&norm);
	return norm;

}
/*}}}*/
/*FUNCTION PetscVec::Scale{{{*/
void PetscVec::Scale(IssmDouble scale_factor){

	_assert_(this->vector);
	VecScale(this->vector,scale_factor); 

}
/*}}}*/
/*FUNCTION PetscVec::Dot{{{*/
IssmDouble PetscVec::Dot(PetscVec* input){

	IssmDouble dot;
	_assert_(this->vector);
	VecDot(this->vector,input->vector,&dot);
	return dot;

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

	_assert_(this->vector);
	VecPointwiseDivide(this->vector,x->vector,y->vector);

}
/*}}}*/
