Index: /issm/trunk/src/c/Reducematrixfromgtofx/Reducematrixfromgton.cpp.petsc2
===================================================================
--- /issm/trunk/src/c/Reducematrixfromgtofx/Reducematrixfromgton.cpp.petsc2	(revision 2062)
+++ /issm/trunk/src/c/Reducematrixfromgtofx/Reducematrixfromgton.cpp.petsc2	(revision 2062)
@@ -0,0 +1,73 @@
+/*!\file: Reducematrixfromgton
+ * \brief
+ */ 
+
+#ifdef HAVE_CONFIG_H
+	#include "config.h"
+#else
+#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
+#endif
+
+#undef __FUNCT__ 
+#define __FUNCT__ "Reducematrixfromgton"
+
+#include "../toolkits/toolkits.h"
+
+
+void  Reducematrixfromgton(Mat* pKnn,Mat Kgg,Mat Gmn,double* pv_m,int msize, double* pv_n,int nsize,int flag){
+
+
+	/*output: */
+	Mat Knn=NULL;
+
+	/*intermediary: */
+	Mat Kmm=NULL;
+	Mat Kmn=NULL;
+	Mat Knm=NULL;
+	Mat KnmGmn=NULL;
+	Mat KmmGmn=NULL;
+	Mat tGmn=NULL;
+	Mat tGmnKmn=NULL;
+	PetscScalar a=1;
+
+	if (msize){
+		
+		/*Partition K_gg*/
+		MatPartition(&Kmm, Kgg, pv_m,msize,pv_m,msize); 
+		MatPartition(&Kmn, Kgg, pv_m,msize,pv_n,nsize); 
+		MatPartition(&Knm, Kgg, pv_n,nsize,pv_m,msize); 
+		MatPartition(&Knn, Kgg, pv_n,nsize,pv_n,nsize);
+	
+		/*Reduce K_gg to K_nn*/
+		//K_nn = K_nn + K_nm * G_mn;
+		MatMatMult(Knm,Gmn,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&KnmGmn);
+		MatAXPY(Knn,a,KnmGmn,DIFFERENT_NONZERO_PATTERN);
+
+		if (flag!=1){
+			//K_mn = K_mn + K_mm * G_mn;
+			MatMatMult(Kmm,Gmn,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&KmmGmn);
+			MatAXPY(Kmn,a,KmmGmn,DIFFERENT_NONZERO_PATTERN);
+			
+			//K_nn = K_nn + G_mn' * K_mn;
+			MatTranspose(Gmn,&tGmn);
+			MatMatMult(tGmn,Kmn,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&tGmnKmn);
+			MatAXPY(Knn,a,tGmnKmn,DIFFERENT_NONZERO_PATTERN);
+		}
+	}
+	else{
+		/*msize=0 just return input*/
+		MatDuplicate(Kgg,MAT_COPY_VALUES,&Knn);
+	}	
+
+	/*Free ressources:*/
+	MatFree(&Kmm);
+	MatFree(&Kmn);
+	MatFree(&Knm);
+	MatFree(&KnmGmn);
+	MatFree(&KmmGmn);
+	MatFree(&tGmn);
+	MatFree(&tGmnKmn);
+	
+	/*Assign output pointers:*/
+	*pKnn=Knn;
+}
Index: /issm/trunk/src/c/Reducematrixfromgtofx/Reducematrixfromgton.cpp.petsc3
===================================================================
--- /issm/trunk/src/c/Reducematrixfromgtofx/Reducematrixfromgton.cpp.petsc3	(revision 2062)
+++ /issm/trunk/src/c/Reducematrixfromgtofx/Reducematrixfromgton.cpp.petsc3	(revision 2062)
@@ -0,0 +1,73 @@
+/*!\file: Reducematrixfromgton
+ * \brief
+ */ 
+
+#ifdef HAVE_CONFIG_H
+	#include "config.h"
+#else
+#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
+#endif
+
+#undef __FUNCT__ 
+#define __FUNCT__ "Reducematrixfromgton"
+
+#include "../toolkits/toolkits.h"
+
+
+void  Reducematrixfromgton(Mat* pKnn,Mat Kgg,Mat Gmn,double* pv_m,int msize, double* pv_n,int nsize,int flag){
+
+
+	/*output: */
+	Mat Knn=NULL;
+
+	/*intermediary: */
+	Mat Kmm=NULL;
+	Mat Kmn=NULL;
+	Mat Knm=NULL;
+	Mat KnmGmn=NULL;
+	Mat KmmGmn=NULL;
+	Mat tGmn=NULL;
+	Mat tGmnKmn=NULL;
+	PetscScalar a=1;
+
+	if (msize){
+		
+		/*Partition K_gg*/
+		MatPartition(&Kmm, Kgg, pv_m,msize,pv_m,msize); 
+		MatPartition(&Kmn, Kgg, pv_m,msize,pv_n,nsize); 
+		MatPartition(&Knm, Kgg, pv_n,nsize,pv_m,msize); 
+		MatPartition(&Knn, Kgg, pv_n,nsize,pv_n,nsize);
+	
+		/*Reduce K_gg to K_nn*/
+		//K_nn = K_nn + K_nm * G_mn;
+		MatMatMult(Knm,Gmn,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&KnmGmn);
+		MatAXPY(Knn,a,KnmGmn,DIFFERENT_NONZERO_PATTERN);
+
+		if (flag!=1){
+			//K_mn = K_mn + K_mm * G_mn;
+			MatMatMult(Kmm,Gmn,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&KmmGmn);
+			MatAXPY(Kmn,a,KmmGmn,DIFFERENT_NONZERO_PATTERN);
+			
+			//K_nn = K_nn + G_mn' * K_mn;
+			MatTranspose(Gmn,MAT_INITIAL_MATRIX,&tGmn);
+			MatMatMult(tGmn,Kmn,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&tGmnKmn);
+			MatAXPY(Knn,a,tGmnKmn,DIFFERENT_NONZERO_PATTERN);
+		}
+	}
+	else{
+		/*msize=0 just return input*/
+		MatDuplicate(Kgg,MAT_COPY_VALUES,&Knn);
+	}	
+
+	/*Free ressources:*/
+	MatFree(&Kmm);
+	MatFree(&Kmn);
+	MatFree(&Knm);
+	MatFree(&KnmGmn);
+	MatFree(&KmmGmn);
+	MatFree(&tGmn);
+	MatFree(&tGmnKmn);
+	
+	/*Assign output pointers:*/
+	*pKnn=Knn;
+}
