Index: /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/sealevelrise_core.cpp
===================================================================
--- /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/sealevelrise_core.cpp	(revision 21316)
+++ /issm/branches/trunk-larour-NatGeoScience2016/src/c/cores/sealevelrise_core.cpp	(revision 21317)
@@ -80,5 +80,5 @@
 		femmodel->parameters->FindParam(&earthid,EarthIdEnum);
 		femmodel->parameters->FindParam(&ismasstransport,TransientIsmasstransportEnum);
-		if((modelid==EarthIdEnum) & !ismasstransport) EarthMassTransport(femmodel);
+		if((modelid==earthid) && (ismasstransport==0)) EarthMassTransport(femmodel);
 	}
 
@@ -96,23 +96,31 @@
 		InputUpdateFromVectorx(femmodel,Sg,SealevelEnum,VertexSIdEnum);
 
-		/*compute other geodetic signatures, such as absolute sea level chagne, components of 3-D crustal motion: */
-		/*Initialize:*/
-		U_radial = new Vector<IssmDouble>(gsize);
-		U_north = new Vector<IssmDouble>(gsize);
-		U_east = new Vector<IssmDouble>(gsize);
-		Sg_absolute = new Vector<IssmDouble>(gsize); 
-		
-		/*call the geodetic main modlule:*/ 
-		femmodel->SealevelriseGeodetic(U_radial,U_north,U_east,Sg,latitude,longitude,radius,xx,yy,zz); 
-
-		/*compute: absolute sea level change = relative sea level change + vertical motion*/
-		Sg->Copy(Sg_absolute); Sg_absolute->AXPY(U_radial,1); 
-		
-		/*get results into elements:*/
-		InputUpdateFromVectorx(femmodel,U_radial,SealevelUmotionEnum,VertexSIdEnum);	// radial displacement 
-		InputUpdateFromVectorx(femmodel,U_north,SealevelNmotionEnum,VertexSIdEnum);	// north motion 
-		InputUpdateFromVectorx(femmodel,U_east,SealevelEmotionEnum,VertexSIdEnum);		// east motion 
-		InputUpdateFromVectorx(femmodel,Sg_absolute,SealevelAbsoluteEnum,VertexSIdEnum);
-		
+		int geodetic=0;
+		if (geodetic){
+
+			/*compute other geodetic signatures, such as absolute sea level chagne, components of 3-D crustal motion: */
+			/*Initialize:*/
+			U_radial = new Vector<IssmDouble>(gsize);
+			U_north = new Vector<IssmDouble>(gsize);
+			U_east = new Vector<IssmDouble>(gsize);
+			Sg_absolute = new Vector<IssmDouble>(gsize); 
+
+			/*call the geodetic main modlule:*/ 
+			femmodel->SealevelriseGeodetic(U_radial,U_north,U_east,Sg,latitude,longitude,radius,xx,yy,zz); 
+
+			/*compute: absolute sea level change = relative sea level change + vertical motion*/
+			Sg->Copy(Sg_absolute); Sg_absolute->AXPY(U_radial,1); 
+
+			/*get results into elements:*/
+			InputUpdateFromVectorx(femmodel,U_radial,SealevelUmotionEnum,VertexSIdEnum);	// radial displacement 
+			InputUpdateFromVectorx(femmodel,U_north,SealevelNmotionEnum,VertexSIdEnum);	// north motion 
+			InputUpdateFromVectorx(femmodel,U_east,SealevelEmotionEnum,VertexSIdEnum);		// east motion 
+			InputUpdateFromVectorx(femmodel,Sg_absolute,SealevelAbsoluteEnum,VertexSIdEnum);
+			delete U_radial;
+			delete U_north;
+			delete U_east;
+			delete Sg_absolute;
+		}	
+
 		if(save_results){
 			if(VerboseSolution()) _printf0_("   saving results\n");
@@ -120,5 +128,5 @@
 			femmodel->RequestedOutputsx(&femmodel->results,requested_outputs,numoutputs);
 		}
-			
+
 		if(solution_type==SealevelriseSolutionEnum)femmodel->RequestedDependentsx();
 
@@ -126,8 +134,4 @@
 		delete Sg;
 		delete Sg_eustatic;
-		delete U_radial;
-		delete U_north;
-		delete U_east;
-		delete Sg_absolute;
 		if(numoutputs){for(int i=0;i<numoutputs;i++){xDelete<char>(requested_outputs[i]);} xDelete<char*>(requested_outputs);}
 	}
@@ -209,5 +213,5 @@
 		 *First, build the global delta thickness vector in the earth model: */
 		nv=femmodel->vertices->NumberOfVertices();
-		forcingglobal= new Vector<IssmDouble>(nv);
+		GetVectorFromInputsx(&forcingglobal,femmodel,forcingenum,VertexSIdEnum);
 
 		/*Retrieve transition vectors, used to plug from each ice cap into the global forcing:*/
@@ -377,4 +381,6 @@
 	Vector<IssmDouble> *deltathickness    = NULL;
 	int nv;
+	
+	if(VerboseSolution()) _printf0_("      computing earth mass transport\n");
 
 	/*No mass transport module was called, so we are just going to retrieve the geometry thickness 
@@ -390,5 +396,5 @@
 	deltathickness = new Vector<IssmDouble>(nv); 
 	newthickness->Copy(deltathickness); deltathickness->AXPY(oldthickness,-1); 
-		
+
 	/*plug into elements:*/
 	InputUpdateFromVectorx(femmodel,deltathickness,SealevelriseDeltathicknessEnum,VertexSIdEnum);