Index: /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 21793)
+++ /issm/trunk-jpl/src/c/analyses/MasstransportAnalysis.cpp	(revision 21794)
@@ -470,5 +470,6 @@
 	/*Intermediaries */
 	IssmDouble  Jdet,dt;
-	IssmDouble  ms,mb,gmb,fmb,thickness,phi;
+	IssmDouble  ms,mb,gmb,fmb,thickness;
+	IssmDouble  gllevelset,phi=1.;
 	IssmDouble* xyz_list = NULL;
 
@@ -485,5 +486,5 @@
 	Input* gmb_input           = element->GetInput(BasalforcingsGroundediceMeltingRateEnum);  _assert_(gmb_input);
 	Input* fmb_input           = element->GetInput(BasalforcingsFloatingiceMeltingRateEnum);  _assert_(fmb_input);
-	Input* groundedice_input   = element->GetInput(MaskGroundediceLevelsetEnum);              _assert_(groundedice_input);
+	Input* gllevelset_input   = element->GetInput(MaskGroundediceLevelsetEnum);              _assert_(gllevelset_input);
 	Input* ms_input            = element->GetInput(SmbMassBalanceEnum);                       _assert_(ms_input);
 	Input* thickness_input     = element->GetInput(ThicknessEnum);                            _assert_(thickness_input);
@@ -500,5 +501,5 @@
 		gmb_input->GetInputValue(&gmb,gauss);
 		fmb_input->GetInputValue(&fmb,gauss);
-		groundedice_input->GetInputValue(&phi,gauss);
+		gllevelset_input->GetInputValue(&phi,gauss);
 		thickness_input->GetInputValue(&thickness,gauss);
 		if(phi>0.) mb=gmb;
@@ -520,7 +521,11 @@
 
 	/*Intermediaries */
-	IssmDouble  Jdet,dt;
-	IssmDouble  ms,mb,gmb,fmb,thickness,phi;
+	int         migration_style,point1;
+	bool        mainlyfloating;
+	IssmDouble  Jdet,fraction1,fraction2;
+	IssmDouble  dt,gllevelset;
+	IssmDouble  ms,mb,gmb,fmb,thickness,phi=1.;
 	IssmDouble* xyz_list = NULL;
+	Gauss*      gauss     = NULL;
 
 	/*Fetch number of nodes and dof for this finite element*/
@@ -537,9 +542,19 @@
 	Input* fmb_input           = element->GetInput(BasalforcingsFloatingiceMeltingRateEnum); _assert_(fmb_input);
 	Input* ms_input            = element->GetInput(SmbMassBalanceEnum);          _assert_(ms_input);
-	Input* groundedice_input   = element->GetInput(MaskGroundediceLevelsetEnum);             _assert_(groundedice_input);
+	Input* gllevelset_input   = element->GetInput(MaskGroundediceLevelsetEnum);             _assert_(gllevelset_input);
 	Input* thickness_input     = element->GetInput(ThicknessEnum);                           _assert_(thickness_input);
 
+	/*Recover portion of element that is grounded*/
+	phi=element->GetGroundedPortion(xyz_list);
+	if(migration_style==SubelementMigration2Enum){
+		gllevelset_input=element->GetInput(MaskGroundediceLevelsetEnum); _assert_(gllevelset_input);
+		element->GetGroundedPart(&point1,&fraction1,&fraction2,&mainlyfloating);
+	   gauss = element->NewGauss(point1,fraction1,fraction2,mainlyfloating,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);
@@ -551,8 +566,22 @@
 		gmb_input->GetInputValue(&gmb,gauss);
 		fmb_input->GetInputValue(&fmb,gauss);
-		groundedice_input->GetInputValue(&phi,gauss);
+		gllevelset_input->GetInputValue(&gllevelset,gauss);
 		thickness_input->GetInputValue(&thickness,gauss);
-		if(phi>0) mb=gmb;
-		else mb=fmb;
+		if(migration_style==SubelementMigrationEnum){
+			if (phi<0.00000001) mb=gmb;
+			else mb=(1-phi)*fmb+phi*gmb; // phi is the fraction of grounded ice so (1-phi) is floating
+		}
+		else if(migration_style==SubelementMigration2Enum){
+			if(gllevelset>0.) mb=gmb;
+			else mb=fmb;
+		}
+		else if(migration_style==SubelementMigration3Enum){
+			if (phi<0.00000001) mb=gmb;
+			else mb=gmb;
+		}
+		else{
+			if(gllevelset>0) mb=gmb;
+			else mb=fmb;
+		}
 
 		for(int i=0;i<numnodes;i++) pe->values[i]+=Jdet*gauss->weight*(thickness+dt*(ms-mb))*basis[i];