Index: ../trunk-jpl/src/m/classes/hydrologydc.py
===================================================================
--- ../trunk-jpl/src/m/classes/hydrologydc.py	(revision 18585)
+++ ../trunk-jpl/src/m/classes/hydrologydc.py	(revision 18586)
@@ -1,3 +1,4 @@
+import numpy
 from fielddisplay import fielddisplay
 from EnumDefinitions import *
 from checkfield import checkfield
@@ -4,14 +5,14 @@
 from WriteData import WriteData
 
 class hydrologydc(object):
-    """
-    Hydrologydc class definition
-    
-    Usage:
-    hydrologydc=hydrologydc();
-    """
+	"""
+	Hydrologydc class definition
 
-    def __init__(self): # {{{
+	Usage:
+		hydrologydc=hydrologydc();
+	"""
+
+	def __init__(self): # {{{
 		self.water_compressibility    = 0
 		self.isefficientlayer         = 0
 		self.penalty_factor           = 0
@@ -22,7 +23,7 @@
 		self.sedimentlimit            = 0
 		self.transfer_flag            = 0
 		self.leakage_factor           = 0
-                self.basal_moulin_input       = float('NaN')
+		self.basal_moulin_input       = float('NaN')
 
 		self.spcsediment_head         = float('NaN')
 		self.sediment_transmitivity   = float('NaN')
@@ -38,153 +39,147 @@
 		self.epl_initial_thickness    = 0
 		self.epl_max_thickness        = 0
 		self.epl_conductivity         = 0
-                
+				 
 		#set defaults
 		self.setdefaultparameters()
+	#}}}
 
-                #}}}
-                
 	def __repr__(self): # {{{
-            string='   hydrology Dual Porous Continuum Equivalent parameters:'
-            string='   - general parameters'
-            string="%s\n%s"%(string,fielddisplay(self,'water_compressibility','compressibility of water [Pa^-1]'))
-            string="%s\n%s"%(string,fielddisplay(self,'isefficientlayer','do we use an efficient drainage system [1: true 0: false]'))
-            string="%s\n%s"%(string,fielddisplay(self,'penalty_factor','exponent of the value used in the penalisation method [dimensionless]'))
-            string="%s\n%s"%(string,fielddisplay(self,'penalty_lock','stabilize unstable constraints that keep zigzagging after n iteration (default is 0, no stabilization)'))
-            string="%s\n%s"%(string,fielddisplay(self,'rel_tol','tolerance of the nonlinear iteration for the transfer between layers [dimensionless]'))
-            string="%s\n%s"%(string,fielddisplay(self,'max_iter','maximum number of nonlinear iteration'))
-            string="%s\n%s"%(string,fielddisplay(self,'sedimentlimit_flag','what kind of upper limit is applied for the inefficient layer'))
-            string='%55s  0: no limit',' '
-            string='%55s  1: user defined: %s',' ','sedimentlimit'
-            string='%55s  2: hydrostatic pressure',' '
-            string='%55s  3: normal stress',' '
+		string='   hydrology Dual Porous Continuum Equivalent parameters:'
+		string='   - general parameters'
+		string="%s\n%s"%(string,fielddisplay(self,'water_compressibility','compressibility of water [Pa^-1]'))
+		string="%s\n%s"%(string,fielddisplay(self,'isefficientlayer','do we use an efficient drainage system [1: true 0: false]'))
+		string="%s\n%s"%(string,fielddisplay(self,'penalty_factor','exponent of the value used in the penalisation method [dimensionless]'))
+		string="%s\n%s"%(string,fielddisplay(self,'penalty_lock','stabilize unstable constraints that keep zigzagging after n iteration (default is 0, no stabilization)'))
+		string="%s\n%s"%(string,fielddisplay(self,'rel_tol','tolerance of the nonlinear iteration for the transfer between layers [dimensionless]'))
+		string="%s\n%s"%(string,fielddisplay(self,'max_iter','maximum number of nonlinear iteration'))
+		string="%s\n%s"%(string,fielddisplay(self,'sedimentlimit_flag','what kind of upper limit is applied for the inefficient layer'))
+		string='%55s  0: no limit',' '
+		string='%55s  1: user defined: %s',' ','sedimentlimit'
+		string='%55s  2: hydrostatic pressure',' '
+		string='%55s  3: normal stress',' '
 
-            if self.sedimentlimit_flag==1:
-                string="%s\n%s"%(string,fielddisplay(self,'sedimentlimit','user defined upper limit for the inefficient layer [m]'))
-                string="%s\n%s"%(string,fielddisplay(self,'basal_moulin_input','water flux at a given point [m3 s-1]'))
-                string="%s\n%s"%(string,fielddisplay(self,'transfer_flag','what kind of transfer method is applied between the layers'))
-                string='%55s  0: no transfer',' '
-                string='%55s  1: constant leakage factor: %s',' ','leakage_factor'
-                
-            if self.transfer_flag is 1:
-                string="%s\n%s"%(string,fielddisplay(self,'leakage_factor','user defined leakage factor [m]'))
-                string='   - for the sediment layer'
-                string="%s\n%s"%(string,fielddisplay(self,'spcsediment_head','sediment water head constraints (NaN means no constraint) [m above MSL]'))
-                string="%s\n%s"%(string,fielddisplay(self,'sediment_compressibility','sediment compressibility [Pa^-1]'))
-                string="%s\n%s"%(string,fielddisplay(self,'sediment_porosity','sediment [dimensionless]'))
-                string="%s\n%s"%(string,fielddisplay(self,'sediment_thickness','sediment thickness [m]'))
-                string="%s\n%s"%(string,fielddisplay(self,'sediment_transmitivity','sediment transmitivity [m^2/s]'))
-            
-            if self.isefficientlayer==1:
-                string='   - for the epl layer'
-                string="%s\n%s"%(string,fielddisplay(self,'spcepl_head','epl water head constraints (NaN means no constraint) [m above MSL]'))
-                string="%s\n%s"%(string,fielddisplay(self,'mask_eplactive_node','active (1) or not (0) EPL'))
-                string="%s\n%s"%(string,fielddisplay(self,'epl_compressibility','epl compressibility [Pa^-1]'))
-                string="%s\n%s"%(string,fielddisplay(self,'epl_porosity','epl [dimensionless]'))
-                string="%s\n%s"%(string,fielddisplay(self,'epl_initial_thickness','epl initial thickness [m]'))
-                string="%s\n%s"%(string,fielddisplay(self,'epl_conductivity','epl conductivity [m^2/s]'))
-            
-            #}}}
-    def setdefaultparameters(self): #{{{ 
-        
-        #Parameters from de Fleurian 2014
-        self.water_compressibility    = 5.04e-10
-        self.isefficientlayer         = 1
-        self.penalty_factor           = 3
-        self.rel_tol                  = 1.0e-06
-        self.max_iter                 = 100
-        self.sedimentlimit_flag       = 0
-        self.sedimentlimit            = 0
-        self.transfer_flag            = 0
-        self.leakage_factor           = 10.0
-        
-        self.sediment_compressibility = 1.0e-08
-        self.sediment_porosity        = 0.4
-        self.sediment_thickness       = 20.0
-        self.sediment_transmitivity   = 8.0e-04
-        
-        self.epl_compressibility      = 1.0e-08
-        self.epl_porosity             = 0.4
-        self.epl_initial_thickness    = 1.0
-        self.epl_max_thickness        = 5.0
-        self.epl_conductivity         = 8.0e-02
-        
-        return self
-    # }}}
-    
-    def initialize(self): # {{{
-        if numpy.all(numpy.isnan(self.basal_moulin_input):
-            self.basal_moulin_input=numpy.zeros(md.mesh.numberofvertices,1)
-            print"      no hydrology.basal_moulin_input specified: values set as zero"
+		if self.sedimentlimit_flag==1:
+			string="%s\n%s"%(string,fielddisplay(self,'sedimentlimit','user defined upper limit for the inefficient layer [m]'))
+			string="%s\n%s"%(string,fielddisplay(self,'basal_moulin_input','water flux at a given point [m3 s-1]'))
+			string="%s\n%s"%(string,fielddisplay(self,'transfer_flag','what kind of transfer method is applied between the layers'))
+			string='%55s  0: no transfer',' '
+			string='%55s  1: constant leakage factor: %s',' ','leakage_factor'
+			 
+		if self.transfer_flag is 1:
+			string="%s\n%s"%(string,fielddisplay(self,'leakage_factor','user defined leakage factor [m]'))
+			string='   - for the sediment layer'
+			string="%s\n%s"%(string,fielddisplay(self,'spcsediment_head','sediment water head constraints (NaN means no constraint) [m above MSL]'))
+			string="%s\n%s"%(string,fielddisplay(self,'sediment_compressibility','sediment compressibility [Pa^-1]'))
+			string="%s\n%s"%(string,fielddisplay(self,'sediment_porosity','sediment [dimensionless]'))
+			string="%s\n%s"%(string,fielddisplay(self,'sediment_thickness','sediment thickness [m]'))
+			string="%s\n%s"%(string,fielddisplay(self,'sediment_transmitivity','sediment transmitivity [m^2/s]'))
 
-        return self
+		if self.isefficientlayer==1:
+			string='   - for the epl layer'
+			string="%s\n%s"%(string,fielddisplay(self,'spcepl_head','epl water head constraints (NaN means no constraint) [m above MSL]'))
+			string="%s\n%s"%(string,fielddisplay(self,'mask_eplactive_node','active (1) or not (0) EPL'))
+			string="%s\n%s"%(string,fielddisplay(self,'epl_compressibility','epl compressibility [Pa^-1]'))
+			string="%s\n%s"%(string,fielddisplay(self,'epl_porosity','epl [dimensionless]'))
+			string="%s\n%s"%(string,fielddisplay(self,'epl_initial_thickness','epl initial thickness [m]'))
+			string="%s\n%s"%(string,fielddisplay(self,'epl_conductivity','epl conductivity [m^2/s]'))
+	#}}}
+	def setdefaultparameters(self): #{{{ 
 
-    # }}}
-    def checkconsistency(self,md,solution,analyses): #{{{ 
-        
-        #Early return
-        if HydrologyDCInefficientAnalysisEnum() not in analyses:
-            return md
+		#Parameters from de Fleurian 2014
+		self.water_compressibility    = 5.04e-10
+		self.isefficientlayer         = 1
+		self.penalty_factor           = 3
+		self.rel_tol                  = 1.0e-06
+		self.max_iter                 = 100
+		self.sedimentlimit_flag       = 0
+		self.sedimentlimit            = 0
+		self.transfer_flag            = 0
+		self.leakage_factor           = 10.0
 
-        md = checkfield(md,'fieldname','hydrology.water_compressibility','>',0,'numel',1)
-        md = checkfield(md,'fieldname','hydrology.isefficientlayer','numel',[1],'values',[0 1])
-        md = checkfield(md,'fieldname','hydrology.penalty_factor','>',0,'numel',1)
-        md = checkfield(md,'fieldname','hydrology.rel_tol','>',0,'numel',1)
-        md = checkfield(md,'fieldname','hydrology.max_iter','>',0,'numel',1)
-        md = checkfield(md,'fieldname','hydrology.sedimentlimit_flag','numel',[1],'values',[0 1 2 3])
-        md = checkfield(md,'fieldname','hydrology.transfer_flag','numel',[1],'values',[0 1])
-        
-        if self.sedimentlimit_flag==1:
-            md = checkfield(md,'fieldname','hydrology.sedimentlimit','>',0,'numel',1)
-    
-        if self.transfer_flag==1:
-            md = checkfield(md,'fieldname','hydrology.leakage_factor','>',0,'numel',1)
+		self.sediment_compressibility = 1.0e-08
+		self.sediment_porosity        = 0.4
+		self.sediment_thickness       = 20.0
+		self.sediment_transmitivity   = 8.0e-04
 
-        md = checkfield(md,'fieldname','hydrology.basal_moulin_input','NaN',1,'forcing',1)
-        md = checkfield(md,'fieldname','hydrology.spcsediment_head','forcing',1)
-        md = checkfield(md,'fieldname','hydrology.sediment_compressibility','>',0,'numel',1)
-        md = checkfield(md,'fieldname','hydrology.sediment_porosity','>',0,'numel',1)
-        md = checkfield(md,'fieldname','hydrology.sediment_thickness','>',0,'numel',1)
-        md = checkfield(md,'fieldname','hydrology.sediment_transmitivity','>=',0,'size',[md.mesh.numberofvertices 1])
-        if self.isefficientlayer==1:
-            md = checkfield(md,'fieldname','hydrology.spcepl_head','forcing',1)
-            md = checkfield(md,'fieldname','hydrology.mask_eplactive_node','size',[md.mesh.numberofvertices 1],'values',[0 1])
-            md = checkfield(md,'fieldname','hydrology.epl_compressibility','>',0,'numel',1)
-            md = checkfield(md,'fieldname','hydrology.epl_porosity','>',0,'numel',1)
-            md = checkfield(md,'fieldname','hydrology.epl_initial_thickness','>',0,'numel',1)
-            md = checkfield(md,'fieldname','hydrology.epl_conductivity','>',0,'numel',1)
+		self.epl_compressibility      = 1.0e-08
+		self.epl_porosity             = 0.4
+		self.epl_initial_thickness    = 1.0
+		self.epl_max_thickness        = 5.0
+		self.epl_conductivity         = 8.0e-02
 
-        # }}}
-    def marshall(self,md,fid): #{{{ 
-        WriteData(fid,'enum',HydrologyModelEnum(),'data',HydrologydcEnum(),'format','Integer')
-        WriteData(fid,'object',self,'fieldname','water_compressibility','format','Double')
-        WriteData(fid,'object',self,'fieldname','isefficientlayer','format','Boolean')
-        WriteData(fid,'object',self,'fieldname','penalty_factor','format','Double')
-        WriteData(fid,'object',self,'fieldname','penalty_lock','format','Integer')
-        WriteData(fid,'object',self,'fieldname','rel_tol','format','Double')
-        WriteData(fid,'object',self,'fieldname','max_iter','format','Integer')
-        WriteData(fid,'object',self,'fieldname','sedimentlimit_flag','format','Integer')
-        WriteData(fid,'object',self,'fieldname','transfer_flag','format','Integer')
+		return self
+	# }}}
 
-        if self.sedimentlimit_flag==1:
-            WriteData(fid,'object',self,'fieldname','sedimentlimit','format','Double')
+	def initialize(self): # {{{
+		if numpy.all(numpy.isnan(self.basal_moulin_input)):
+			self.basal_moulin_input=numpy.zeros(md.mesh.numberofvertices,1)
+			print"      no hydrology.basal_moulin_input specified: values set as zero"
 
-        if self.transfer_flag==1:
-            WriteData(fid,'object',self,'fieldname','leakage_factor','format','Double')
+		return self
+	# }}}
+	def checkconsistency(self,md,solution,analyses): #{{{ 
 
-	WriteData(fid,'object',self,'fieldname','basal_moulin_input','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1)
-        WriteData(fid,'object',self,'fieldname','spcsediment_head','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1)
-        WriteData(fid,'object',self,'fieldname','sediment_compressibility','format','Double')
-        WriteData(fid,'object',self,'fieldname','sediment_porosity','format','Double')			
-        WriteData(fid,'object',self,'fieldname','sediment_thickness','format','Double')
-        WriteData(fid,'object',self,'fieldname','sediment_transmitivity','format','DoubleMat','mattype',1)		
+		#Early return
+		if HydrologyDCInefficientAnalysisEnum() not in analyses:
+			return md
 
-        if self.isefficientlayer==1:	
-            WriteData(fid,'object',self,'fieldname','spcepl_head','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1)	
-            WriteData(fid,'object',self,'fieldname','mask_eplactive_node','format','DoubleMat','mattype',1)
-            WriteData(fid,'object',self,'fieldname','epl_compressibility','format','Double')			
-            WriteData(fid,'object',self,'fieldname','epl_porosity','format','Double')			
-            WriteData(fid,'object',self,'fieldname','epl_initial_thickness','format','Double')
-            WriteData(fid,'object',self,'fieldname','epl_conductivity','format','Double')
-# }}}
+		md = checkfield(md,'fieldname','hydrology.water_compressibility','>',0,'numel',1)
+		md = checkfield(md,'fieldname','hydrology.isefficientlayer','numel',[1],'values',[0,1])
+		md = checkfield(md,'fieldname','hydrology.penalty_factor','>',0,'numel',1)
+		md = checkfield(md,'fieldname','hydrology.rel_tol','>',0,'numel',1)
+		md = checkfield(md,'fieldname','hydrology.max_iter','>',0,'numel',1)
+		md = checkfield(md,'fieldname','hydrology.sedimentlimit_flag','numel',[1],'values',[0,1,2,3])
+		md = checkfield(md,'fieldname','hydrology.transfer_flag','numel',[1],'values',[0,1])
 
+		if self.sedimentlimit_flag==1:
+			md = checkfield(md,'fieldname','hydrology.sedimentlimit','>',0,'numel',1)
+
+		if self.transfer_flag==1:
+			md = checkfield(md,'fieldname','hydrology.leakage_factor','>',0,'numel',1)
+
+		md = checkfield(md,'fieldname','hydrology.basal_moulin_input','NaN',1,'forcing',1)
+		md = checkfield(md,'fieldname','hydrology.spcsediment_head','forcing',1)
+		md = checkfield(md,'fieldname','hydrology.sediment_compressibility','>',0,'numel',1)
+		md = checkfield(md,'fieldname','hydrology.sediment_porosity','>',0,'numel',1)
+		md = checkfield(md,'fieldname','hydrology.sediment_thickness','>',0,'numel',1)
+		md = checkfield(md,'fieldname','hydrology.sediment_transmitivity','>=',0,'size',[md.mesh.numberofvertices,1])
+		if self.isefficientlayer==1:
+			md = checkfield(md,'fieldname','hydrology.spcepl_head','forcing',1)
+			md = checkfield(md,'fieldname','hydrology.mask_eplactive_node','size',[md.mesh.numberofvertices,1],'values',[0,1])
+			md = checkfield(md,'fieldname','hydrology.epl_compressibility','>',0,'numel',1)
+			md = checkfield(md,'fieldname','hydrology.epl_porosity','>',0,'numel',1)
+			md = checkfield(md,'fieldname','hydrology.epl_initial_thickness','>',0,'numel',1)
+			md = checkfield(md,'fieldname','hydrology.epl_conductivity','>',0,'numel',1)
+	# }}}
+	def marshall(self,md,fid): #{{{ 
+		WriteData(fid,'enum',HydrologyModelEnum(),'data',HydrologydcEnum(),'format','Integer')
+		WriteData(fid,'object',self,'fieldname','water_compressibility','format','Double')
+		WriteData(fid,'object',self,'fieldname','isefficientlayer','format','Boolean')
+		WriteData(fid,'object',self,'fieldname','penalty_factor','format','Double')
+		WriteData(fid,'object',self,'fieldname','penalty_lock','format','Integer')
+		WriteData(fid,'object',self,'fieldname','rel_tol','format','Double')
+		WriteData(fid,'object',self,'fieldname','max_iter','format','Integer')
+		WriteData(fid,'object',self,'fieldname','sedimentlimit_flag','format','Integer')
+		WriteData(fid,'object',self,'fieldname','transfer_flag','format','Integer')
+		if self.sedimentlimit_flag==1:
+			WriteData(fid,'object',self,'fieldname','sedimentlimit','format','Double')
+
+		if self.transfer_flag==1:
+			WriteData(fid,'object',self,'fieldname','leakage_factor','format','Double')
+
+		WriteData(fid,'object',self,'fieldname','basal_moulin_input','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1)
+		WriteData(fid,'object',self,'fieldname','spcsediment_head','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1)
+		WriteData(fid,'object',self,'fieldname','sediment_compressibility','format','Double')
+		WriteData(fid,'object',self,'fieldname','sediment_porosity','format','Double')			
+		WriteData(fid,'object',self,'fieldname','sediment_thickness','format','Double')
+		WriteData(fid,'object',self,'fieldname','sediment_transmitivity','format','DoubleMat','mattype',1)		
+
+		if self.isefficientlayer==1:	
+			WriteData(fid,'object',self,'fieldname','spcepl_head','format','DoubleMat','mattype',1,'forcinglength',md.mesh.numberofvertices+1)	
+			WriteData(fid,'object',self,'fieldname','mask_eplactive_node','format','DoubleMat','mattype',1)
+			WriteData(fid,'object',self,'fieldname','epl_compressibility','format','Double')			
+			WriteData(fid,'object',self,'fieldname','epl_porosity','format','Double')			
+			WriteData(fid,'object',self,'fieldname','epl_initial_thickness','format','Double')
+			WriteData(fid,'object',self,'fieldname','epl_conductivity','format','Double')
+	# }}}
