Index: /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp	(revision 22280)
+++ /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp	(revision 22281)
@@ -315,29 +315,57 @@
 void TransientInput::GetInputAverageOnTimes(IssmDouble** pvalues, IssmDouble init_time){/*{{{*/
 
-	IssmDouble  time;
-	IssmDouble  preceding_time;
+	int					numnodes;
+	IssmDouble  subtimestep;
 	IssmDouble* values= NULL;
-
-	values=xNew<IssmDouble>(3);
-	
+	Input* input=NULL;
+
+	/*Initialize numnode from transientinput out of time loop*/
 	for(int i=0;i<this->numtimesteps;i++){
-		time						 = this->timesteps[i];
-		TriaInput*	input		 = (TriaInput*)this->inputs->GetObjectByOffset(i); _assert_(input);
-		int					numnodes = input->NumberofNodes(input->interpolation_type);
+		/*First compute the lengt of the sub-timestep*/
 		if(i==0){
-			preceding_time = init_time;
-			for(int j=0;j<numnodes;j++){
-				values[j]=(input->values[j]*(time-preceding_time));
-			}
+			subtimestep = this->timesteps[i]-init_time;
 		}
 		else{
-			preceding_time = this->timesteps[i-1];
+			subtimestep = this->timesteps[i]-this->timesteps[i-1];
+		}	
+		/*Figure out type of input and get it*/
+		input = xDynamicCast<Input*>(this->inputs->GetObjectByOffset(i)); _assert_(input);
+		switch(input->ObjectEnum()){
+		case TriaInputEnum:{
+			TriaInput* triainput = xDynamicCast<TriaInput*>(this->inputs->GetObjectByOffset(i)); _assert_(triainput);
+			numnodes= triainput->NumberofNodes(triainput->interpolation_type);
+			values=xNewZeroInit<IssmDouble>(numnodes);
+			/*Sum the values weighted by their respective sub-timestep*/
 			for(int j=0;j<numnodes;j++){
-				values[j]+=(input->values[j]*(time-preceding_time));
-				if(i==this->numtimesteps-1){
-					values[j]=values[j]/(this->timesteps[this->numtimesteps-1]-init_time);
-				}
+				values[j]+=(triainput->values[j]*subtimestep);
 			}
+			break;
 		}
+		case PentaInputEnum:{
+			PentaInput*	pentainput = xDynamicCast<PentaInput*>(this->inputs->GetObjectByOffset(i)); _assert_(pentainput);
+			numnodes= pentainput->NumberofNodes(pentainput->interpolation_type);
+			values=xNewZeroInit<IssmDouble>(numnodes);
+			/*Sum the values weighted by their respective sub-timestep*/
+			for(int j=0;j<numnodes;j++){
+				values[j]+=(pentainput->values[j]*subtimestep);
+			}
+			break;
+		}
+		case TetraInputEnum:{
+			TetraInput*	tetrainput = xDynamicCast<TetraInput*>(this->inputs->GetObjectByOffset(i)); _assert_(tetrainput);
+			numnodes= tetrainput->NumberofNodes(tetrainput->interpolation_type);
+			values=xNewZeroInit<IssmDouble>(numnodes);
+			/*Sum the values weighted by their respective sub-timestep*/
+			for(int j=0;j<numnodes;j++){
+				values[j]+=(tetrainput->values[j]*subtimestep);
+			}
+			break;
+		}
+		default: _error_("not implemented yet");
+		}
+	}
+	/*Compute the average based on the length of the full timestep*/
+	for(int j=0;j<numnodes;j++){
+		values[j]=values[j]/(this->timesteps[this->numtimesteps-1]-init_time);
 	}
 	*pvalues=values;
Index: /issm/trunk-jpl/src/c/cores/hydrology_core.cpp
===================================================================
--- /issm/trunk-jpl/src/c/cores/hydrology_core.cpp	(revision 22280)
+++ /issm/trunk-jpl/src/c/cores/hydrology_core.cpp	(revision 22281)
@@ -24,14 +24,11 @@
 	femmodel->parameters->FindParam(&solution_type,SolutionTypeEnum);	
 
-	if(VerboseSolution()) _printf0_("   computing water heads\n");
-			
-	/*first compute slopes: */
+	/*Using the Shreve based Model*/
 	if (hydrology_model==HydrologyshreveEnum){
+		if(VerboseSolution()) _printf0_("   computing water heads\n");
+		/*first compute slopes: */
 		surfaceslope_core(femmodel);
 		bedslope_core(femmodel);
-	}
-
-	/*Using the Shreve based Model*/
-	if (hydrology_model==HydrologyshreveEnum){
+		/*and then go to water column*/
 		if(VerboseSolution()) _printf0_("   computing water column\n");
 		femmodel->SetCurrentConfiguration(HydrologyShreveAnalysisEnum);
@@ -45,8 +42,4 @@
 			int outputs[3] = {WatercolumnEnum,HydrologyWaterVxEnum,HydrologyWaterVyEnum};
 			femmodel->RequestedOutputsx(&femmodel->results,&outputs[0],3);
-			
-			/*unload results*/
-			if(VerboseSolution()) _printf0_("   saving temporary results\n");
-			OutputResultsx(femmodel);
 		}
 	}
@@ -56,5 +49,5 @@
 		/*intermediary: */
 		int        step,hydroslices;
-		int        numoutputs;
+		int        numoutputs=0;
 		char       **requested_outputs = NULL;
 		IssmDouble time,init_time,hydrotime,yts;
@@ -111,7 +104,4 @@
 			if(VerboseSolution()) _printf0_("   saving results \n");
 			femmodel->RequestedOutputsx(&femmodel->results,requested_outputs,numoutputs);
-			/*unload results removed 23/11 needs checking*/
-			/* if(VerboseSolution()) _printf0_("   saving temporary results\n"); */
-			/* OutputResultsx(femmodel); */
 		}
 		/*Free ressources:*/	
@@ -119,5 +109,5 @@
 			for(int i=0;i<numoutputs;i++){
 				xDelete<char>(requested_outputs[i]);
-			} 
+			}
 			xDelete<char*>(requested_outputs);
 		}
@@ -138,8 +128,4 @@
 			int outputs[5] = {HydrologyHeadEnum,HydrologyGapHeightEnum,EffectivePressureEnum,HydrologyBasalFluxEnum,DegreeOfChannelizationEnum};
 			femmodel->RequestedOutputsx(&femmodel->results,&outputs[0],5);
-			
-			/*unload results*/
-			if(VerboseSolution()) _printf0_("   saving temporary results\n");
-			OutputResultsx(femmodel);
 		}
 	}