/*! \file VecTranspose.cpp
 *  \brief: transpose of a petsc vector
 */

#include "./petscpatches.h"
#include "../../../shared/shared.h"

int VecTranspose(Vec* ptvector,Vec vector){

	int i;
	int size;
	int lower_row,upper_row,range;
	
	int* idxm=NULL;
	double* values=NULL;

	int* tidxm=NULL;
	double* tvalues=NULL;
	
	/*output: */
	Vec tvector=NULL;

	/*Get size of input vector: */
	VecGetSize(vector,&size);
	
	/*Create new vector of same size: */
	tvector=NewVec(size);

	/*Extract values locally from input vector: */
	VecGetOwnershipRange(vector,&lower_row,&upper_row);
	upper_row--;
	range=upper_row-lower_row+1;    

	if (range){
		idxm=xNew<int>(range); 
		tidxm=xNew<int>(range); 
		for (i=0;i<range;i++){
			*(idxm+i)=lower_row+i;
		} 
		values=xNew<double>(range);
		tvalues=xNew<double>(range);
		
		VecGetValues(vector,range,idxm,values);
		
		/*Transfer values into tidxm, and idxm into tvalues: */
		for (i=0;i<range;i++){
			tidxm[i]=(int)values[i];
			tvalues[i]=(double)idxm[i];
		}
		VecSetValues(tvector,range,tidxm,tvalues,INSERT_VALUES);
	}

	/*Assemble: */
	VecAssemblyBegin(tvector);
	VecAssemblyEnd(tvector);

	/*Free ressources: */
	xDelete<int>(idxm);
	xDelete<double>(values);
	xDelete<int>(tidxm);
	xDelete<double>(tvalues);

	/*Assign output pointers: */
	*ptvector=tvector;
	return 1;
}
