Index: /issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.cpp
===================================================================
--- /issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.cpp	(revision 11733)
+++ /issm/trunk-jpl/src/c/objects/Numerics/ElementMatrix.cpp	(revision 11734)
@@ -260,45 +260,41 @@
 	this->CheckConsistency();
 
-	#ifdef _HAVE_PETSC_
-		if(this->dofsymmetrical){
-			/*only use row dofs to add values into global matrices: */
-			
-			if(this->row_fsize){
-				/*first, retrieve values that are in the f-set from the g-set values matrix: */
-				localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
-				for(i=0;i<this->row_fsize;i++){
-					for(j=0;j<this->row_fsize;j++){
-						*(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
-					}
+	if(this->dofsymmetrical){
+		/*only use row dofs to add values into global matrices: */
+
+		if(this->row_fsize){
+			/*first, retrieve values that are in the f-set from the g-set values matrix: */
+			localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
+			for(i=0;i<this->row_fsize;i++){
+				for(j=0;j<this->row_fsize;j++){
+					*(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
 				}
-				/*add local values into global  matrix, using the fglobaldoflist: */
-				MatSetValues(Kff->matrix,this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,(const double*)localvalues,ADD_VALUES);
-
-				/*Free ressources:*/
-				xfree((void**)&localvalues);
-			}
-
-
-			if((this->row_ssize!=0) && (this->row_fsize!=0)){
-				/*first, retrieve values that are in the f and s-set from the g-set values matrix: */
-				localvalues=(double*)xmalloc(this->row_fsize*this->row_ssize*sizeof(double));
-				for(i=0;i<this->row_fsize;i++){
-					for(j=0;j<this->row_ssize;j++){
-						*(localvalues+this->row_ssize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_slocaldoflist[j]);
-					}
+			}
+			/*add local values into global  matrix, using the fglobaldoflist: */
+			Kff->SetValues(this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,localvalues,ADD_VAL);
+
+			/*Free ressources:*/
+			xfree((void**)&localvalues);
+		}
+
+
+		if((this->row_ssize!=0) && (this->row_fsize!=0)){
+			/*first, retrieve values that are in the f and s-set from the g-set values matrix: */
+			localvalues=(double*)xmalloc(this->row_fsize*this->row_ssize*sizeof(double));
+			for(i=0;i<this->row_fsize;i++){
+				for(j=0;j<this->row_ssize;j++){
+					*(localvalues+this->row_ssize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_slocaldoflist[j]);
 				}
-				/*add local values into global  matrix, using the fglobaldoflist: */
-				MatSetValues(Kfs->matrix,this->row_fsize,this->row_fglobaldoflist,this->row_ssize,this->row_sglobaldoflist,(const double*)localvalues,ADD_VALUES);
-
-				/*Free ressources:*/
-				xfree((void**)&localvalues);
-			}
-		}
-		else{
-			_error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
-		}
-	#else
-		_error_("not supported yet!");
-	#endif
+			}
+			/*add local values into global  matrix, using the fglobaldoflist: */
+			Kfs->SetValues(this->row_fsize,this->row_fglobaldoflist,this->row_ssize,this->row_sglobaldoflist,localvalues,ADD_VAL);
+
+			/*Free ressources:*/
+			xfree((void**)&localvalues);
+		}
+	}
+	else{
+		_error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
+	}
 
 }
@@ -316,30 +312,26 @@
 	this->CheckConsistency();
 
-	#ifdef _HAVE_PETSC_
-		if(this->dofsymmetrical){
-			/*only use row dofs to add values into global matrices: */
-
-			if(this->row_fsize){
-				/*first, retrieve values that are in the f-set from the g-set values matrix: */
-				localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
-				for(i=0;i<this->row_fsize;i++){
-					for(j=0;j<this->row_fsize;j++){
-						*(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
-					}
+	if(this->dofsymmetrical){
+		/*only use row dofs to add values into global matrices: */
+
+		if(this->row_fsize){
+			/*first, retrieve values that are in the f-set from the g-set values matrix: */
+			localvalues=(double*)xmalloc(this->row_fsize*this->row_fsize*sizeof(double));
+			for(i=0;i<this->row_fsize;i++){
+				for(j=0;j<this->row_fsize;j++){
+					*(localvalues+this->row_fsize*i+j)=*(this->values+this->ncols*this->row_flocaldoflist[i]+this->row_flocaldoflist[j]);
 				}
-				/*add local values into global  matrix, using the fglobaldoflist: */
-				MatSetValues(Jff->matrix,this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,(const double*)localvalues,ADD_VALUES);
-
-				/*Free ressources:*/
-				xfree((void**)&localvalues);
-			}
-
-		}
-		else{
-			_error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
-		}
-	#else
-		_error_("not supported yet!");
-	#endif
+			}
+			/*add local values into global  matrix, using the fglobaldoflist: */
+			Jff->SetValues(this->row_fsize,this->row_fglobaldoflist,this->row_fsize,this->row_fglobaldoflist,localvalues,ADD_VAL);
+
+			/*Free ressources:*/
+			xfree((void**)&localvalues);
+		}
+
+	}
+	else{
+		_error_(" non dofsymmetrical matrix AddToGlobal routine not support yet!");
+	}
 
 }
Index: /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp
===================================================================
--- /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp	(revision 11733)
+++ /issm/trunk-jpl/src/c/objects/Numerics/ElementVector.cpp	(revision 11734)
@@ -166,20 +166,16 @@
 	double* localvalues=NULL;
 
-	#ifdef _HAVE_PETSC_
-		if(this->fsize){
-			/*first, retrieve values that are in the f-set from the g-set values vector: */
-			localvalues=(double*)xmalloc(this->fsize*sizeof(double));
-			for(i=0;i<this->fsize;i++){
-				localvalues[i]=this->values[this->flocaldoflist[i]];
-			}
-			/*add local values into global  vector, using the fglobaldoflist: */
-			VecSetValues(pf->vector,this->fsize,this->fglobaldoflist,(const double*)localvalues,ADD_VALUES);
-
-			/*Free ressources:*/
-			xfree((void**)&localvalues);
-		}
-	#else
-		_error_("not supported yet!");
-	#endif
+	if(this->fsize){
+		/*first, retrieve values that are in the f-set from the g-set values vector: */
+		localvalues=(double*)xmalloc(this->fsize*sizeof(double));
+		for(i=0;i<this->fsize;i++){
+			localvalues[i]=this->values[this->flocaldoflist[i]];
+		}
+		/*add local values into global  vector, using the fglobaldoflist: */
+		pf->SetValues(this->fsize,this->fglobaldoflist,localvalues,ADD_VAL);
+
+		/*Free ressources:*/
+		xfree((void**)&localvalues);
+	}
 	
 }
@@ -191,20 +187,16 @@
 	double* localvalues=NULL;
 
-	#ifdef _HAVE_PETSC_
-		if(this->fsize){
-			/*first, retrieve values that are in the f-set from the g-set values vector: */
-			localvalues=(double*)xmalloc(this->fsize*sizeof(double));
-			for(i=0;i<this->fsize;i++){
-				localvalues[i]=this->values[this->flocaldoflist[i]];
-			}
-			/*add local values into global  vector, using the fglobaldoflist: */
-			VecSetValues(pf->vector,this->fsize,this->fglobaldoflist,(const double*)localvalues,INSERT_VALUES);
-
-			/*Free ressources:*/
-			xfree((void**)&localvalues);
-		}
-	#else
-		_error_("not supported yet!");
-	#endif
+	if(this->fsize){
+		/*first, retrieve values that are in the f-set from the g-set values vector: */
+		localvalues=(double*)xmalloc(this->fsize*sizeof(double));
+		for(i=0;i<this->fsize;i++){
+			localvalues[i]=this->values[this->flocaldoflist[i]];
+		}
+		/*add local values into global  vector, using the fglobaldoflist: */
+		pf->SetValues(this->fsize,this->fglobaldoflist,localvalues,INS_VAL);
+
+		/*Free ressources:*/
+		xfree((void**)&localvalues);
+	}
 
 }
Index: /issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp
===================================================================
--- /issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp	(revision 11733)
+++ /issm/trunk-jpl/src/c/objects/Numerics/Matrix.cpp	(revision 11734)
@@ -48,5 +48,5 @@
 	this->matrix=NewMat(pM,pN);
 	#else
-	this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
+	this->matrix=new SeqMat(pM,pN);
 	#endif
 	#ifdef _HAVE_ADOLC_
@@ -64,5 +64,5 @@
 	this->matrix=NewMat(pM,pN,sparsity);
 	#else
-	this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
+	this->matrix=new SeqMat(pM,pN,sparsity);
 	#endif
 	#ifdef _HAVE_ADOLC_
@@ -97,5 +97,5 @@
 	xfree((void**)&idxn);
 	#else
-	this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
+	this->matrix=new SeqMat(serial_mat,pM,pN,sparsity);
 	#endif
 	#ifdef _HAVE_ADOLC_
@@ -113,5 +113,5 @@
 	this->matrix=NewMat(pM,pN,connectivity,numberofdofspernode);
 	#else
-	this->matrix=(double*)xcalloc(pM*pN,sizeof(double));
+	this->matrix=new SeqMat(pM,pN,connectivity,numberofdofspernode);
 	#endif
 	#ifdef _HAVE_ADOLC_
@@ -126,5 +126,5 @@
 	MatFree(&this->matrix);
 	#else
-	xfree((void**)&this->matrix);
+	delete this->matrix;
 	#endif
 	#ifdef _HAVE_ADOLC_
@@ -143,11 +143,5 @@
 	MatView(this->matrix,PETSC_VIEWER_STDOUT_WORLD);
 	#else
-	printf("Matrix size: %i-%i\n",M,N);
-	for(i=0;i<M;i++){
-		for(j=0;j<N;j++){
-			printf("%g ",*(matrix+N*i+j));
-		}
-		printf("\n");
-	}
+	this->matrix->Echo();
 	#endif
 
@@ -173,8 +167,26 @@
 	PetscMatrixToMatlabMatrix(&dataref,this->matrix);
 	#else
-	_error_("not implemented yet!");
+	dataref=this->matrix->ToMatlabMatrix();
 	#endif
 	return dataref;
 
+}
+/*}}}*/
+/*FUNCTION MatlabMatrixToMatrix{{{1*/
+Matrix* MatlabMatrixToMatrix(const mxArray* mxmatrix){
+
+	int dummy;
+	Matrix* matrix=NULL;
+
+	/*allocate matrix object: */
+	matrix=new Matrix();
+
+	#ifdef _HAVE_PETSC_
+	MatlabMatrixToPetscMatrix(&matrix->matrix,&matrix->M,&matrix->N,mxmatrix);
+	#else
+	matrix->matrix=MatlabMatrixToSeqMat(mxmatrix);
+	#endif
+	
+	return matrix;
 }
 /*}}}*/
@@ -190,5 +202,5 @@
 		#endif
 	#else
-		/*do nothing:*/
+		this->matrix->Assemble();
 	#endif
 
@@ -203,5 +215,5 @@
 		MatNorm(this->matrix,ISSMToPetscNormMode(norm_type),&norm);
 	#else
-		_error_("not implemented yet!");
+		norm=this->matrix->Norm(norm_type);
 	#endif
 	return norm;
@@ -215,5 +227,5 @@
 		MatGetSize(this->matrix,pM,pN);
 	#else
-		_error_("not implemented yet!");
+		this->matrix->GetSize(pM,pN);
 	#endif
 }
@@ -226,5 +238,5 @@
 		MatGetLocalSize(this->matrix,pM,pN);
 	#else
-		_error_("not implemented yet!");
+		this->matrix->GetLocalSize(pM,pN);
 	#endif
 }
@@ -238,5 +250,5 @@
 		MatMultPatch(this->matrix,X->vector,AX->vector);
 	#else
-		_error_("not implemented yet!");
+		this->matrix->MatMult(X->vector,AX->vector);
 	#endif
 }
@@ -249,9 +261,9 @@
 	output=new Matrix();
 
-	#ifdef _HAVE_PETSC
+	#ifdef _HAVE_PETSC_
 		_assert_(this->matrix);
 		MatDuplicate(this->matrix,MAT_COPY_VALUES,&output->matrix);
 	#else
-		_error_("not implemented yet!");
+		output->matrix=this->matrix->Duplicate();
 	#endif
 }
@@ -262,10 +274,20 @@
 	double* output=NULL;
 
-	#ifdef _HAVE_PETSC
+	#ifdef _HAVE_PETSC_
 		MatToSerial(&output,this->matrix);
 	#else
-		_error_("not implemented yet!");
+		output=this->matrix->ToSerial();
 	#endif
 	return output;
 }
 /*}}}*/
+/*FUNCTION Matrix::SetValues{{{1*/
+void Matrix::SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode){
+
+	#ifdef _HAVE_PETSC_
+		MatSetValues(this->matrix,m,idxm,n,idxn,values,ISSMToPetscInsertMode(mode));
+	#else
+		this->matrix->SetValues(m,idxm,n,idxn,values,mode);
+	#endif
+}
+/*}}}*/
Index: /issm/trunk-jpl/src/c/objects/Numerics/Matrix.h
===================================================================
--- /issm/trunk-jpl/src/c/objects/Numerics/Matrix.h	(revision 11733)
+++ /issm/trunk-jpl/src/c/objects/Numerics/Matrix.h	(revision 11734)
@@ -37,5 +37,5 @@
 		Mat matrix; 
 		#else
-		double* matrix; 
+		SeqMat* matrix; 
 		#endif
 		#ifdef _HAVE_ADOLC_
@@ -63,6 +63,12 @@
 		Matrix* Duplicate(void);
 		double* ToSerial(void);
+		void SetValues(int m,int* idxm,int n,int* idxn,double* values,InsMode mode);
 		/*}}}*/
 
 };
+/*API: */
+#ifdef _SERIAL_
+Matrix* MatlabMatrixToMatrix(const mxArray* mxmatrix);
+#endif
+
 #endif //#ifndef _MATRIX_H_
Index: /issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp
===================================================================
--- /issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp	(revision 11733)
+++ /issm/trunk-jpl/src/c/objects/Numerics/Vector.cpp	(revision 11734)
@@ -25,5 +25,4 @@
 Vector::Vector(){
 
-	this->M=0;
 	#ifdef _HAVE_PETSC_
 	this->vector=NULL;
@@ -42,6 +41,5 @@
 	this->vector=NewVec(pM,fromlocalsize);
 	#else
-	this->M=pM;
-	this->vector=(double*)xcalloc(pM,sizeof(double));
+	this->vector=new SeqVec(pM,fromlocalsize);
 	#endif
 	#ifdef _HAVE_ADOLC_
@@ -51,5 +49,5 @@
 /*}}}*/
 /*FUNCTION Vector::Vector(double* serial_vec,int M){{{1*/
-Vector::Vector(double* serial_vec,int pM){
+Vector::Vector(double* serial_vec,int M){
 
 	int i,j;
@@ -58,5 +56,5 @@
 		int* idxm=NULL;
 
-		this->vector=NewVec(pM);
+		this->vector=NewVec(M);
 		
 		idxm=(int*)xmalloc(M*sizeof(int));
@@ -70,9 +68,8 @@
 
 	#else
-		this->M=pM;
-		this->vector=(double*)xcalloc(pM,sizeof(double));
-	#endif
-	#ifdef _HAVE_ADOLC_
-		this->avector=(adouble*)xmalloc(pM*sizeof(adouble));
+		this->vector=new SeqVec(serial_vec,M);
+	#endif
+	#ifdef _HAVE_ADOLC_
+		this->avector=(adouble*)xmalloc(M*sizeof(adouble));
 	#endif
 }
@@ -82,7 +79,4 @@
 Vector::Vector(Vec petsc_vec){
 
-	/*Get Vector size*/
-	VecGetSize(petsc_vec,&this->M);
-
 	/*copy vector*/
 	VecDuplicate(petsc_vec,&this->vector);
@@ -98,5 +92,5 @@
 	VecFree(&this->vector);
 	#else
-	xfree((void**)&this->vector);
+	delete this->vector;
 	#endif
 	#ifdef _HAVE_ADOLC_
@@ -115,16 +109,9 @@
 	VecView(this->vector,PETSC_VIEWER_STDOUT_WORLD);
 	#else
-	printf("Vector size: %i\n",M);
-	for(i=0;i<M;i++){
-		printf("%g\n ",*(vector+i));
-	}
-	#endif
-
-	#ifdef _HAVE_ADOLC_
-	/*Not sure about that one. Should we use the overloaded operator >>?*/
-	printf("ADOLC Vector equivalent:" );
-	for(i=0;i<M;i++){
-		printf("%g\n ",*(avector+i));
-	}
+	this->vector->Echo();
+	#endif
+
+	#ifdef _HAVE_ADOLC_
+	/*do nothing for now: */
 	#endif
 }
@@ -139,8 +126,26 @@
 	PetscVectorToMatlabVector(&dataref,this->vector);
 	#else
-	_error_("not implemented yet!");
+	dataref=this->vector->ToMatlabVector();
 	#endif
 	return dataref;
 
+}
+/*}}}*/
+/*FUNCTION MatlabVectorToVector{{{1*/
+Vector* MatlabVectorToVector(const mxArray* mxvector){
+
+	int dummy;
+	Vector* vector=NULL;
+
+	/*allocate vector object: */
+	vector=new Vector();
+
+	#ifdef _HAVE_PETSC_
+	MatlabVectorToPetscVector(&vector->vector,&dummy,mxvector);
+	#else
+	vector->vector=MatlabVectorToSeqVec(mxvector);
+	#endif
+	
+	return vector;
 }
 /*}}}*/
@@ -154,5 +159,5 @@
 		VecAssemblyEnd(this->vector);
 	#else
-		/*do nothing*/
+		this->vector->Assemble();
 	#endif
 
@@ -167,5 +172,5 @@
 		VecSetValues(this->vector,ssize,list,values,ISSMToPetscInsertMode(mode));
 	#else
-		_error_("not implemented yet!");
+		this->vector->SetValues(ssize,list,values,mode);
 	#endif
 
@@ -179,5 +184,5 @@
 		VecSetValues(this->vector,1,&dof,&value,ISSMToPetscInsertMode(mode));
 	#else
-		_error_("not implemented yet!");
+		this->vector->SetValue(dof,value,mode);
 	#endif
 
@@ -191,5 +196,5 @@
 		VecGetValues(this->vector,1,&dof,pvalue);
 	#else
-		_error_("not implemented yet!");
+	this->vector->GetValue(pvalue,dof);
 	#endif
 }
@@ -202,5 +207,5 @@
 		VecGetSize(this->vector,pM);
 	#else
-		_error_("not implemented yet!");
+		this->vector->GetSize(pM);
 	#endif
 
@@ -214,5 +219,5 @@
 		VecGetLocalSize(this->vector,pM);
 	#else
-		_error_("not implemented yet!");
+		this->vector->GetLocalSize(pM);
 	#endif
 
@@ -229,7 +234,6 @@
 		VecDuplicate(this->vector,&vec_output);
 		output->vector=vec_output;
-		VecGetSize(output->vector,&output->M);
-	#else
-		_error_("not implemented yet!");
+	#else
+		output->vector=this->vector->Duplicate();
 	#endif
 	return output;
@@ -244,5 +248,5 @@
 		VecSet(this->vector,value);
 	#else
-		_error_("not implemented yet!");
+		this->vector->Set(value);
 	#endif
 
@@ -256,5 +260,5 @@
 		VecAXPY(this->vector,a,X->vector);
 	#else
-		_error_("not implemented yet!");
+		this->vector->AXPY(X->vector,a);
 	#endif
 }
@@ -267,5 +271,5 @@
 		VecAYPX(this->vector,a,X->vector);
 	#else
-		_error_("not implemented yet!");
+		this->vector->AYPX(X->vector,a);
 	#endif
 
@@ -280,5 +284,5 @@
 		VecToMPISerial(&vec_serial, this->vector);
 	#else
-		_error_("not implemented yet!");
+		vec_serial=this->vector->ToMPISerial();
 	#endif
 
@@ -293,5 +297,5 @@
 		VecCopy(this->vector,to->vector);
 	#else
-		_error_("not implemented yet!");
+		this->vector->Copy(to->vector);
 	#endif
 
@@ -306,5 +310,5 @@
 		VecNorm(this->vector,ISSMToPetscNormMode(norm_type),&norm);
 	#else
-		_error_("not implemented yet!");
+		norm=this->vector->Norm(norm_type);
 	#endif
 	return norm;
@@ -318,5 +322,5 @@
 		VecScale(this->vector,scale_factor); 
 	#else
-		_error_("not implemented yet!");
+		this->vector->Scale(scale_factor);
 	#endif
 }
@@ -330,5 +334,5 @@
 		VecDot(this->vector,vector->vector,&dot);
 	#else
-		_error_("not implemented yet!");
+		dot=this->vector->Dot(vector->vector);
 	#endif
 	return dot;
@@ -342,6 +346,6 @@
 		VecPointwiseDivide(this->vector,x->vector,y->vector);
 	#else
-		_error_("not implemented yet!");
-	#endif
-}
-/*}}}*/
+		this->vector->PointwiseDivide(x->vector,y->vector);
+	#endif
+}
+/*}}}*/
Index: /issm/trunk-jpl/src/c/objects/Numerics/Vector.h
===================================================================
--- /issm/trunk-jpl/src/c/objects/Numerics/Vector.h	(revision 11733)
+++ /issm/trunk-jpl/src/c/objects/Numerics/Vector.h	(revision 11734)
@@ -31,10 +31,8 @@
 	public:
 	
-		int M;
-
 		#ifdef _HAVE_PETSC_
 		Vec vector; 
 		#else
-		double* vector; 
+		SeqVec* vector;
 		#endif
 		#ifdef _HAVE_ADOLC_
@@ -74,3 +72,9 @@
 		/*}}}*/
 };
+
+/*API: */
+#ifdef _SERIAL_
+Vector* MatlabVectorToVector(const mxArray* mxvector);
+#endif
+
 #endif //#ifndef _VECTOR_H_