Index: /issm/trunk/src/c/objects/Pengrid.cpp
===================================================================
--- /issm/trunk/src/c/objects/Pengrid.cpp	(revision 925)
+++ /issm/trunk/src/c/objects/Pengrid.cpp	(revision 926)
@@ -309,4 +309,5 @@
 void  Pengrid::PenaltyCreateKMatrixMelting(Mat Kgg,void* vinputs,double kmax,int analysis_type,int sub_analysis_type){
 
+
 	int found=0;
 	const int numgrids=1;
@@ -325,4 +326,7 @@
 	ParameterInputs* inputs=NULL;
 
+	/*check that pengrid is not a clone (penalty to be added only once)*/
+	if (node->IsClone()) return;
+
 	/*recover pointers: */
 	inputs=(ParameterInputs*)vinputs;
@@ -345,4 +349,9 @@
 	if (temperature<t_pmp){ //If T<Tpmp, there must be no melting. Therefore, melting should be  constrained to 0 when T<Tpmp, instead of using spcs, use penalties
 		Ke[0][0]=kmax*pow(10,penalty_offset);
+		if (id==65){
+			printf("Ke_terms[0] = %10.10lf\n",Ke[0]);
+			printf("penalty_offset(K) = %g\n",penalty_offset);
+			printf("kmax = %10.10lf\n",kmax);
+		}
 	}
 	
@@ -458,4 +467,7 @@
 	ParameterInputs* inputs=NULL;
 
+	/*check that pengrid is not a clone (penalty to be added only once)*/
+	if (node->IsClone()) return;
+
 	/*recover pointers: */
 	inputs=(ParameterInputs*)vinputs;
@@ -495,4 +507,9 @@
 		if (sub_analysis_type==SteadyAnalysisEnum()){
 			P_terms[0]=melting_offset*pow(10,penalty_offset)*(temperature-t_pmp);
+			if (id==65){
+				printf("penalty_offset(P) = %g\n",penalty_offset);
+				printf("melting_offset(P) = %10.10lf\n",melting_offset);
+				printf("P_terms[0] = %10.10lf\n",P_terms[0]);
+			}
 		}
 		else{
@@ -553,7 +570,13 @@
 	ParameterInputs* inputs=NULL;
 
+	/*check that pengrid is not a clone (penalty to be added only once)*/
+	if (node->IsClone()){
+		unstable=0;
+		*punstable=unstable;
+		return;
+	}
+
 	/*recover pointers: */
 	inputs=(ParameterInputs*)vinputs;
-
 
 	//First recover beta, pressure and temperature vectors;