Index: /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp	(revision 18033)
+++ /issm/trunk-jpl/src/c/analyses/DamageEvolutionAnalysis.cpp	(revision 18034)
@@ -88,7 +88,5 @@
 void DamageEvolutionAnalysis::CreateConstraints(Constraints* constraints,IoModel* iomodel){/*{{{*/
 
-	int stabilization;
 	int finiteelement;
-	iomodel->Constant(&stabilization,DamageStabilizationEnum);
 	iomodel->Constant(&finiteelement,DamageElementinterpEnum);
 
@@ -450,6 +448,6 @@
 
 	/*Fetch number of vertices and allocate output*/
-	int numvertices = element->GetNumberOfVertices();
-	IssmDouble* f   = xNew<IssmDouble>(numvertices);
+	int numnodes = element->GetNumberOfNodes();
+	IssmDouble* f   = xNew<IssmDouble>(numnodes);
 
 	/*retrieve parameters:*/
@@ -481,6 +479,6 @@
 	/*Calculate damage evolution source term: */
 	Gauss* gauss=element->NewGauss();
-	for (int iv=0;iv<numvertices;iv++){
-		gauss->GaussVertex(iv);
+	for (int i=0;i<numnodes;i++){
+		gauss->GaussNode(element->GetElementType(),i);
 		
 		damage_input->GetInputValue(&damage,gauss);
@@ -509,5 +507,5 @@
 		NegPsi=max(-Chi,0.); /* healing only for compressive stresses */
 
-		f[iv]= c1*(pow(PosPsi,c2) - healing*pow(NegPsi,c2))*pow((1.-damage),-c3);
+		f[i]= c1*(pow(PosPsi,c2) - healing*pow(NegPsi,c2))*pow((1.-damage),-c3);
 	}