diff --git a/src/c/analyses/EnthalpyAnalysis.cpp b/src/c/analyses/EnthalpyAnalysis.cpp
index ff8086672..31126e1a9 100644
--- a/src/c/analyses/EnthalpyAnalysis.cpp
+++ b/src/c/analyses/EnthalpyAnalysis.cpp
@@ -1350,35 +1350,6 @@ void           EnthalpyAnalysis::GetBasalConstraintsTransient(Vector<IssmDouble>
 	delete gauss;
 	delete gaussup;
 }/*}}}*/
-void           EnthalpyAnalysis::GetBConduct(IssmDouble* B,Element* element,IssmDouble* xyz_list,Gauss* gauss){/*{{{*/
-	/*Compute B  matrix. B=[B1 B2 B3 B4 B5 B6] where Bi is of size 5*1.
-	 * For node i, Bi' can be expressed in the actual coordinate system
-	 * by:
-	 *       Bi_conduct=[ dh/dx ]
-	 *                  [ dh/dy ]
-	 *                  [ dh/dz ]
-	 * where h is the interpolation function for node i.
-	 *
-	 * We assume B has been allocated already, of size: 3x(1*numnodes)
-	 */
-
-	/*Fetch number of nodes for this finite element*/
-	int numnodes = element->GetNumberOfNodes();
-
-	/*Get nodal functions derivatives*/
-	IssmDouble* dbasis=xNew<IssmDouble>(3*numnodes);
-	element->NodalFunctionsDerivatives(dbasis,xyz_list,gauss);
-
-	/*Build B: */
-	for(int i=0;i<numnodes;i++){
-		B[numnodes*0+i] = dbasis[0*numnodes+i];
-		B[numnodes*1+i] = dbasis[1*numnodes+i];
-		B[numnodes*2+i] = dbasis[2*numnodes+i];
-	}
-
-	/*Clean-up*/
-	xDelete<IssmDouble>(dbasis);
-}/*}}}*/
 void           EnthalpyAnalysis::GetSolutionFromInputs(Vector<IssmDouble>* solution,Element* element){/*{{{*/
 	element->GetSolutionFromInputsOneDof(solution,EnthalpyEnum);
 }/*}}}*/
diff --git a/src/c/analyses/EnthalpyAnalysis.h b/src/c/analyses/EnthalpyAnalysis.h
index dabe76e76..98930e5a8 100644
--- a/src/c/analyses/EnthalpyAnalysis.h
+++ b/src/c/analyses/EnthalpyAnalysis.h
@@ -44,7 +44,6 @@ class EnthalpyAnalysis: public Analysis{
 		static void       GetBasalConstraints(Vector<IssmDouble>* vec_spc,Element* element);
 		static void       GetBasalConstraintsSteadystate(Vector<IssmDouble>* vec_spc,Element* element);
 		static void       GetBasalConstraintsTransient(Vector<IssmDouble>* vec_spc,Element* element);
-		void              GetBConduct(IssmDouble* B,Element* element,IssmDouble* xyz_list,Gauss* gauss);
 		void              GetSolutionFromInputs(Vector<IssmDouble>* solution,Element* element);
 		static int        GetThermalBasalCondition(Element* element, IssmDouble enthalpy, IssmDouble enthalpy_up, IssmDouble pressure, IssmDouble pressure_up, IssmDouble watercolumn, IssmDouble meltingrate);
 		static IssmDouble GetWetIceConductivity(Element* element, IssmDouble enthalpy, IssmDouble pressure);