Index: /issm/trunk-jpl/src/c/Makefile.am
===================================================================
--- /issm/trunk-jpl/src/c/Makefile.am	(revision 19748)
+++ /issm/trunk-jpl/src/c/Makefile.am	(revision 19749)
@@ -71,4 +71,5 @@
 					./classes/Loads/Moulin.cpp\
 					./classes/Loads/Numericalflux.cpp\
+					./classes/Loads/Neumannflux.cpp\
 					./classes/matrix/ElementMatrix.cpp\
 					./classes/matrix/ElementVector.cpp\
Index: /issm/trunk-jpl/src/c/analyses/HydrologySommersAnalysis.cpp
===================================================================
--- /issm/trunk-jpl/src/c/analyses/HydrologySommersAnalysis.cpp	(revision 19748)
+++ /issm/trunk-jpl/src/c/analyses/HydrologySommersAnalysis.cpp	(revision 19749)
@@ -30,5 +30,5 @@
 	if(hydrology_model!=HydrologysommersEnum) return;
 
-	//create penalties for nodes: no node can have water above the max
+	/*Create discrete loads for Moulins*/
 	CreateSingleNodeToElementConnectivity(iomodel);
 	for(int i=0;i<iomodel->numberofvertices;i++){
@@ -46,4 +46,19 @@
 	}
 	iomodel->DeleteData(1,MeshVertexonbaseEnum);
+
+	/*Deal with Neumann BC*/
+	int M,N;
+	int *segments = NULL;
+	iomodel->FetchData(&segments,&M,&N,MeshSegmentsEnum);
+
+	/*Check that the size seem right*/
+	_assert_(N==3); _assert_(M>=3);
+	for(int i=0;i<M;i++){
+		if(iomodel->my_elements[segments[i*3+2]-1]){
+			loads->AddObject(new Neumannflux(iomodel->loadcounter+i+1,i,iomodel,segments,HydrologySommersAnalysisEnum));
+		}
+	}
+
+	xDelete<int>(segments);
 
 }/*}}}*/
@@ -100,4 +115,5 @@
 	iomodel->FetchDataToInput(elements,HydrologyBumpHeightEnum);
 	iomodel->FetchDataToInput(elements,HydrologyReynoldsEnum);
+	iomodel->FetchDataToInput(elements,HydrologyNeumannfluxEnum);
 	iomodel->FetchDataToInput(elements,VxEnum);
 	iomodel->FetchDataToInput(elements,VyEnum);
@@ -260,5 +276,5 @@
 		IssmDouble pressure_ice   = rho_ice*g*thickness;    _assert_(pressure_ice>0.); 
 		IssmDouble pressure_water = rho_water*g*(head-bed);
-		_assert_(pressure_water<=pressure_ice);
+		if(pressure_water>pressure_ice) pressure_water = pressure_ice;
 
 		meltrate = 1/latentheat*(G+frictionheat+rho_water*g*conductivity*(dh[0]*dh[0]+dh[1]*dh[1]));
@@ -457,5 +473,5 @@
 		IssmDouble pressure_ice   = rho_ice*g*thickness;    _assert_(pressure_ice>0.); 
 		IssmDouble pressure_water = rho_water*g*(head-bed);
-		_assert_(pressure_water<=pressure_ice);
+		if(pressure_water>pressure_ice) pressure_water = pressure_ice;
 
 		meltrate = 1/latentheat*(G+frictionheat+rho_water*g*conductivity*(dh[0]*dh[0]+dh[1]*dh[1]));
Index: /issm/trunk-jpl/src/c/classes/Loads/Neumannflux.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Loads/Neumannflux.cpp	(revision 19749)
+++ /issm/trunk-jpl/src/c/classes/Loads/Neumannflux.cpp	(revision 19749)
@@ -0,0 +1,391 @@
+/*!\file Neumannflux.c
+ * \brief: implementation of the Neumannflux object
+ */
+
+/*Headers:*/
+/*{{{*/
+#ifdef HAVE_CONFIG_H
+#include <config.h>
+#else
+#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
+#endif
+
+#include "shared/shared.h"
+#include "../classes.h"
+/*}}}*/	
+
+/*Load macros*/
+#define NUMVERTICES 2
+#define NUMNODES_BOUNDARY 2
+
+/*Neumannflux constructors and destructor*/
+Neumannflux::Neumannflux(){/*{{{*/
+	this->parameters = NULL;
+	this->helement   = NULL;
+	this->element    = NULL;
+	this->hnodes     = NULL;
+	this->hvertices  = NULL;
+	this->nodes      = NULL;
+}
+/*}}}*/
+Neumannflux::Neumannflux(int neumannflux_id,int i,IoModel* iomodel,int* segments,int in_analysis_type){/*{{{*/
+
+
+	/*Some sanity checks*/
+	_assert_(segments);
+
+	/*neumannflux constructor data: */
+	int neumannflux_elem_id;
+	int neumannflux_vertex_ids[2];
+	int neumannflux_node_ids[2];
+
+	/*1: Get vertices ids*/
+	neumannflux_vertex_ids[0]=segments[3*i+0];
+	neumannflux_vertex_ids[1]=segments[3*i+1];
+
+	/*2: Get the ids of the nodes*/
+	neumannflux_node_ids[0]=iomodel->nodecounter+neumannflux_vertex_ids[0];
+	neumannflux_node_ids[1]=iomodel->nodecounter+neumannflux_vertex_ids[1];
+
+	/*Get element id*/
+	neumannflux_elem_id = segments[3*i+2];
+
+	/*Ok, we have everything to build the object: */
+	this->id=neumannflux_id;
+	this->analysis_type=in_analysis_type;
+
+	/*Hooks: */
+	this->hnodes    =new Hook(&neumannflux_node_ids[0],2);
+	this->hvertices =new Hook(&neumannflux_vertex_ids[0],2);
+	this->helement  =new Hook(&neumannflux_elem_id,1);
+
+	//this->parameters: we still can't point to it, it may not even exist. Configure will handle this.
+	this->parameters=NULL;
+	this->element=NULL;
+	this->nodes=NULL;
+}
+/*}}}*/
+Neumannflux::~Neumannflux(){/*{{{*/
+	this->parameters=NULL;
+	delete helement;
+	delete hnodes;
+	delete hvertices;
+}
+/*}}}*/
+
+/*Object virtual functions definitions:*/
+Object* Neumannflux::copy() {/*{{{*/
+
+	Neumannflux* neumannflux=NULL;
+
+	neumannflux=new Neumannflux();
+
+	/*copy fields: */
+	neumannflux->id=this->id;
+	neumannflux->analysis_type=this->analysis_type;
+
+	/*point parameters: */
+	neumannflux->parameters=this->parameters;
+
+	/*now deal with hooks and objects: */
+	neumannflux->hnodes    = (Hook*)this->hnodes->copy();
+	neumannflux->hvertices = (Hook*)this->hvertices->copy();
+	neumannflux->helement  = (Hook*)this->helement->copy();
+
+	/*corresponding fields*/
+	neumannflux->nodes    = (Node**)neumannflux->hnodes->deliverp();
+	neumannflux->vertices = (Vertex**)neumannflux->hvertices->deliverp();
+	neumannflux->element  = (Element*)neumannflux->helement->delivers();
+
+	return neumannflux;
+}
+/*}}}*/
+void    Neumannflux::Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){ /*{{{*/
+
+	_assert_(this);
+
+	/*ok, marshall operations: */
+	MARSHALLING_ENUM(NeumannfluxEnum);
+	MARSHALLING(id);
+	MARSHALLING(analysis_type);
+
+	if(marshall_direction==MARSHALLING_BACKWARD){
+		this->hnodes      = new Hook();
+		this->hvertices   = new Hook();
+		this->helement    = new Hook();
+	}
+
+	this->hnodes->Marshall(pmarshalled_data,pmarshalled_data_size,marshall_direction);
+	this->helement->Marshall(pmarshalled_data,pmarshalled_data_size,marshall_direction);
+	this->hvertices->Marshall(pmarshalled_data,pmarshalled_data_size,marshall_direction);
+
+	/*corresponding fields*/
+	nodes    =(Node**)this->hnodes->deliverp();
+	vertices =(Vertex**)this->hvertices->deliverp();
+	element  =(Element*)this->helement->delivers();
+
+}
+/*}}}*/
+void    Neumannflux::DeepEcho(void){/*{{{*/
+
+	_printf_("Neumannflux:\n");
+	_printf_("   id: " << id << "\n");
+	_printf_("   analysis_type: " << EnumToStringx(analysis_type) << "\n");
+	hnodes->DeepEcho();
+	hvertices->DeepEcho();
+	helement->DeepEcho();
+	_printf_("   parameters\n");
+	if(parameters)
+	 parameters->DeepEcho();
+	else
+	 _printf_("      NULL\n");
+}		
+/*}}}*/
+void    Neumannflux::Echo(void){/*{{{*/
+	_printf_("Neumannflux:\n");
+	_printf_("   id: " << id << "\n");
+	_printf_("   analysis_type: " << EnumToStringx(analysis_type) << "\n");
+	hnodes->Echo();
+	hvertices->Echo();
+	helement->Echo();
+	_printf_("   parameters: " << parameters << "\n");
+}
+/*}}}*/
+int     Neumannflux::Id(void){/*{{{*/
+	return id;
+}
+/*}}}*/
+int     Neumannflux::ObjectEnum(void){/*{{{*/
+
+	return NeumannfluxEnum;
+
+}
+/*}}}*/
+
+/*Load virtual functions definitions:*/
+void  Neumannflux::Configure(Elements* elementsin,Loads* loadsin,Nodes* nodesin,Vertices* verticesin,Materials* materialsin,Parameters* parametersin){/*{{{*/
+
+	/*Take care of hooking up all objects for this element, ie links the objects in the hooks to their respective 
+	 * datasets, using internal ids and offsets hidden in hooks: */
+	hnodes->configure((DataSet*)nodesin);
+	hvertices->configure((DataSet*)verticesin);
+	helement->configure((DataSet*)elementsin);
+
+	/*Initialize hooked fields*/
+	this->nodes    = (Node**)hnodes->deliverp();
+	this->vertices = (Vertex**)hvertices->deliverp();
+	this->element  = (Element*)helement->delivers();
+
+	/*point parameters to real dataset: */
+	this->parameters=parametersin;
+}
+/*}}}*/
+void  Neumannflux::CreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* Kfs){/*{{{*/
+
+	/*recover some parameters*/
+	ElementMatrix* Ke=NULL;
+	int analysis_type;
+	this->parameters->FindParam(&analysis_type,AnalysisTypeEnum);
+
+	/*Just branch to the correct element stiffness matrix generator, according to the type of analysis we are carrying out: */
+	switch(analysis_type){
+		case HydrologySommersAnalysisEnum:
+			/*Nothing!*/
+			break;
+		default:
+			_error_("analysis " << analysis_type << " (" << EnumToStringx(analysis_type) << ") not supported yet");
+	}
+
+	/*Add to global matrix*/
+	if(Ke){
+		Ke->AddToGlobal(Kff,Kfs);
+		delete Ke;
+	}
+
+}
+/*}}}*/
+void  Neumannflux::CreatePVector(Vector<IssmDouble>* pf){/*{{{*/
+
+	/*recover some parameters*/
+	ElementVector* pe=NULL;
+	int analysis_type;
+	this->parameters->FindParam(&analysis_type,AnalysisTypeEnum);
+
+	switch(analysis_type){
+		case HydrologySommersAnalysisEnum:
+			pe=CreatePVectorHydrologySommers();
+			break;
+		default:
+			_error_("analysis " << analysis_type << " (" << EnumToStringx(analysis_type) << ") not supported yet");
+	}
+
+	/*Add to global matrix*/
+	if(pe){
+		pe->AddToGlobal(pf);
+		delete pe;
+	}
+
+}
+/*}}}*/
+void  Neumannflux::GetNodesLidList(int* lidlist){/*{{{*/
+
+	_assert_(lidlist);
+	_assert_(nodes);
+
+	for(int i=0;i<NUMNODES_BOUNDARY;i++) lidlist[i]=nodes[i]->Lid();
+}
+/*}}}*/
+void  Neumannflux::GetNodesSidList(int* sidlist){/*{{{*/
+
+	_assert_(sidlist);
+	_assert_(nodes);
+
+	for(int i=0;i<NUMNODES_BOUNDARY;i++) sidlist[i]=nodes[i]->Sid();
+}
+/*}}}*/
+int   Neumannflux::GetNumberOfNodes(void){/*{{{*/
+
+	return NUMNODES_BOUNDARY;
+}
+/*}}}*/
+bool  Neumannflux::InAnalysis(int in_analysis_type){/*{{{*/
+	if (in_analysis_type==this->analysis_type) return true;
+	else return false;
+}
+/*}}}*/
+bool  Neumannflux::IsPenalty(void){/*{{{*/
+	return false;
+}
+/*}}}*/
+void  Neumannflux::PenaltyCreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* Kfs,IssmDouble kmax){/*{{{*/
+
+	/*No stiffness loads applied, do nothing: */
+	return;
+
+}
+/*}}}*/
+void  Neumannflux::PenaltyCreatePVector(Vector<IssmDouble>* pf,IssmDouble kmax){/*{{{*/
+
+	/*No penalty loads applied, do nothing: */
+	return;
+
+}
+/*}}}*/
+void  Neumannflux::ResetHooks(){/*{{{*/
+
+	this->nodes=NULL;
+	this->vertices=NULL;
+	this->element=NULL;
+	this->parameters=NULL;
+
+	/*Get Element type*/
+	this->hnodes->reset();
+	this->hvertices->reset();
+	this->helement->reset();
+
+}
+/*}}}*/
+void  Neumannflux::SetCurrentConfiguration(Elements* elementsin,Loads* loadsin,Nodes* nodesin,Vertices* verticesin,Materials* materialsin,Parameters* parametersin){/*{{{*/
+
+}
+/*}}}*/
+void  Neumannflux::SetwiseNodeConnectivity(int* pd_nz,int* po_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum){/*{{{*/
+
+	/*Output */
+	int d_nz = 0;
+	int o_nz = 0;
+
+	/*Loop over all nodes*/
+	for(int i=0;i<this->GetNumberOfNodes();i++){
+
+		if(!flags[this->nodes[i]->Lid()]){
+
+			/*flag current node so that no other element processes it*/
+			flags[this->nodes[i]->Lid()]=true;
+
+			int counter=0;
+			while(flagsindices[counter]>=0) counter++;
+			flagsindices[counter]=this->nodes[i]->Lid();
+
+			/*if node is clone, we have an off-diagonal non-zero, else it is a diagonal non-zero*/
+			switch(set2_enum){
+				case FsetEnum:
+					if(nodes[i]->indexing.fsize){
+						if(this->nodes[i]->IsClone())
+						 o_nz += 1;
+						else
+						 d_nz += 1;
+					}
+					break;
+				case GsetEnum:
+					if(nodes[i]->indexing.gsize){
+						if(this->nodes[i]->IsClone())
+						 o_nz += 1;
+						else
+						 d_nz += 1;
+					}
+					break;
+				case SsetEnum:
+					if(nodes[i]->indexing.ssize){
+						if(this->nodes[i]->IsClone())
+						 o_nz += 1;
+						else
+						 d_nz += 1;
+					}
+					break;
+				default: _error_("not supported");
+			}
+		}
+	}
+
+	/*Assign output pointers: */
+	*pd_nz=d_nz;
+	*po_nz=o_nz;
+}
+/*}}}*/
+
+/*Neumannflux management*/
+ElementVector* Neumannflux::CreatePVectorHydrologySommers(void){/*{{{*/
+
+	/* constants*/
+	const int numdof=2;
+
+	/* Intermediaries*/
+	IssmDouble Jdet,flux;
+	IssmDouble xyz_list[NUMVERTICES][3];
+	IssmDouble basis[numdof];
+
+	/*Initialize Load Vector and return if necessary*/
+	Tria*  tria=(Tria*)element;
+	_assert_(tria->FiniteElement()==P1Enum); 
+	if(!tria->IsIceInElement() || tria->IsFloating()) return NULL;
+
+	/*Initialize Element vector and other vectors*/
+	ElementVector* pe=new ElementVector(nodes,NUMNODES_BOUNDARY,this->parameters);
+
+	/*Retrieve all inputs and parameters*/
+	GetVerticesCoordinates(&xyz_list[0][0],vertices,NUMVERTICES);
+	Input* flux_input = tria->inputs->GetInput(HydrologyNeumannfluxEnum);  _assert_(flux_input); 
+
+	/*Check wether it is an inflow or outflow BC (0 is the middle of the segment)*/
+	int index1=tria->GetNodeIndex(nodes[0]);
+	int index2=tria->GetNodeIndex(nodes[1]);
+
+	/* Start  looping on the number of gaussian points: */
+	GaussTria* gauss=new GaussTria(index1,index2,2);
+	for(int ig=gauss->begin();ig<gauss->end();ig++){
+
+		gauss->GaussPoint(ig);
+
+		tria->GetSegmentJacobianDeterminant(&Jdet,&xyz_list[0][0],gauss);
+		tria->GetSegmentNodalFunctions(&basis[0],gauss,index1,index2,tria->FiniteElement());
+		flux_input->GetInputValue(&flux,gauss);
+
+		for(int i=0;i<numdof;i++) pe->values[i] += gauss->weight*Jdet*flux*basis[i];
+	}
+
+	/*Clean up and return*/
+	delete gauss;
+	return pe;
+}
+/*}}}*/
Index: /issm/trunk-jpl/src/c/classes/Loads/Neumannflux.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/Loads/Neumannflux.h	(revision 19749)
+++ /issm/trunk-jpl/src/c/classes/Loads/Neumannflux.h	(revision 19749)
@@ -0,0 +1,80 @@
+/*!\file Neumannflux.h
+ * \brief: header file for icefront object
+ */
+
+#ifndef _NEUMANNFLUX_H_
+#define _NEUMANNFLUX_H_
+
+/*Headers:*/
+#include "./Load.h"
+class Hook;
+class Parameters;
+class IoModel;
+class Element;
+class Vertex;
+class ElementMatrix;
+class ElementVector;
+
+class Neumannflux: public Load {
+
+	public: 
+		int id;
+		int analysis_type;
+
+		/*Hooks*/
+		Hook *helement;
+		Hook *hnodes;
+		Hook *hvertices;
+
+		/*Corresponding fields*/
+		Element     *element;
+		Vertex     **vertices;
+		Node       **nodes;
+		Parameters  *parameters;
+
+		/*Neumannflux constructors,destructors {{{*/
+		Neumannflux();
+		Neumannflux(int numericalflux_id,int i,IoModel* iomodel,int* segments,int analysis_type);
+		~Neumannflux();
+		/*}}}*/
+		/*Object virtual functions definitions:{{{ */
+		Object *copy();
+		void    DeepEcho();
+		void    Echo();
+		int     Id();
+		int     ObjectEnum();
+		void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction);
+		/*}}}*/
+		/*Update virtual functions resolution: {{{*/
+		void InputUpdateFromConstant(IssmDouble constant, int name){/*Do nothing*/};
+		void InputUpdateFromConstant(int constant, int name){/*Do nothing*/};
+		void InputUpdateFromConstant(bool constant, int name){/*Do nothing*/};
+		void InputUpdateFromIoModel(int index, IoModel* iomodel){_error_("not implemented yet");};
+		void InputUpdateFromMatrixDakota(IssmDouble* matrix, int nrows, int ncols, int name, int type){/*Do nothing*/}
+		void InputUpdateFromVector(IssmDouble* vector, int name, int type){/*Do nothing*/}
+		void InputUpdateFromVectorDakota(IssmDouble* vector, int name, int type){/*Do nothing*/}
+		/*}}}*/
+		/*Load virtual functions definitions: {{{*/
+		void Configure(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
+		void CreateJacobianMatrix(Matrix<IssmDouble>* Jff){_error_("Not implemented yet");};
+		void CreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* Kfs);
+		void CreatePVector(Vector<IssmDouble>* pf);
+		void GetNodesLidList(int* lidlist);
+		void GetNodesSidList(int* sidlist);
+		int  GetNumberOfNodes(void);
+		bool InAnalysis(int analysis_type);
+		bool IsPenalty(void);
+		void PenaltyCreateJacobianMatrix(Matrix<IssmDouble>* Jff,IssmDouble kmax){_error_("Not implemented yet");};
+		void PenaltyCreateKMatrix(Matrix<IssmDouble>* Kff, Matrix<IssmDouble>* kfs, IssmDouble kmax);
+		void PenaltyCreatePVector(Vector<IssmDouble>* pf, IssmDouble kmax);
+		void ResetHooks();
+		void SetwiseNodeConnectivity(int* d_nz,int* o_nz,Node* node,bool* flags,int* flagsindices,int set1_enum,int set2_enum);
+		void SetCurrentConfiguration(Elements* elements,Loads* loads,Nodes* nodes,Vertices* vertices,Materials* materials,Parameters* parameters);
+		/*}}}*/
+		/*Neumannflux management:{{{*/
+		ElementVector* CreatePVectorHydrologySommers(void);
+		/*}}}*/
+
+};
+
+#endif  /* _NEUMANNFLUX_H_ */
Index: /issm/trunk-jpl/src/c/classes/classes.h
===================================================================
--- /issm/trunk-jpl/src/c/classes/classes.h	(revision 19748)
+++ /issm/trunk-jpl/src/c/classes/classes.h	(revision 19749)
@@ -34,4 +34,5 @@
 #include "./Loads/Friction.h"
 #include "./Loads/Numericalflux.h"
+#include "./Loads/Neumannflux.h"
 #include "./Loads/Riftfront.h"
 #include "./Loads/Penpair.h"
Index: /issm/trunk-jpl/src/c/cores/hydrology_core.cpp
===================================================================
--- /issm/trunk-jpl/src/c/cores/hydrology_core.cpp	(revision 19748)
+++ /issm/trunk-jpl/src/c/cores/hydrology_core.cpp	(revision 19749)
@@ -79,5 +79,4 @@
 
 	else if (hydrology_model==HydrologysommersEnum){
-		if(VerboseSolution()) _printf0_("   computing water head\n");
 		femmodel->SetCurrentConfiguration(HydrologySommersAnalysisEnum);
 		solutionsequence_nonlinear(femmodel,modify_loads);
Index: /issm/trunk-jpl/src/c/cores/transient_core.cpp
===================================================================
--- /issm/trunk-jpl/src/c/cores/transient_core.cpp	(revision 19748)
+++ /issm/trunk-jpl/src/c/cores/transient_core.cpp	(revision 19749)
@@ -79,5 +79,5 @@
 		femmodel->parameters->SetParam(step,StepEnum);
 
-		if(VerboseSolution()) _printf0_("iteration " << step << "/" << floor((finaltime-time)/dt)+step << "  time [yr]: " << setprecision(8) << time/yts << " (time step: " << dt/yts << ")\n");
+		if(VerboseSolution()) _printf0_("iteration " << step << "/" << floor((finaltime-time)/dt)+step << "  time [yr]: " << setprecision(4) << time/yts << " (time step: " << dt/yts << ")\n");
 		if(step%output_frequency==0 || (time >= finaltime - (yts*DBL_EPSILON)) || step==1)
 		 save_results=true;
Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/ModelProcessorx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/ModelProcessorx/ModelProcessorx.cpp	(revision 19748)
+++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/ModelProcessorx.cpp	(revision 19749)
@@ -59,7 +59,10 @@
 		 * constraints, and ids for objects created in next call to CreateDataSets
 		 * will need to start at the end of the updated counters: */
-		iomodel->nodecounter       = nodes->MaximumId();
+		if(nodes->Size()) iomodel->nodecounter = nodes->MaximumId();
 		iomodel->loadcounter       = loads->NumberOfLoads();
 		iomodel->constraintcounter = constraints->NumberOfConstraints();
+
+		/*Make sure nodecounter is at least 0 (if no node exists, maxid will be -1*/
+		_assert_(iomodel->nodecounter>=0);
 	}
 
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 19748)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h	(revision 19749)
@@ -57,5 +57,4 @@
 	BaseEnum,
 	ConstantsGEnum,
-	ConstantsOmegaEnum,
 	ConstantsReferencetemperatureEnum,
 	ConstantsYtsEnum,
@@ -164,4 +163,5 @@
 	HydrologyMoulinInputEnum,
 	HydrologyReynoldsEnum,
+	HydrologyNeumannfluxEnum,
 	HydrologySpcheadEnum,
 	HydrologyConductivityEnum,
@@ -277,4 +277,5 @@
 	MeshZEnum,
 	MeshElementtypeEnum,
+	MeshSegmentsEnum,
 	DomainTypeEnum,
 	DomainDimensionEnum,
@@ -579,4 +580,5 @@
 	NumericalfluxEnum,
 	NumericalfluxTypeEnum,
+	NeumannfluxEnum,
 	ParamEnum,
 	MoulinEnum,
Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 19748)
+++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp	(revision 19749)
@@ -63,5 +63,4 @@
 		case BaseEnum : return "Base";
 		case ConstantsGEnum : return "ConstantsG";
-		case ConstantsOmegaEnum : return "ConstantsOmega";
 		case ConstantsReferencetemperatureEnum : return "ConstantsReferencetemperature";
 		case ConstantsYtsEnum : return "ConstantsYts";
@@ -170,4 +169,5 @@
 		case HydrologyMoulinInputEnum : return "HydrologyMoulinInput";
 		case HydrologyReynoldsEnum : return "HydrologyReynolds";
+		case HydrologyNeumannfluxEnum : return "HydrologyNeumannflux";
 		case HydrologySpcheadEnum : return "HydrologySpchead";
 		case HydrologyConductivityEnum : return "HydrologyConductivity";
@@ -283,4 +283,5 @@
 		case MeshZEnum : return "MeshZ";
 		case MeshElementtypeEnum : return "MeshElementtype";
+		case MeshSegmentsEnum : return "MeshSegments";
 		case DomainTypeEnum : return "DomainType";
 		case DomainDimensionEnum : return "DomainDimension";
@@ -571,4 +572,5 @@
 		case NumericalfluxEnum : return "Numericalflux";
 		case NumericalfluxTypeEnum : return "NumericalfluxType";
+		case NeumannfluxEnum : return "Neumannflux";
 		case ParamEnum : return "Param";
 		case MoulinEnum : return "Moulin";
Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 19748)
+++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp	(revision 19749)
@@ -63,5 +63,4 @@
 	      else if (strcmp(name,"Base")==0) return BaseEnum;
 	      else if (strcmp(name,"ConstantsG")==0) return ConstantsGEnum;
-	      else if (strcmp(name,"ConstantsOmega")==0) return ConstantsOmegaEnum;
 	      else if (strcmp(name,"ConstantsReferencetemperature")==0) return ConstantsReferencetemperatureEnum;
 	      else if (strcmp(name,"ConstantsYts")==0) return ConstantsYtsEnum;
@@ -137,9 +136,9 @@
 	      else if (strcmp(name,"HydrologydcSedimentTransmitivity")==0) return HydrologydcSedimentTransmitivityEnum;
 	      else if (strcmp(name,"HydrologydcWaterCompressibility")==0) return HydrologydcWaterCompressibilityEnum;
+	      else if (strcmp(name,"HydrologydcSpceplHead")==0) return HydrologydcSpceplHeadEnum;
          else stage=2;
    }
    if(stage==2){
-	      if (strcmp(name,"HydrologydcSpceplHead")==0) return HydrologydcSpceplHeadEnum;
-	      else if (strcmp(name,"HydrologydcMaskEplactiveNode")==0) return HydrologydcMaskEplactiveNodeEnum;
+	      if (strcmp(name,"HydrologydcMaskEplactiveNode")==0) return HydrologydcMaskEplactiveNodeEnum;
 	      else if (strcmp(name,"HydrologydcMaskEplactiveElt")==0) return HydrologydcMaskEplactiveEltEnum;
 	      else if (strcmp(name,"HydrologydcEplCompressibility")==0) return HydrologydcEplCompressibilityEnum;
@@ -173,4 +172,5 @@
 	      else if (strcmp(name,"HydrologyMoulinInput")==0) return HydrologyMoulinInputEnum;
 	      else if (strcmp(name,"HydrologyReynolds")==0) return HydrologyReynoldsEnum;
+	      else if (strcmp(name,"HydrologyNeumannflux")==0) return HydrologyNeumannfluxEnum;
 	      else if (strcmp(name,"HydrologySpchead")==0) return HydrologySpcheadEnum;
 	      else if (strcmp(name,"HydrologyConductivity")==0) return HydrologyConductivityEnum;
@@ -289,4 +289,5 @@
 	      else if (strcmp(name,"MeshZ")==0) return MeshZEnum;
 	      else if (strcmp(name,"MeshElementtype")==0) return MeshElementtypeEnum;
+	      else if (strcmp(name,"MeshSegments")==0) return MeshSegmentsEnum;
 	      else if (strcmp(name,"DomainType")==0) return DomainTypeEnum;
 	      else if (strcmp(name,"DomainDimension")==0) return DomainDimensionEnum;
@@ -382,9 +383,9 @@
 	      else if (strcmp(name,"SmbV")==0) return SmbVEnum;
 	      else if (strcmp(name,"SmbDswrf")==0) return SmbDswrfEnum;
-	      else if (strcmp(name,"SmbDlwrf")==0) return SmbDlwrfEnum;
          else stage=4;
    }
    if(stage==4){
-	      if (strcmp(name,"SmbP")==0) return SmbPEnum;
+	      if (strcmp(name,"SmbDlwrf")==0) return SmbDlwrfEnum;
+	      else if (strcmp(name,"SmbP")==0) return SmbPEnum;
 	      else if (strcmp(name,"SmbSwf")==0) return SmbSwfEnum;
 	      else if (strcmp(name,"SmbEAir")==0) return SmbEAirEnum;
@@ -505,9 +506,9 @@
 	      else if (strcmp(name,"FreeSurfaceBaseAnalysis")==0) return FreeSurfaceBaseAnalysisEnum;
 	      else if (strcmp(name,"FreeSurfaceTopAnalysis")==0) return FreeSurfaceTopAnalysisEnum;
-	      else if (strcmp(name,"SurfaceNormalVelocity")==0) return SurfaceNormalVelocityEnum;
          else stage=5;
    }
    if(stage==5){
-	      if (strcmp(name,"ExtrudeFromBaseAnalysis")==0) return ExtrudeFromBaseAnalysisEnum;
+	      if (strcmp(name,"SurfaceNormalVelocity")==0) return SurfaceNormalVelocityEnum;
+	      else if (strcmp(name,"ExtrudeFromBaseAnalysis")==0) return ExtrudeFromBaseAnalysisEnum;
 	      else if (strcmp(name,"ExtrudeFromTopAnalysis")==0) return ExtrudeFromTopAnalysisEnum;
 	      else if (strcmp(name,"DepthAverageAnalysis")==0) return DepthAverageAnalysisEnum;
@@ -583,4 +584,5 @@
 	      else if (strcmp(name,"Numericalflux")==0) return NumericalfluxEnum;
 	      else if (strcmp(name,"NumericalfluxType")==0) return NumericalfluxTypeEnum;
+	      else if (strcmp(name,"Neumannflux")==0) return NeumannfluxEnum;
 	      else if (strcmp(name,"Param")==0) return ParamEnum;
 	      else if (strcmp(name,"Moulin")==0) return MoulinEnum;
@@ -627,10 +629,10 @@
 	      else if (strcmp(name,"Water")==0) return WaterEnum;
 	      else if (strcmp(name,"Closed")==0) return ClosedEnum;
-	      else if (strcmp(name,"Free")==0) return FreeEnum;
-	      else if (strcmp(name,"Open")==0) return OpenEnum;
          else stage=6;
    }
    if(stage==6){
-	      if (strcmp(name,"Adjointp")==0) return AdjointpEnum;
+	      if (strcmp(name,"Free")==0) return FreeEnum;
+	      else if (strcmp(name,"Open")==0) return OpenEnum;
+	      else if (strcmp(name,"Adjointp")==0) return AdjointpEnum;
 	      else if (strcmp(name,"Adjointx")==0) return AdjointxEnum;
 	      else if (strcmp(name,"Adjointy")==0) return AdjointyEnum;
@@ -750,10 +752,10 @@
 	      else if (strcmp(name,"OneLayerP4z")==0) return OneLayerP4zEnum;
 	      else if (strcmp(name,"CrouzeixRaviart")==0) return CrouzeixRaviartEnum;
-	      else if (strcmp(name,"LACrouzeixRaviart")==0) return LACrouzeixRaviartEnum;
-	      else if (strcmp(name,"SaveResults")==0) return SaveResultsEnum;
          else stage=7;
    }
    if(stage==7){
-	      if (strcmp(name,"BoolExternalResult")==0) return BoolExternalResultEnum;
+	      if (strcmp(name,"LACrouzeixRaviart")==0) return LACrouzeixRaviartEnum;
+	      else if (strcmp(name,"SaveResults")==0) return SaveResultsEnum;
+	      else if (strcmp(name,"BoolExternalResult")==0) return BoolExternalResultEnum;
 	      else if (strcmp(name,"DoubleExternalResult")==0) return DoubleExternalResultEnum;
 	      else if (strcmp(name,"DoubleMatExternalResult")==0) return DoubleMatExternalResultEnum;
@@ -873,10 +875,10 @@
 	      else if (strcmp(name,"MisfitDefinitionenum")==0) return MisfitDefinitionenumEnum;
 	      else if (strcmp(name,"MisfitModelEnum")==0) return MisfitModelEnumEnum;
-	      else if (strcmp(name,"MisfitObservation")==0) return MisfitObservationEnum;
-	      else if (strcmp(name,"MisfitObservationEnum")==0) return MisfitObservationEnumEnum;
          else stage=8;
    }
    if(stage==8){
-	      if (strcmp(name,"MisfitLocal")==0) return MisfitLocalEnum;
+	      if (strcmp(name,"MisfitObservation")==0) return MisfitObservationEnum;
+	      else if (strcmp(name,"MisfitObservationEnum")==0) return MisfitObservationEnumEnum;
+	      else if (strcmp(name,"MisfitLocal")==0) return MisfitLocalEnum;
 	      else if (strcmp(name,"MisfitTimeinterpolation")==0) return MisfitTimeinterpolationEnum;
 	      else if (strcmp(name,"MisfitWeights")==0) return MisfitWeightsEnum;
Index: /issm/trunk-jpl/src/m/classes/constants.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/constants.m	(revision 19748)
+++ /issm/trunk-jpl/src/m/classes/constants.m	(revision 19749)
@@ -61,5 +61,4 @@
 		function marshall(self,md,fid) % {{{
 			WriteData(fid,'object',self,'fieldname','g','format','Double');
-			WriteData(fid,'object',self,'fieldname','omega','format','Double');
 			WriteData(fid,'object',self,'fieldname','yts','format','Double');
 			WriteData(fid,'object',self,'fieldname','referencetemperature','format','Double');
Index: /issm/trunk-jpl/src/m/classes/hydrologysommers.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/hydrologysommers.m	(revision 19748)
+++ /issm/trunk-jpl/src/m/classes/hydrologysommers.m	(revision 19749)
@@ -14,7 +14,8 @@
 		reynolds        = NaN;
 		spchead         = NaN;
+		neumannflux     = NaN;
 	end
 	methods
-		function self = extrude(self,md) % {{{
+		function self = extrud,e(self,md) % {{{
 		end % }}}
 		function self = hydrologysommers(varargin) % {{{
@@ -44,4 +45,5 @@
 			md = checkfield(md,'fieldname','hydrology.moulin_input','>=',0,'NaN',1,'timeseries',1);
 			md = checkfield(md,'fieldname','hydrology.reynolds','>',0,'size',[md.mesh.numberofelements 1],'NaN',1);
+			md = checkfield(md,'fieldname','hydrology.neumannflux','timeseries',1,'NaN',1);
 			md = checkfield(md,'fieldname','hydrology.spchead','size',[md.mesh.numberofvertices 1]);
 		end % }}}
@@ -55,4 +57,5 @@
 			fielddisplay(self,'moulin_input','liquid water input from moulins (at the vertices) to subglacial system (m^3/s)');
 			fielddisplay(self,'reynolds','Reynolds'' number');
+			fielddisplay(self,'neumannflux','water flux applied along the model boundary (m^2/s)');
 			fielddisplay(self,'spchead','water head constraints (NaN means no constraint) (m)');
 		end % }}}
@@ -69,4 +72,5 @@
 			WriteData(fid,'object',self,'class','hydrology','fieldname','moulin_input','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1);
 			WriteData(fid,'object',self,'class','hydrology','fieldname','reynolds','format','DoubleMat','mattype',2);
+			WriteData(fid,'object',self,'class','hydrology','fieldname','neumannflux','format','DoubleMat','mattype',2,'timeserieslength',md.mesh.numberofelements+1);
 			WriteData(fid,'object',self,'class','hydrology','fieldname','spchead','format','DoubleMat','mattype',1);
 		end % }}}
Index: /issm/trunk-jpl/src/m/classes/mesh2d.m
===================================================================
--- /issm/trunk-jpl/src/m/classes/mesh2d.m	(revision 19748)
+++ /issm/trunk-jpl/src/m/classes/mesh2d.m	(revision 19749)
@@ -91,4 +91,5 @@
 			md = checkfield(md,'fieldname','mesh.numberofvertices','>',0);
 			md = checkfield(md,'fieldname','mesh.average_vertex_connectivity','>=',9,'message','''mesh.average_vertex_connectivity'' should be at least 9 in 2d');
+			md = checkfield(md,'fieldname','mesh.segments','NaN',1,'>',0,'size',[NaN 3]);
 
 			switch(solution),
@@ -126,46 +127,4 @@
 			fielddisplay(self,'epsg','EPSG code (ex: 3413 for UPS Greenland, 3031 for UPS Antarctica)');
 		end % }}}
-		function createxml(self,fid) % {{{
-			fprintf(fid, '<!-- 2D tria Mesh (horizontal) -->\n');
-
-			%elements and vertices
-			fprintf(fid,'%s\n%s\n%s\n','<frame key="1" label="Elements and vertices">','<section name="mesh" />');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n','  <parameter key ="numberofelements" type="',class(self.numberofelements),'" default="',convert2str(self.numberofelements),'">','     <section name="mesh" />','     <help> number of elements </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n','  <parameter key ="numberofvertices" type="',class(self.numberofvertices),'" default="',convert2str(self.numberofvertices),'">','     <section name="mesh" />','     <help> number of vertices </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n','  <parameter key ="elements" type="',class(self.elements),'" default="',convert2str(self.elements),'">','     <section name="mesh" />','     <help> vertex indices of the mesh elements </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="x" type="',class(self.x),'" default="',convert2str(self.x),'">','     <section name="mesh" />','     <help> vertices x coordinate [m] </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="y" type="',class(self.y),'" default="',convert2str(self.y),'">','     <section name="mesh" />','     <help> vertices y coordinate [m] </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="edges" type="',class(self.edges),'" default="',convert2str(self.edges),'">','     <section name="mesh" />','     <help> edges of the 2d mesh (vertex1 vertex2 element1 element2) </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="numberofedges" type="',class(self.numberofedges),'" default="',convert2str(self.numberofedges),'">','     <section name="mesh" />','     <help> number of edges of the 2d mesh </help>','  </parameter>');
-			fprintf(fid,'%s\n%s\n','</frame>');
-
-			% properties
-			fprintf(fid,'%s\n%s\n%s\n','<frame key="2" label="Properties">','<section name="mesh" />');             
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="vertexonboundary" type="',class(self.vertexonboundary),'" default="',convert2str(self.vertexonboundary),'">','     <section name="mesh" />','     <help> vertices on the boundary of the domain flag list </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="segments" type="',class(self.segments),'" default="',convert2str(self.segments),'">','     <section name="mesh" />','     <help> edges on domain boundary (vertex1 vertex2 element) </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="segmentmarkers" type="',class(self.segmentmarkers),'" default="',convert2str(self.segmentmarkers),'">','     <section name="mesh" />','     <help> number associated to each segment </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="vertexconnectivity" type="',class(self.vertexconnectivity),'" default="',convert2str(self.vertexconnectivity),'">','     <section name="mesh" />','     <help> list of vertices connected to vertex_i </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="elementconnectivity" type="',class(self.elementconnectivity),'" default="',convert2str(self.elementconnectivity),'">','     <section name="mesh" />','     <help> list of vertices connected to element_i </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="average_vertex_connectivity" type="',class(self.average_vertex_connectivity),'" default="',convert2str(self.average_vertex_connectivity),'">','     <section name="mesh" />','     <help> average number of vertices connected to one vertex </help>','  </parameter>');
-			fprintf(fid,'%s\n%s\n','</frame>');
-
-			%extracted model
-			fprintf(fid,'%s\n%s\n%s\n','<frame key="3" label="Extracted Model">','<section name="mesh" />'); 
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="extractedvertices" type="',class(self.extractedvertices),'" default="',convert2str(self.extractedvertices),'">','     <section name="mesh" />','     <help> vertices extracted from the model </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="extractedelements" type="',class(self.extractedelements),'" default="',convert2str(self.extractedelements),'">','     <section name="mesh" />','     <help> elements extracted from the model </help>','  </parameter>');
-			fprintf(fid,'%s\n%s\n','</frame>');
-
-			%projection
-			fprintf(fid,'%s\n%s\n%s\n','<frame key="4" label="Projection">','<section name="mesh" />'); 
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="lat" type="',class(self.lat),'" default="',convert2str(self.lat),'">','     <section name="mesh" />','     <help> vertices latitude [degrees] </help>','  </parameter>');
-			fprintf(fid,'%s%s%s%s%s\n%s\n%s\n%s\n','  <parameter key ="long" type="',class(self.long),'" default="',convert2str(self.long),'">','     <section name="mesh" />','     <help> verticies longitude [degrees] </help>','  </parameter>');
-			% choice (epsg) 'n' or 's'
-			fprintf(fid,'%s\n%s\n%s\n','  <parameter key ="epsg" type="alternative" optional="false">','     <section name="mesh" />','     <help> Indicate epsg ''n'' or ''s'' </help>');
-			fprintf(fid,'%s\n','       <option value="n" type="string" default="true"> </option>');
-			fprintf(fid,'%s\n','       <option value="s" type="string" default="false"> </option>');
-			fprintf(fid,'%s\n','  </parameter>');
-			fprintf(fid,'%s\n%s\n','</frame>');
-
-		end % }}}
 		function marshall(self,md,fid) % {{{
 			WriteData(fid,'enum',DomainTypeEnum(),'data',StringToEnum(['Domain' domaintype(self)]),'format','Integer');
@@ -180,4 +139,5 @@
 			WriteData(fid,'object',self,'class','mesh','fieldname','average_vertex_connectivity','format','Integer');
 			WriteData(fid,'object',self,'class','mesh','fieldname','vertexonboundary','format','DoubleMat','mattype',1);
+			WriteData(fid,'object',self,'class','mesh','fieldname','segments','format','DoubleMat','mattype',3);
 		end % }}}
 		function t = domaintype(self) % {{{
Index: sm/trunk-jpl/src/m/enum/ConstantsOmegaEnum.m
===================================================================
--- /issm/trunk-jpl/src/m/enum/ConstantsOmegaEnum.m	(revision 19748)
+++ 	(revision )
@@ -1,11 +1,0 @@
-function macro=ConstantsOmegaEnum()
-%CONSTANTSOMEGAENUM - Enum of ConstantsOmega
-%
-%   WARNING: DO NOT MODIFY THIS FILE
-%            this file has been automatically generated by src/c/shared/Enum/Synchronize.sh
-%            Please read src/c/shared/Enum/README for more information
-%
-%   Usage:
-%      macro=ConstantsOmegaEnum()
-
-macro=StringToEnum('ConstantsOmega');
Index: /issm/trunk-jpl/src/m/enum/EnumDefinitions.py
===================================================================
--- /issm/trunk-jpl/src/m/enum/EnumDefinitions.py	(revision 19748)
+++ /issm/trunk-jpl/src/m/enum/EnumDefinitions.py	(revision 19749)
@@ -55,5 +55,4 @@
 def BaseEnum(): return StringToEnum("Base")[0]
 def ConstantsGEnum(): return StringToEnum("ConstantsG")[0]
-def ConstantsOmegaEnum(): return StringToEnum("ConstantsOmega")[0]
 def ConstantsReferencetemperatureEnum(): return StringToEnum("ConstantsReferencetemperature")[0]
 def ConstantsYtsEnum(): return StringToEnum("ConstantsYts")[0]
@@ -162,4 +161,5 @@
 def HydrologyMoulinInputEnum(): return StringToEnum("HydrologyMoulinInput")[0]
 def HydrologyReynoldsEnum(): return StringToEnum("HydrologyReynolds")[0]
+def HydrologyNeumannfluxEnum(): return StringToEnum("HydrologyNeumannflux")[0]
 def HydrologySpcheadEnum(): return StringToEnum("HydrologySpchead")[0]
 def HydrologyConductivityEnum(): return StringToEnum("HydrologyConductivity")[0]
@@ -275,4 +275,5 @@
 def MeshZEnum(): return StringToEnum("MeshZ")[0]
 def MeshElementtypeEnum(): return StringToEnum("MeshElementtype")[0]
+def MeshSegmentsEnum(): return StringToEnum("MeshSegments")[0]
 def DomainTypeEnum(): return StringToEnum("DomainType")[0]
 def DomainDimensionEnum(): return StringToEnum("DomainDimension")[0]
@@ -563,4 +564,5 @@
 def NumericalfluxEnum(): return StringToEnum("Numericalflux")[0]
 def NumericalfluxTypeEnum(): return StringToEnum("NumericalfluxType")[0]
+def NeumannfluxEnum(): return StringToEnum("Neumannflux")[0]
 def ParamEnum(): return StringToEnum("Param")[0]
 def MoulinEnum(): return StringToEnum("Moulin")[0]
Index: /issm/trunk-jpl/src/m/enum/HydrologyNeumannfluxEnum.m
===================================================================
--- /issm/trunk-jpl/src/m/enum/HydrologyNeumannfluxEnum.m	(revision 19749)
+++ /issm/trunk-jpl/src/m/enum/HydrologyNeumannfluxEnum.m	(revision 19749)
@@ -0,0 +1,11 @@
+function macro=HydrologyNeumannfluxEnum()
+%HYDROLOGYNEUMANNFLUXENUM - Enum of HydrologyNeumannflux
+%
+%   WARNING: DO NOT MODIFY THIS FILE
+%            this file has been automatically generated by src/c/shared/Enum/Synchronize.sh
+%            Please read src/c/shared/Enum/README for more information
+%
+%   Usage:
+%      macro=HydrologyNeumannfluxEnum()
+
+macro=StringToEnum('HydrologyNeumannflux');
Index: /issm/trunk-jpl/src/m/enum/MeshSegmentsEnum.m
===================================================================
--- /issm/trunk-jpl/src/m/enum/MeshSegmentsEnum.m	(revision 19749)
+++ /issm/trunk-jpl/src/m/enum/MeshSegmentsEnum.m	(revision 19749)
@@ -0,0 +1,11 @@
+function macro=MeshSegmentsEnum()
+%MESHSEGMENTSENUM - Enum of MeshSegments
+%
+%   WARNING: DO NOT MODIFY THIS FILE
+%            this file has been automatically generated by src/c/shared/Enum/Synchronize.sh
+%            Please read src/c/shared/Enum/README for more information
+%
+%   Usage:
+%      macro=MeshSegmentsEnum()
+
+macro=StringToEnum('MeshSegments');
Index: /issm/trunk-jpl/src/m/enum/NeumannfluxEnum.m
===================================================================
--- /issm/trunk-jpl/src/m/enum/NeumannfluxEnum.m	(revision 19749)
+++ /issm/trunk-jpl/src/m/enum/NeumannfluxEnum.m	(revision 19749)
@@ -0,0 +1,11 @@
+function macro=NeumannfluxEnum()
+%NEUMANNFLUXENUM - Enum of Neumannflux
+%
+%   WARNING: DO NOT MODIFY THIS FILE
+%            this file has been automatically generated by src/c/shared/Enum/Synchronize.sh
+%            Please read src/c/shared/Enum/README for more information
+%
+%   Usage:
+%      macro=NeumannfluxEnum()
+
+macro=StringToEnum('Neumannflux');
