Index: /issm/trunk-jpl/src/m/classes/matenhancedice.py =================================================================== --- /issm/trunk-jpl/src/m/classes/matenhancedice.py (revision 25170) +++ /issm/trunk-jpl/src/m/classes/matenhancedice.py (revision 25171) @@ -9,9 +9,9 @@ MATICE class definition - Usage: - matenhancedice = matenhancedice() + Usage: + matenhancedice = matenhancedice() """ - def __init__(self): # {{{ + def __init__(self): #{{{ self.rho_ice = 0. self.rho_water = 0. @@ -32,5 +32,5 @@ self.rheology_law = '' - #giaivins: + #GIA self.lithosphere_shear_modulus = 0. self.lithosphere_density = 0. @@ -38,37 +38,38 @@ self.mantle_density = 0. - #SLR - self.earth_density = 0 # average density of the Earth, (kg / m^3) + #SLR + self.earth_density = 0 # average density of the Earth, (kg/m^3) self.setdefaultparameters() #}}} - def __repr__(self): # {{{ - string = " Materials:" - string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "water density [kg / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]")) - string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "effectiveconductivity_averaging", "computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)")) - string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K")) - string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]")) - string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]")) - string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]")) - string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]")) - string = "%s\n%s" % (string, fielddisplay(self, "rheology_E", "enhancement factor")) - string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]")) - string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent")) - string = "%s\n%s" % (string, fielddisplay(self, "rheology_law", "law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius' or 'LliboutryDuval'")) - string = "%s\n%s" % (string, fielddisplay(self, "lithosphere_shear_modulus", "Lithosphere shear modulus [Pa]")) - string = "%s\n%s" % (string, fielddisplay(self, "lithosphere_density", "Lithosphere density [g / cm^ - 3]")) - string = "%s\n%s" % (string, fielddisplay(self, "mantle_shear_modulus", "Mantle shear modulus [Pa]")) - string = "%s\n%s" % (string, fielddisplay(self, "mantle_density", "Mantle density [g / cm^ - 3]")) - string = "%s\n%s" % (string, fielddisplay(self, "earth_density", "Mantle density [kg / m^ - 3]")) - return string + def __repr__(self): #{{{ + s = " Materials:" + s = "%s\n%s" % (s, fielddisplay(self, "rho_ice", "ice density [kg/m^3]")) + s = "%s\n%s" % (s, fielddisplay(self, "rho_water", "water density [kg/m^3]")) + s = "%s\n%s" % (s, fielddisplay(self, "rho_freshwater", "fresh water density [kg/m^3]")) + s = "%s\n%s" % (s, fielddisplay(self, "mu_water", "water viscosity [N s/m^2]")) + s = "%s\n%s" % (s, fielddisplay(self, "heatcapacity", "heat capacity [J/kg/K]")) + s = "%s\n%s" % (s, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W/m/K]")) + s = "%s\n%s" % (s, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W/m/K]")) + s = "%s\n%s" % (s, fielddisplay(self, "effectiveconductivity_averaging", "computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)")) + s = "%s\n%s" % (s, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K")) + s = "%s\n%s" % (s, fielddisplay(self, "latentheat", "latent heat of fusion [J/m^3]")) + s = "%s\n%s" % (s, fielddisplay(self, "beta", "rate of change of melting point with pressure [K/Pa]")) + s = "%s\n%s" % (s, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W/kg/K]")) + s = "%s\n%s" % (s, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m/s]")) + s = "%s\n%s" % (s, fielddisplay(self, "rheology_E", "enhancement factor")) + s = "%s\n%s" % (s, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1/n)]")) + s = "%s\n%s" % (s, fielddisplay(self, "rheology_n", "Glen's flow law exponent")) + s = "%s\n%s" % (s, fielddisplay(self, "rheology_law", "law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius' or 'LliboutryDuval'")) + s = "%s\n%s" % (s, fielddisplay(self, "lithosphere_shear_modulus", "Lithosphere shear modulus [Pa]")) + s = "%s\n%s" % (s, fielddisplay(self, "lithosphere_density", "Lithosphere density [g/cm^-3]")) + s = "%s\n%s" % (s, fielddisplay(self, "mantle_shear_modulus", "Mantle shear modulus [Pa]")) + s = "%s\n%s" % (s, fielddisplay(self, "mantle_density", "Mantle density [g/cm^-3]")) + s = "%s\n%s" % (s, fielddisplay(self, "earth_density", "Mantle density [kg/m^-3]")) + + return s #}}} - def extrude(self, md): # {{{ + def extrude(self, md): #{{{ self.rheology_E = project3d(md, 'vector', self.rheology_E, 'type', 'node') self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node') @@ -77,5 +78,5 @@ #}}} - def setdefaultparameters(self): # {{{ + def setdefaultparameters(self): #{{{ #ice density (kg / m^3) self.rho_ice = 917. @@ -120,5 +121,5 @@ #}}} - def checkconsistency(self, md, solution, analyses): # {{{ + def checkconsistency(self, md, solution, analyses): #{{{ md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0) md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0) @@ -136,4 +137,5 @@ md = checkfield(md, 'fieldname', 'materials.mantle_shear_modulus', '>', 0, 'numel', 1) md = checkfield(md, 'fieldname', 'materials.mantle_density', '>', 0, 'numel', 1) + if 'SealevelriseAnalysis' in analyses: md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', 1) @@ -141,5 +143,5 @@ # }}} - def marshall(self, prefix, md, fid): # {{{ + def marshall(self, prefix, md, fid): #{{{ WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 4, 'format', 'Integer') WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_ice', 'format', 'Double') Index: /issm/trunk-jpl/src/m/classes/matestar.py =================================================================== --- /issm/trunk-jpl/src/m/classes/matestar.py (revision 25170) +++ /issm/trunk-jpl/src/m/classes/matestar.py (revision 25171) @@ -1,17 +1,19 @@ +import numpy as np + +from checkfield import checkfield from fielddisplay import fielddisplay from project3d import project3d -from checkfield import checkfield from WriteData import WriteData class matestar(object): - """ - matestar class definition + ''' + MATESTAR class definition - Usage: - matestar = matestar() - """ + Usage: + matestar = matestar() + ''' - def __init__(self): # {{{ + def __init__(self): #{{{ self.rho_ice = 0. self.rho_water = 0. @@ -27,10 +29,10 @@ self.mixed_layer_capacity = 0. self.thermal_exchange_velocity = 0. - self.rheology_B = float('NaN') - self.rheology_Ec = float('NaN') - self.rheology_Es = float('NaN') + self.rheology_B = np.nan + self.rheology_Ec = np.nan + self.rheology_Es = np.nan self.rheology_law = '' - #giaivins: + #GIA self.lithosphere_shear_modulus = 0. self.lithosphere_density = 0. @@ -38,46 +40,48 @@ self.mantle_density = 0. - #slr + #SLR self.earth_density = 0 - #set default parameters: + #set default parameters self.setdefaultparameters() #}}} - def __repr__(self): # {{{ - string = " Materials:" - string = "%s\n%s" % (string, fielddisplay(self, 'rho_ice', 'ice density [kg / m^3]')) - string = "%s\n%s" % (string, fielddisplay(self, 'rho_water', 'ocean water density [kg / m^3]')) - string = "%s\n%s" % (string, fielddisplay(self, 'rho_freshwater', 'fresh water density [kg / m^3]')) - string = "%s\n%s" % (string, fielddisplay(self, 'mu_water', 'water viscosity [N s / m^2]')) - string = "%s\n%s" % (string, fielddisplay(self, 'heatcapacity', 'heat capacity [J / kg / K]')) - string = "%s\n%s" % (string, fielddisplay(self, 'thermalconductivity', ['ice thermal conductivity [W / m / K]'])) - string = "%s\n%s" % (string, fielddisplay(self, 'temperateiceconductivity', 'temperate ice thermal conductivity [W / m / K]')) - string = "%s\n%s" % (string, fielddisplay(self, "effectiveconductivity_averaging", "computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)")) - string = "%s\n%s" % (string, fielddisplay(self, 'meltingpoint', 'melting point of ice at 1atm in K')) - string = "%s\n%s" % (string, fielddisplay(self, 'latentheat', 'latent heat of fusion [J / kg]')) - string = "%s\n%s" % (string, fielddisplay(self, 'beta', 'rate of change of melting point with pressure [K / Pa]')) - string = "%s\n%s" % (string, fielddisplay(self, 'mixed_layer_capacity', 'mixed layer capacity [W / kg / K]')) - string = "%s\n%s" % (string, fielddisplay(self, 'thermal_exchange_velocity', 'thermal exchange velocity [m / s]')) - string = "%s\n%s" % (string, fielddisplay(self, 'rheology_B', 'flow law parameter [Pa s^(1 / 3)]')) - string = "%s\n%s" % (string, fielddisplay(self, 'rheology_Ec', 'compressive enhancement factor')) - string = "%s\n%s" % (string, fielddisplay(self, 'rheology_Es', 'shear enhancement factor')) - string = "%s\n%s" % (string, fielddisplay(self, 'rheology_law', ['law for the temperature dependance of the rheology: ''None'', ''BuddJacka'', ''Cuffey'', ''CuffeyTemperate'', ''Paterson'', ''Arrhenius'' or ''LliboutryDuval'''])) - string = "%s\n%s" % (string, fielddisplay(self, 'lithosphere_shear_modulus', 'Lithosphere shear modulus [Pa]')) - string = "%s\n%s" % (string, fielddisplay(self, 'lithosphere_density', 'Lithosphere density [g / cm^ - 3]')) - string = "%s\n%s" % (string, fielddisplay(self, 'mantle_shear_modulus', 'Mantle shear modulus [Pa]')) - string = "%s\n%s" % (string, fielddisplay(self, 'mantle_density', 'Mantle density [g / cm^ - 3]')) - string = "%s\n%s" % (string, fielddisplay(self, 'earth_density', 'Mantle density [kg / m^ - 3]')) - return string + def __repr__(self): #{{{ + s = " Materials:" + s = "%s\n%s" % (s, fielddisplay(self, 'rho_ice', 'ice density [kg/m^3]')) + s = "%s\n%s" % (s, fielddisplay(self, 'rho_water', 'ocean water density [kg/m^3]')) + s = "%s\n%s" % (s, fielddisplay(self, 'rho_freshwater', 'fresh water density [kg/m^3]')) + s = "%s\n%s" % (s, fielddisplay(self, 'mu_water', 'water viscosity [N s/m^2]')) + s = "%s\n%s" % (s, fielddisplay(self, 'heatcapacity', 'heat capacity [J/kg/K]')) + s = "%s\n%s" % (s, fielddisplay(self, 'thermalconductivity', ['ice thermal conductivity [W/m/K]'])) + s = "%s\n%s" % (s, fielddisplay(self, 'temperateiceconductivity', 'temperate ice thermal conductivity [W/m/K]')) + s = "%s\n%s" % (s, fielddisplay(self, "effectiveconductivity_averaging", "computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)")) + s = "%s\n%s" % (s, fielddisplay(self, 'meltingpoint', 'melting point of ice at 1atm in K')) + s = "%s\n%s" % (s, fielddisplay(self, 'latentheat', 'latent heat of fusion [J / kg]')) + s = "%s\n%s" % (s, fielddisplay(self, 'beta', 'rate of change of melting point with pressure [K/Pa]')) + s = "%s\n%s" % (s, fielddisplay(self, 'mixed_layer_capacity', 'mixed layer capacity [W/kg/K]')) + s = "%s\n%s" % (s, fielddisplay(self, 'thermal_exchange_velocity', 'thermal exchange velocity [m/s]')) + s = "%s\n%s" % (s, fielddisplay(self, 'rheology_B', 'flow law parameter [Pa s^(1/3)]')) + s = "%s\n%s" % (s, fielddisplay(self, 'rheology_Ec', 'compressive enhancement factor')) + s = "%s\n%s" % (s, fielddisplay(self, 'rheology_Es', 'shear enhancement factor')) + s = "%s\n%s" % (s, fielddisplay(self, 'rheology_law', ['law for the temperature dependance of the rheology: ''None'', ''BuddJacka'', ''Cuffey'', ''CuffeyTemperate'', ''Paterson'', ''Arrhenius'' or ''LliboutryDuval'''])) + s = "%s\n%s" % (s, fielddisplay(self, 'lithosphere_shear_modulus', 'Lithosphere shear modulus [Pa]')) + s = "%s\n%s" % (s, fielddisplay(self, 'lithosphere_density', 'Lithosphere density [g/cm^-3]')) + s = "%s\n%s" % (s, fielddisplay(self, 'mantle_shear_modulus', 'Mantle shear modulus [Pa]')) + s = "%s\n%s" % (s, fielddisplay(self, 'mantle_density', 'Mantle density [g/cm^-3]')) + s = "%s\n%s" % (s, fielddisplay(self, 'earth_density', 'Mantle density [kg/m^-3]')) + + return s #}}} - def extrude(self, md): # {{{ + def extrude(self, md): #{{{ self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node') self.rheology_Ec = project3d(md, 'vector', self.rheology_Ec, 'type', 'node') self.rheology_Es = project3d(md, 'vector', self.rheology_Es, 'type', 'node') + return self #}}} - def setdefaultparameters(self): # {{{ + def setdefaultparameters(self): #{{{ #ice density (kg / m^3) self.rho_ice = 917. @@ -109,10 +113,10 @@ #available: none, paterson and arrhenius self.rheology_law = 'Paterson' - # GIA: + # GIA self.lithosphere_shear_modulus = 6.7e10 # (Pa) self.lithosphere_density = 3.32 # (g / cm^ - 3) self.mantle_shear_modulus = 1.45e11 # (Pa) self.mantle_density = 3.34 # (g / cm^ - 3) - #SLR + #SLR self.earth_density = 5512 # average density of the Earth, (kg / m^3) @@ -120,5 +124,5 @@ #}}} - def checkconsistency(self, md, solution, analyses): # {{{ + def checkconsistency(self, md, solution, analyses): #{{{ md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0) md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0) @@ -143,5 +147,5 @@ # }}} - def marshall(self, prefix, md, fid): # {{{ + def marshall(self, prefix, md, fid): #{{{ WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 2, 'format', 'Integer') WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_ice', 'format', 'Double') Index: /issm/trunk-jpl/src/m/classes/rotational.py =================================================================== --- /issm/trunk-jpl/src/m/classes/rotational.py (revision 25170) +++ /issm/trunk-jpl/src/m/classes/rotational.py (revision 25171) @@ -46,5 +46,5 @@ def checkconsistency(self, md, solution, analyses): # {{{ - if ('SealevelriseAnalysis' not in analyses) or (solution == 'TransientSolution' and md.transient.isrotational == 0): + if ('SealevelriseAnalysis' not in analyses) or (solution == 'TransientSolution' and md.transient.isslr == 0): return md Index: /issm/trunk-jpl/src/m/classes/slr.py =================================================================== --- /issm/trunk-jpl/src/m/classes/slr.py (revision 25170) +++ /issm/trunk-jpl/src/m/classes/slr.py (revision 25171) @@ -122,5 +122,5 @@ def checkconsistency(self, md, solution, analyses): # {{{ #Early return - if (solution != 'SealevelriseAnalysis') or (solution == 'TransientSolution' and not md.transient.isslr): + if (solution != 'SealevelriseAnalysis') or (solution == 'TransientSolution' and md.transient.isslr == 0): return md @@ -147,6 +147,6 @@ #check that love numbers are provided at the same level of accuracy: - if (size(self.love_h, 0) != size(self.love_k, 0) | size(self.love_h, 0) != size(self.love_l, 0)): - error('slr error message: love numbers should be provided at the same level of accuracy') + if (size(self.love_h, 0) != size(self.love_k, 0)) or (size(self.love_h, 0) != size(self.love_l, 0)): + raise Exception('slr error message: love numbers should be provided at the same level of accuracy') #cross check that whereever we have an ice load, the mask is < 0 on each vertex: @@ -160,5 +160,5 @@ #a coupler to a planet model is provided. if self.geodetic and not md.transient.iscoupler and domaintype(md.mesh) != 'mesh3dsurface': - error('model is requesting geodetic computations without being a mesh3dsurface, or being coupled to one!') + raise Exception('model is requesting geodetic computations without being a mesh3dsurface, or being coupled to one!') return md # }}} Index: /issm/trunk-jpl/src/m/classes/surfaceload.py =================================================================== --- /issm/trunk-jpl/src/m/classes/surfaceload.py (revision 25170) +++ /issm/trunk-jpl/src/m/classes/surfaceload.py (revision 25171) @@ -41,14 +41,14 @@ def checkconsistency(self, md, solution, analyses): # {{{ - if ('SealevelriseAnalysis' not in analyses) or (solution == 'TransientSolution' and md.transient.issurfaceload == 0): + if ('SealevelriseAnalysis' not in analyses) or (solution == 'TransientSolution' and md.transient.isslr == 0): return md - if len(self.icethicknesschange) > 0: + if len(self.icethicknesschange): md = checkfield(md,'fieldname', 'solidearth.surfaceload.icethicknesschange', 'timeseries', 1, 'NaN', 1, 'Inf', 1) - if len(self.waterheightchange) > 0: + if len(self.waterheightchange): md = checkfield(md,'fieldname', 'solidearth.surfaceload.waterheightchange', 'timeseries', 1, 'NaN', 1, 'Inf', 1) - if len(self.other) > 0: + if len(self.other): md = checkfield(md,'fieldname', 'solidearth.surfaceload.other', 'timeseries', 1, 'NaN', 1, 'Inf', 1) Index: /issm/trunk-jpl/test/NightlyRun/test2002.py =================================================================== --- /issm/trunk-jpl/test/NightlyRun/test2002.py (revision 25170) +++ /issm/trunk-jpl/test/NightlyRun/test2002.py (revision 25171) @@ -20,5 +20,5 @@ md.solidearth.surfaceload.icethicknesschange = np.zeros((md.mesh.numberofelements, 1)) md.solidearth.sealevel = np.zeros((md.mesh.numberofvertices, 1)) -md.dsl.global_average_thermosteric_sea_level_change = np.zeros((1, 1)) +md.dsl.global_average_thermosteric_sea_level_change = np.zeros((2, 1)) md.dsl.sea_surface_height_change_above_geoid = np.zeros((md.mesh.numberofvertices + 1, 1)) md.dsl.sea_water_pressure_change_at_sea_floor = np.zeros((md.mesh.numberofvertices + 1, 1))