Index: /issm/trunk-jpl/src/c/classes/Loads/Penpair.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Loads/Penpair.cpp	(revision 15793)
+++ /issm/trunk-jpl/src/c/classes/Loads/Penpair.cpp	(revision 15794)
@@ -303,5 +303,5 @@
 				case SSAApproximationEnum: return PenaltyCreateKMatrixStressbalanceSSAHO(kmax); 
 				case HOApproximationEnum:   return PenaltyCreateKMatrixStressbalanceSSAHO(kmax); 
-				default: _error_("not supported yet");
+				default: _error_("Approximation "<<EnumToStringx(approximation1)<<" not supported yet");
 			}
 		case HOApproximationEnum:
@@ -309,19 +309,18 @@
 				case SSAApproximationEnum: return PenaltyCreateKMatrixStressbalanceSSAHO(kmax); 
 				case HOApproximationEnum:   return PenaltyCreateKMatrixStressbalanceSSAHO(kmax); 
-				default: _error_("not supported yet");
+				default: _error_("Approximation "<<EnumToStringx(approximation1)<<" not supported yet");
 			}
-		case FSApproximationEnum:
+		case FSvelocityEnum:
 			switch(approximation1){
-				case FSApproximationEnum: return PenaltyCreateKMatrixStressbalanceFS(kmax); 
+				case FSvelocityEnum: return PenaltyCreateKMatrixStressbalanceFS(kmax); 
 				case NoneApproximationEnum: return   PenaltyCreateKMatrixStressbalanceFS(kmax); 
-				default: _error_("not supported yet");
+				default: _error_("Approximation "<<EnumToStringx(approximation1)<<" not supported yet");
 			}
 		case NoneApproximationEnum:
 			switch(approximation1){
-				case FSApproximationEnum: return PenaltyCreateKMatrixStressbalanceFS(kmax); 
+				case FSvelocityEnum: return PenaltyCreateKMatrixStressbalanceFS(kmax); 
 				case NoneApproximationEnum: return   PenaltyCreateKMatrixStressbalanceFS(kmax); 
-				default: _error_("not supported yet");
 			}
-		default: _error_("not supported yet");
+		default: _error_("Approximation "<<EnumToStringx(approximation0)<<" not supported yet");
 	}
 }
@@ -357,5 +356,5 @@
 ElementMatrix* Penpair::PenaltyCreateKMatrixStressbalanceFS(IssmDouble kmax){
 
-	const int numdof=NUMVERTICES*NDOF4;
+	const int  numdof=NUMVERTICES*NDOF3;
 	IssmDouble penalty_offset;
 
@@ -368,22 +367,17 @@
 	//Create elementary matrix: add penalty to 
 	Ke->values[0*numdof+0]=+kmax*pow(10.,penalty_offset);
-	Ke->values[0*numdof+4]=-kmax*pow(10.,penalty_offset);
-	Ke->values[4*numdof+0]=-kmax*pow(10.,penalty_offset);
+	Ke->values[0*numdof+3]=-kmax*pow(10.,penalty_offset);
+	Ke->values[3*numdof+0]=-kmax*pow(10.,penalty_offset);
+	Ke->values[3*numdof+3]=+kmax*pow(10.,penalty_offset);
+
+	Ke->values[1*numdof+1]=+kmax*pow(10.,penalty_offset);
+	Ke->values[1*numdof+4]=-kmax*pow(10.,penalty_offset);
+	Ke->values[4*numdof+1]=-kmax*pow(10.,penalty_offset);
 	Ke->values[4*numdof+4]=+kmax*pow(10.,penalty_offset);
 
-	Ke->values[1*numdof+1]=+kmax*pow(10.,penalty_offset);
-	Ke->values[1*numdof+5]=-kmax*pow(10.,penalty_offset);
-	Ke->values[5*numdof+1]=-kmax*pow(10.,penalty_offset);
+	Ke->values[2*numdof+2]=+kmax*pow(10.,penalty_offset);
+	Ke->values[2*numdof+5]=-kmax*pow(10.,penalty_offset);
+	Ke->values[5*numdof+2]=-kmax*pow(10.,penalty_offset);
 	Ke->values[5*numdof+5]=+kmax*pow(10.,penalty_offset);
-
-	Ke->values[2*numdof+2]=+kmax*pow(10.,penalty_offset);
-	Ke->values[2*numdof+6]=-kmax*pow(10.,penalty_offset);
-	Ke->values[6*numdof+2]=-kmax*pow(10.,penalty_offset);
-	Ke->values[6*numdof+6]=+kmax*pow(10.,penalty_offset);
-
-	Ke->values[3*numdof+3]=+kmax*pow(10.,penalty_offset);
-	Ke->values[3*numdof+7]=-kmax*pow(10.,penalty_offset);
-	Ke->values[7*numdof+3]=-kmax*pow(10.,penalty_offset);
-	Ke->values[7*numdof+7]=+kmax*pow(10.,penalty_offset);
 
 	/*Clean up and return*/