Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/model/solve.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/model/solve.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/model/solve.py (revision 12949) @@ -15,16 +15,16 @@ where varargin is a list of paired arguments of string OR enums solution types available comprise: - - DiagnosticSolutionEnum - - PrognosticSolutionEnum - - ThermalSolutionEnum - - SteadystateSolutionEnum - - TransientSolutionEnum... - - BalancethicknessSolutionEnum - - BedSlopeSolutionEnum - - SurfaceSlopeSolutionEnum - - HydrologySolutionEnum - - FlaimSolutionEnum + - DiagnosticSolutionEnum + - PrognosticSolutionEnum + - ThermalSolutionEnum + - SteadystateSolutionEnum + - TransientSolutionEnum... + - BalancethicknessSolutionEnum + - BedSlopeSolutionEnum + - SurfaceSlopeSolutionEnum + - HydrologySolutionEnum + - FlaimSolutionEnum extra options: - loadonly : does not solve. only load results Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/model/setflowequation.m =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/model/setflowequation.m (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/model/setflowequation.m (revision 12949) @@ -34,11 +34,11 @@ %Flag the elements that have not been flagged as filltype if strcmpi(filltype,'hutter'), - hutterflag(find(~macayealflag & ~pattynflag))=1; + hutterflag(find(~(macayealflag | pattynflag)))=1; elseif strcmpi(filltype,'macayeal'), - macayealflag(find(~hutterflag & ~pattynflag & ~stokesflag))=1; + macayealflag(find(~(hutterflag | pattynflag | ~stokesflag)))=1; elseif strcmpi(filltype,'pattyn'), - pattynflag(find(~hutterflag & ~macayealflag & ~stokesflag))=1; + pattynflag(find(~(hutterflag | ~macayealflag | ~stokesflag)))=1; end %check that each element has at least one flag Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/initialization.m =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/initialization.m (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/initialization.m (revision 12949) @@ -28,7 +28,7 @@ end % }}} function md = checkconsistency(obj,md,solution,analyses) % {{{ if ismember(DiagnosticHorizAnalysisEnum,analyses) - if ~isnan(md.initialization.vx) & ~isnan(md.initialization.vy), + if ~(isnan(md.initialization.vx) | isnan(md.initialization.vy)), md = checkfield(md,'initialization.vx','NaN',1,'size',[md.mesh.numberofvertices 1]); md = checkfield(md,'initialization.vy','NaN',1,'size',[md.mesh.numberofvertices 1]); end Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/surfaceforcings.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/surfaceforcings.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/surfaceforcings.py (revision 12949) @@ -1,5 +1,8 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class surfaceforcings: #properties @@ -7,6 +10,22 @@ # {{{ Properties self.precipitation = float('NaN') self.mass_balance = float('NaN') + self.ispdd = 0 + self.issmbgradients = 0 + self.isdelta18o = 0 + self.hc = float('NaN') + self.smb_pos_max = float('NaN') + self.smb_pos_min = float('NaN') + self.a_pos = float('NaN') + self.b_pos = float('NaN') + self.a_neg = float('NaN') + self.b_neg = float('NaN') + self.monthlytemperatures = float('NaN') + self.delta18o = float('NaN') + self.delta18o_surface = float('NaN') + self.temperatures_presentday = float('NaN') + self.temperatures_lgm = float('NaN') + self.precipitations_presentday = float('NaN') #set defaults self.setdefaultparameters() @@ -18,12 +37,96 @@ string="%s\n\n%s"%(string,fielddisplay(obj,'precipitation','surface precipitation [m/yr water eq]')) string="%s\n%s"%(string,fielddisplay(obj,'mass_balance','surface mass balance [m/yr ice eq]')) + string="%s\n%s"%(string,fielddisplay(obj,'ispdd','is pdd activated (0 or 1, default is 0)')) + string="%s\n%s"%(string,fielddisplay(obj,'isdelta18o','is temperature and precipitation delta18o parametrisation activated (0 or 1, default is 0)')) + string="%s\n%s"%(string,fielddisplay(obj,'monthlytemperatures','monthly surface temperatures [Kelvin], required if pdd is activated and delta18o not activated')) + string="%s\n%s"%(string,fielddisplay(obj,'precipitation','surface precipitation [m/yr water eq]')) + string="%s\n%s"%(string,fielddisplay(obj,'temperatures_presentday','monthly present day surface temperatures [Kelvin], required if pdd is activated and delta18o activated')) + string="%s\n%s"%(string,fielddisplay(obj,'temperatures_lgm','monthly LGM surface temperatures [Kelvin], required if pdd is activated and delta18o activated')) + string="%s\n%s"%(string,fielddisplay(obj,'precipitations_presentday','monthly surface precipitation [m/yr water eq], required if pdd is activated and delta18o activated')) + string="%s\n%s"%(string,fielddisplay(obj,'delta18o','delta18o, required if pdd is activated and delta18o activated')) + string="%s\n%s"%(string,fielddisplay(obj,'delta18o_surface','surface elevation of the delta18o site, required if pdd is activated and delta18o activated')) + string="%s\n%s"%(string,fielddisplay(obj,'issmbgradients','is smb gradients method activated (0 or 1, default is 0)')) + string="%s\n%s"%(string,fielddisplay(obj,'hc',' elevation of intersection between accumulation and ablation regime required if smb gradients is activated')) + string="%s\n%s"%(string,fielddisplay(obj,'smb_pos_max',' maximum value of positive smb required if smb gradients is activated')) + string="%s\n%s"%(string,fielddisplay(obj,'smb_pos_min',' minimum value of positive smb required if smb gradients is activated')) + string="%s\n%s"%(string,fielddisplay(obj,'a_pos',' intercept of hs - smb regression line for accumulation regime required if smb gradients is activated')) + string="%s\n%s"%(string,fielddisplay(obj,'b_pos',' slope of hs - smb regression line for accumulation regime required if smb gradients is activated')) + string="%s\n%s"%(string,fielddisplay(obj,'a_neg',' intercept of hs - smb regression line for ablation regime required if smb gradients is activated')) + string="%s\n%s"%(string,fielddisplay(obj,'b_neg',' slope of hs - smb regression line for ablation regime required if smb gradients is activated')) return string #}}} def setdefaultparameters(obj): # {{{setdefaultparameters + + #pdd method not used in default mode + obj.ispdd=0 + obj.issmbgradients=0 + obj.isdelta18o=0 + return obj #}}} + def checkconsistency(self,md,solution,analyses): # {{{ + + if PrognosticAnalysisEnum in analyses: + md = checkfield(md,'surfaceforcings.ispdd','numel',1,'values',[0,1]) + md = checkfield(md,'surfaceforcings.issmbgradients','numel',1,'values',[0,1]) + if self.ispdd: + if not self.isdelta18o: + md = checkfield(md,'surfaceforcings.monthlytemperatures','forcing',1,'NaN',1) + md = checkfield(md,'surfaceforcings.precipitation','forcing',1,'NaN',1) + else: + md = checkfield(md,'surfaceforcings.delta18o','NaN',1) + md = checkfield(md,'surfaceforcings.delta18o_surface','NaN',1) + md = checkfield(md,'surfaceforcings.temperatures_presentday','size',[md.mesh.numberofvertices+1,12],'NaN',1) + md = checkfield(md,'surfaceforcings.temperatures_lgm','size',[md.mesh.numberofvertices+1,12],'NaN',1) + md = checkfield(md,'surfaceforcings.precipitations_presentday','size',[md.mesh.numberofvertices+1,12],'NaN',1) + elif self.issmbgradients: + md = checkfield(md,'surfaceforcings.hc','forcing',1,'NaN',1) + md = checkfield(md,'surfaceforcings.smb_pos_max','forcing',1,'NaN',1) + md = checkfield(md,'surfaceforcings.smb_pos_min','forcing',1,'NaN',1) + md = checkfield(md,'surfaceforcings.a_pos','forcing',1,'NaN',1) + md = checkfield(md,'surfaceforcings.b_pos','forcing',1,'NaN',1) + md = checkfield(md,'surfaceforcings.a_neg','forcing',1,'NaN',1) + md = checkfield(md,'surfaceforcings.b_neg','forcing',1,'NaN',1) + else: + md = checkfield(md,'surfaceforcings.mass_balance','forcing',1,'NaN',1) + + if BalancethicknessAnalysisEnum in analyses: + md = checkfield(md,'surfaceforcings.mass_balance','size',[md.mesh.numberofvertices],'NaN',1) + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','precipitation','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','mass_balance','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','ispdd','format','Boolean') + WriteData(fid,'object',self,'fieldname','isdelta18o','format','Boolean') + + if self.ispdd: + if self.isdelta18o: + WriteData(fid,'object',self,'fieldname','temperatures_presentday','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','temperatures_lgm','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','precipitations_presentday','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','delta18o_surface','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','delta18o','format','DoubleMat','mattype',1) + else: + WriteData(fid,'object',self,'fieldname','monthlytemperatures','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','precipitation','format','DoubleMat','mattype',1) + + WriteData(fid,'object',self,'fieldname','issmbgradients','format','Boolean') + + if self.issmbgradients: + WriteData(fid,'object',self,'fieldname','hc','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','smb_pos_max','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','smb_pos_min','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','a_pos','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','b_pos','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','a_neg','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','b_neg','format','DoubleMat','mattype',1) + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/mesh.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/mesh.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/mesh.py (revision 12949) @@ -1,8 +1,8 @@ #module imports import numpy from fielddisplay import fielddisplay +from EnumDefinitions import * from checkfield import * -from EnumDefinitions import * from MatlabFuncs import * class mesh: @@ -130,7 +130,7 @@ else: md = checkfield(md,'mesh.elements','size',[md.mesh.numberofelements,6]) if any(numpy.logical_not(ismember(range(1,md.mesh.numberofvertices+1),md.mesh.elements))): - md = checkmessage(md,"orphan nodes have been found. Check the mesh outline") + 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) Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/timestepping.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/timestepping.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/timestepping.py (revision 12949) @@ -1,12 +1,16 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class timestepping: #properties def __init__(self): # {{{ Properties + self.start_time = 0 + self.final_time = 0 self.time_step = 0 - self.final_time = 0 self.time_adapt = 0 self.cfl_coefficient = 0 @@ -17,8 +21,9 @@ def __repr__(obj): # {{{ Display string=" timestepping parameters:" - string="%s\n\n%s"%(string,fielddisplay(obj,"time_step","length of time steps [yrs]")) + string="%s\n\n%s"%(string,fielddisplay(obj,"start_time","simulation starting time [yrs]")) string="%s\n%s"%(string,fielddisplay(obj,"final_time","final time to stop the simulation [yrs]")) + string="%s\n%s"%(string,fielddisplay(obj,"time_step","length of time steps [yrs]")) string="%s\n%s"%(string,fielddisplay(obj,"time_adapt","use cfl condition to define time step ? (0 or 1) ")) string="%s\n%s"%(string,fielddisplay(obj,"cfl_coefficient","coefficient applied to cfl condition")) return string @@ -39,3 +44,25 @@ return obj #}}} + + def checkconsistency(self,md,solution,analyses): # {{{ + + md = checkfield(md,'timestepping.start_time','numel',1,'NaN',1) + md = checkfield(md,'timestepping.final_time','numel',1,'NaN',1) + md = checkfield(md,'timestepping.time_step','numel',1,'>=',0,'NaN',1) + md = checkfield(md,'timestepping.time_adapt','numel',1,'values',[0,1]) + md = checkfield(md,'timestepping.cfl_coefficient','numel',1,'>',0,'<=',1) + if self.final_time-self.start_time<0: + md.checkmessage("timestepping.final_time should be larger than timestepping.start_time") + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','start_time','format','Double') + WriteData(fid,'object',self,'fieldname','final_time','format','Double') + WriteData(fid,'object',self,'fieldname','time_step','format','Double') + WriteData(fid,'object',self,'fieldname','time_adapt','format','Boolean') + WriteData(fid,'object',self,'fieldname','cfl_coefficient','format','Double') + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/constants.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/constants.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/constants.py (revision 12949) @@ -1,5 +1,8 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class constants: #properties @@ -39,3 +42,18 @@ return obj #}}} + def checkconsistency(self,md,solution,analyses): # {{{ + + md = checkfield(md,'constants.g','>',0,'size',[1]) + md = checkfield(md,'constants.yts','>',0,'size',[1]) + md = checkfield(md,'constants.referencetemperature','size',[1]) + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','g','format','Double') + WriteData(fid,'object',self,'fieldname','yts','format','Double') + WriteData(fid,'object',self,'fieldname','referencetemperature','format','Double') + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/friction.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/friction.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/friction.py (revision 12949) @@ -1,5 +1,8 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class friction: #properties @@ -27,3 +30,22 @@ return obj #}}} + def checkconsistency(self,md,solution,analyses): # {{{ + + #Early return + if not DiagnosticHorizAnalysisEnum in analyses and not ThermalAnalysisEnum in analyses: + return md + + md = checkfield(md,'friction.coefficient','NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'friction.q','NaN',1,'size',[md.mesh.numberofelements]) + md = checkfield(md,'friction.p','NaN',1,'size',[md.mesh.numberofelements]) + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','coefficient','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','p','format','DoubleMat','mattype',2) + WriteData(fid,'object',self,'fieldname','q','format','DoubleMat','mattype',2) + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/materials.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/materials.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/materials.py (revision 12949) @@ -1,23 +1,26 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class materials: #properties def __init__(self): # {{{ Properties - self.rho_ice = 0; - self.rho_water = 0; - self.mu_water = 0; - self.heatcapacity = 0; - self.latentheat = 0; - self.thermalconductivity = 0; - self.meltingpoint = 0; - self.beta = 0; - self.mixed_layer_capacity = 0; - self.thermal_exchange_velocity = 0; + self.rho_ice = 0 + self.rho_water = 0 + self.mu_water = 0 + self.heatcapacity = 0 + self.latentheat = 0 + self.thermalconductivity = 0 + self.meltingpoint = 0 + self.beta = 0 + self.mixed_layer_capacity = 0 + self.thermal_exchange_velocity = 0 self.rheology_B = float('NaN') self.rheology_n = float('NaN') - self.rheology_law = ""; + self.rheology_law = "" self.setdefaultparameters() #}}} @@ -78,3 +81,32 @@ obj.rheology_law='Paterson' return obj #}}} + + def checkconsistency(self,md,solution,analyses): # {{{ + md = checkfield(md,'materials.rho_ice','>',0) + md = checkfield(md,'materials.rho_water','>',0) + md = checkfield(md,'materials.rho_freshwater','>',0) + md = checkfield(md,'materials.mu_water','>',0) + md = checkfield(md,'materials.rheology_B','>',0,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'materials.rheology_n','>',0,'size',[md.mesh.numberofelements]) + md = checkfield(md,'materials.rheology_law','values',['None','Paterson','Arrhenius']) + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','rho_ice','format','Double') + WriteData(fid,'object',self,'fieldname','rho_water','format','Double') + WriteData(fid,'object',self,'fieldname','rho_freshwater','format','Double') + WriteData(fid,'object',self,'fieldname','mu_water','format','Double') + WriteData(fid,'object',self,'fieldname','heatcapacity','format','Double') + WriteData(fid,'object',self,'fieldname','latentheat','format','Double') + WriteData(fid,'object',self,'fieldname','thermalconductivity','format','Double') + WriteData(fid,'object',self,'fieldname','meltingpoint','format','Double') + WriteData(fid,'object',self,'fieldname','beta','format','Double') + WriteData(fid,'object',self,'fieldname','mixed_layer_capacity','format','Double') + WriteData(fid,'object',self,'fieldname','thermal_exchange_velocity','format','Double') + WriteData(fid,'object',self,'fieldname','rheology_B','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','rheology_n','format','DoubleMat','mattype',2) + WriteData(fid,'data',StringToEnum(self.rheology_law),'enum',MaterialsRheologyLawEnum,'format','Integer') + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/initialization.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/initialization.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/initialization.py (revision 12949) @@ -1,5 +1,9 @@ #module imports +import numpy from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class initialization: #properties @@ -40,3 +44,41 @@ return obj #}}} + def checkconsistency(obj,md,solution,analyses): # {{{ + if DiagnosticHorizAnalysisEnum in analyses: + if not any(numpy.isnan(md.initialization.vx) or numpy.isnan(md.initialization.vy)): + md = checkfield(md,'initialization.vx','NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'initialization.vy','NaN',1,'size',[md.mesh.numberofvertices]) + if PrognosticAnalysisEnum in analyses: + md = checkfield(md,'initialization.vx','NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'initialization.vy','NaN',1,'size',[md.mesh.numberofvertices]) + if HydrologyAnalysisEnum in analyses: + md = checkfield(md,'initialization.watercolumn','NaN',1,'size',[md.mesh.numberofvertices]) + if BalancethicknessAnalysisEnum in analyses: + md = checkfield(md,'initialization.vx','NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'initialization.vy','NaN',1,'size',[md.mesh.numberofvertices]) + #Triangle with zero velocity + if any(numpy.logical_and(numpy.sum(numpy.abs(md.initialization.vx[md.mesh.elements.astype(int)-1]),axis=1)==0,\ + numpy.sum(numpy.abs(md.initialization.vy[md.mesh.elements.astype(int)-1]),axis=1)==0)): + md.checkmessage("at least one triangle has all its vertices with a zero velocity") + if ThermalAnalysisEnum in analyses: + md = checkfield(md,'initialization.vx','NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'initialization.vy','NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'initialization.vz','NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'initialization.pressure','NaN',1,'size',[md.mesh.numberofvertices]) + if (EnthalpyAnalysisEnum in analyses and md.thermal.isenthalpy) or solution==EnthalpySolutionEnum: + md = checkfield(md,'initialization.waterfraction','>=',0,'size',[md.mesh.numberofvertices]) + + return md + # }}} + + def marshall(obj,fid): # {{{ + WriteData(fid,'data',obj.vx,'format','DoubleMat','mattype',1,'enum',VxEnum) + WriteData(fid,'data',obj.vy,'format','DoubleMat','mattype',1,'enum',VyEnum) + WriteData(fid,'data',obj.vz,'format','DoubleMat','mattype',1,'enum',VzEnum) + WriteData(fid,'data',obj.pressure,'format','DoubleMat','mattype',1,'enum',PressureEnum) + WriteData(fid,'data',obj.temperature,'format','DoubleMat','mattype',1,'enum',TemperatureEnum) + WriteData(fid,'data',obj.watercolumn,'format','DoubleMat','mattype',1,'enum',WatercolumnEnum) + WriteData(fid,'data',obj.waterfraction,'format','DoubleMat','mattype',1,'enum',WaterfractionEnum) + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/mask.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/mask.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/mask.py (revision 12949) @@ -1,5 +1,8 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class mask: #properties @@ -34,3 +37,24 @@ return obj #}}} + def checkconsistency(self,md,solution,analyses): # {{{ + + md = checkfield(md,'mask.elementonfloatingice','size',[md.mesh.numberofelements],'values',[0,1]) + md = checkfield(md,'mask.elementongroundedice','size',[md.mesh.numberofelements],'values',[0,1]) + md = checkfield(md,'mask.elementonwater' ,'size',[md.mesh.numberofelements],'values',[0,1]) + md = checkfield(md,'mask.vertexonfloatingice' ,'size',[md.mesh.numberofvertices],'values',[0,1]) + md = checkfield(md,'mask.vertexongroundedice' ,'size',[md.mesh.numberofvertices],'values',[0,1]) + md = checkfield(md,'mask.vertexonwater' ,'size',[md.mesh.numberofvertices],'values',[0,1]) + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','elementonfloatingice','format','BooleanMat','mattype',2) + WriteData(fid,'object',self,'fieldname','elementongroundedice','format','BooleanMat','mattype',2) + WriteData(fid,'object',self,'fieldname','elementonwater','format','BooleanMat','mattype',2) + WriteData(fid,'object',self,'fieldname','vertexonfloatingice','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','vertexongroundedice','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','vertexonwater','format','DoubleMat','mattype',1) + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/surfaceforcings.m =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/surfaceforcings.m (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/surfaceforcings.m (revision 12949) @@ -10,7 +10,7 @@ ispdd = 0; issmbgradients = 0; isdelta18o = 0; - hc = NaN; + hc = NaN; smb_pos_max = NaN; smb_pos_min = NaN; a_pos = NaN; @@ -45,26 +45,26 @@ if ismember(PrognosticAnalysisEnum,analyses), md = checkfield(md,'surfaceforcings.ispdd','numel',1,'values',[0 1]); - checkfield(md,'surfaceforcings.issmbgradients','numel',1,'values',[0 1]); + md = checkfield(md,'surfaceforcings.issmbgradients','numel',1,'values',[0 1]); if(obj.ispdd) - if(obj.isdelta18o==0) - md = checkfield(md,'surfaceforcings.monthlytemperatures','forcing',1,'NaN',1); - md = checkfield(md,'surfaceforcings.precipitation','forcing',1,'NaN',1); - else - md = checkfield(md,'surfaceforcings.delta18o','NaN',1); - md = checkfield(md,'surfaceforcings.delta18o_surface','NaN',1); - md = checkfield(md,'surfaceforcings.temperatures_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1); - md = checkfield(md,'surfaceforcings.temperatures_lgm','size',[md.mesh.numberofvertices+1 12],'NaN',1); - md = checkfield(md,'surfaceforcings.precipitations_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1); - end + if(obj.isdelta18o==0) + md = checkfield(md,'surfaceforcings.monthlytemperatures','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.precipitation','forcing',1,'NaN',1); + else + md = checkfield(md,'surfaceforcings.delta18o','NaN',1); + md = checkfield(md,'surfaceforcings.delta18o_surface','NaN',1); + md = checkfield(md,'surfaceforcings.temperatures_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1); + md = checkfield(md,'surfaceforcings.temperatures_lgm','size',[md.mesh.numberofvertices+1 12],'NaN',1); + md = checkfield(md,'surfaceforcings.precipitations_presentday','size',[md.mesh.numberofvertices+1 12],'NaN',1); + end elseif(obj.issmbgradients) - checkfield(md,'surfaceforcings.hc','forcing',1,'NaN',1); - checkfield(md,'surfaceforcings.smb_pos_max','forcing',1,'NaN',1); - checkfield(md,'surfaceforcings.smb_pos_min','forcing',1,'NaN',1); - checkfield(md,'surfaceforcings.a_pos','forcing',1,'NaN',1); - checkfield(md,'surfaceforcings.b_pos','forcing',1,'NaN',1); - checkfield(md,'surfaceforcings.a_neg','forcing',1,'NaN',1); - checkfield(md,'surfaceforcings.b_neg','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.hc','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.smb_pos_max','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.smb_pos_min','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.a_pos','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.b_pos','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.a_neg','forcing',1,'NaN',1); + md = checkfield(md,'surfaceforcings.b_neg','forcing',1,'NaN',1); else md = checkfield(md,'surfaceforcings.mass_balance','forcing',1,'NaN',1); end @@ -102,16 +102,16 @@ WriteData(fid,'object',obj,'fieldname','ispdd','format','Boolean'); WriteData(fid,'object',obj,'fieldname','isdelta18o','format','Boolean'); if obj.ispdd, - if obj.isdelta18o - WriteData(fid,'object',obj,'fieldname','temperatures_presentday','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','temperatures_lgm','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','precipitations_presentday','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','delta18o_surface','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','delta18o','format','DoubleMat','mattype',1); - else - WriteData(fid,'object',obj,'fieldname','monthlytemperatures','format','DoubleMat','mattype',1); - WriteData(fid,'object',obj,'fieldname','precipitation','format','DoubleMat','mattype',1); - end + if obj.isdelta18o + WriteData(fid,'object',obj,'fieldname','temperatures_presentday','format','DoubleMat','mattype',1); + WriteData(fid,'object',obj,'fieldname','temperatures_lgm','format','DoubleMat','mattype',1); + WriteData(fid,'object',obj,'fieldname','precipitations_presentday','format','DoubleMat','mattype',1); + WriteData(fid,'object',obj,'fieldname','delta18o_surface','format','DoubleMat','mattype',1); + WriteData(fid,'object',obj,'fieldname','delta18o','format','DoubleMat','mattype',1); + else + WriteData(fid,'object',obj,'fieldname','monthlytemperatures','format','DoubleMat','mattype',1); + WriteData(fid,'object',obj,'fieldname','precipitation','format','DoubleMat','mattype',1); + end end WriteData(fid,'object',obj,'fieldname','issmbgradients','format','Boolean'); if obj.issmbgradients, Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/geometry.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/geometry.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/geometry.py (revision 12949) @@ -1,5 +1,8 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class geometry: #properties @@ -33,3 +36,25 @@ return obj #}}} + def checkconsistency(self,md,solution,analyses): # {{{ + + md = checkfield(md,'geometry.surface' ,'NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'geometry.bed' ,'NaN',1,'size',[md.mesh.numberofvertices]) + md = checkfield(md,'geometry.thickness','NaN',1,'size',[md.mesh.numberofvertices],'>',0) + if any((self.thickness-self.surface+self.bed)>10**-9): + md.checkmessage("equality thickness=surface-bed violated") + end + if solution==TransientSolutionEnum and md.transient.isgroundingline: + md = checkfield(md,'geometry.bathymetry','NaN',1,'size',[md.mesh.numberofvertices]) + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'data',self.surface,'format','DoubleMat','mattype',1,'enum',SurfaceEnum) + WriteData(fid,'data',self.thickness,'format','DoubleMat','mattype',1,'enum',ThicknessEnum) + WriteData(fid,'data',self.bed,'format','DoubleMat','mattype',1,'enum',BedEnum) + WriteData(fid,'data',self.bathymetry,'format','DoubleMat','mattype',1,'enum',BathymetryEnum) + WriteData(fid,'object',self,'fieldname','hydrostatic_ratio','format','DoubleMat','mattype',1) + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/private.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/private.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/private.py (revision 12949) @@ -5,10 +5,10 @@ #properties def __init__(self): # {{{ Properties - self.isconsistent = True; + self.isconsistent = True self.runtimename = '' self.bamg = {} - self.solution = ''; + self.solution = '' #set defaults self.setdefaultparameters() Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/basalforcings.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/basalforcings.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/basalforcings.py (revision 12949) @@ -1,5 +1,8 @@ #module imports from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class basalforcings: #properties @@ -27,3 +30,24 @@ return obj #}}} + def checkconsistency(self,md,solution,analyses): # {{{ + + if PrognosticAnalysisEnum in analyses and not (solution==TransientSolutionEnum and not md.transient.isprognostic): + md = checkfield(md,'basalforcings.melting_rate','NaN',1,'forcing',1) + + if BalancethicknessAnalysisEnum in analyses: + md = checkfield(md,'basalforcings.melting_rate','NaN',1,'size',[md.mesh.numberofvertices]) + + if ThermalAnalysisEnum in analyses and not (solution==TransientSolutionEnum and not md.transient.isthermal): + md = checkfield(md,'basalforcings.melting_rate','NaN',1,'forcing',1) + md = checkfield(md,'basalforcings.geothermalflux','NaN',1,'forcing',1,'>=',0) + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','melting_rate','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','melting_rate_correction','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','geothermalflux','format','DoubleMat','mattype',1) + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/flowequation.py =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/flowequation.py (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/flowequation.py (revision 12949) @@ -1,5 +1,9 @@ #module imports +import numpy from fielddisplay import fielddisplay +from EnumDefinitions import * +from checkfield import * +from WriteData import * class flowequation: #properties @@ -39,3 +43,59 @@ return obj #}}} + def checkconsistency(self,md,solution,analyses): # {{{ + + if DiagnosticHorizAnalysisEnum in analyses: + md = checkfield(md,'flowequation.ismacayealpattyn','numel',1,'values',[0,1]) + md = checkfield(md,'flowequation.ishutter','numel',1,'values',[0,1]) + md = checkfield(md,'flowequation.isstokes','numel',1,'values',[0,1]) + md = checkfield(md,'flowequation.bordermacayeal','size',[md.mesh.numberofvertices],'values',[0,1]) + md = checkfield(md,'flowequation.borderpattyn','size',[md.mesh.numberofvertices],'values',[0,1]) + md = checkfield(md,'flowequation.borderstokes','size',[md.mesh.numberofvertices],'values',[0,1]) + if md.mesh.dimension==2: + md = checkfield(md,'flowequation.vertex_equation','size',[md.mesh.numberofvertices],'values',[1,2]) + md = checkfield(md,'flowequation.element_equation','size',[md.mesh.numberofelements],'values',[1,2]) + else: + md = checkfield(md,'flowequation.vertex_equation','size',[md.mesh.numberofvertices],'values',range(0,7+1)) + md = checkfield(md,'flowequation.element_equation','size',[md.mesh.numberofelements],'values',range(0,7+1)) + if not (md.flowequation.ismacayealpattyn or md.flowequation.ishutter or md.flowequation.isstokes): + md.checkmessage("no element types set for this model. At least one of ismacayealpattyn, ishutter or isstokes need to be =1") + + if DiagnosticHutterAnalysisEnum in analyses: + if any(md.flowequation.element_equation==1): + if any(numpy.logical_and(md.flowequation.element_equation,md.mask.elementonfloatingice)): + print "\n !!! Warning: Hutter's model is not consistent on ice shelves !!!\n" + + return md + # }}} + + def marshall(self,fid): # {{{ + WriteData(fid,'object',self,'fieldname','ismacayealpattyn','format','Boolean') + WriteData(fid,'object',self,'fieldname','ishutter','format','Boolean') + WriteData(fid,'object',self,'fieldname','isstokes','format','Boolean') + WriteData(fid,'object',self,'fieldname','bordermacayeal','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','borderpattyn','format','DoubleMat','mattype',1) + WriteData(fid,'object',self,'fieldname','borderstokes','format','DoubleMat','mattype',1) + #convert approximations to enums + data=self.vertex_equation + data[[i for i,item in enumerate(data) if item==0]]=NoneApproximationEnum + data[[i for i,item in enumerate(data) if item==1]]=HutterApproximationEnum + data[[i for i,item in enumerate(data) if item==2]]=MacAyealApproximationEnum + data[[i for i,item in enumerate(data) if item==3]]=PattynApproximationEnum + data[[i for i,item in enumerate(data) if item==4]]=StokesApproximationEnum + data[[i for i,item in enumerate(data) if item==5]]=MacAyealPattynApproximationEnum + data[[i for i,item in enumerate(data) if item==6]]=MacAyealStokesApproximationEnum + data[[i for i,item in enumerate(data) if item==7]]=PattynStokesApproximationEnum + WriteData(fid,'data',data,'enum',FlowequationVertexEquationEnum,'format','DoubleMat','mattype',1) + data=self.element_equation + data[[i for i,item in enumerate(data) if item==0]]=NoneApproximationEnum + data[[i for i,item in enumerate(data) if item==1]]=HutterApproximationEnum + data[[i for i,item in enumerate(data) if item==2]]=MacAyealApproximationEnum + data[[i for i,item in enumerate(data) if item==3]]=PattynApproximationEnum + data[[i for i,item in enumerate(data) if item==4]]=StokesApproximationEnum + data[[i for i,item in enumerate(data) if item==5]]=MacAyealPattynApproximationEnum + data[[i for i,item in enumerate(data) if item==6]]=MacAyealStokesApproximationEnum + data[[i for i,item in enumerate(data) if item==7]]=PattynStokesApproximationEnum + WriteData(fid,'data',data,'enum',FlowequationElementEquationEnum,'format','DoubleMat','mattype',2) + # }}} + Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/flowequation.m =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/flowequation.m (revision 12948) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/m/classes/flowequation.m (revision 12949) @@ -43,7 +43,7 @@ md = checkfield(md,'flowequation.vertex_equation','size',[md.mesh.numberofvertices 1],'values',[0:7]); md = checkfield(md,'flowequation.element_equation','size',[md.mesh.numberofelements 1],'values',[0:7]); end - if (md.flowequation.ismacayealpattyn==0 && md.flowequation.ishutter==0 && md.flowequation.isstokes==0), + if ~(md.flowequation.ismacayealpattyn || md.flowequation.ishutter || md.flowequation.isstokes), md = checkmessage(md,['no element types set for this model. At least one of ismacayealpattyn, ishutter or isstokes need to be =1']); end end