source: issm/oecreview/Archive/24684-25833/ISSM-24755-24756.diff@ 25834

Last change on this file since 25834 was 25834, checked in by Mathieu Morlighem, 4 years ago

CHG: added 24684-25833
File size: 29.1 KB
RevLine 
[25834]1Index: ../trunk-jpl/src/m/classes/materials.py
2===================================================================
3--- ../trunk-jpl/src/m/classes/materials.py (revision 24755)
4+++ ../trunk-jpl/src/m/classes/materials.py (revision 24756)
5@@ -4,30 +4,6 @@
6 from checkfield import checkfield
7 from WriteData import WriteData
8
9-
10-def naturetointeger(strnat): #{{{
11-
12- intnat = np.zeros(len(strnat))
13- for i in range(len(intnat)):
14- if strnat[i] == 'damageice':
15- intnat[i] = 1
16- elif strnat[i] == 'estar':
17- intnat[i] = 2
18- elif strnat[i] == 'ice':
19- intnat[i] = 3
20- elif strnat[i] == 'enhancedice':
21- intnat[i] = 4
22- elif strnat[i] == 'litho':
23- intnat[i] = 6
24- elif strnat[i] == 'hydro':
25- intnat[i] = 7
26- else:
27- raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')")
28-
29- return intnat
30- # }}}
31-
32-
33 class materials(object):
34 """
35 MATERIALS class definition
36@@ -44,8 +20,8 @@
37 self.nature = args
38
39 for i in range(len(self.nature)):
40- if not(self.nature[i] == 'litho' or self.nature[i] == 'ice'):
41- raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
42+ if not(self.nature[i] == 'litho' or self.nature[i] == 'ice' or self.nature[i] == 'hydro'):
43+ raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
44
45 #start filling in the dynamic fields:
46 for i in range(len(self.nature)):
47@@ -52,7 +28,6 @@
48 nat = self.nature[i]
49 if nat == 'ice':
50 setattr(self, 'rho_ice', 0)
51- setattr(self, 'rho_ice', 0)
52 setattr(self, 'rho_water', 0)
53 setattr(self, 'rho_freshwater', 0)
54 setattr(self, 'mu_water', 0)
55@@ -78,61 +53,17 @@
56 setattr(self, 'isburgers', 0)
57 setattr(self, 'density', 0)
58 setattr(self, 'issolid', 0)
59+ elif nat == 'hydro':
60+ setattr(self, 'rho_ice', 0)
61+ setattr(self, 'rho_water', 0)
62+ setattr(self, 'earth_density', 0)
63 else:
64- raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
65- #set default parameters:
66+ raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
67+
68+ #set default parameters:
69 self.setdefaultparameters()
70 #}}}
71
72- def __repr__(self): # {{{
73- string = " Materials:"
74- for i in range(len(self.nature)):
75- nat = self.nature[i]
76- if nat == 'ice':
77- string = "%s\n%s" % (string, 'Ice:')
78- string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]"))
79- string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "water density [kg / m^3]"))
80- string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]"))
81- string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]"))
82- string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]"))
83- string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]"))
84- string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]"))
85- string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K"))
86- string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]"))
87- string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]"))
88- string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]"))
89- string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]"))
90- string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]"))
91- string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent"))
92- string = "%s\n%s" % (string, fielddisplay(self, "rheology_law", "law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', or 'NyeH2O'"))
93- elif nat == 'litho':
94- string = "%s\n%s" % (string, 'Litho:')
95- string = "%s\n%s" % (string, fielddisplay(self, 'numlayers', 'number of layers (default 2)'))
96- string = "%s\n%s" % (string, fielddisplay(self, 'radius', 'array describing the radius for each interface (numlayers + 1) [m]'))
97- string = "%s\n%s" % (string, fielddisplay(self, 'viscosity', 'array describing each layer''s viscosity (numlayers) [Pa.s]'))
98- string = "%s\n%s" % (string, fielddisplay(self, 'lame_lambda', 'array describing the lame lambda parameter (numlayers) [Pa]'))
99- string = "%s\n%s" % (string, fielddisplay(self, 'lame_mu', 'array describing the shear modulus for each layers (numlayers) [Pa]'))
100- string = "%s\n%s" % (string, fielddisplay(self, 'burgers_viscosity', 'array describing each layer''s transient viscosity, only for Burgers rheologies (numlayers) [Pa.s]'))
101- string = "%s\n%s" % (string, fielddisplay(self, 'burgers_mu', 'array describing each layer''s transient shear modulus, only for Burgers rheologies (numlayers) [Pa]'))
102- string = "%s\n%s" % (string, fielddisplay(self, 'isburgers', 'array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)'))
103- string = "%s\n%s" % (string, fielddisplay(self, 'density', 'array describing each layer''s density (numlayers) [kg / m^3]'))
104- string = "%s\n%s" % (string, fielddisplay(self, 'issolid', 'array describing whether the layer is solid or liquid (default 1) (numlayers)'))
105-
106- else:
107- raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho')")
108-
109- return string
110- #}}}
111-
112- def extrude(self, md): # {{{
113- for i in range(len(self.nature)):
114- nat = self.nature[i]
115- if nat == 'ice':
116- self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node')
117- self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element')
118- return self
119- #}}}
120-
121 def setdefaultparameters(self): # {{{
122 for i in range(len(self.nature)):
123 nat = self.nature[i]
124@@ -164,7 +95,6 @@
125 #Rheology law: what is the temperature dependence of B with T
126 #available: none, paterson and arrhenius
127 self.rheology_law = 'Paterson'
128-
129 elif nat == 'litho':
130 #we default to a configuration that enables running GIA solutions using giacaron and / or giaivins.
131 self.numlayers = 2
132@@ -179,13 +109,63 @@
133 self.isburgers = [0, 0]
134 self.density = [5.51 * 1e3, 5.50 * 1e3] # (Pa) #mantle and lithosphere density [kg / m^3]
135 self.issolid = [1, 1] # is layer solid or liquid.
136-
137+ elif nat == 'hydro':
138+ #ice density (kg / m^3)
139+ self.rho_ice = 917.
140+ #ocean water density (kg / m^3)
141+ self.rho_water = 1023.
142+ #average density of the Earth (kg / m^3)
143+ self.earth_density = 5512
144 else:
145- raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho')")
146+ raise RuntimeError("materials setdefaultparameters error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
147
148 return self
149 #}}}
150
151+ def __repr__(self): # {{{
152+ string = " Materials:"
153+ for i in range(len(self.nature)):
154+ nat = self.nature[i]
155+ if nat == 'ice':
156+ string = "%s\n%s" % (string, 'Ice:')
157+ string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]"))
158+ string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "ocean water density [kg / m^3]"))
159+ string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]"))
160+ string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]"))
161+ string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]"))
162+ string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]"))
163+ string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]"))
164+ string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K"))
165+ string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]"))
166+ string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]"))
167+ string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]"))
168+ string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]"))
169+ string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]"))
170+ string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent"))
171+ string = "%s\n%s" % (string, fielddisplay(self, "rheology_law", "law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', or 'NyeH2O'"))
172+ elif nat == 'litho':
173+ string = "%s\n%s" % (string, 'Litho:')
174+ string = "%s\n%s" % (string, fielddisplay(self, 'numlayers', 'number of layers (default 2)'))
175+ string = "%s\n%s" % (string, fielddisplay(self, 'radius', 'array describing the radius for each interface (numlayers + 1) [m]'))
176+ string = "%s\n%s" % (string, fielddisplay(self, 'viscosity', 'array describing each layer''s viscosity (numlayers) [Pa.s]'))
177+ string = "%s\n%s" % (string, fielddisplay(self, 'lame_lambda', 'array describing the lame lambda parameter (numlayers) [Pa]'))
178+ string = "%s\n%s" % (string, fielddisplay(self, 'lame_mu', 'array describing the shear modulus for each layers (numlayers) [Pa]'))
179+ string = "%s\n%s" % (string, fielddisplay(self, 'burgers_viscosity', 'array describing each layer''s transient viscosity, only for Burgers rheologies (numlayers) [Pa.s]'))
180+ string = "%s\n%s" % (string, fielddisplay(self, 'burgers_mu', 'array describing each layer''s transient shear modulus, only for Burgers rheologies (numlayers) [Pa]'))
181+ string = "%s\n%s" % (string, fielddisplay(self, 'isburgers', 'array describing whether we adopt a MaxWell (0) or Burgers (1) rheology (default 0)'))
182+ string = "%s\n%s" % (string, fielddisplay(self, 'density', 'array describing each layer''s density (numlayers) [kg / m^3]'))
183+ string = "%s\n%s" % (string, fielddisplay(self, 'issolid', 'array describing whether the layer is solid or liquid (default 1) (numlayers)'))
184+ elif nat == 'hydro':
185+ string = "%s\n%s" % (string, 'Hydro:')
186+ string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]"))
187+ string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "ocean water density [kg / m^3]"))
188+ string = "%s\n%s" % (string, fielddisplay(self, "earth_density", "mantle density [kg / m^3]"))
189+ else:
190+ raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
191+
192+ return string
193+ #}}}
194+
195 def checkconsistency(self, md, solution, analyses): # {{{
196 for i in range(len(self.nature)):
197 nat = self.nature[i]
198@@ -197,7 +177,6 @@
199 md = checkfield(md, 'fieldname', 'materials.rheology_B', '>', 0, 'timeseries', 1, 'NaN', 1, 'Inf', 1)
200 md = checkfield(md, 'fieldname', 'materials.rheology_n', '>', 0, 'size', [md.mesh.numberofelements])
201 md = checkfield(md, 'fieldname', 'materials.rheology_law', 'values', ['None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', 'NyeH2O'])
202-
203 elif nat == 'litho':
204 if 'LoveAnalysis' not in analyses:
205 return md
206@@ -223,9 +202,12 @@
207 for i in range(md.materials.numlayers - 1):
208 if (not md.materials.issolid[i]) and (not md.materials.issolid[i + 1]): #if there are at least two consecutive indices that contain issolid = 0
209 raise RuntimeError("%s%i%s" % ('2 or more adjacent fluid layers detected starting at layer ', i, '. This is not supported yet. Consider merging them.'))
210-
211+ elif nat == 'hydro':
212+ md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0)
213+ md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0)
214+ md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', [1])
215 else:
216- raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho')")
217+ raise RuntimeError("materials checkconsistency error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
218
219 return md
220 # }}}
221@@ -232,8 +214,8 @@
222
223 def marshall(self, prefix, md, fid): # {{{
224 #1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum
225- WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 6, 'format', 'Integer')
226 WriteData(fid, prefix, 'name', 'md.materials.nature', 'data', naturetointeger(self.nature), 'format', 'IntMat', 'mattype', 3)
227+ WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 5, 'format', 'Integer') #DANGER, this can evolve if you have classes
228
229 for i in range(len(self.nature)):
230 nat = self.nature[i]
231@@ -254,7 +236,6 @@
232 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_B', 'format', 'DoubleMat', 'mattype', 1, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts)
233 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_n', 'format', 'DoubleMat', 'mattype', 2)
234 WriteData(fid, prefix, 'data', self.rheology_law, 'name', 'md.materials.rheology_law', 'format', 'String')
235-
236 elif nat == 'litho':
237 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'numlayers', 'format', 'Integer')
238 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'radius', 'format', 'DoubleMat', 'mattype', 3)
239@@ -266,7 +247,43 @@
240 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'isburgers', 'format', 'DoubleMat', 'mattype', 3)
241 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'burgers_viscosity', 'format', 'DoubleMat', 'mattype', 3)
242 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'burgers_mu', 'format', 'DoubleMat', 'mattype', 3)
243-
244+ elif nat == 'hydro':
245+ WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_ice', 'format', 'Double')
246+ WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_water', 'format', 'Double')
247+ WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'earth_density', 'format', 'Double')
248 else:
249- raise RuntimeError("materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'')")
250+ raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
251 # }}}
252+
253+ def extrude(self, md): # {{{
254+ for i in range(len(self.nature)):
255+ nat = self.nature[i]
256+ if nat == 'ice':
257+ self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node')
258+ self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element')
259+ return self
260+ #}}}
261+
262+def naturetointeger(strnat): #{{{
263+ intnat = np.zeros(len(strnat))
264+
265+ for i in range(len(intnat)):
266+ if strnat[i] == 'damageice':
267+ intnat[i] = 1
268+ elif strnat[i] == 'estar':
269+ intnat[i] = 2
270+ elif strnat[i] == 'ice':
271+ intnat[i] = 3
272+ elif strnat[i] == 'enhancedice':
273+ intnat[i] = 4
274+ # elif strnat[i] == 'materials':
275+ # intnat[i] = 5 #this case will never happen, kept to ensure equivalent of codes between IoCodeToMeterialsEnum and IoCodeToNatureEnum
276+ elif strnat[i] == 'litho':
277+ intnat[i] = 6
278+ elif strnat[i] == 'hydro':
279+ intnat[i] = 7
280+ else:
281+ raise RuntimeError("materials constructor error message: nature of the material not supported yet! ('ice' or 'litho' or 'hydro')")
282+
283+ return intnat
284+# }}}
285Index: ../trunk-jpl/src/m/classes/materials.m
286===================================================================
287--- ../trunk-jpl/src/m/classes/materials.m (revision 24755)
288+++ ../trunk-jpl/src/m/classes/materials.m (revision 24756)
289@@ -4,7 +4,7 @@
290 % materials=materials();
291
292 classdef materials < dynamicprops
293- properties (SetAccess=public)
294+ properties (SetAccess=public)
295 nature={};
296 %all properties are dynamic.
297 end
298@@ -12,20 +12,20 @@
299 function self = materials(varargin) % {{{
300 if nargin==0
301 self.nature={'ice'};
302- else
303+ else
304 self.nature=varargin;
305 end
306-
307- %check this is acceptable:
308+
309+ %check this is acceptable:
310 for i=1:length(self.nature),
311- if ~(strcmpi(self.nature{i},'litho') | strcmpi(self.nature{i},'ice') | strcmpi(self.nature{i},'hydro')),
312+ if ~(strcmpi(self.nature{i},'litho') | strcmpi(self.nature{i},'ice') | strcmpi(self.nature{i},'hydro')),
313 error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
314 end
315 end
316-
317- %start filling in the dynamic fields:
318+
319+ %start filling in the dynamic fields:
320 for i=1:length(self.nature),
321- nat=self.nature{i};
322+ nat=self.nature{i};
323 switch nat
324 case 'ice'
325 self.addprop('rho_ice');
326@@ -36,7 +36,7 @@
327 self.addprop('latentheat');
328 self.addprop('thermalconductivity');
329 self.addprop('temperateiceconductivity');
330- self.addprop('meltingpoint');
331+ self.addprop('meltingpoint');
332 self.addprop('beta');
333 self.addprop('mixed_layer_capacity');
334 self.addprop('thermal_exchange_velocity');
335@@ -45,7 +45,7 @@
336 self.addprop('rheology_law');
337 case 'litho'
338 self.addprop('numlayers');
339- self.addprop('radius');
340+ self.addprop('radius');
341 self.addprop('viscosity');
342 self.addprop('lame_lambda');
343 self.addprop('lame_mu');
344@@ -69,7 +69,7 @@
345 function self = setdefaultparameters(self) % {{{
346
347 for i=1:length(self.nature),
348- nat=self.nature{i};
349+ nat=self.nature{i};
350 switch nat
351 case 'ice'
352 %ice density (kg/m^3)
353@@ -82,7 +82,7 @@
354 self.rho_freshwater=1000.;
355
356 %water viscosity (N.s/m^2)
357- self.mu_water=0.001787;
358+ self.mu_water=0.001787;
359
360 %ice heat capacity cp (J/kg/K)
361 self.heatcapacity=2093.;
362@@ -92,7 +92,7 @@
363
364 %ice thermal conductivity (W/m/K)
365 self.thermalconductivity=2.4;
366-
367+
368 %wet ice thermal conductivity (W/m/K)
369 self.temperateiceconductivity=.24;
370
371@@ -113,11 +113,11 @@
372 self.rheology_law='Paterson';
373
374 case 'litho'
375- %we default to a configuration that enables running GIA solutions using giacaron and/or giaivins.
376+ %we default to a configuration that enables running GIA solutions using giacaron and/or giaivins.
377 self.numlayers=2;
378
379 %center of the earth (approximation, must not be 0), then the lab (lithosphere/asthenosphere boundary) then the surface
380- %(with 1d3 to avoid numerical singularities)
381+ %(with 1d3 to avoid numerical singularities)
382 self.radius=[1e3;6278*1e3;6378*1e3];
383
384 self.viscosity=[1e21;1e40]; %mantle and lithosphere viscosity (respectively) [Pa.s]
385@@ -135,8 +135,8 @@
386
387 %ocean water density (kg/m^3)
388 self.rho_water=1023.;
389- %SLR
390- self.earth_density= 5512; % average density of the Earth, (kg/m^3)
391+ % average density of the Earth (kg/m^3)
392+ self.earth_density=5512;
393
394 otherwise
395 error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
396@@ -147,7 +147,7 @@
397 disp(sprintf(' Materials:'));
398
399 for i=1:length(self.nature),
400- nat=self.nature{i};
401+ nat=self.nature{i};
402 switch nat
403 case 'ice'
404 disp(sprintf(' \nIce:'));
405@@ -182,7 +182,7 @@
406 disp(sprintf(' \nHydro:'));
407 fielddisplay(self,'rho_ice','ice density [kg/m^3]');
408 fielddisplay(self,'rho_water','ocean water density [kg/m^3]');
409- fielddisplay(self,'earth_density','Mantle density [kg/m^-3]');
410+ fielddisplay(self,'earth_density','mantle density [kg/m^3]');
411
412 otherwise
413 error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
414@@ -192,7 +192,7 @@
415 function md = checkconsistency(self,md,solution,analyses) % {{{
416
417 for i=1:length(self.nature),
418- nat=self.nature{i};
419+ nat=self.nature{i};
420 switch nat
421 case 'ice'
422 md = checkfield(md,'fieldname','materials.rho_ice','>',0);
423@@ -216,11 +216,11 @@
424 md = checkfield(md,'fieldname','materials.burgers_mu','Inf',1,'size',[md.materials.numlayers 1],'>=',0);
425
426 for i=1:md.materials.numlayers,
427- if md.materials.isburgers(i) & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))),
428+ if md.materials.isburgers(i) & (isnan(md.materials.burgers_viscosity(i) | isnan(md.materials.burgers_mu(i)))),
429 error('materials checkconsistency error message: Litho burgers_viscosity or burgers_mu has NaN values, inconsistent with isburgers choice');
430 end
431 end
432- if md.materials.issolid(1)==0 | md.materials.lame_mu(1)==0
433+ if md.materials.issolid(1)==0 | md.materials.lame_mu(1)==0
434 error('First layer must be solid (issolid(1) > 0 AND lame_mu(1) > 0). Add a weak inner core if necessary.');
435 end
436 ind=find(md.materials.issolid==0);
437@@ -241,12 +241,12 @@
438 end % }}}
439 function marshall(self,prefix,md,fid) % {{{
440
441- %1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum
442+ %1: MatdamageiceEnum 2: MatestarEnum 3: MaticeEnum 4: MatenhancediceEnum 5: MaterialsEnum
443 WriteData(fid,prefix,'name','md.materials.nature','data',naturetointeger(self.nature),'format','IntMat','mattype',3);
444- WriteData(fid,prefix,'name','md.materials.type','data',5,'format','Integer'); %DANGER, this can evolve if you have classes.
445+ WriteData(fid,prefix,'name','md.materials.type','data',5,'format','Integer'); %DANGER, this can evolve if you have classes.
446
447 for i=1:length(self.nature),
448- nat=self.nature{i};
449+ nat=self.nature{i};
450 switch nat
451 case 'ice'
452 WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_ice','format','Double');
453@@ -270,11 +270,11 @@
454 WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_mu','format','DoubleMat','mattype',3);
455 WriteData(fid,prefix,'object',self,'class','materials','fieldname','lame_lambda','format','DoubleMat','mattype',3);
456 WriteData(fid,prefix,'object',self,'class','materials','fieldname','issolid','format','DoubleMat','mattype',3);
457- WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3);
458- WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3);
459- WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3);
460- WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3);
461- WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3);
462+ WriteData(fid,prefix,'object',self,'class','materials','fieldname','density','format','DoubleMat','mattype',3);
463+ WriteData(fid,prefix,'object',self,'class','materials','fieldname','viscosity','format','DoubleMat','mattype',3);
464+ WriteData(fid,prefix,'object',self,'class','materials','fieldname','isburgers','format','DoubleMat','mattype',3);
465+ WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_viscosity','format','DoubleMat','mattype',3);
466+ WriteData(fid,prefix,'object',self,'class','materials','fieldname','burgers_mu','format','DoubleMat','mattype',3);
467 case 'hydro'
468 WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_ice','format','Double');
469 WriteData(fid,prefix,'object',self,'class','materials','fieldname','rho_water','format','Double');
470@@ -287,7 +287,7 @@
471 end % }}}
472 function self = extrude(self,md) % {{{
473 for i=1:length(self.nature),
474- nat=self.nature{i};
475+ nat=self.nature{i};
476 switch nat
477 case 'ice'
478 self.rheology_B=project3d(md,'vector',self.rheology_B,'type','node');
479@@ -296,9 +296,9 @@
480 end
481 end % }}}
482 function savemodeljs(self,fid,modelname) % {{{
483-
484+
485 for i=1:length(self.nature),
486- nat=self.nature{i};
487+ nat=self.nature{i};
488 switch nat
489 case 'ice'
490 writejsdouble(fid,[modelname '.materials.rho_ice'],self.rho_ice);
491@@ -323,11 +323,11 @@
492 writejsdouble(fid,[modelname '.materials.lame_mu'],self.lame_mu);
493 writejsdouble(fid,[modelname '.materials.lame_lambda'],self.lame_lambda);
494 writejsdouble(fid,[modelname '.materials.issolid'],self.issolid);
495- writejsdouble(fid,[modelname '.materials.density'],self.density);
496- writejsdouble(fid,[modelname '.materials.viscosity'],self.viscosity);
497- writejsdouble(fid,[modelname '.materials.isburgers'],self.isburgers);
498- writejsdouble(fid,[modelname '.materials.burgers_viscosity'],self.burgers_viscosity);
499- writejsdouble(fid,[modelname '.materials.burgers_mu'],self.burgers_mu);
500+ writejsdouble(fid,[modelname '.materials.density'],self.density);
501+ writejsdouble(fid,[modelname '.materials.viscosity'],self.viscosity);
502+ writejsdouble(fid,[modelname '.materials.isburgers'],self.isburgers);
503+ writejsdouble(fid,[modelname '.materials.burgers_viscosity'],self.burgers_viscosity);
504+ writejsdouble(fid,[modelname '.materials.burgers_mu'],self.burgers_mu);
505
506 case 'hydro'
507 writejsdouble(fid,[modelname '.materials.rho_ice'],self.rho_ice);
508@@ -344,26 +344,27 @@
509 end % }}}
510 end
511 end
512- function intnat = naturetointeger(strnat) % {{{
513- intnat=zeros(length(strnat),1);
514- for i=1:length(strnat),
515- switch strnat{i},
516- case 'damageice'
517- intnat(i)=1;
518- case 'estar'
519- intnat(i)=2;
520- case 'ice'
521- intnat(i)=3;
522- case 'enhancedice'
523- intnat(i)=4;
524- %case 'materials' %this case will never happen, kept to ensure equivalent of codes between IoCodeToMaterialsEnumm and IoCodeToNatureEnum
525- % intnat(i)=5;
526- case 'litho'
527- intnat(i)=6;
528- case 'hydro'
529- intnat(i)=7;
530- otherwise
531- error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
532- end
533- end
534- end % }}}
535+
536+function intnat = naturetointeger(strnat) % {{{
537+ intnat=zeros(length(strnat),1);
538+ for i=1:length(strnat),
539+ switch strnat{i},
540+ case 'damageice'
541+ intnat(i)=1;
542+ case 'estar'
543+ intnat(i)=2;
544+ case 'ice'
545+ intnat(i)=3;
546+ case 'enhancedice'
547+ intnat(i)=4;
548+ %case 'materials' %this case will never happen, kept to ensure equivalent of codes between IoCodeToMaterialsEnum and IoCodeToNatureEnum
549+ % intnat(i)=5;
550+ case 'litho'
551+ intnat(i)=6;
552+ case 'hydro'
553+ intnat(i)=7;
554+ otherwise
555+ error('materials constructor error message: nature of the material not supported yet! (''ice'' or ''litho'' or ''hydro'')');
556+ end
557+ end
558+end % }}}
Note: See TracBrowser for help on using the repository browser.