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

class thermal(object):
	"""
	THERMAL class definition

	   Usage:
	      thermal=thermal();
	"""

	#properties
	def __init__(self):
		# {{{ Properties
		self.spctemperature    = float('NaN')
		self.penalty_threshold = 0
		self.stabilization     = 0
		self.maxiter           = 0
		self.penalty_lock      = 0
		self.penalty_factor    = 0
		self.isenthalpy        = 0

		#set defaults
		self.setdefaultparameters()

		#}}}
	def __repr__(self):
		# {{{ Display
		string='   Thermal solution parameters:'
		string="%s\n%s"%(string,fielddisplay(self,'spctemperature','temperature constraints (NaN means no constraint)'))
		string="%s\n%s"%(string,fielddisplay(self,'stabilization','0->no, 1->artificial_diffusivity, 2->SUPG'))
		string="%s\n%s"%(string,fielddisplay(self,'maxiter','maximum number of non linear iterations'))
		string="%s\n%s"%(string,fielddisplay(self,'penalty_lock','stabilize unstable thermal constraints that keep zigzagging after n iteration (default is 0, no stabilization)'))
		string="%s\n%s"%(string,fielddisplay(self,'penalty_threshold','threshold to declare convergence of thermal solution (default is 0)'))
		string="%s\n%s"%(string,fielddisplay(self,'isenthalpy','use an enthalpy formulation to include temperate ice (default is 0)'))
		return string
		#}}}
		
	def setdefaultparameters(self):
		# {{{setdefaultparameters
		
		#Number of unstable constraints acceptable
		self.penalty_threshold=0

		#Type of stabilization used
		self.stabilization=1

		#Maximum number of iterations
		self.maxiter=100

		#factor used to compute the values of the penalties: kappa=max(stiffness matrix)*10^penalty_factor
		self.penalty_factor=3

		#Should we use cold ice (default) or enthalpy formulation
		self.isenthalpy=0

		return self
	#}}}

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

		#Early return
		if (ThermalAnalysisEnum() not in analyses and EnthalpyAnalysisEnum() not in analyses) or (solution==TransientSolutionEnum() and not md.transient.isthermal):
			return md

		md = checkfield(md,'thermal.stabilization','numel',[1],'values',[0,1,2])
		md = checkfield(md,'thermal.spctemperature','forcing',1)
		if EnthalpyAnalysisEnum() in analyses and (md.thermal.isenthalpy or solution==EnthalpySolutionEnum()) and md.mesh.dimension==3:
			pos=numpy.nonzero(numpy.logical_not(numpy.isnan(md.thermal.spctemperature[0:md.mesh.numberofvertices])))
			replicate=numpy.tile(md.geometry.surface-md.mesh.z,(1,numpy.size(md.thermal.spctemperature,axis=1)))
			md = checkfield(md,'thermal.spctemperature[numpy.nonzero(numpy.logical_not(numpy.isnan(md.thermal.spctemperature[0:md.mesh.numberofvertices,:])))]','<',md.materials.meltingpoint-md.materials.beta*md.materials.rho_ice*md.constants.g*replicate[pos],'message',"spctemperature should be below the adjusted melting point")
			md = checkfield(md,'thermal.isenthalpy','numel',[1],'values',[0,1])

		return md
	# }}}

	def marshall(self,fid):    # {{{
		WriteData(fid,'object',self,'fieldname','spctemperature','format','DoubleMat','mattype',1)
		WriteData(fid,'object',self,'fieldname','penalty_threshold','format','Integer')
		WriteData(fid,'object',self,'fieldname','stabilization','format','Integer')
		WriteData(fid,'object',self,'fieldname','maxiter','format','Integer')
		WriteData(fid,'object',self,'fieldname','penalty_lock','format','Integer')
		WriteData(fid,'object',self,'fieldname','penalty_factor','format','Double')
		WriteData(fid,'object',self,'fieldname','isenthalpy','format','Boolean')
	# }}}

