source: issm/trunk-jpl/src/m/classes/matenhancedice.py

Last change on this file was 26358, checked in by jdquinn, 4 years ago

CHG: Completed MATLAB -> Python updates for SE; archive updates now that GMSH can be used on macOS and Linux; various minor bug fixes; formatting; cleanup
File size: 8.3 KB
Line 
1import numpy as np
2
3from checkfield import checkfield
4from fielddisplay import fielddisplay
5from project3d import project3d
6from WriteData import WriteData
7
8
9class matenhancedice(object):
10 """MATICE class definition
11
12 Usage:
13 matenhancedice = matenhancedice()
14 """
15
16 def __init__(self): # {{{
17 self.rho_ice = 0
18 self.rho_water = 0
19 self.rho_freshwater = 0
20 self.mu_water = 0
21 self.heatcapacity = 0
22 self.latentheat = 0
23 self.thermalconductivity = 0
24 self.temperateiceconductivity = 0
25 self.effectiveconductivity_averaging = 0
26 self.meltingpoint = 0
27 self.beta = 0
28 self.mixed_layer_capacity = 0
29 self.thermal_exchange_velocity = 0
30 self.rheology_E = np.nan
31 self.rheology_B = np.nan
32 self.rheology_n = np.nan
33 self.rheology_law = ''
34
35 #SLC
36 self.earth_density = 0
37
38 self.setdefaultparameters()
39 # }}}
40
41 def __repr__(self): # {{{
42 s = ' Materials:\n'
43 s += '{}\n'.format(fielddisplay(self, 'rho_ice', 'ice density [kg/m^3]'))
44 s += '{}\n'.format(fielddisplay(self, 'rho_water', 'water density [kg/m^3]'))
45 s += '{}\n'.format(fielddisplay(self, 'rho_freshwater', 'fresh water density [kg/m^3]'))
46 s += '{}\n'.format(fielddisplay(self, 'mu_water', 'water viscosity [N s/m^2]'))
47 s += '{}\n'.format(fielddisplay(self, 'heatcapacity', 'heat capacity [J/kg/K]'))
48 s += '{}\n'.format(fielddisplay(self, 'thermalconductivity', 'ice thermal conductivity [W/m/K]'))
49 s += '{}\n'.format(fielddisplay(self, 'temperateiceconductivity', 'temperate ice thermal conductivity [W/m/K]'))
50 s += '{}\n'.format(fielddisplay(self, 'effectiveconductivity_averaging', 'computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)'))
51 s += '{}\n'.format(fielddisplay(self, 'meltingpoint', 'melting point of ice at 1atm in K'))
52 s += '{}\n'.format(fielddisplay(self, 'latentheat', 'latent heat of fusion [J/m^3]'))
53 s += '{}\n'.format(fielddisplay(self, 'beta', 'rate of change of melting point with pressure [K/Pa]'))
54 s += '{}\n'.format(fielddisplay(self, 'mixed_layer_capacity', 'mixed layer capacity [W/kg/K]'))
55 s += '{}\n'.format(fielddisplay(self, 'thermal_exchange_velocity', 'thermal exchange velocity [m/s]'))
56 s += '{}\n'.format(fielddisplay(self, 'rheology_E', 'enhancement factor'))
57 s += '{}\n'.format(fielddisplay(self, 'rheology_B', 'flow law parameter [Pa s^(1/n)]'))
58 s += '{}\n'.format(fielddisplay(self, 'rheology_n', 'Glen\'s flow law exponent'))
59 s += '{}\n'.format(fielddisplay(self, 'rheology_law', 'law for the temperature dependance of the rheology: \'None\', \'BuddJacka\', \'Cuffey\', \'CuffeyTemperate\', \'Paterson\', \'Arrhenius\' or \'LliboutryDuval\''))
60 s += '{}\n'.format(fielddisplay(self, 'earth_density', 'Mantle density [kg/m^-3]'))
61
62 return s
63 # }}}
64
65 def extrude(self, md): # {{{
66 self.rheology_E = project3d(md, 'vector', self.rheology_E, 'type', 'node')
67 self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node')
68 self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element')
69 return self
70 # }}}
71
72 def setdefaultparameters(self): # {{{
73 # Ice density (kg/m^3)
74 self.rho_ice = 917.0
75 # Ocean water density (kg/m^3)
76 self.rho_water = 1023.0
77 # Fresh water density (kg/m^3)
78 self.rho_freshwater = 1000.0
79 # Water viscosity (N.s/m^2)
80 self.mu_water = 0.001787
81 # Ice heat capacity cp (J/kg/K)
82 self.heatcapacity = 2093.0
83 # Ice latent heat of fusion L (J/kg)
84 self.latentheat = 3.34 * pow(10, 5)
85 # Ice thermal conductivity (W/m/K)
86 self.thermalconductivity = 2.4
87 # Temperate ice thermal conductivity (W/m/K)
88 self.temperateiceconductivity = 0.24
89 # Computation of effective conductivity
90 self.effectiveconductivity_averaging = 1
91 # The melting point of ice at 1 atmosphere of pressure in K
92 self.meltingpoint = 273.15
93 # Rate of change of melting point with pressure (K/Pa)
94 self.beta = 9.8 * pow(10, -8)
95 # Mixed layer (ice-water interface) heat capacity (J/kg/K)
96 self.mixed_layer_capacity = 3974.0
97 # Thermal exchange velocity (ice-water interface) (m/s)
98 self.thermal_exchange_velocity = 1.00 * pow(10, -4)
99 # Rheology law: what is the temperature dependence of B with T
100 # available: none, paterson and arrhenius
101 self.rheology_law = 'Paterson'
102
103 # SLC
104 self.earth_density = 5512 # average density of the Earth, (kg/m^3)
105
106 return self
107 # }}}
108
109 def checkconsistency(self, md, solution, analyses): # {{{
110 md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0)
111 md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0)
112 md = checkfield(md, 'fieldname', 'materials.rho_freshwater', '>', 0)
113 md = checkfield(md, 'fieldname', 'materials.mu_water', '>', 0)
114 md = checkfield(md, 'fieldname', 'materials.rheology_E', '>', 0, 'timeseries', 1, 'NaN', 1, 'Inf', 1)
115 md = checkfield(md, 'fieldname', 'materials.rheology_B', '>', 0, 'timeseries', 1, 'NaN', 1, 'Inf', 1)
116 md = checkfield(md, 'fieldname', 'materials.rheology_n', '>', 0, 'size', [md.mesh.numberofelements])
117 md = checkfield(md, 'fieldname', 'materials.rheology_law', 'values', ['None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval'])
118 md = checkfield(md, 'fieldname', 'materials.effectiveconductivity_averaging', 'numel', [1], 'values', [0, 1, 2])
119
120 if 'SealevelchangeAnalysis' in analyses:
121 md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', 1)
122 return md
123 # }}}
124
125 def marshall(self, prefix, md, fid): # {{{
126 WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 4, 'format', 'Integer')
127 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_ice', 'format', 'Double')
128 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_water', 'format', 'Double')
129 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_freshwater', 'format', 'Double')
130 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'mu_water', 'format', 'Double')
131 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'heatcapacity', 'format', 'Double')
132 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'latentheat', 'format', 'Double')
133 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'thermalconductivity', 'format', 'Double')
134 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'temperateiceconductivity', 'format', 'Double')
135 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'effectiveconductivity_averaging', 'format', 'Integer')
136 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'meltingpoint', 'format', 'Double')
137 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'beta', 'format', 'Double')
138 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'mixed_layer_capacity', 'format', 'Double')
139 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'thermal_exchange_velocity', 'format', 'Double')
140 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_E', 'format', 'DoubleMat', 'mattype', 1, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts)
141 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_B', 'format', 'DoubleMat', 'mattype', 1, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts)
142 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_n', 'format', 'DoubleMat', 'mattype', 2)
143 WriteData(fid, prefix, 'data', self.rheology_law, 'name', 'md.materials.rheology_law', 'format', 'String')
144 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'earth_density', 'format', 'Double')
145 # }}}
Note: See TracBrowser for help on using the repository browser.