Index: /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 24116)
+++ /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 24117)
@@ -694,4 +694,7 @@
 
 	/*Intermediaries */
+	int melt_style, point1;
+	bool mainlyfloating;
+	IssmDouble  fraction1,fraction2,gllevelset;
 	IssmDouble  Jdet,dt;
 	IssmDouble  ms,mb,gmb,fmb,thickness,phi=1.;
@@ -708,4 +711,5 @@
 	element->GetVerticesCoordinates(&xyz_list);
 	element->FindParam(&dt,TimesteppingTimeStepEnum);
+	element->FindParam(&melt_style,GroundinglineMeltInterpolationEnum);
 	Input* gmb_input        = element->GetInput(BasalforcingsGroundediceMeltingRateEnum); _assert_(gmb_input);
 	Input* fmb_input        = element->GetInput(BasalforcingsFloatingiceMeltingRateEnum); _assert_(fmb_input);
@@ -714,6 +718,16 @@
 	Input* thickness_input  = element->GetInput(ThicknessEnum);                           _assert_(thickness_input);
 
+   /*Recover portion of element that is grounded*/
+   Gauss* gauss=NULL;
+   phi=element->GetGroundedPortion(xyz_list);
+   if(melt_style==SubelementMelt2Enum){
+      element->GetGroundedPart(&point1,&fraction1,&fraction2,&mainlyfloating);
+      gauss = element->NewGauss(point1,fraction1,fraction2,3);
+   }
+   else{
+      gauss = element->NewGauss(3);
+   }
+
 	/* Start  looping on the number of gaussian points: */
-	Gauss* gauss=element->NewGauss(2);
 	for(int ig=gauss->begin();ig<gauss->end();ig++){
 		gauss->GaussPoint(ig);
@@ -725,9 +739,24 @@
 		gmb_input->GetInputValue(&gmb,gauss);
 		fmb_input->GetInputValue(&fmb,gauss);
-		gllevelset_input->GetInputValue(&phi,gauss);
+		gllevelset_input->GetInputValue(&gllevelset,gauss);
 		thickness_input->GetInputValue(&thickness,gauss);
 
-		if(phi>0.) mb=gmb;
-		else mb=fmb;
+ 		if(melt_style==SubelementMelt1Enum){
+         if (phi>0.999999999) mb=gmb;
+         else mb=(1-phi)*fmb+phi*gmb; // phi is the fraction of grounded ice so (1-phi) is floating
+      }
+      else if(melt_style==SubelementMelt2Enum){
+         if(gllevelset>0.) mb=gmb;
+         else mb=fmb;
+      }
+      else if(melt_style==NoMeltOnPartiallyFloatingEnum){
+         if (phi<0.00000001) mb=fmb;
+         else mb=gmb;
+      }
+      else if(melt_style==FullMeltOnPartiallyFloatingEnum){
+         if (phi<0.99999999) mb=fmb;
+         else mb=gmb;
+      }
+      else  _error_("melt interpolation "<<EnumToStringx(melt_style)<<" not implemented yet");
 
 		for(int i=0;i<numnodes;i++) pe->values[i]+=Jdet*gauss->weight*(thickness+dt*(ms-mb))*basis[i];
