#module imports
import numpy
from fielddisplay import fielddisplay
from EnumDefinitions import *
from checkfield import *
from MatlabFuncs import *

class mesh(object):
	"""
	MESH class definition

	   Usage:
	      mesh=mesh();
	"""

	#properties
	def __init__(self):
		# {{{ Properties
		self.x                           = float('NaN');
		self.y                           = float('NaN');
		self.z                           = float('NaN');
		self.elements                    = float('NaN');
		self.dimension                   = 0;
		self.numberoflayers              = 0;
		self.numberofelements            = 0;
		self.numberofvertices            = 0;
		self.numberofedges               = 0;
		
		self.lat                         = float('NaN');
		self.long                        = float('NaN');
		self.hemisphere                  = float('NaN');

		self.elementonbed                = float('NaN');
		self.elementonsurface            = float('NaN');
		self.vertexonbed                 = float('NaN');
		self.vertexonsurface             = float('NaN');
		self.lowerelements               = float('NaN');
		self.lowervertex                 = float('NaN');
		self.upperelements               = float('NaN');
		self.uppervertex                 = float('NaN');
		self.vertexonboundary            = float('NaN');

		self.edges                       = float('NaN');
		self.segments                    = float('NaN');
		self.segmentmarkers              = float('NaN');
		self.vertexconnectivity          = float('NaN');
		self.elementconnectivity         = float('NaN');
		self.average_vertex_connectivity = 0;

		self.x2d                         = float('NaN');
		self.y2d                         = float('NaN');
		self.elements2d                  = float('NaN');
		self.numberofvertices2d          = 0;
		self.numberofelements2d          = 0;

		self.extractedvertices           = float('NaN');
		self.extractedelements           = float('NaN');

		#set defaults
		self.setdefaultparameters()

		#}}}
	def __repr__(self):
		# {{{ Display
		string="   Mesh:" 


		if self.dimension==3:
			string="%s\n%s"%(string,"\n      Elements and vertices of the original 2d mesh:")
			
			string="%s\n%s"%(string,fielddisplay(self,"numberofelements2d","number of elements"))
			string="%s\n%s"%(string,fielddisplay(self,"numberofvertices2d","number of vertices"))
			string="%s\n%s"%(string,fielddisplay(self,"elements2d","index into (x,y,z), coordinates of the vertices"))
			string="%s\n%s"%(string,fielddisplay(self,"x2d","vertices x coordinate"))
			string="%s\n%s"%(string,fielddisplay(self,"y2d","vertices y coordinate"))

			string="%s\n%s"%(string,"\n\n      Elements and vertices of the extruded 3d mesh:")
		else:
			string="%s\n%s"%(string,"\n      Elements and vertices:")
		string="%s\n%s"%(string,fielddisplay(self,"numberofelements","number of elements"))
		string="%s\n%s"%(string,fielddisplay(self,"numberofvertices","number of vertices"))
		string="%s\n%s"%(string,fielddisplay(self,"elements","index into (x,y,z), coordinates of the vertices"))
		string="%s\n%s"%(string,fielddisplay(self,"x","vertices x coordinate"))
		string="%s\n%s"%(string,fielddisplay(self,"y","vertices y coordinate"))
		string="%s\n%s"%(string,fielddisplay(self,"z","vertices z coordinate"))
		string="%s\n%s"%(string,fielddisplay(self,"edges","edges of the 2d mesh (vertex1 vertex2 element1 element2)"))
		string="%s\n%s"%(string,fielddisplay(self,"numberofedges","number of edges of the 2d mesh"))

		string="%s%s"%(string,"\n\n      Properties:")
		string="%s\n%s"%(string,fielddisplay(self,"dimension","mesh dimension (2d or 3d)"))
		string="%s\n%s"%(string,fielddisplay(self,"numberoflayers","number of extrusion layers"))
		string="%s\n%s"%(string,fielddisplay(self,"vertexonbed","lower vertices flags list"))
		string="%s\n%s"%(string,fielddisplay(self,"elementonbed","lower elements flags list"))
		string="%s\n%s"%(string,fielddisplay(self,"vertexonsurface","upper vertices flags list"))
		string="%s\n%s"%(string,fielddisplay(self,"elementonsurface","upper elements flags list"))
		string="%s\n%s"%(string,fielddisplay(self,"uppervertex","upper vertex list (-1 for vertex on the upper surface)"))
		string="%s\n%s"%(string,fielddisplay(self,"upperelements","upper element list (-1 for element on the upper layer)"))
		string="%s\n%s"%(string,fielddisplay(self,"lowervertex","lower vertex list (-1 for vertex on the lower surface)"))
		string="%s\n%s"%(string,fielddisplay(self,"lowerelements","lower element list (-1 for element on the lower layer)"))
		string="%s\n%s"%(string,fielddisplay(self,"vertexonboundary","vertices on the boundary of the domain flag list"))
		string="%s\n%s"%(string,fielddisplay(self,"segments","edges on domain boundary (vertex1 vertex2 element)"))
		string="%s\n%s"%(string,fielddisplay(self,"segmentmarkers","number associated to each segment"))
		string="%s\n%s"%(string,fielddisplay(self,"vertexconnectivity","list of vertices connected to vertex_i"))
		string="%s\n%s"%(string,fielddisplay(self,"elementconnectivity","list of vertices connected to element_i"))
		string="%s\n%s"%(string,fielddisplay(self,"average_vertex_connectivity","average number of vertices connected to one vertex"))

		string="%s%s"%(string,"\n\n      Extracted model:")
		string="%s\n%s"%(string,fielddisplay(self,"extractedvertices","vertices extracted from the model"))
		string="%s\n%s"%(string,fielddisplay(self,"extractedelements","elements extracted from the model"))

		string="%s%s"%(string,"\n\n      Projection:")
		string="%s\n%s"%(string,fielddisplay(self,"lat","vertices latitude"))
		string="%s\n%s"%(string,fielddisplay(self,"long","vertices longitude"))
		string="%s\n%s"%(string,fielddisplay(self,"hemisphere","Indicate hemisphere 'n' or 's'"))
		return string
		#}}}
		
	def setdefaultparameters(self):
		# {{{setdefaultparameters
		
		#the connectivity is the averaged number of nodes linked to a
		#given node through an edge. This connectivity is used to initially
		#allocate memory to the stiffness matrix. A value of 16 seems to
		#give a good memory/time ration. This value can be checked in
		#trunk/test/Miscellaneous/runme.m
		self.average_vertex_connectivity=25

		return self
	#}}}

	def checkconsistency(self,md,solution,analyses):    # {{{

		md = checkfield(md,'mesh.x','NaN',1,'size',[md.mesh.numberofvertices])
		md = checkfield(md,'mesh.y','NaN',1,'size',[md.mesh.numberofvertices])
		md = checkfield(md,'mesh.z','NaN',1,'size',[md.mesh.numberofvertices])
		md = checkfield(md,'mesh.elements','NaN',1,'>',0,'values',numpy.arange(1,md.mesh.numberofvertices+1))
		if md.mesh.dimension==2:
			md = checkfield(md,'mesh.elements','size',[md.mesh.numberofelements,3])
		else:
			md = checkfield(md,'mesh.elements','size',[md.mesh.numberofelements,6])
		if numpy.any(numpy.logical_not(ismember(numpy.arange(1,md.mesh.numberofvertices+1),md.mesh.elements))):
			md.checkmessage("orphan nodes have been found. Check the mesh outline")
		md = checkfield(md,'mesh.dimension','values',[2,3])
		md = checkfield(md,'mesh.numberoflayers','>=',0)
		md = checkfield(md,'mesh.numberofelements','>',0)
		md = checkfield(md,'mesh.numberofvertices','>',0)
		#no checks for numberofedges lat long and hemisphere
		md = checkfield(md,'mesh.elementonbed','size',[md.mesh.numberofelements],'values',[0,1])
		md = checkfield(md,'mesh.elementonsurface','size',[md.mesh.numberofelements],'values',[0,1])
		md = checkfield(md,'mesh.vertexonbed','size',[md.mesh.numberofvertices],'values',[0,1])
		md = checkfield(md,'mesh.vertexonsurface','size',[md.mesh.numberofvertices],'values',[0,1])
		if md.mesh.dimension==2:
			md = checkfield(md,'mesh.average_vertex_connectivity','>=',9,'message',"'mesh.average_vertex_connectivity' should be at least 9 in 2d")
		else:
			md = checkfield(md,'mesh.average_vertex_connectivity','>=',24,'message',"'mesh.average_vertex_connectivity' should be at least 24 in 3d")
		md = checkfield(md,'mesh.elementconnectivity','size',[md.mesh.numberofelements,3],'NaN',1)

		#Solution specific checks
		if   solution==PrognosticSolutionEnum():
			if md.prognostic.stabilization==3:
				md = checkfield(md,'mesh.dimension','values',2,'message',"Discontinuous Galerkin only supported for 2d meshes")
				md = checkfield(md,'mesh.edges','size',[float('NaN'),4])
				md = checkfield(md,'mesh.edges[:,1:3]','>',0)
		elif solution==BalancethicknessSolutionEnum():
			if md.balancethickness.stabilization==3:
				md = checkfield(md,'mesh.dimension','values',2,'message',"Discontinuous Galerkin only supported for 2d meshes")
				md = checkfield(md,'mesh.edges','size',[float('NaN'),4])
				md = checkfield(md,'mesh.edges[:,1:3]','>',0)
		elif solution==TransientSolutionEnum():
			if md.transient.isprognostic and md.prognostic.stabilization==3:
				md = checkfield(md,'mesh.dimension','values',2,'message',"Discontinuous Galerkin only supported for 2d meshes")
				md = checkfield(md,'mesh.edges','size',[float('NaN'),4])
				md = checkfield(md,'mesh.edges[:,1:3]','>',0)
		elif solution==ThermalSolutionEnum():
			md = checkfield(md,'mesh.dimension','values',3,'message','thermal solution only supported on 3d meshes')

		return md
	# }}}

	def marshall(self,fid):    # {{{
		WriteData(fid,'object',self,'fieldname','x','format','DoubleMat','mattype',1)
		WriteData(fid,'object',self,'fieldname','y','format','DoubleMat','mattype',1)
		WriteData(fid,'object',self,'fieldname','z','format','DoubleMat','mattype',1)
		WriteData(fid,'object',self,'fieldname','elements','format','DoubleMat','mattype',2)
		WriteData(fid,'object',self,'fieldname','dimension','format','Integer')
		WriteData(fid,'object',self,'fieldname','numberoflayers','format','Integer')
		WriteData(fid,'object',self,'fieldname','numberofelements','format','Integer')
		WriteData(fid,'object',self,'fieldname','numberofvertices','format','Integer')
		WriteData(fid,'object',self,'fieldname','numberofedges','format','Integer')
		WriteData(fid,'object',self,'fieldname','elementonbed','format','BooleanMat','mattype',2)
		WriteData(fid,'object',self,'fieldname','elementonsurface','format','BooleanMat','mattype',2)
		WriteData(fid,'object',self,'fieldname','vertexonbed','format','BooleanMat','mattype',1)
		WriteData(fid,'object',self,'fieldname','vertexonsurface','format','BooleanMat','mattype',1)
		WriteData(fid,'object',self,'fieldname','lowerelements','format','DoubleMat','mattype',2)
		WriteData(fid,'object',self,'fieldname','upperelements','format','DoubleMat','mattype',2)
		WriteData(fid,'object',self,'fieldname','edges','format','DoubleMat','mattype',3)
		WriteData(fid,'object',self,'fieldname','elementconnectivity','format','DoubleMat','mattype',3)
		WriteData(fid,'object',self,'fieldname','average_vertex_connectivity','format','Integer')
		WriteData(fid,'object',self,'fieldname','elements2d','format','DoubleMat','mattype',3)
		WriteData(fid,'object',self,'fieldname','numberofvertices2d','format','Integer')
		WriteData(fid,'object',self,'fieldname','numberofelements2d','format','Integer')
	# }}}

