Index: /issm/trunk/src/c/objects/Elements/Penta.cpp
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 5768)
+++ /issm/trunk/src/c/objects/Elements/Penta.cpp	(revision 5769)
@@ -894,7 +894,11 @@
 			GetSolutionFromInputsDiagnosticStokes(solution);
 		}
-		else{
+		else if (approximation==MacAyealApproximationEnum || approximation==PattynApproximationEnum ||approximation==HutterApproximationEnum){
 			GetSolutionFromInputsDiagnosticHoriz(solution);
 		}
+		else if(approximation==MacAyealPattynApproximationEnum || approximation==PattynStokesApproximationEnum){
+			return; //do nothing the elements around with update the solution
+		}
+		else ISSMERROR("approximation type : %i (%s) not implemented yet",approximation,EnumToString(approximation));
 	}
 	else if(analysis_type==DiagnosticHutterAnalysisEnum){
@@ -1888,4 +1892,7 @@
 			this->inputs->AddInput(new IntInput(ApproximationEnum,StokesApproximationEnum));
 		}
+		else if (*(iomodel->elements_type+index)==PattynStokesApproximationEnum){
+			this->inputs->AddInput(new IntInput(ApproximationEnum,PattynStokesApproximationEnum));
+		}
 		else if (*(iomodel->elements_type+index)==NoneApproximationEnum){
 			this->inputs->AddInput(new IntInput(ApproximationEnum,NoneApproximationEnum));
@@ -3988,10 +3995,5 @@
 	/*If the element is a coupling, do nothing: every grid is also on an other elements 
 	 * (as coupling is between MacAyeal and Pattyn) so the other element will take care of it*/
-	if(approximation==MacAyealPattynApproximationEnum){
-	 return;
-	}
-	else{
-		GetDofList(&doflist,approximation);
-	}
+	GetDofList(&doflist,approximation);
 
 	/*Ok, we have vx and vy in values, fill in vx and vy arrays: */
@@ -4445,4 +4447,7 @@
 	else if (approximation==MacAyealPattynApproximationEnum){
 		InputUpdateFromSolutionDiagnosticMacAyealPattyn(solution);
+	}
+	else if (approximation==PattynStokesApproximationEnum){
+		InputUpdateFromSolutionDiagnosticPattynStokes(solution);
 	}
 }
@@ -4707,4 +4712,91 @@
 }
 
+/*}}}*/
+/*FUNCTION Penta::InputUpdateFromSolutionDiagnosticPattynStokes {{{1*/
+void  Penta::InputUpdateFromSolutionDiagnosticPattynStokes(double* solution){
+
+	int i;
+
+	const int    numdofpervertexp=2;
+	const int    numdofpervertexs=4;
+	const int    numdofp=numdofpervertexp*NUMVERTICES;
+	const int    numdofs=numdofpervertexs*NUMVERTICES;
+	int*         doflistp=NULL;
+	int*         doflists=NULL;
+	double       pattyn_values[numdofp];
+	double       stokes_values[numdofs];
+	double       vx[NUMVERTICES];
+	double       vy[NUMVERTICES];
+	double       vz[NUMVERTICES];
+	double       vel[NUMVERTICES];
+	int          dummy;
+	double       pressure[NUMVERTICES];
+	double       xyz_list[NUMVERTICES][3];
+
+	double *vz_ptr         = NULL;
+	Penta  *penta          = NULL;
+
+	/*OK, we have to add results of this element for pattyn 
+	 * and results from the penta at base for macayeal. Now recover results*/
+
+	/*Get dof listof this element (pattyn dofs and stokes dofs): */
+	GetDofList(&doflistp,PattynApproximationEnum);
+	GetDofList(&doflists,StokesApproximationEnum);
+
+	/*Get node data: */
+	GetVerticesCoordinates(&xyz_list[0][0], nodes, NUMVERTICES);
+
+	/*Use the dof list to index into the solution vector: */
+	for(i=0;i<numdofp;i++){
+		pattyn_values[i]=solution[doflistp[i]];
+	}
+	for(i=0;i<numdofs;i++){
+		stokes_values[i]=solution[doflists[i]];
+	}
+
+	/*Ok, we have vx and vy in values, fill in vx and vy arrays: */
+	for(i=0;i<NUMVERTICES;i++){
+		vx[i]=stokes_values[i*numdofpervertexs+0]+pattyn_values[i*numdofpervertexp+0];
+		vy[i]=stokes_values[i*numdofpervertexs+1]+pattyn_values[i*numdofpervertexp+1];
+		vz[i]=stokes_values[i*numdofpervertexs+2];
+		pressure[i]=stokes_values[i*numdofpervertexs+3];
+	}
+
+	/*Get Vz*/
+	Input* vz_input=inputs->GetInput(VzEnum);
+	if (vz_input){
+		if (vz_input->Enum()!=PentaVertexInputEnum){
+			ISSMERROR("Cannot compute Vel as Vz is of type %s",EnumToString(vz_input->Enum()));
+		}
+		vz_input->GetValuesPtr(&vz_ptr,&dummy);
+		for(i=0;i<NUMVERTICES;i++) vz[i]=vz[i]+vz_ptr[i];
+	}
+	else{
+		ISSMERROR("Cannot compute Vel there is no Vz input");
+	}
+
+	/*Now Compute vel*/
+	for(i=0;i<NUMVERTICES;i++) {
+		vel[i]=pow( pow(vx[i],2.0) + pow(vy[i],2.0) + pow(vz[i],2.0) , 0.5);
+	}
+
+	/*Now, we have to move the previous Vx and Vy inputs  to old 
+	 * status, otherwise, we'll wipe them off: */
+	this->inputs->ChangeEnum(VxEnum,VxOldEnum);
+	this->inputs->ChangeEnum(VyEnum,VyOldEnum);
+	this->inputs->ChangeEnum(VzEnum,VzOldEnum);
+	this->inputs->ChangeEnum(PressureEnum,PressureOldEnum);
+
+	/*Add vx and vy as inputs to the tria element: */
+	this->inputs->AddInput(new PentaVertexInput(VxEnum,vx));
+	this->inputs->AddInput(new PentaVertexInput(VyEnum,vy));
+	this->inputs->AddInput(new PentaVertexInput(VzEnum,vz));
+	this->inputs->AddInput(new PentaVertexInput(VelEnum,vel));
+	this->inputs->AddInput(new PentaVertexInput(PressureEnum,pressure));
+
+	/*Free ressources:*/
+	xfree((void**)&doflistp);
+	xfree((void**)&doflists);
+}
 /*}}}*/
 /*FUNCTION Penta::InputUpdateFromSolutionDiagnosticHutter {{{1*/
Index: /issm/trunk/src/c/objects/Elements/Penta.h
===================================================================
--- /issm/trunk/src/c/objects/Elements/Penta.h	(revision 5768)
+++ /issm/trunk/src/c/objects/Elements/Penta.h	(revision 5769)
@@ -168,4 +168,5 @@
 		void    InputUpdateFromSolutionDiagnosticMacAyealPattyn( double* solutiong);
 		void    InputUpdateFromSolutionDiagnosticPattyn( double* solutiong);
+		void    InputUpdateFromSolutionDiagnosticPattynStokes( double* solutiong);
 		void    InputUpdateFromSolutionDiagnosticHutter( double* solutiong);
 		void    InputUpdateFromSolutionDiagnosticVert( double* solutiong);
