Index: /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp	(revision 18815)
+++ /issm/trunk-jpl/src/c/analyses/LevelsetAnalysis.cpp	(revision 18816)
@@ -48,7 +48,9 @@
 			case DefaultCalvingEnum:
 				iomodel->FetchDataToInput(elements,CalvingCalvingrateEnum);
+				iomodel->FetchDataToInput(elements,CalvingMeltingrateEnum);
 				break;
 			case CalvingLevermannEnum:
 				iomodel->FetchDataToInput(elements,CalvinglevermannCoeffEnum);
+				iomodel->FetchDataToInput(elements,CalvinglevermannMeltingrateEnum);
 				break;
 			default:
@@ -125,5 +127,5 @@
 	IssmDouble h,hx,hy,hz;
 	IssmDouble vel;
-	IssmDouble norm_dlsf, calvingrate;
+	IssmDouble norm_dlsf, norm_calving, calvingrate, meltingrate;
 	IssmDouble* xyz_list = NULL;
 
@@ -151,4 +153,5 @@
 	IssmDouble*    w        = xNew<IssmDouble>(dim);
 	IssmDouble*    c        = xNewZeroInit<IssmDouble>(dim);
+	IssmDouble*    m        = xNewZeroInit<IssmDouble>(dim);
 	IssmDouble*    dlsf     = xNew<IssmDouble>(dim);
 
@@ -163,4 +166,5 @@
 	Input* lsf_slopey_input   = NULL;
 	Input* calvingrate_input  = NULL;
+	Input* meltingrate_input  = NULL;
 
 	/*Load velocities*/
@@ -186,4 +190,5 @@
 				if(dim==2) lsf_slopey_input  = basalelement->GetInput(LevelsetfunctionSlopeYEnum); _assert_(lsf_slopey_input);
 				calvingrate_input = basalelement->GetInput(CalvingCalvingrateEnum);     _assert_(calvingrate_input);
+				meltingrate_input = basalelement->GetInput(CalvingMeltingrateEnum);     _assert_(meltingrate_input);
 				break;
 			case CalvingLevermannEnum:
@@ -201,4 +206,5 @@
 					}
 				}
+				meltingrate_input = basalelement->GetInput(CalvinglevermannMeltingrateEnum);     _assert_(meltingrate_input);
 				break;
 			default:
@@ -238,4 +244,5 @@
 					if(dim==2) lsf_slopey_input->GetInputValue(&dlsf[1],gauss);
 					calvingrate_input->GetInputValue(&calvingrate,gauss);
+					meltingrate_input->GetInputValue(&meltingrate,gauss);
 
 					norm_dlsf=0.;
@@ -244,11 +251,26 @@
 
 					if(norm_dlsf>1.e-10)
-					 for(i=0;i<dim;i++) c[i]=calvingrate*dlsf[i]/norm_dlsf;
+					 for(i=0;i<dim;i++){ 
+					  c[i]=calvingrate*dlsf[i]/norm_dlsf;
+					  m[i]=meltingrate*dlsf[i]/norm_dlsf;
+					 }
 					else
-					 for(i=0;i<dim;i++) c[i]=0.;
+					 for(i=0;i<dim;i++){
+						 c[i]=0.;
+						 m[i]=0.;
+					 }
+
 					break;
 				case CalvingLevermannEnum:
 					calvingratex_input->GetInputValue(&c[0],gauss); // in 3D case, add mesh velocity 
 					if(dim==2) calvingratey_input->GetInputValue(&c[1],gauss);
+					meltingrate_input->GetInputValue(&meltingrate,gauss);
+
+					norm_calving=0.;
+					for(i=0;i<dim;i++) norm_calving+=pow(c[i],2);
+					norm_calving=sqrt(norm_calving)+1.e-14;
+
+					for(i=0;i<dim;i++) m[i]=meltingrate*c[i]/norm_calving;
+					
 					break;
 				default:
@@ -258,5 +280,5 @@
 
 		/*Levelset speed is ice velocity - calving rate*/
-		for(i=0;i<dim;i++) w[i]=v[i]-c[i];
+		for(i=0;i<dim;i++) w[i]=v[i]-c[i]-m[i];
 
 		/*Compute D*/
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 18815)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 18816)
@@ -213,8 +213,10 @@
 	CalvingLawEnum,
 	CalvingCalvingrateEnum,
+	CalvingMeltingrateEnum,
 	CalvingLevermannEnum,
 	DefaultCalvingEnum,
 	CalvingRequestedOutputsEnum,
 	CalvinglevermannCoeffEnum,
+	CalvinglevermannMeltingrateEnum,
 	CalvingratexEnum,
 	CalvingrateyEnum,
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 18815)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 18816)
@@ -221,8 +221,10 @@
 		case CalvingLawEnum : return "CalvingLaw";
 		case CalvingCalvingrateEnum : return "CalvingCalvingrate";
+		case CalvingMeltingrateEnum : return "CalvingMeltingrate";
 		case CalvingLevermannEnum : return "CalvingLevermann";
 		case DefaultCalvingEnum : return "DefaultCalving";
 		case CalvingRequestedOutputsEnum : return "CalvingRequestedOutputs";
 		case CalvinglevermannCoeffEnum : return "CalvinglevermannCoeff";
+		case CalvinglevermannMeltingrateEnum : return "CalvinglevermannMeltingrate";
 		case CalvingratexEnum : return "Calvingratex";
 		case CalvingrateyEnum : return "Calvingratey";
Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 18815)
+++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 18816)
@@ -224,8 +224,10 @@
 	      else if (strcmp(name,"CalvingLaw")==0) return CalvingLawEnum;
 	      else if (strcmp(name,"CalvingCalvingrate")==0) return CalvingCalvingrateEnum;
+	      else if (strcmp(name,"CalvingMeltingrate")==0) return CalvingMeltingrateEnum;
 	      else if (strcmp(name,"CalvingLevermann")==0) return CalvingLevermannEnum;
 	      else if (strcmp(name,"DefaultCalving")==0) return DefaultCalvingEnum;
 	      else if (strcmp(name,"CalvingRequestedOutputs")==0) return CalvingRequestedOutputsEnum;
 	      else if (strcmp(name,"CalvinglevermannCoeff")==0) return CalvinglevermannCoeffEnum;
+	      else if (strcmp(name,"CalvinglevermannMeltingrate")==0) return CalvinglevermannMeltingrateEnum;
 	      else if (strcmp(name,"Calvingratex")==0) return CalvingratexEnum;
 	      else if (strcmp(name,"Calvingratey")==0) return CalvingrateyEnum;
@@ -258,10 +260,10 @@
 	      else if (strcmp(name,"MeshVertexonsurface")==0) return MeshVertexonsurfaceEnum;
 	      else if (strcmp(name,"MeshVertexonboundary")==0) return MeshVertexonboundaryEnum;
-	      else if (strcmp(name,"MeshX")==0) return MeshXEnum;
-	      else if (strcmp(name,"MeshY")==0) return MeshYEnum;
          else stage=3;
    }
    if(stage==3){
-	      if (strcmp(name,"MeshZ")==0) return MeshZEnum;
+	      if (strcmp(name,"MeshX")==0) return MeshXEnum;
+	      else if (strcmp(name,"MeshY")==0) return MeshYEnum;
+	      else if (strcmp(name,"MeshZ")==0) return MeshZEnum;
 	      else if (strcmp(name,"MeshElementtype")==0) return MeshElementtypeEnum;
 	      else if (strcmp(name,"DomainType")==0) return DomainTypeEnum;
@@ -381,10 +383,10 @@
 	      else if (strcmp(name,"SolutionType")==0) return SolutionTypeEnum;
 	      else if (strcmp(name,"AnalysisType")==0) return AnalysisTypeEnum;
-	      else if (strcmp(name,"ConfigurationType")==0) return ConfigurationTypeEnum;
-	      else if (strcmp(name,"AdjointBalancethicknessAnalysis")==0) return AdjointBalancethicknessAnalysisEnum;
          else stage=4;
    }
    if(stage==4){
-	      if (strcmp(name,"AdjointBalancethickness2Analysis")==0) return AdjointBalancethickness2AnalysisEnum;
+	      if (strcmp(name,"ConfigurationType")==0) return ConfigurationTypeEnum;
+	      else if (strcmp(name,"AdjointBalancethicknessAnalysis")==0) return AdjointBalancethicknessAnalysisEnum;
+	      else if (strcmp(name,"AdjointBalancethickness2Analysis")==0) return AdjointBalancethickness2AnalysisEnum;
 	      else if (strcmp(name,"AdjointHorizAnalysis")==0) return AdjointHorizAnalysisEnum;
 	      else if (strcmp(name,"AnalysisCounter")==0) return AnalysisCounterEnum;
@@ -504,10 +506,10 @@
 	      else if (strcmp(name,"MassconaxpbyNamex")==0) return MassconaxpbyNamexEnum;
 	      else if (strcmp(name,"MassconaxpbyNamey")==0) return MassconaxpbyNameyEnum;
-	      else if (strcmp(name,"MassconaxpbyAlpha")==0) return MassconaxpbyAlphaEnum;
-	      else if (strcmp(name,"MassconaxpbyBeta")==0) return MassconaxpbyBetaEnum;
          else stage=5;
    }
    if(stage==5){
-	      if (strcmp(name,"NodeSId")==0) return NodeSIdEnum;
+	      if (strcmp(name,"MassconaxpbyAlpha")==0) return MassconaxpbyAlphaEnum;
+	      else if (strcmp(name,"MassconaxpbyBeta")==0) return MassconaxpbyBetaEnum;
+	      else if (strcmp(name,"NodeSId")==0) return NodeSIdEnum;
 	      else if (strcmp(name,"VectorParam")==0) return VectorParamEnum;
 	      else if (strcmp(name,"Riftfront")==0) return RiftfrontEnum;
@@ -627,10 +629,10 @@
 	      else if (strcmp(name,"StrainRateyz")==0) return StrainRateyzEnum;
 	      else if (strcmp(name,"StrainRatezz")==0) return StrainRatezzEnum;
-	      else if (strcmp(name,"Divergence")==0) return DivergenceEnum;
-	      else if (strcmp(name,"MaxDivergence")==0) return MaxDivergenceEnum;
          else stage=6;
    }
    if(stage==6){
-	      if (strcmp(name,"GiaCrossSectionShape")==0) return GiaCrossSectionShapeEnum;
+	      if (strcmp(name,"Divergence")==0) return DivergenceEnum;
+	      else if (strcmp(name,"MaxDivergence")==0) return MaxDivergenceEnum;
+	      else if (strcmp(name,"GiaCrossSectionShape")==0) return GiaCrossSectionShapeEnum;
 	      else if (strcmp(name,"GiadWdt")==0) return GiadWdtEnum;
 	      else if (strcmp(name,"GiaW")==0) return GiaWEnum;
@@ -750,10 +752,10 @@
 	      else if (strcmp(name,"MpiDense")==0) return MpiDenseEnum;
 	      else if (strcmp(name,"MpiSparse")==0) return MpiSparseEnum;
-	      else if (strcmp(name,"Seq")==0) return SeqEnum;
-	      else if (strcmp(name,"Mpi")==0) return MpiEnum;
          else stage=7;
    }
    if(stage==7){
-	      if (strcmp(name,"Mumps")==0) return MumpsEnum;
+	      if (strcmp(name,"Seq")==0) return SeqEnum;
+	      else if (strcmp(name,"Mpi")==0) return MpiEnum;
+	      else if (strcmp(name,"Mumps")==0) return MumpsEnum;
 	      else if (strcmp(name,"Gsl")==0) return GslEnum;
 	      else if (strcmp(name,"Option")==0) return OptionEnum;
Index: /issm/trunk-jpl/src/m/classes/calving.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/calving.m	(revision 18815)
+++ /issm/trunk-jpl/src/m/classes/calving.m	(revision 18816)
@@ -7,4 +7,5 @@
 	properties (SetAccess=public) 
 		 calvingrate            = NaN;
+		 meltingrate            = NaN;
 	end
 	methods (Static)
@@ -60,8 +61,10 @@
 			if (solution~=TransientSolutionEnum() | md.transient.iscalving==0), return; end
 			md = checkfield(md,'fieldname','calving.calvingrate','NaN',1,'size',[md.mesh.numberofvertices 1],'>=',0);
+			md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices 1],'>=',0);
 		end % }}}
 		function disp(obj) % {{{
 			disp(sprintf('   Calving parameters:'));
 			fielddisplay(obj,'calvingrate','calving rate at given location [m/a]');
+			fielddisplay(obj,'meltingrate','melting rate at given location [m/a]');
 		end % }}}
 		function marshall(obj,md,fid) % {{{
@@ -69,4 +72,5 @@
 			WriteData(fid,'enum',CalvingLawEnum(),'data',DefaultCalvingEnum(),'format','Integer');
 			WriteData(fid,'object',obj,'fieldname','calvingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
+			WriteData(fid,'object',obj,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
 		end % }}}
 	end
Index: /issm/trunk-jpl/src/m/classes/calving.py
===================================================================
--- /issm/trunk-jpl/src/m/classes/calving.py	(revision 18815)
+++ /issm/trunk-jpl/src/m/classes/calving.py	(revision 18816)
@@ -15,4 +15,5 @@
 	def __init__(self): # {{{
 		self.calvingrate            = float('NaN')
+		self.meltingrate            = float('NaN')
 
 		#set defaults
@@ -23,4 +24,5 @@
 		string='   Calving parameters:'
 		string="%s\n%s"%(string,fielddisplay(self,'calvingrate','calving rate at given location [m/a]'))
+		string="%s\n%s"%(string,fielddisplay(self,'meltingrate','melting rate at given location [m/a]'))
 
 		return string
@@ -37,4 +39,5 @@
 
 		md = checkfield(md,'fieldname','calving.calvingrate','NaN',1,'size',[md.mesh.numberofvertices],'>=',0)
+		md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices],'>=',0)
 		return md
 	# }}}
@@ -45,3 +48,4 @@
 		WriteData(fid,'enum',CalvingLawEnum(),'data',DefaultCalvingEnum(),'format','Integer');
 		WriteData(fid,'object',self,'fieldname','calvingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts)
+		WriteData(fid,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts)
 	# }}}
Index: /issm/trunk-jpl/src/m/classes/calvinglevermann.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/calvinglevermann.m	(revision 18815)
+++ /issm/trunk-jpl/src/m/classes/calvinglevermann.m	(revision 18816)
@@ -6,5 +6,6 @@
 classdef calvinglevermann
 	properties (SetAccess=public) 
-		 coeff = NaN;
+		 coeff       = NaN;
+		 meltingrate = NaN;
 	end
 	methods
@@ -35,13 +36,17 @@
 			if (solution~=TransientSolutionEnum() | md.transient.iscalving==0), return; end
 			md = checkfield(md,'fieldname','calving.coeff','>',0,'size',[md.mesh.numberofvertices 1]);
+			md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices 1],'>=',0);
 		end % }}}
 		function disp(obj) % {{{
 			disp(sprintf('   Calving Levermann parameters:'));
 			fielddisplay(obj,'coeff','proportionality coefficient in Levermann model');
+			fielddisplay(obj,'meltingrate','melting rate at given location [m/a]');
 
 		end % }}}
 		function marshall(obj,md,fid) % {{{
+			yts=365.0*24.0*3600.0;
 			WriteData(fid,'enum',CalvingLawEnum(),'data',CalvingLevermannEnum(),'format','Integer');
 			WriteData(fid,'enum',CalvinglevermannCoeffEnum(),'data',obj.coeff,'format','DoubleMat','mattype',1);
+			WriteData(fid,'object',obj,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts);
 		end % }}}
 	end
Index: /issm/trunk-jpl/src/m/classes/calvinglevermann.py
===================================================================
--- /issm/trunk-jpl/src/m/classes/calvinglevermann.py	(revision 18815)
+++ /issm/trunk-jpl/src/m/classes/calvinglevermann.py	(revision 18816)
@@ -14,5 +14,6 @@
 
 	def __init__(self): # {{{
-		self.coeff = float('NaN')
+		self.coeff       = float('NaN')
+		self.meltingrate = float('NaN')
 
 		#set defaults
@@ -23,4 +24,5 @@
 		string='   Calving Levermann parameters:'
 		string="%s\n%s"%(string,fielddisplay(self,'coeff','proportionality coefficient in Levermann model'))
+		string="%s\n%s"%(string,fielddisplay(self,'meltingrate','melting rate at given location [m/a]'))
 
 		return string
@@ -38,8 +40,11 @@
 
 		md = checkfield(md,'fieldname','calving.coeff','size',[md.mesh.numberofvertices],'>',0)
+		md = checkfield(md,'fieldname','calving.meltingrate','NaN',1,'size',[md.mesh.numberofvertices],'>=',0)
 		return md
 	# }}}
 	def marshall(self,md,fid):    # {{{
+		yts=365.*24.*3600.
 		WriteData(fid,'enum',CalvingLawEnum(),'data',CalvingLevermannEnum(),'format','Integer');
 		WriteData(fid,'enum',CalvinglevermannCoeffEnum(),'data',self.coeff,'format','DoubleMat','mattype',1)
+		WriteData(fid,'object',self,'fieldname','meltingrate','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1,'scale',1./yts)
 	# }}}
Index: /issm/trunk-jpl/test/NightlyRun/IdToName.m
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/IdToName.m	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/IdToName.m	(revision 18816)
@@ -195,5 +195,6 @@
 	case 804, name='ValleyGlacierLevelsetCalvingSIA2d';
 	case 805, name='ValleyGlacierLevelsetEnthCalvingHO3d';
-	case 806, name='SquareShelfLevelsetCalvingSIA2dLevermann';
+	case 806, name='SquareShelfLevelsetCalvingSSA2dLevermann';
+	case 807, name='SquareShelfLevelsetMeltingSSA2d';
 	case 1101, name='ISMIPAHO';
 	case 1102, name='ISMIPAFS';
Index: /issm/trunk-jpl/test/NightlyRun/IdToName.py
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/IdToName.py	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/IdToName.py	(revision 18816)
@@ -180,5 +180,6 @@
 	804  : 'ValleyGlacierLevelsetCalvingSIA2d',	
 	805  : 'ValleyGlacierLevelsetEnthCalvingHO3d',	
-	806  : 'SquareShelfLevelsetCalvingSIA2dLevermann',	
+	806  : 'SquareShelfLevelsetCalvingSSA2dLevermann',	
+	807  : 'SquareShelfLevelsetMeltingSSA2d',
 	1101 : 'ISMIPAHO',
 	1102 : 'ISMIPAFS',
Index: /issm/trunk-jpl/test/NightlyRun/test804.m
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test804.m	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/test804.m	(revision 18816)
@@ -15,4 +15,5 @@
 
 md.calving.calvingrate=1000.*ones(md.mesh.numberofvertices,1);
+md.calving.meltingrate=zeros(md.mesh.numberofvertices,1);
 
 md=solve(md,TransientSolutionEnum());
Index: /issm/trunk-jpl/test/NightlyRun/test804.py
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test804.py	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/test804.py	(revision 18816)
@@ -25,4 +25,5 @@
 
 md.calving.calvingrate=1000.*numpy.ones((md.mesh.numberofvertices,1))
+md.calving.meltingrate=numpy.zeros((md.mesh.numberofvertices,1))
 
 md=solve(md,TransientSolutionEnum())
Index: /issm/trunk-jpl/test/NightlyRun/test805.m
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test805.m	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/test805.m	(revision 18816)
@@ -22,4 +22,5 @@
 
 md.calving.calvingrate=1000.*ones(md.mesh.numberofvertices,1);
+md.calving.meltingrate=zeros(md.mesh.numberofvertices,1);
 
 md=solve(md,TransientSolutionEnum());
Index: /issm/trunk-jpl/test/NightlyRun/test805.py
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test805.py	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/test805.py	(revision 18816)
@@ -32,4 +32,5 @@
 
 md.calving.calvingrate=1000.*numpy.ones((md.mesh.numberofvertices,1))
+md.calving.meltingrate=numpy.zeros((md.mesh.numberofvertices,1))
 
 md=solve(md,TransientSolutionEnum())
Index: /issm/trunk-jpl/test/NightlyRun/test806.m
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test806.m	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/test806.m	(revision 18816)
@@ -26,4 +26,5 @@
 md.calving=calvinglevermann();
 md.calving.coeff=4.89e13*ones(md.mesh.numberofvertices,1);
+md.calving.meltingrate=zeros(md.mesh.numberofvertices,1);
 
 md.transient.requested_outputs={'default','StrainRateparallel','StrainRateperpendicular','Calvingratex','Calvingratey','CalvingCalvingrate'};
Index: /issm/trunk-jpl/test/NightlyRun/test806.py
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test806.py	(revision 18815)
+++ /issm/trunk-jpl/test/NightlyRun/test806.py	(revision 18816)
@@ -37,6 +37,7 @@
 md.transient.iscalving=True;
 
-md.calving=calvinglevermann();
-md.calving.coeff=4.89e13*numpy.ones((md.mesh.numberofvertices,1));
+md.calving=calvinglevermann()
+md.calving.coeff=4.89e13*numpy.ones((md.mesh.numberofvertices,1))
+md.calving.meltingrate=numpy.zeros((md.mesh.numberofvertices,1))
 
 md.transient.requested_outputs=['default','StrainRateparallel','StrainRateperpendicular','Calvingratex','Calvingratey','CalvingCalvingrate']