Index: /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp	(revision 16733)
+++ /issm/trunk-jpl/src/c/analyses/EnthalpyAnalysis.cpp	(revision 16734)
@@ -188,4 +188,90 @@
 }/*}}}*/
 void EnthalpyAnalysis::InputUpdateFromSolution(IssmDouble* solution,Element* element){/*{{{*/
-	_error_("not implemented yet");
-}/*}}}*/
+
+	bool        converged;
+	int         i,rheology_law;
+	IssmDouble  B_average,s_average,T_average=0.,P_average=0.;
+	int        *doflist   = NULL;
+	IssmDouble *xyz_list  = NULL;
+
+	/*Fetch number of nodes and dof for this finite element*/
+	int numnodes    = element->GetNumberOfNodes();
+
+	/*Fetch dof list and allocate solution vector*/
+	element->GetDofList(&doflist,NoneApproximationEnum,GsetEnum);
+	IssmDouble* values        = xNew<IssmDouble>(numnodes);
+	IssmDouble* pressure      = xNew<IssmDouble>(numnodes);
+	IssmDouble* temperature   = xNew<IssmDouble>(numnodes);
+	IssmDouble* waterfraction = xNew<IssmDouble>(numnodes);
+
+	/*Use the dof list to index into the solution vector: */
+	for(i=0;i<numnodes;i++){
+		values[i]=solution[doflist[i]];
+
+		/*Check solution*/
+		if(xIsNan<IssmDouble>(values[i])) _error_("NaN found in solution vector");
+	}
+
+	/*Get all inputs and parameters*/
+	element->GetInputValue(&converged,ConvergedEnum);
+	element->GetInputListOnNodes(pressure,PressureEnum);
+	if(converged){
+		for(i=0;i<numnodes;i++){
+			element->EnthalpyToThermal(&temperature[i],&waterfraction[i],values[i],pressure[i]);
+			if(waterfraction[i]<0.) _error_("Negative water fraction found in solution vector");
+			//if(waterfraction[i]>1.) _error_("Water fraction >1 found in solution vector");
+		}
+		element->AddInput(EnthalpyEnum,values,P1Enum);
+		element->AddInput(WaterfractionEnum,waterfraction,P1Enum);
+		element->AddInput(TemperatureEnum,temperature,P1Enum);
+
+		/*Update Rheology only if converged (we must make sure that the temperature is below melting point
+		 * otherwise the rheology could be negative*/
+		element->FindParam(&rheology_law,MaterialsRheologyLawEnum);
+		switch(rheology_law){
+			case NoneEnum:
+				/*Do nothing: B is not temperature dependent*/
+				break;
+			case PatersonEnum:
+				for(i=0;i<numnodes;i++) T_average+=values[i]/reCast<IssmDouble>(numnodes);
+				B_average=Paterson(T_average);
+				element->AddMaterialInput(MaterialsRheologyBEnum,&B_average,P0Enum);
+				break;
+			case ArrheniusEnum:{
+				Input* surface_input=element->GetInput(SurfaceEnum); _assert_(surface_input);
+				surface_input->GetInputAverage(&s_average);
+				element->GetVerticesCoordinates(&xyz_list);
+				for(i=0;i<numnodes;i++) T_average+=values[i]/reCast<IssmDouble>(numnodes);
+				//B_average=Arrhenius(T_average,
+							//s_average-((xyz_list[0][2]+xyz_list[1][2]+xyz_list[2][2]+xyz_list[3][2]+xyz_list[4][2]+xyz_list[5][2])/6.0),
+							//element->GetMaticeParameter(MaterialsRheologyNEnum));
+				element->AddMaterialInput(MaterialsRheologyBEnum,&B_average,P0Enum);
+				break;
+				}
+			case LliboutryDuvalEnum:
+				for(i=0;i<numnodes;i++) T_average+=values[i]/reCast<IssmDouble>(numnodes);
+				for(i=0;i<numnodes;i++) P_average+=pressure[i]/reCast<IssmDouble>(numnodes);
+				B_average=LliboutryDuval(T_average,P_average,
+							element->GetMaterialParameter(MaterialsRheologyNEnum),
+							element->GetMaterialParameter(MaterialsBetaEnum),
+							element->GetMaterialParameter(ConstantsReferencetemperatureEnum),
+							element->GetMaterialParameter(MaterialsHeatcapacityEnum),
+							element->GetMaterialParameter(MaterialsLatentheatEnum));
+				element->AddMaterialInput(MaterialsRheologyBEnum,&B_average,P0Enum);
+				break;
+			default:
+				_error_("Rheology law " << EnumToStringx(rheology_law) << " not supported yet");
+		}
+	}
+	else{
+		element->AddInput(EnthalpyPicardEnum,values,P1Enum);
+	}
+
+	/*Free ressources:*/
+	xDelete<IssmDouble>(values);
+	xDelete<IssmDouble>(pressure);
+	xDelete<IssmDouble>(temperature);
+	xDelete<IssmDouble>(waterfraction);
+	xDelete<IssmDouble>(xyz_list);
+	xDelete<int>(doflist);
+}/*}}}*/
Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 16733)
+++ /issm/trunk-jpl/src/c/classes/Elements/Element.h	(revision 16734)
@@ -46,4 +46,5 @@
 		virtual void   CreatePVector(Vector<IssmDouble>* pf)=0;
 		virtual void   CreateJacobianMatrix(Matrix<IssmDouble>* Jff)=0;
+		virtual void   EnthalpyToThermal(IssmDouble* ptemperature,IssmDouble* pwaterfraction,IssmDouble enthalpy,IssmDouble pressure)=0;
 		virtual void   FindParam(int* pvalue,int paramenum)=0;
 		virtual void   FindParam(IssmDouble* pvalue,int paramenum)=0;
Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 16733)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 16734)
@@ -963,4 +963,9 @@
 void Penta::Echo(void){
 	this->DeepEcho();
+}
+/*}}}*/
+/*FUNCTION Penta::EnthalpyToThermal{{{*/
+void Penta::EnthalpyToThermal(IssmDouble* ptemperature,IssmDouble* pwaterfraction,IssmDouble enthalpy,IssmDouble pressure){
+	matpar->EnthalpyToThermal(ptemperature,pwaterfraction,enthalpy,pressure);
 }
 /*}}}*/
Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 16733)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.h	(revision 16734)
@@ -84,4 +84,5 @@
 		void   CreateJacobianMatrix(Matrix<IssmDouble>* Jff);
 		void   Delta18oParameterization(void);
+		void   EnthalpyToThermal(IssmDouble* ptemperature,IssmDouble* pwaterfraction,IssmDouble enthalpy,IssmDouble pressure);
 		void	 GetDofList(int** pdoflist,int approximation_enum,int setenum);
 		void	 GetDofListVelocity(int** pdoflist,int setenum);
Index: /issm/trunk-jpl/src/c/classes/Elements/Seg.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Seg.h	(revision 16733)
+++ /issm/trunk-jpl/src/c/classes/Elements/Seg.h	(revision 16734)
@@ -83,4 +83,5 @@
 		void        CreateJacobianMatrix(Matrix<IssmDouble>* Jff){_error_("not implemented yet");};
 		void        Delta18oParameterization(void){_error_("not implemented yet");};
+		void        EnthalpyToThermal(IssmDouble* ptemperature,IssmDouble* pwaterfraction,IssmDouble enthalpy,IssmDouble pressure){_error_("not implemented yet");};
 		void        FindParam(int* pvalue,int paramenum);
 		void        FindParam(IssmDouble* pvalue,int paramenum);
Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 16733)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.h	(revision 16734)
@@ -79,4 +79,5 @@
 		void        CreateJacobianMatrix(Matrix<IssmDouble>* Jff);
 		void        Delta18oParameterization(void);
+		void        EnthalpyToThermal(IssmDouble* ptemperature,IssmDouble* pwaterfraction,IssmDouble enthalpy,IssmDouble pressure){_error_("not implemented yet");};
 		void        FindParam(int* pvalue,int paramenum);
 		void        FindParam(IssmDouble* pvalue,int paramenum);