source: issm/oecreview/Archive/24684-25833/ISSM-25168-25169.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: 67.7 KB
RevLine 
[25834]1Index: ../trunk-jpl/src/m/classes/matdamageice.py
2===================================================================
3--- ../trunk-jpl/src/m/classes/matdamageice.py (revision 25168)
4+++ ../trunk-jpl/src/m/classes/matdamageice.py (revision 25169)
5@@ -129,6 +129,16 @@
6 md = checkfield(md, 'fieldname', 'materials.mantle_shear_modulus', '>', 0, 'numel', [1])
7 md = checkfield(md, 'fieldname', 'materials.mantle_density', '>', 0, 'numel', [1])
8 md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', [1])
9+
10+ if 'GiaAnalysis' in analyses:
11+ md = checkfield(md, 'fieldname', 'materials.lithosphere_shear_modulus', '>', 0, 'numel', 1)
12+ md = checkfield(md, 'fieldname', 'materials.lithosphere_density', '>', 0, 'numel', 1)
13+ md = checkfield(md,'fieldname', 'materials.mantle_shear_modulus', '>', 0, 'numel', 1)
14+ md = checkfield(md,'fieldname', 'materials.mantle_density', '>', 0, 'numel', 1)
15+
16+ if 'SealevelriseAnalysis' in analyses:
17+ md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', 1)
18+
19 return md
20 # }}}
21
22Index: ../trunk-jpl/src/m/classes/toolkits.m
23===================================================================
24--- ../trunk-jpl/src/m/classes/toolkits.m (revision 25168)
25+++ ../trunk-jpl/src/m/classes/toolkits.m (revision 25169)
26@@ -4,170 +4,176 @@
27 % self=toolkits();
28
29 classdef toolkits < dynamicprops
30- properties (SetAccess=public)
31- DefaultAnalysis = struct();
32- RecoveryAnalysis = struct();
33- %The other properties are dynamic
34- end
35- methods (Static)
36- function self = loadobj(self) % {{{
37- % This function is directly called by matlab when a model object is
38- % loaded. Update old properties here
39+ properties (SetAccess=public)
40+ DefaultAnalysis = struct();
41+ RecoveryAnalysis = struct();
42+ %The other properties are dynamic
43+ end
44+ methods (Static)
45+ function self = loadobj(self) % {{{
46+ % This function is directly called by matlab when a model object is
47+ % loaded. Update old properties here
48
49- if isempty(fieldnames(self.RecoveryAnalysis));
50- disp('WARNING: updating toolkits (RecoveryAnalysis not set)');
51- self.RecoveryAnalysis = self.DefaultAnalysis;
52- end
53- end% }}}
54- end
55- methods
56- function self = toolkits(varargin) % {{{
57- switch nargin
58- case 0
59- self=setdefaultparameters(self);
60- case 1
61- self=structtoobj(self,varargin{1});
62- otherwise
63- error('constructor not supported');
64- end
65- end % }}}
66- function self = addoptions(self,analysis,varargin) % {{{
67- % Usage example:
68- % md.toolkits=addoptions(md.toolkits,'StressbalanceAnalysis',FSoptions());
69- % md.toolkits=addoptions(md.toolkits,'StressbalanceAnalysis');
70+ if isempty(fieldnames(self.RecoveryAnalysis));
71+ disp('WARNING: updating toolkits (RecoveryAnalysis not set)');
72+ self.RecoveryAnalysis = self.DefaultAnalysis;
73+ end
74+ end% }}}
75+ end
76+ methods
77+ function self = toolkits(varargin) % {{{
78+ switch nargin
79+ case 0
80+ self=setdefaultparameters(self);
81+ case 1
82+ self=structtoobj(self,varargin{1});
83+ otherwise
84+ error('constructor not supported');
85+ end
86+ end % }}}
87+ function self = addoptions(self,analysis,varargin) % {{{
88+ %ADDOPTIONS - add analysis to md.toolkits.analysis
89+ %
90+ % Optional third parameter adds toolkits options to analysis.
91+ %
92+ % Usage:
93+ % md.toolkits=addoptions(md.toolkits,'StressbalanceAnalysis',FSoptions());
94+ % md.toolkits=addoptions(md.toolkits,'StressbalanceAnalysis');
95
96- %Create dynamic property if property does not exist yet
97- if ~ismember(analysis,properties(self)),
98- self.addprop(analysis);
99- end
100+ %Create dynamic property if property does not exist yet
101+ if ~ismember(analysis,properties(self)),
102+ self.addprop(analysis);
103+ end
104
105- %Add toolkits options to analysis
106- if nargin==3, self.(analysis) = varargin{1}; end
107- end
108- %}}}
109- function self = setdefaultparameters(self) % {{{
110+ %Add toolkits options to analysis
111+ if nargin==3,
112+ self.(analysis) = varargin{1};
113+ end
114+ end
115+ %}}}
116+ function self = setdefaultparameters(self) % {{{
117
118- %default toolkits:
119- if IssmConfig('_HAVE_PETSC_'),
120- %MUMPS is the default toolkits
121- if IssmConfig('_HAVE_MUMPS_'),
122- self.DefaultAnalysis = mumpsoptions();
123- else
124- self.DefaultAnalysis = iluasmoptions();
125- end
126- else
127- if IssmConfig('_HAVE_MUMPS_'),
128- self.DefaultAnalysis = issmmumpssolver();
129- elseif IssmConfig('_HAVE_GSL_'),
130- self.DefaultAnalysis = issmgslsolver();
131- else
132- disp('WARNING: Need at least Mumps or Gsl to define an issm solver type, no default solver assigned');
133- end
134- end
135+ %default toolkits:
136+ if IssmConfig('_HAVE_PETSC_'),
137+ %MUMPS is the default toolkits
138+ if IssmConfig('_HAVE_MUMPS_'),
139+ self.DefaultAnalysis = mumpsoptions();
140+ else
141+ self.DefaultAnalysis = iluasmoptions();
142+ end
143+ else
144+ if IssmConfig('_HAVE_MUMPS_'),
145+ self.DefaultAnalysis = issmmumpssolver();
146+ elseif IssmConfig('_HAVE_GSL_'),
147+ self.DefaultAnalysis = issmgslsolver();
148+ else
149+ disp('WARNING: Need at least Mumps or Gsl to define an issm solver type, no default solver assigned');
150+ end
151+ end
152
153- %Use same solver for Recovery mode
154- self.RecoveryAnalysis = self.DefaultAnalysis;
155+ %Use same solver for Recovery mode
156+ self.RecoveryAnalysis = self.DefaultAnalysis;
157
158
159- end % }}}
160- function disp(self) % {{{
161- analyses=properties(self);
162- disp(sprintf('List of toolkits options per analysis:\n'));
163- for i=1:numel(analyses),
164- analysis=analyses{i};
165- disp([analysis ':']);
166- disp(self.(analysis));
167- end
168- end % }}}
169- function md = checkconsistency(self,md,solution,analyses) % {{{
170- analyses=properties(self);
171- for i=1:numel(analyses),
172- switch analyses{i}
173- case 'DefaultAnalysis'
174- case 'RecoveryAnalysis'
175- case 'StressbalanceAnalysis'
176- case 'GLheightadvectionAnalysis'
177- case 'MasstransportAnalysis'
178- case 'ThermalAnalysis'
179- case 'EnthalpyAnalysis'
180- case 'AdjointBalancethicknessAnalysis'
181- case 'BalancethicknessAnalysis'
182- case 'Balancethickness2Analysis'
183- case 'BalancethicknessSoftAnalysis'
184- case 'BalancevelocityAnalysis'
185- case 'DamageEvolutionAnalysis'
186- case 'LoveAnalysis'
187- case 'EsaAnalysis'
188- case 'SealevelriseAnalysis'
189- otherwise
190+ end % }}}
191+ function disp(self) % {{{
192+ analyses=properties(self);
193+ disp(sprintf('List of toolkits options per analysis:\n'));
194+ for i=1:numel(analyses),
195+ analysis=analyses{i};
196+ disp([analysis ':']);
197+ disp(self.(analysis));
198+ end
199+ end % }}}
200+ function md = checkconsistency(self,md,solution,analyses) % {{{
201+ analyses=properties(self);
202+ for i=1:numel(analyses),
203+ switch analyses{i}
204+ case 'DefaultAnalysis'
205+ case 'RecoveryAnalysis'
206+ case 'StressbalanceAnalysis'
207+ case 'GLheightadvectionAnalysis'
208+ case 'MasstransportAnalysis'
209+ case 'ThermalAnalysis'
210+ case 'EnthalpyAnalysis'
211+ case 'AdjointBalancethicknessAnalysis'
212+ case 'BalancethicknessAnalysis'
213+ case 'Balancethickness2Analysis'
214+ case 'BalancethicknessSoftAnalysis'
215+ case 'BalancevelocityAnalysis'
216+ case 'DamageEvolutionAnalysis'
217+ case 'LoveAnalysis'
218+ case 'EsaAnalysis'
219+ case 'SealevelriseAnalysis'
220+ otherwise
221 md = checkmessage(md,['md.toolkits.' analyses{i} ' not supported yet']);
222- end
223- if isempty(fieldnames(self.(analyses{i})))
224- md = checkmessage(md,['md.toolkits.' analyses{i} ' is empty']);
225- end
226- end
227- end % }}}
228- function ToolkitsFile(toolkits,filename) % {{{
229- %TOOLKITSFILE - build toolkits file
230- %
231- % Build a Petsc compatible options file, from the toolkits model field + return options string.
232- % This file will also be used when the toolkit used is 'issm' instead of 'petsc'
233- %
234- % Usage: ToolkitsFile(toolkits,filename);
235+ end
236+ if isempty(fieldnames(self.(analyses{i})))
237+ md = checkmessage(md,['md.toolkits.' analyses{i} ' is empty']);
238+ end
239+ end
240+ end % }}}
241+ function ToolkitsFile(toolkits,filename) % {{{
242+ %TOOLKITSFILE - build toolkits file
243+ %
244+ % Build a Petsc compatible options file, from the toolkits model field + return options string.
245+ % This file will also be used when the toolkit used is 'issm' instead of 'petsc'
246+ %
247+ % Usage: ToolkitsFile(toolkits,filename);
248
249- %open file for writing
250- fid=fopen(filename,'w');
251- if fid==-1,
252- error(['ToolkitsFile error: could not open ' filename ' for writing']);
253- end
254+ %open file for writing
255+ fid=fopen(filename,'w');
256+ if fid==-1,
257+ error(['ToolkitsFile error: could not open ' filename ' for writing']);
258+ end
259
260- %write header
261- fprintf(fid,'%s%s%s\n','%Toolkits options file: ',filename,' written from Matlab toolkits array');
262+ %write header
263+ fprintf(fid,'%s%s%s\n','%Toolkits options file: ',filename,' written from Matlab toolkits array');
264
265- %start writing options
266- analyses=properties(toolkits);
267- for i=1:numel(analyses),
268- analysis=analyses{i};
269- options=toolkits.(analysis);
270+ %start writing options
271+ analyses=properties(toolkits);
272+ for i=1:numel(analyses),
273+ analysis=analyses{i};
274+ options=toolkits.(analysis);
275
276- %first write analysis:
277- fprintf(fid,'\n+%s\n',analysis); %append a + to recognize it's an analysis enum
278+ %first write analysis:
279+ fprintf(fid,'\n+%s\n',analysis); %append a + to recognize it's an analysis enum
280
281- %now, write options
282- optionslist=fieldnames(options);
283- for j=1:numel(optionslist),
284- optionname=optionslist{j};
285- optionvalue=options.(optionname);
286+ %now, write options
287+ optionslist=fieldnames(options);
288+ for j=1:numel(optionslist),
289+ optionname=optionslist{j};
290+ optionvalue=options.(optionname);
291
292- if isempty(optionvalue),
293- %this option has only one argument
294- fprintf(fid,'-%s\n',optionname);
295- else
296- %option with value. value can be string or scalar
297- if isnumeric(optionvalue),
298- fprintf(fid,'-%s %g\n',optionname,optionvalue);
299- elseif ischar(optionvalue),
300- fprintf(fid,'-%s %s\n',optionname,optionvalue);
301- else
302- error(['ToolkitsFile error: option ' optionname ' is not well formatted']);
303- end
304- end
305- end
306- end
307+ if isempty(optionvalue),
308+ %this option has only one argument
309+ fprintf(fid,'-%s\n',optionname);
310+ else
311+ %option with value. value can be string or scalar
312+ if isnumeric(optionvalue),
313+ fprintf(fid,'-%s %g\n',optionname,optionvalue);
314+ elseif ischar(optionvalue),
315+ fprintf(fid,'-%s %s\n',optionname,optionvalue);
316+ else
317+ error(['ToolkitsFile error: option ' optionname ' is not well formatted']);
318+ end
319+ end
320+ end
321+ end
322
323- fclose(fid);
324- end %}}}
325- function savemodeljs(self,fid,modelname) % {{{
326+ fclose(fid);
327+ end %}}}
328+ function savemodeljs(self,fid,modelname) % {{{
329
330- analyses=properties(self);
331- for i=1:numel(analyses),
332- if isempty(fieldnames(self.(analyses{i})))
333- error(['md.toolkits.' analyses{i} ' is empty']);
334- else
335- writejsstruct(fid,[modelname '.toolkits.' analyses{i}],self.(analyses{i}));
336- end
337- end
338+ analyses=properties(self);
339+ for i=1:numel(analyses),
340+ if isempty(fieldnames(self.(analyses{i})))
341+ error(['md.toolkits.' analyses{i} ' is empty']);
342+ else
343+ writejsstruct(fid,[modelname '.toolkits.' analyses{i}],self.(analyses{i}));
344+ end
345+ end
346
347- end % }}}
348- end
349+ end % }}}
350+ end
351 end
352Index: ../trunk-jpl/src/m/classes/giamme.py
353===================================================================
354--- ../trunk-jpl/src/m/classes/giamme.py (revision 25168)
355+++ ../trunk-jpl/src/m/classes/giamme.py (revision 25169)
356@@ -55,8 +55,8 @@
357 WriteData(fid, prefix, 'name', 'md.gia.model', 'data', 3, 'format', 'Integer')
358
359 if isinstance(self.Ngia, list) or self.Ngia.ndim == 1: #list or 1D numpy.ndarray
360- WriteData(fid, prefix, 'data', np.zeros(md.mesh.numberofvertices), 'format', 'DoubleMat', 'mattype', 1, 'name', 'md.gia.Ngia')
361- WriteData(fid, prefix, 'data', np.zeros(md.mesh.numberofvertices), 'format', 'DoubleMat', 'mattype', 1, 'name', 'md.gia.Ugia')
362+ WriteData(fid, prefix, 'data', np.zeros((md.mesh.numberofvertices, 1)), 'format', 'DoubleMat', 'mattype', 1, 'name', 'md.gia.Ngia')
363+ WriteData(fid, prefix, 'data', np.zeros((md.mesh.numberofvertices, 1)), 'format', 'DoubleMat', 'mattype', 1, 'name', 'md.gia.Ugia')
364 WriteData(fid, prefix, 'data', 1, 'name', 'md.gia.modelid', 'format', 'Double')
365 WriteData(fid, prefix, 'name', 'md.gia.nummodels', 'data', 1, 'format', 'Integer')
366 else: #2D numpy.ndarray
367Index: ../trunk-jpl/src/m/classes/dslmme.py
368===================================================================
369--- ../trunk-jpl/src/m/classes/dslmme.py (nonexistent)
370+++ ../trunk-jpl/src/m/classes/dslmme.py (revision 25169)
371@@ -0,0 +1,80 @@
372+import numpy as np
373+
374+from checkfield import checkfield
375+from fielddisplay import fielddisplay
376+from project3d import project3d
377+from WriteData import WriteData
378+
379+
380+class dslmme(object):
381+ '''
382+ DSLMME class definition
383+
384+ Usage:
385+ dsl = dslmme() #dynamic sea level class based on a multi-model ensemble of CMIP5 outputs
386+ '''
387+
388+ def __init__(self, *args): #{{{
389+ self.modelid = 0 #index into the multi-model ensemble
390+ self.global_average_thermosteric_sea_level_change = [] #corresponds to zostoga fields in CMIP5 archives. Specified as a temporally variable global rate (mm/yr) for each ensemble.
391+ self.sea_surface_height_change_above_geoid = [] #corresponds to zos fields in CMIP5 archives. Spatial average is 0. Specified as a spatio-temporally variable rate (mm/yr) for each ensemble
392+ self.sea_water_pressure_change_at_sea_floor = [] #corresponds to bpo fields in CMIP5 archives. Specified as a spatio-temporally variable rate (in mm/yr equivalent, not in Pa/yr!) for each ensemble
393+ self.compute_fingerprints = 0 #corresponds to zos fields in CMIP5 archives. Spatial average is 0. Specified as a spatio-temporally variable rate (mm/yr) for each ensemble
394+
395+ nargin = len(args)
396+
397+ if nargin == 0:
398+ self.setdefaultparameters()
399+ else:
400+ raise Exception('constructor not supported')
401+ #}}}
402+
403+ def __repr__(self): # {{{
404+ s = ' gia mme parameters:\n'
405+ s += '{}\n'.format(fielddisplay(self, 'modelid', 'index into the multi-model ensemble, determines which field will be used.'))
406+ s += '{}\n'.format(fielddisplay(self, 'Ngia', 'geoid (mm/yr).'))
407+ s += '{}\n'.format(fielddisplay(self, 'Ugia', 'uplift (mm/yr).'))
408+
409+ return s
410+ #}}}
411+
412+ def setdefaultparameters(self): # {{{
413+ return
414+ #}}}
415+
416+ def checkconsistency(self, md, solution, analyses): # {{{
417+ if ('SealevelriseAnalysis' not in analyses) or (solution == 'TransientSolution' and md.transient.isslr == 0):
418+ return md
419+
420+ for i in range(len(self.global_average_thermosteric_sea_level_change)):
421+ md = checkfield(md, 'field', self.global_average_thermosteric_sea_level_change[i], 'NaN', 1, 'Inf', 1)
422+ md = checkfield(md, 'field', self.sea_surface_height_change_above_geoid[i], 'NaN', 1, 'Inf', 1, 'timeseries', 1)
423+ md = checkfield(md, 'field', self.sea_water_pressure_change_at_sea_floor[i], 'NaN', 1, 'Inf', 1, 'timeseries', 1)
424+
425+ md = checkfield(md, 'field', self.modelid, 'NaN', 1, 'Inf', 1, '>=', 1, '<=',len(self.global_average_thermosteric_sea_level_change))
426+
427+ if self.compute_fingerprints:
428+ #check geodetic flag of slr is on
429+ if md.slr.geodetic==0,
430+ raise Exception('DSL checkconsistency error message: if bottom pressure fingerprints computations are requested, slr class should have geodetic flag on')
431+
432+ return md
433+ #}}}
434+
435+ def marshall(self, prefix, md, fid): #{{{
436+ WriteData(fid, prefix, 'name', 'md.dsl.model', 'data', 2, 'format', 'Integer')
437+ WriteData(fid, prefix, 'object', self, 'fieldname', 'compute_fingerprints', 'format', 'Integer')
438+ WriteData(fid, prefix, 'object', self, 'fieldname', 'modelid', 'format', 'Integer')
439+ WriteData(fid, prefix, 'name', 'md.dsl.nummodels', 'data', len(self.global_average_thermosteric_sea_level_change), 'format', 'Integer')
440+ WriteData(fid, prefix, 'object', self, 'fieldname', 'global_average_thermosteric_sea_level_change', 'format', 'MatArray', 'timeseries', 1, 'timeserieslength', 2, 'yts', md.constants.yts, 'scale', 1e-3 / md.constants.yts)
441+ WriteData(fid, prefix, 'object', self, 'fieldname', 'sea_water_pressure_change_at_sea_floor', 'format', 'MatArray', 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts, 'scale', 1e-3 / md.constants.yts)
442+ WriteData(fid, prefix, 'object', self, 'fieldname', 'sea_surface_height_change_above_geoid', 'format', 'MatArray', 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts, 'scale', 1e-3 / md.constants.yts)
443+ #}}}
444+
445+ def extrude(self, md): #{{{
446+ for i in range(len(self.global_average_thermosteric_sea_level_change)):
447+ self.sea_surface_height_change_above_geoid[i] = project3d(md, 'vector', self.self.sea_surface_height_change_above_geoid[i], 'type', 'node', 'layer', 1)
448+ self.sea_water_pressure_change_at_sea_floor[i] = project3d(md, 'vector', self.sea_water_pressure_change_at_sea_floor[i], 'type', 'node', 'layer', 1)
449+
450+ return self
451+ #}}}
452Index: ../trunk-jpl/src/m/classes/geometry.py
453===================================================================
454--- ../trunk-jpl/src/m/classes/geometry.py (revision 25168)
455+++ ../trunk-jpl/src/m/classes/geometry.py (revision 25169)
456@@ -1,41 +1,41 @@
457 import numpy as np
458+
459+from checkfield import checkfield
460+from fielddisplay import fielddisplay
461 from project3d import project3d
462-from fielddisplay import fielddisplay
463-from checkfield import checkfield
464 from WriteData import WriteData
465
466
467 class geometry(object):
468- """
469+ '''
470 GEOMETRY class definition
471
472- Usage:
473- geometry = geometry()
474- """
475+ Usage:
476+ geometry = geometry()
477+ '''
478
479- def __init__(self): # {{{
480- self.surface = float('NaN')
481- self.thickness = float('NaN')
482- self.base = float('NaN')
483- self.bed = float('NaN')
484- self.hydrostatic_ratio = float('NaN')
485+ def __init__(self): #{{{
486+ self.surface = np.nan
487+ self.thickness = np.nan
488+ self.base = np.nan
489+ self.bed = np.nan
490+ self.hydrostatic_ratio = np.nan
491
492- #set defaults
493+ #set defaults
494 self.setdefaultparameters()
495-
496 #}}}
497
498- def __repr__(self): # {{{
499-
500- string = " geometry parameters:"
501- string = "%s\n%s" % (string, fielddisplay(self, 'surface', 'ice upper surface elevation [m]'))
502- string = "%s\n%s" % (string, fielddisplay(self, 'thickness', 'ice thickness [m]'))
503- string = "%s\n%s" % (string, fielddisplay(self, 'base', 'ice base elevation [m]'))
504- string = "%s\n%s" % (string, fielddisplay(self, 'bed', 'bed elevation [m]'))
505- return string
506+ def __repr__(self): #{{{
507+ s = " geometry parameters:"
508+ s = "%s\n%s" % (s, fielddisplay(self, 'surface', 'ice upper surface elevation [m]'))
509+ s = "%s\n%s" % (s, fielddisplay(self, 'thickness', 'ice thickness [m]'))
510+ s = "%s\n%s" % (s, fielddisplay(self, 'base', 'ice base elevation [m]'))
511+ s = "%s\n%s" % (s, fielddisplay(self, 'bed', 'bed elevation [m]'))
512+
513+ return s
514 #}}}
515
516- def extrude(self, md): # {{{
517+ def extrude(self, md): #{{{
518 self.surface = project3d(md, 'vector', self.surface, 'type', 'node')
519 self.thickness = project3d(md, 'vector', self.thickness, 'type', 'node')
520 self.hydrostatic_ratio = project3d(md, 'vector', self.hydrostatic_ratio, 'type', 'node')
521@@ -44,13 +44,15 @@
522 return self
523 #}}}
524
525- def setdefaultparameters(self): # {{{
526+ def setdefaultparameters(self): #{{{
527 return self
528 #}}}
529
530- def checkconsistency(self, md, solution, analyses): # {{{
531+ def checkconsistency(self, md, solution, analyses): #{{{
532 if (solution == 'TransientSolution' and md.transient.isgia) or (solution == 'GiaSolution'):
533- md = checkfield(md, 'fieldname', 'geometry.thickness', 'NaN', 1, 'Inf', 1, 'timeseries', 1)
534+ md = checkfield(md, 'fieldname', 'geometry.thickness', 'timeseries', 1, 'NaN', 1, 'Inf', 1)
535+ elif solution == 'SealevelriseSolution':
536+ md = checkfield(md, 'fieldname', 'geometry.bed', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
537 elif solution == 'LoveSolution':
538 return
539 else:
540@@ -71,7 +73,7 @@
541 return md
542 # }}}
543
544- def marshall(self, prefix, md, fid): # {{{
545+ def marshall(self, prefix, md, fid): #{{{
546 WriteData(fid, prefix, 'object', self, 'fieldname', 'surface', 'format', 'DoubleMat', 'mattype', 1)
547 WriteData(fid, prefix, 'object', self, 'fieldname', 'thickness', 'format', 'DoubleMat', 'mattype', 1, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts)
548 WriteData(fid, prefix, 'object', self, 'fieldname', 'base', 'format', 'DoubleMat', 'mattype', 1)
549Index: ../trunk-jpl/src/m/classes/geometry.m
550===================================================================
551--- ../trunk-jpl/src/m/classes/geometry.m (revision 25168)
552+++ ../trunk-jpl/src/m/classes/geometry.m (revision 25169)
553@@ -60,7 +60,7 @@
554 elseif strcmpi(solution,'LoveSolution'),
555 return;
556 else
557- md = checkfield(md,'fieldname','geometry.surface' ,'NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]);
558+ md = checkfield(md,'fieldname','geometry.surface' ,'NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]);
559 md = checkfield(md,'fieldname','geometry.base' ,'NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1]);
560 md = checkfield(md,'fieldname','geometry.thickness','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1],'>',0);
561 if any(abs(self.thickness-self.surface+self.base)>10^-9),
562Index: ../trunk-jpl/src/m/classes/matice.m
563===================================================================
564--- ../trunk-jpl/src/m/classes/matice.m (revision 25168)
565+++ ../trunk-jpl/src/m/classes/matice.m (revision 25169)
566@@ -67,7 +67,7 @@
567 self.rho_freshwater=1000.;
568
569 %water viscosity (N.s/m^2)
570- self.mu_water=0.001787;
571+ self.mu_water=0.001787;
572
573 %ice heat capacity cp (J/kg/K)
574 self.heatcapacity=2093.;
575@@ -177,7 +177,6 @@
576 WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_B','format','DoubleMat','mattype',1,'timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts);
577 WriteData(fid,prefix,'object',self,'class','materials','fieldname','rheology_n','format','DoubleMat','mattype',2);
578 WriteData(fid,prefix,'data',self.rheology_law,'name','md.materials.rheology_law','format','String');
579-
580 WriteData(fid,prefix,'object',self,'class','materials','fieldname','lithosphere_shear_modulus','format','Double');
581 WriteData(fid,prefix,'object',self,'class','materials','fieldname','lithosphere_density','format','Double','scale',10^3);
582 WriteData(fid,prefix,'object',self,'class','materials','fieldname','mantle_shear_modulus','format','Double');
583Index: ../trunk-jpl/src/m/classes/dsl.py
584===================================================================
585--- ../trunk-jpl/src/m/classes/dsl.py (revision 25168)
586+++ ../trunk-jpl/src/m/classes/dsl.py (revision 25169)
587@@ -1,19 +1,20 @@
588 import numpy as np
589+
590+from checkfield import checkfield
591 from fielddisplay import fielddisplay
592-from checkfield import checkfield
593+from project3d import project3d
594 from WriteData import WriteData
595-from project3d import project3d
596
597
598 class dsl(object):
599- """
600- dsl Class definition
601+ '''
602+ DSL - slass definition
603
604- Usage:
605- dsl = dsl()
606- """
607+ Usage:
608+ dsl = dsl()
609+ '''
610
611- def __init__(self): # {{{
612+ def __init__(self): #{{{
613 self.global_average_thermosteric_sea_level_change = 0 #corresponds to zostoga field in CMIP5 archives. Specified as a temporally variable global rate (mm/yr)
614 self.sea_surface_height_change_above_geoid = float('NaN') #corresponds to zos field in CMIP5 archives. Spatial average is 0. Specified as a spatio-temporally variable rate (mm/yr)
615 self.sea_water_pressure_change_at_sea_floor = float('NaN') #corresponds to bpo field in CMIP5 archives. Specified as a spatio-temporally variable rate (in Pa/yr)
616@@ -20,28 +21,29 @@
617 self.compute_fingerprints = 0; #do we use the sea water pressure change to compute fingerprints and correct sea_surface_height_change_above_geoid
618 #}}}
619
620- def __repr__(self): # {{{
621- string = " dsl parameters:"
622- string = "%s\n%s" % (string, fielddisplay(self, 'global_average_thermosteric_sea_level_change', 'corresponds to zostoga field in CMIP5 archives. Specified as a temporally variable global rate (mm/yr)'))
623- string = "%s\n%s" % (string, fielddisplay(self, 'sea_surface_height_change_above_geoid', 'corresponds to zos field in CMIP5 archives. Spatial average is 0. Specified as a spatio-temporally variable rate (mm/yr)'))
624- string = "%s\n%s" % (string, fielddisplay(self, 'sea_water_pressure_change_at_sea_floor', 'corresponds to bpo field in CMIP5 archives. Specified as a spatio-temporally variable rate (in Pa/yr)'))
625- string = "%s\n%s" % (string, fielddisplay(self, 'compute_fingerprints', 'do we use the sea water pressure change to compute fingerprints and correct sea_surface_height_change_above_geoid'))
626- return string
627+ def __repr__(self): #{{{
628+ s = " dsl parameters:"
629+ s = "%s\n%s" % (s, fielddisplay(self, 'global_average_thermosteric_sea_level_change', 'corresponds to zostoga field in CMIP5 archives. Specified as a temporally variable global rate (mm/yr)'))
630+ s = "%s\n%s" % (s, fielddisplay(self, 'sea_surface_height_change_above_geoid', 'corresponds to zos field in CMIP5 archives. Spatial average is 0. Specified as a spatio-temporally variable rate (mm/yr)'))
631+ s = "%s\n%s" % (s, fielddisplay(self, 'sea_water_pressure_change_at_sea_floor', 'corresponds to bpo field in CMIP5 archives. Specified as a spatio-temporally variable rate (in Pa/yr)'))
632+ s = "%s\n%s" % (s, fielddisplay(self, 'compute_fingerprints', 'do we use the sea water pressure change to compute fingerprints and correct sea_surface_height_change_above_geoid'))
633+
634+ return s
635 #}}}
636
637- def extrude(self, md): # {{{
638+ def extrude(self, md): #{{{
639 self.sea_surface_height_change_above_geoid = project3d(md, 'vector', self.sea_surface_height_change_above_geoid, 'type', 'node')
640 self.sea_water_pressure_change_at_sea_floor = project3d(md, 'vector', self.sea_water_pressure_change_at_sea_floor, 'type', 'node')
641 return self
642 #}}}
643
644- def defaultoutputs(self, md): # {{{
645+ def defaultoutputs(self, md): #{{{
646 return []
647 #}}}
648
649- def checkconsistency(self, md, solution, analyses): # {{{
650+ def checkconsistency(self, md, solution, analyses): #{{{
651 # Early retun
652- if not 'SealevelriseAnalysis' in analyses:
653+ if 'SealevelriseAnalysis' not in analyses:
654 return md
655
656 if solution == 'TransientSolution' and md.transient.isslr == 0:
657@@ -59,7 +61,7 @@
658 return md
659 # }}}
660
661- def marshall(self, prefix, md, fid): # {{{
662+ def marshall(self, prefix, md, fid): #{{{
663 yts = md.constants.yts
664 WriteData(fid, prefix, 'name', 'md.dsl.model', 'data', 1, 'format', 'Integer')
665 WriteData(fid, prefix, 'object', self, 'class', 'dsl', 'fieldname', 'compute_fingerprints', 'format', 'Integer')
666Index: ../trunk-jpl/src/m/classes/nodalvalue.py
667===================================================================
668--- ../trunk-jpl/src/m/classes/nodalvalue.py (nonexistent)
669+++ ../trunk-jpl/src/m/classes/nodalvalue.py (revision 25169)
670@@ -0,0 +1,69 @@
671+import numpy as np
672+
673+from checkfield import checkfield
674+from fielddisplay import fielddisplay
675+from WriteData import WriteData
676+
677+
678+class dslmme(object):
679+ '''
680+ NODALVALUE class definition
681+
682+ Usage:
683+ nodalvalue=nodalvalue()
684+ nodalvalue=nodalvalue(
685+ 'name', 'SealevelriseSNodalValue',
686+ 'definitionstring', 'Outputdefinition1',
687+ 'model_string', 'SealevelriseS',
688+ 'node', 1
689+ )
690+ '''
691+
692+ def __init__(self, *args): #{{{
693+ self.name = ''
694+ self.definitionstring = '' #string that identifies this output definition uniquely, from 'Outputdefinition[1-10]'
695+ self.model_string = '' #string for field that is being retrieved
696+ self.node = np.nan #for which node are we retrieving the value?
697+
698+ #use provided options to change fields
699+ options = pairoptions(*args)
700+
701+ #get name
702+ self.name = options.getfieldvalue(options, 'name', '')
703+ self.definitionstring = options.getfieldvalue(options, 'definitionstring', '')
704+ self.model_string = options.getfieldvalue(options, 'model_string', '')
705+ self.node = options.getfieldvalue(options, 'node', '')
706+ #}}}
707+
708+ def __repr__(self): # {{{
709+ s = ' Nodalvalue:\n'
710+ s += '{}\n'.format(fielddisplay(self, 'name', 'identifier for this nodalvalue response'))
711+ s += '{}\n'.format(fielddisplay(self, 'definitionstring', 'string that identifies this output definition uniquely, from \'Outputdefinition[1-10]\''))
712+ s += '{}\n'.format(fielddisplay(self, 'model_string', 'string for field that is being retrieved'))
713+ s += '{}\n'.format(fielddisplay(self, 'node', 'vertex index at which we retrieve the value'))
714+
715+ return s
716+ #}}}
717+
718+ def setdefaultparameters(self): # {{{
719+ return
720+ #}}}
721+
722+ def checkconsistency(self, md, solution, analyses): # {{{
723+ if not isinstance(self.name, basestring):
724+ raise Exception('nodalvalue error message: \'name\' field should be a string!')
725+ OutputdefinitionStringArray = []
726+ for i in range(100):
727+ OutputdefinitionStringArray[i] = 'Outputdefinition{}'.format(i)
728+ md = checkfield(md, 'fieldname', 'self.definitionstring', 'field', self.definitionstring, 'values', OutputdefinitionStringArray)
729+ md = checkfield(md, 'fieldname', 'self.node', 'field', self.node, 'values', range(md.mesh.numberofvertices))
730+
731+ return md
732+ #}}}
733+
734+ def marshall(self, prefix, md, fid): #{{{
735+ WriteData(fid, prefix, 'data', self.name, 'name', 'md.nodalvalue.name', 'format', 'String')
736+ WriteData(fid, prefix, 'data', self.definitionstring, 'name', 'md.nodalvalue.definitionenum', 'format', 'String')
737+ WriteData(fid, prefix, 'data', self.model_string, 'name', 'md.nodalvalue.model_enum', 'format', 'String')
738+ WriteData(fid, prefix, 'data', self.node, 'name', 'md.nodalvalue.node', 'format', 'Integer')
739+ #}}}
740Index: ../trunk-jpl/src/m/classes/solidearth.py
741===================================================================
742--- ../trunk-jpl/src/m/classes/solidearth.py (revision 25168)
743+++ ../trunk-jpl/src/m/classes/solidearth.py (revision 25169)
744@@ -23,7 +23,7 @@
745 self.sealevel = np.nan
746 self.settings = solidearthsettings()
747 self.surfaceload = surfaceload()
748- self.love = lovenumbers()
749+ self.lovenumbers = lovenumbers()
750 self.rotational = rotational()
751 self.planetradius = planetradius('earth')
752 self.requested_outputs = ['default']
753@@ -64,7 +64,7 @@
754
755 self.settings.checkconsistency(md, solution, analyses)
756 self.surfaceload.checkconsistency(md, solution, analyses)
757- self.love.checkconsistency(md, solution, analyses)
758+ self.lovenumbers.checkconsistency(md, solution, analyses)
759 self.rotational.checkconsistency(md, solution, analyses)
760
761 return md
762@@ -81,7 +81,7 @@
763
764 self.settings.marshall(prefix, md, fid)
765 self.surfaceload.marshall(prefix, md, fid)
766- self.love.marshall(prefix, md, fid)
767+ self.lovenumbers.marshall(prefix, md, fid)
768 self.rotational.marshall(prefix, md, fid)
769
770 #process requested outputs
771Index: ../trunk-jpl/src/m/classes/toolkits.py
772===================================================================
773--- ../trunk-jpl/src/m/classes/toolkits.py (revision 25168)
774+++ ../trunk-jpl/src/m/classes/toolkits.py (revision 25169)
775@@ -7,14 +7,14 @@
776
777
778 class toolkits(object):
779- """
780+ '''
781 TOOLKITS class definition
782
783- Usage:
784- self = toolkits()
785- """
786+ Usage:
787+ self = toolkits()
788+ '''
789
790- def __init__(self): # {{{
791+ def __init__(self): #{{{
792 #default toolkits
793 if IssmConfig('_HAVE_PETSC_')[0]:
794 #MUMPS is the default toolkits
795@@ -33,18 +33,18 @@
796 #Use same solver for Recovery mode
797 self.RecoveryAnalysis = self.DefaultAnalysis
798
799- #The other properties are dynamic
800- # }}}
801+ #The other properties are dynamic
802+ #}}}
803
804- def __repr__(self): # {{{
805+ def __repr__(self): #{{{
806 s = "List of toolkits options per analysis:\n\n"
807 for analysis in list(vars(self).keys()):
808- s += "%s\n" % fielddisplay(self, analysis, '')
809+ s += "{}\n".format(fielddisplay(self, analysis, ''))
810
811- return s
812- # }}}
813+ return s
814+ #}}}
815
816- def addoptions(self, analysis, *args): # {{{
817+ def addoptions(self, analysis, *args): #{{{
818 # Usage example:
819 # md.toolkits = addoptions(md.toolkits, 'StressbalanceAnalysis', FSoptions())
820 # md.toolkits = addoptions(md.toolkits, 'StressbalanceAnalysis')
821@@ -58,43 +58,49 @@
822 setattr(self, analysis, args[0])
823
824 return self
825- # }}}
826+ #}}}
827
828- def checkconsistency(self, md, solution, analyses): # {{{
829+ def checkconsistency(self, md, solution, analyses): #{{{
830+ # TODO
831+ # - Implement something closer to a switch as in
832+ # src/m/classes/toolkits.m?
833+ #
834 for analysis in list(vars(self).keys()):
835 if not getattr(self, analysis):
836- md.checkmessage("md.toolkits.%s is empty" % analysis)
837+ md.checkmessage("md.toolkits.{} is empty".format(analysis))
838
839 return md
840- # }}}
841+ #}}}
842
843- def ToolkitsFile(self, filename): # {{{
844- """
845+ def ToolkitsFile(self, filename): #{{{
846+ '''
847 TOOLKITSFILE - build toolkits file
848
849- Build a Petsc compatible options file, from the toolkits model field + return options string
850- This file will also be used when the toolkit used is 'issm' instead of 'petsc'
851+ Build a Petsc compatible options file, from the toolkits model
852+ field + return options string.
853+ This file will also be used when the toolkit used is 'issm' instead
854+ of 'petsc'.s
855
856+ Usage:
857+ ToolkitsFile(toolkits, filename)
858+ '''
859
860- Usage: ToolkitsFile(toolkits, filename)
861- """
862-
863- #open file for writing
864+ #open file for writing
865 try:
866 fid = open(filename, 'w')
867 except IOError as e:
868 raise IOError("ToolkitsFile error: could not open {}' for writing due to".format(filename), e)
869
870- #write header
871+ #write header
872 fid.write("%s%s%s\n" % ('%Toolkits options file: ', filename, ' written from Python toolkits array'))
873
874- #start writing options
875+ #start writing options
876 for analysis in list(vars(self).keys()):
877 options = getattr(self, analysis)
878
879- #first write analysis:
880+ #first write analysis:
881 fid.write("\n+{}\n".format(analysis)) #append a + to recognize it's an analysis enum
882- #now, write options
883+ #now, write options
884 for optionname, optionvalue in list(options.items()):
885
886 if not optionvalue:
887@@ -110,4 +116,4 @@
888 raise TypeError("ToolkitsFile error: option '{}' is not well formatted.".format(optionname))
889
890 fid.close()
891- # }}}
892+ #}}}
893Index: ../trunk-jpl/src/m/classes/matdamageice.m
894===================================================================
895--- ../trunk-jpl/src/m/classes/matdamageice.m (revision 25168)
896+++ ../trunk-jpl/src/m/classes/matdamageice.m (revision 25169)
897@@ -119,7 +119,7 @@
898 md = checkfield(md,'fieldname','materials.rheology_n','>',0,'size',[md.mesh.numberofelements 1]);
899 md = checkfield(md,'fieldname','materials.rheology_law','values',{'None' 'BuddJacka' 'Cuffey' 'CuffeyTemperate' 'Paterson' 'Arrhenius' 'LliboutryDuval'});
900 md = checkfield(md,'fieldname','materials.effectiveconductivity_averaging','numel',[1],'values',[0 1 2]);
901-
902+
903 if ismember('GiaAnalysis',analyses),
904 md = checkfield(md,'fieldname','materials.lithosphere_shear_modulus','>',0,'numel',1);
905 md = checkfield(md,'fieldname','materials.lithosphere_density','>',0,'numel',1);
906Index: ../trunk-jpl/src/m/classes/dslmme.m
907===================================================================
908--- ../trunk-jpl/src/m/classes/dslmme.m (revision 25168)
909+++ ../trunk-jpl/src/m/classes/dslmme.m (revision 25169)
910@@ -14,12 +14,6 @@
911
912 end
913 methods
914- function self = extrude(self,md) % {{{
915- for i=1:length(self.global_average_thermosteric_sea_level_change),
916- self.sea_surface_height_change_above_geoid{i}=project3d(md,'vector',self.sea_surface_height_change_above_geoid{i},'type','node','layer',1);
917- self.sea_water_pressure_change_at_sea_floor{i}=project3d(md,'vector',self.sea_water_pressure_change_at_sea_floor{i},'type','node','layer',1);
918- end
919- end % }}}
920 function self = dslmme(varargin) % {{{
921 switch nargin
922 case 0
923@@ -79,6 +73,12 @@
924 WriteData(fid,prefix,'object',self,'fieldname','sea_surface_height_change_above_geoid','format','MatArray','timeserieslength',md.mesh.numberofvertices+1,'yts',md.constants.yts,'scale',1e-3/md.constants.yts);
925
926 end % }}}
927+ function self = extrude(self,md) % {{{
928+ for i=1:length(self.global_average_thermosteric_sea_level_change),
929+ self.sea_surface_height_change_above_geoid{i}=project3d(md,'vector',self.sea_surface_height_change_above_geoid{i},'type','node','layer',1);
930+ self.sea_water_pressure_change_at_sea_floor{i}=project3d(md,'vector',self.sea_water_pressure_change_at_sea_floor{i},'type','node','layer',1);
931+ end
932+ end % }}}
933 function savemodeljs(self,fid,modelname) % {{{
934
935 writejsdouble(fid,[modelname '.dsl.modelid'],self.modelid);
936Index: ../trunk-jpl/src/m/solve/solve.py
937===================================================================
938--- ../trunk-jpl/src/m/solve/solve.py (revision 25168)
939+++ ../trunk-jpl/src/m/solve/solve.py (revision 25169)
940@@ -1,13 +1,13 @@
941 import datetime
942 import os
943 import shutil
944-from pairoptions import pairoptions
945+
946 from ismodelselfconsistent import ismodelselfconsistent
947+from loadresultsfromcluster import loadresultsfromcluster
948 from marshall import marshall
949+from pairoptions import pairoptions
950+from preqmu import *
951 from waitonlock import waitonlock
952-from loadresultsfromcluster import loadresultsfromcluster
953-from preqmu import *
954-#from MatlabFuncs import *
955
956
957 def solve(md, solutionstring, *args):
958@@ -14,35 +14,37 @@
959 '''
960 SOLVE - apply solution sequence for this model
961
962- Usage:
963- md = solve(md, solutionstring, varargin)
964- where varargin is a list of paired arguments of string OR enums
965+ Usage:
966+ md = solve(md, solutionstring, varargin)
967+ where varargin is a list of paired arguments of string OR enums
968
969 solution types available comprise:
970- - 'Stressbalance' or 'sb'
971- - 'Masstransport' or 'mt'
972- - 'Thermal' or 'th'
973- - 'Steadystate' or 'ss'
974- - 'Transient' or 'tr'
975- - 'Balancethickness' or 'mc'
976- - 'Balancevelocity' or 'bv'
977- - 'BedSlope' or 'bsl'
978- - 'SurfaceSlope' or 'ssl'
979- - 'Hydrology' or 'hy'
980- - 'DamageEvolution' or 'da'
981- - 'Gia' or 'gia'
982- - 'Esa' or 'esa'
983- - 'Sealevelrise' or 'slr'
984- - 'Love' or 'lv'
985+ - 'Stressbalance' or 'sb'
986+ - 'Masstransport' or 'mt'
987+ - 'Thermal' or 'th'
988+ - 'Steadystate' or 'ss'
989+ - 'Transient' or 'tr'
990+ - 'Balancethickness' or 'mc'
991+ - 'Balancevelocity' or 'bv'
992+ - 'BedSlope' or 'bsl'
993+ - 'SurfaceSlope' or 'ssl'
994+ - 'Hydrology' or 'hy'
995+ - 'DamageEvolution' or 'da'
996+ - 'Gia' or 'gia'
997+ - 'Esa' or 'esa'
998+ - 'Sealevelrise' or 'slr'
999+ - 'Love' or 'lv'
1000
1001- extra options:
1002- - loadonly : does not solve. only load results
1003- - checkconsistency : 'yes' or 'no' (default is 'yes'), ensures checks on consistency of model
1004- - restart: 'directory name (relative to the execution directory) where the restart file is located.
1005+ extra options:
1006+ - loadonly : does not solve. only load results
1007+ - checkconsistency : 'yes' or 'no' (default is 'yes'), ensures
1008+ checks on consistency of model
1009+ - restart : 'directory name (relative to the execution
1010+ directory) where the restart file is located
1011
1012- Examples:
1013- md = solve(md, 'Stressbalance')
1014- md = solve(md, 'sb')
1015+ Examples:
1016+ md = solve(md, 'Stressbalance')
1017+ md = solve(md, 'sb')
1018 '''
1019
1020 #recover and process solve options
1021Index: ../trunk-jpl/src/m/consistency/checkfield.py
1022===================================================================
1023--- ../trunk-jpl/src/m/consistency/checkfield.py (revision 25168)
1024+++ ../trunk-jpl/src/m/consistency/checkfield.py (revision 25169)
1025@@ -7,32 +7,35 @@
1026
1027
1028 def checkfield(md, *args):
1029- """
1030+ '''
1031 CHECKFIELD - check field consistency
1032
1033- Used to check model consistency.,
1034- Requires:
1035- 'field' or 'fieldname' option. If 'fieldname' is provided, it will retrieve it from the model md. (md.(fieldname))
1036- If 'field' is provided, it will assume the argument following 'field' is a numeric array.
1037+ Used to check model consistency
1038+
1039+ Requires:
1040+ 'field' or 'fieldname' option. If 'fieldname' is provided, it will retrieve it from the model md. (md.(fieldname))
1041+ If 'field' is provided, it will assume the argument following
1042+ 'field' is a numeric array.
1043
1044- Available options:
1045- - NaN: 1 if check that there is no NaN
1046- - size: [lines cols], NaN for non checked dimensions, or 'universal' for any input type (nodal, element, time series, etc)
1047- -> : greater than provided value
1048- ->= : greater or equal to provided value
1049- - < : smallerthan provided value
1050- - <=: smaller or equal to provided value
1051- - < vec: smallerthan provided values on each vertex
1052- - timeseries: 1 if check time series consistency (size and time)
1053- - values: cell of strings or vector of acceptable values
1054- - numel: list of acceptable number of elements
1055- - cell: 1 if check that is cell
1056- - empty: 1 if check that non empty
1057- - message: overloaded error message
1058+ Available options:
1059+ - NaN: 1 if check that there is no NaN
1060+ - size: [lines cols], NaN for non checked dimensions, or 'universal'
1061+ for any input type (nodal, element, time series, etc)
1062+ - > : greater than provided value
1063+ - >= : greater or equal to provided value
1064+ - < : smallerthan provided value
1065+ - <=: smaller or equal to provided value
1066+ - < vec: smallerthan provided values on each vertex
1067+ - timeseries: 1 if check time series consistency (size and time)
1068+ - values: list of acceptable values
1069+ - numel: list of acceptable number of elements
1070+ - cell: 1 if check that is cell
1071+ - empty: 1 if check that non empty
1072+ - message: overloaded error message
1073
1074- Usage:
1075- md = checkfield(md, fieldname, options)
1076- """
1077+ Usage:
1078+ md = checkfield(md, fieldname, options)
1079+ '''
1080
1081 #get options
1082 options = pairoptions(*args)
1083@@ -142,7 +145,7 @@
1084 #check cell
1085 if options.getfieldvalue('cell', 0):
1086 if not isinstance(field, (tuple, list, dict)):
1087- md = md.checkmessage(options.getfieldvalue('message', "field '{}' should be a cell".format(fieldname)))
1088+ md = md.checkmessage(options.getfieldvalue('message', "field '{}' should be a tuple, list, or dict".format(fieldname)))
1089
1090 #check values
1091 if options.exist('values'):
1092Index: ../trunk-jpl/test/NightlyRun/test2002.py
1093===================================================================
1094--- ../trunk-jpl/test/NightlyRun/test2002.py (revision 25168)
1095+++ ../trunk-jpl/test/NightlyRun/test2002.py (revision 25169)
1096@@ -38,7 +38,7 @@
1097
1098 #mask: {{{
1099 mask = gmtmask(md.mesh.lat, md.mesh.long)
1100-icemask = np.ones(md.mesh.numberofvertices)
1101+icemask = np.ones((md.mesh.numberofvertices, 1))
1102 pos = np.where(mask == 0)[0]
1103 icemask[pos] = -1
1104 pos = np.where(np.sum(mask[md.mesh.elements - 1], axis=1) < 3)
1105@@ -47,8 +47,8 @@
1106 md.mask.ocean_levelset = -icemask
1107
1108 #make sure that the ice level set is all inclusive:
1109-md.mask.land_levelset = np.zeros((md.mesh.numberofvertices))
1110-md.mask.ocean_levelset = -np.ones((md.mesh.numberofvertices))
1111+md.mask.land_levelset = np.zeros((md.mesh.numberofvertices, 1))
1112+md.mask.ocean_levelset = -np.ones((md.mesh.numberofvertices, 1))
1113
1114 #make sure that the elements that have loads are fully grounded
1115 pos = np.nonzero(md.solidearth.surfaceload.icethicknesschange)[0]
1116Index: ../trunk-jpl/test/NightlyRun/test2002.m
1117===================================================================
1118--- ../trunk-jpl/test/NightlyRun/test2002.m (revision 25168)
1119+++ ../trunk-jpl/test/NightlyRun/test2002.m (revision 25169)
1120@@ -64,7 +64,9 @@
1121 md.solidearth.settings.maxiter=10;
1122
1123 %eustatic run:
1124-md.solidearth.settings.rigid=0; md.solidearth.settings.elastic=0;md.solidearth.settings.rotation=0;
1125+md.solidearth.settings.rigid=0;
1126+md.solidearth.settings.elastic=0;
1127+md.solidearth.settings.rotation=0;
1128 md=solve(md,'Sealevelrise');
1129 Seustatic=md.results.SealevelriseSolution.Sealevel;
1130
1131Index: ../trunk-jpl/src/m/classes/matice.py
1132===================================================================
1133--- ../trunk-jpl/src/m/classes/matice.py (revision 25168)
1134+++ ../trunk-jpl/src/m/classes/matice.py (revision 25169)
1135@@ -1,3 +1,5 @@
1136+import numpy as np
1137+
1138 from fielddisplay import fielddisplay
1139 from project3d import project3d
1140 from checkfield import checkfield
1141@@ -5,14 +7,14 @@
1142
1143
1144 class matice(object):
1145- """
1146+ '''
1147 MATICE class definition
1148
1149- Usage:
1150- matice = matice()
1151- """
1152+ Usage:
1153+ matice = matice()
1154+ '''
1155
1156- def __init__(self): # {{{
1157+ def __init__(self): #{{{
1158 self.rho_ice = 0.
1159 self.rho_water = 0.
1160 self.rho_freshwater = 0.
1161@@ -26,80 +28,80 @@
1162 self.beta = 0.
1163 self.mixed_layer_capacity = 0.
1164 self.thermal_exchange_velocity = 0.
1165- self.rheology_B = float('NaN')
1166- self.rheology_n = float('NaN')
1167+ self.rheology_B = np.nan
1168+ self.rheology_n = np.nan
1169 self.rheology_law = ''
1170
1171- #giaivins:
1172+ #giaivins
1173 self.lithosphere_shear_modulus = 0.
1174 self.lithosphere_density = 0.
1175 self.mantle_shear_modulus = 0.
1176 self.mantle_density = 0.
1177
1178- #SLR
1179- self.earth_density = 5512
1180+ #slr
1181+ self.earth_density = 0
1182 self.setdefaultparameters()
1183 #}}}
1184
1185- def __repr__(self): # {{{
1186- string = " Materials:"
1187+ def __repr__(self): #{{{
1188+ s = " Materials:"
1189+ s = "%s\n%s" % (s, fielddisplay(self, "rho_ice", "ice density [kg/m^3]"))
1190+ s = "%s\n%s" % (s, fielddisplay(self, "rho_water", "water density [kg/m^3]"))
1191+ s = "%s\n%s" % (s, fielddisplay(self, "rho_freshwater", "fresh water density [kg/m^3]"))
1192+ s = "%s\n%s" % (s, fielddisplay(self, "mu_water", "water viscosity [Ns/m^2]"))
1193+ s = "%s\n%s" % (s, fielddisplay(self, "heatcapacity", "heat capacity [J/kg/K]"))
1194+ s = "%s\n%s" % (s, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W/m/K]"))
1195+ s = "%s\n%s" % (s, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W/m/K]"))
1196+ s = "%s\n%s" % (s, fielddisplay(self, "effectiveconductivity_averaging", "computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)"))
1197+ s = "%s\n%s" % (s, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K"))
1198+ s = "%s\n%s" % (s, fielddisplay(self, "latentheat", "latent heat of fusion [J/m^3]"))
1199+ s = "%s\n%s" % (s, fielddisplay(self, "beta", "rate of change of melting point with pressure [K/Pa]"))
1200+ s = "%s\n%s" % (s, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W/kg/K]"))
1201+ s = "%s\n%s" % (s, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m/s]"))
1202+ s = "%s\n%s" % (s, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1/n)]"))
1203+ s = "%s\n%s" % (s, fielddisplay(self, "rheology_n", "Glen's flow law exponent"))
1204+ s = "%s\n%s" % (s, fielddisplay(self, "rheology_law", "law for the temperature dependance of the rheology: 'None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', or 'NyeH2O'"))
1205+ s = "%s\n%s" % (s, fielddisplay(self, "lithosphere_shear_modulus", "Lithosphere shear modulus [Pa]"))
1206+ s = "%s\n%s" % (s, fielddisplay(self, "lithosphere_density", "Lithosphere density [g/cm^-3]"))
1207+ s = "%s\n%s" % (s, fielddisplay(self, "mantle_shear_modulus", "Mantle shear modulus [Pa]"))
1208+ s = "%s\n%s" % (s, fielddisplay(self, "mantle_density", "Mantle density [g/cm^-3]"))
1209+ s = "%s\n%s" % (s, fielddisplay(self, "earth_density", "Mantle density [kg/m^-3]"))
1210
1211- string = "%s\n%s" % (string, fielddisplay(self, "rho_ice", "ice density [kg / m^3]"))
1212- string = "%s\n%s" % (string, fielddisplay(self, "rho_water", "water density [kg / m^3]"))
1213- string = "%s\n%s" % (string, fielddisplay(self, "rho_freshwater", "fresh water density [kg / m^3]"))
1214- string = "%s\n%s" % (string, fielddisplay(self, "mu_water", "water viscosity [N s / m^2]"))
1215- string = "%s\n%s" % (string, fielddisplay(self, "heatcapacity", "heat capacity [J / kg / K]"))
1216- string = "%s\n%s" % (string, fielddisplay(self, "thermalconductivity", "ice thermal conductivity [W / m / K]"))
1217- string = "%s\n%s" % (string, fielddisplay(self, "temperateiceconductivity", "temperate ice thermal conductivity [W / m / K]"))
1218- string = "%s\n%s" % (string, fielddisplay(self, "effectiveconductivity_averaging", "computation of effectiveconductivity: (0) arithmetic mean, (1) harmonic mean, (2) geometric mean (default)"))
1219- string = "%s\n%s" % (string, fielddisplay(self, "meltingpoint", "melting point of ice at 1atm in K"))
1220- string = "%s\n%s" % (string, fielddisplay(self, "latentheat", "latent heat of fusion [J / m^3]"))
1221- string = "%s\n%s" % (string, fielddisplay(self, "beta", "rate of change of melting point with pressure [K / Pa]"))
1222- string = "%s\n%s" % (string, fielddisplay(self, "mixed_layer_capacity", "mixed layer capacity [W / kg / K]"))
1223- string = "%s\n%s" % (string, fielddisplay(self, "thermal_exchange_velocity", "thermal exchange velocity [m / s]"))
1224- string = "%s\n%s" % (string, fielddisplay(self, "rheology_B", "flow law parameter [Pa s^(1 / n)]"))
1225- string = "%s\n%s" % (string, fielddisplay(self, "rheology_n", "Glen's flow law exponent"))
1226- 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'"))
1227- string = "%s\n%s" % (string, fielddisplay(self, "lithosphere_shear_modulus", "Lithosphere shear modulus [Pa]"))
1228- string = "%s\n%s" % (string, fielddisplay(self, "lithosphere_density", "Lithosphere density [g / cm^ - 3]"))
1229- string = "%s\n%s" % (string, fielddisplay(self, "mantle_shear_modulus", "Mantle shear modulus [Pa]"))
1230- string = "%s\n%s" % (string, fielddisplay(self, "mantle_density", "Mantle density [g / cm^ - 3]"))
1231- string = "%s\n%s" % (string, fielddisplay(self, "earth_density", "Mantle density [kg / m^ - 3]"))
1232- return string
1233+ return s
1234 #}}}
1235
1236- def extrude(self, md): # {{{
1237+ def extrude(self, md): #{{{
1238 self.rheology_B = project3d(md, 'vector', self.rheology_B, 'type', 'node')
1239 self.rheology_n = project3d(md, 'vector', self.rheology_n, 'type', 'element')
1240 return self
1241 #}}}
1242
1243- def setdefaultparameters(self): # {{{
1244- #ice density (kg / m^3)
1245+ def setdefaultparameters(self): #{{{
1246+ #ice density (kg/m^3)
1247 self.rho_ice = 917.
1248- #ocean water density (kg / m^3)
1249+ #ocean water density (kg/m^3)
1250 self.rho_water = 1023.
1251- #fresh water density (kg / m^3)
1252+ #fresh water density (kg/m^3)
1253 self.rho_freshwater = 1000.
1254- #water viscosity (N.s / m^2)
1255+ #water viscosity (N.s/m^2)
1256 self.mu_water = 0.001787
1257- #ice heat capacity cp (J / kg / K)
1258+ #ice heat capacity cp (J/kg/K)
1259 self.heatcapacity = 2093.
1260- #ice latent heat of fusion L (J / kg)
1261+ #ice latent heat of fusion L (J/kg)
1262 self.latentheat = 3.34e5
1263- #ice thermal conductivity (W / m / K)
1264+ #ice thermal conductivity (W/m/K)
1265 self.thermalconductivity = 2.4
1266+ #temperate ice thermal conductivity (W/m/K)
1267+ self.temperateiceconductivity = 0.24
1268 #computation of effective conductivity
1269 self.effectiveconductivity_averaging = 1
1270- #temperate ice thermal conductivity (W / m / K)
1271- self.temperateiceconductivity = 0.24
1272 #the melting point of ice at 1 atmosphere of pressure in K
1273 self.meltingpoint = 273.15
1274- #rate of change of melting point with pressure (K / Pa)
1275+ #rate of change of melting point with pressure (K/Pa)
1276 self.beta = 9.8e-8
1277- #mixed layer (ice-water interface) heat capacity (J / kg / K)
1278+ #mixed layer (ice-water interface) heat capacity (J/kg/K)
1279 self.mixed_layer_capacity = 3974.
1280- #thermal exchange velocity (ice-water interface) (m / s)
1281+ #thermal exchange velocity (ice-water interface) (m/s)
1282 self.thermal_exchange_velocity = 1.00e-4
1283 #Rheology law: what is the temperature dependence of B with T
1284 #available: none, paterson and arrhenius
1285@@ -106,34 +108,39 @@
1286 self.rheology_law = 'Paterson'
1287
1288 # GIA:
1289- self.lithosphere_shear_modulus = 6.7e10 # (Pa)
1290- self.lithosphere_density = 3.32 # (g / cm^ - 3)
1291+ self.lithosphere_shear_modulus = 6.7e10 # (Pa)
1292+ self.lithosphere_density = 3.32 # (g/cm^-3)
1293 self.mantle_shear_modulus = 1.45e11 # (Pa)
1294- self.mantle_density = 3.34 # (g / cm^ - 3)
1295+ self.mantle_density = 3.34 # (g/cm^-3)
1296
1297- #SLR
1298- self.earth_density = 5512 # average density of the Earth, (kg / m^3)
1299- return self
1300+ # SLR
1301+ self.earth_density = 5512 # average density of the Earth, (kg/m^3)
1302 #}}}
1303
1304- def checkconsistency(self, md, solution, analyses): # {{{
1305- md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0)
1306- md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0)
1307- md = checkfield(md, 'fieldname', 'materials.rho_freshwater', '>', 0)
1308- md = checkfield(md, 'fieldname', 'materials.mu_water', '>', 0)
1309- md = checkfield(md, 'fieldname', 'materials.rheology_B', '>', 0, 'timeseries', 1, 'NaN', 1, 'Inf', 1)
1310- md = checkfield(md, 'fieldname', 'materials.rheology_n', '>', 0, 'size', [md.mesh.numberofelements])
1311- md = checkfield(md, 'fieldname', 'materials.rheology_law', 'values', ['None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', 'NyeH2O'])
1312- md = checkfield(md, 'fieldname', 'materials.effectiveconductivity_averaging', 'numel', [1], 'values', [0, 1, 2])
1313- md = checkfield(md, 'fieldname', 'materials.lithosphere_shear_modulus', '>', 0, 'numel', [1])
1314- md = checkfield(md, 'fieldname', 'materials.lithosphere_density', '>', 0, 'numel', [1])
1315- md = checkfield(md, 'fieldname', 'materials.mantle_shear_modulus', '>', 0, 'numel', [1])
1316- md = checkfield(md, 'fieldname', 'materials.mantle_density', '>', 0, 'numel', [1])
1317- md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', [1])
1318+ def checkconsistency(self, md, solution, analyses): #{{{
1319+ if solution != 'SealevelriseSolution':
1320+ md = checkfield(md, 'fieldname', 'materials.rho_ice', '>', 0)
1321+ md = checkfield(md, 'fieldname', 'materials.rho_water', '>', 0)
1322+ md = checkfield(md, 'fieldname', 'materials.rho_freshwater', '>', 0)
1323+ md = checkfield(md, 'fieldname', 'materials.mu_water', '>', 0)
1324+ md = checkfield(md, 'fieldname', 'materials.rheology_B', '>', 0, 'timeseries', 1, 'NaN', 1, 'Inf', 1)
1325+ md = checkfield(md, 'fieldname', 'materials.rheology_n', '>', 0, 'size', [md.mesh.numberofelements])
1326+ md = checkfield(md, 'fieldname', 'materials.rheology_law', 'values', ['None', 'BuddJacka', 'Cuffey', 'CuffeyTemperate', 'Paterson', 'Arrhenius', 'LliboutryDuval', 'NyeCO2', 'NyeH2O'])
1327+ md = checkfield(md, 'fieldname', 'materials.effectiveconductivity_averaging', 'numel', [1], 'values', [0, 1, 2])
1328+
1329+ if 'GiaAnalysis' in analyses:
1330+ md = checkfield(md, 'fieldname', 'materials.lithosphere_shear_modulus', '>', 0, 'numel', [1])
1331+ md = checkfield(md, 'fieldname', 'materials.lithosphere_density', '>', 0, 'numel', [1])
1332+ md = checkfield(md, 'fieldname', 'materials.mantle_shear_modulus', '>', 0, 'numel', [1])
1333+ md = checkfield(md, 'fieldname', 'materials.mantle_density', '>', 0, 'numel', [1])
1334+
1335+ if 'SealevelriseAnalysis' in analyses:
1336+ md = checkfield(md, 'fieldname', 'materials.earth_density', '>', 0, 'numel', [1])
1337+
1338 return md
1339- # }}}
1340+ #}}}
1341
1342- def marshall(self, prefix, md, fid): # {{{
1343+ def marshall(self, prefix, md, fid): #{{{
1344 WriteData(fid, prefix, 'name', 'md.materials.type', 'data', 3, 'format', 'Integer')
1345 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_ice', 'format', 'Double')
1346 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rho_water', 'format', 'Double')
1347@@ -150,12 +157,10 @@
1348 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'thermal_exchange_velocity', 'format', 'Double')
1349 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_B', 'format', 'DoubleMat', 'mattype', 1, 'timeserieslength', md.mesh.numberofvertices + 1, 'yts', md.constants.yts)
1350 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'rheology_n', 'format', 'DoubleMat', 'mattype', 2)
1351- WriteData(fid, prefix, 'data', self.rheology_law, 'name', 'md.materials.rheology_law', 'format', 'String')
1352-
1353+ WriteData(fid, prefix, 'data', self.rheology_law, 'name', 'md.materials.rheology_law', 'format', 's')
1354 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'lithosphere_shear_modulus', 'format', 'Double')
1355 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'lithosphere_density', 'format', 'Double', 'scale', 10.**3.)
1356 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'mantle_shear_modulus', 'format', 'Double')
1357 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'mantle_density', 'format', 'Double', 'scale', 10.**3.)
1358 WriteData(fid, prefix, 'object', self, 'class', 'materials', 'fieldname', 'earth_density', 'format', 'Double')
1359-
1360- # }}}
1361+ #}}}
1362Index: ../trunk-jpl/src/m/classes/slr.py
1363===================================================================
1364--- ../trunk-jpl/src/m/classes/slr.py (revision 25168)
1365+++ ../trunk-jpl/src/m/classes/slr.py (revision 25169)
1366@@ -49,34 +49,34 @@
1367 #}}}
1368
1369 def __repr__(self): # {{{
1370- string = ' slr parameters:'
1371- string = "%s\n%s" % (string, fielddisplay(self, 'deltathickness', 'thickness change: ice height equivalent [m]'))
1372- string = "%s\n%s" % (string, fielddisplay(self, 'sealevel', 'current sea level (prior to computation) [m]'))
1373- string = "%s\n%s" % (string, fielddisplay(self, 'spcthickness', 'thickness constraints (NaN means no constraint) [m]'))
1374- string = "%s\n%s" % (string, fielddisplay(self, 'reltol', 'sea level rise relative convergence criterion, (NaN: not applied)'))
1375- string = "%s\n%s" % (string, fielddisplay(self, 'abstol', 'sea level rise absolute convergence criterion, (default, NaN: not applied)'))
1376- string = "%s\n%s" % (string, fielddisplay(self, 'maxiter', 'maximum number of nonlinear iterations'))
1377- string = "%s\n%s" % (string, fielddisplay(self, 'love_h', 'load Love number for radial displacement'))
1378- string = "%s\n%s" % (string, fielddisplay(self, 'love_k', 'load Love number for gravitational potential perturbation'))
1379- string = "%s\n%s" % (string, fielddisplay(self, 'love_l', 'load Love number for horizontal displaements'))
1380- string = "%s\n%s" % (string, fielddisplay(self, 'tide_love_k', 'tidal load Love number (degree 2)'))
1381- string = "%s\n%s" % (string, fielddisplay(self, 'tide_love_h', 'tidal load Love number (degree 2)'))
1382- string = "%s\n%s" % (string, fielddisplay(self, 'fluid_love', 'secular fluid Love number'))
1383- string = "%s\n%s" % (string, fielddisplay(self, 'equatorial_moi', 'mean equatorial moment of inertia [kg m^2]'))
1384- string = "%s\n%s" % (string, fielddisplay(self, 'polar_moi', 'polar moment of inertia [kg m^2]'))
1385- string = "%s\n%s" % (string, fielddisplay(self, 'angular_velocity', 'mean rotational velocity of earth [per second]'))
1386- string = "%s\n%s" % (string, fielddisplay(self, 'ocean_area_scaling', 'correction for model representation of ocean area [default: No correction]'))
1387- string = "%s\n%s" % (string, fielddisplay(self, 'hydro_rate', 'rate of hydrological expansion [mm / yr]'))
1388- string = "%s\n%s" % (string, fielddisplay(self, 'geodetic', 'compute geodetic SLR? (in addition to steric?) default 0'))
1389- string = "%s\n%s" % (string, fielddisplay(self, 'geodetic_run_frequency', 'how many time steps we skip before we run SLR solver during transient (default: 1)'))
1390- string = "%s\n%s" % (string, fielddisplay(self, 'rigid', 'rigid earth graviational potential perturbation'))
1391- string = "%s\n%s" % (string, fielddisplay(self, 'elastic', 'elastic earth graviational potential perturbation'))
1392- string = "%s\n%s" % (string, fielddisplay(self, 'rotation', 'earth rotational potential perturbation'))
1393- string = "%s\n%s" % (string, fielddisplay(self, 'degacc', 'accuracy (default .01 deg) for numerical discretization of the Green''s functions'))
1394- string = "%s\n%s" % (string, fielddisplay(self, 'transitions', 'indices into parts of the mesh that will be icecaps'))
1395- string = "%s\n%s" % (string, fielddisplay(self, 'requested_outputs', 'additional outputs requested'))
1396+ s = ' slr parameters:'
1397+ s = "%s\n%s" % (s, fielddisplay(self, 'deltathickness', 'thickness change: ice height equivalent [m]'))
1398+ s = "%s\n%s" % (s, fielddisplay(self, 'sealevel', 'current sea level (prior to computation) [m]'))
1399+ s = "%s\n%s" % (s, fielddisplay(self, 'spcthickness', 'thickness constraints (NaN means no constraint) [m]'))
1400+ s = "%s\n%s" % (s, fielddisplay(self, 'reltol', 'sea level rise relative convergence criterion, (NaN: not applied)'))
1401+ s = "%s\n%s" % (s, fielddisplay(self, 'abstol', 'sea level rise absolute convergence criterion, (default, NaN: not applied)'))
1402+ s = "%s\n%s" % (s, fielddisplay(self, 'maxiter', 'maximum number of nonlinear iterations'))
1403+ s = "%s\n%s" % (s, fielddisplay(self, 'love_h', 'load Love number for radial displacement'))
1404+ s = "%s\n%s" % (s, fielddisplay(self, 'love_k', 'load Love number for gravitational potential perturbation'))
1405+ s = "%s\n%s" % (s, fielddisplay(self, 'love_l', 'load Love number for horizontal displaements'))
1406+ s = "%s\n%s" % (s, fielddisplay(self, 'tide_love_k', 'tidal load Love number (degree 2)'))
1407+ s = "%s\n%s" % (s, fielddisplay(self, 'tide_love_h', 'tidal load Love number (degree 2)'))
1408+ s = "%s\n%s" % (s, fielddisplay(self, 'fluid_love', 'secular fluid Love number'))
1409+ s = "%s\n%s" % (s, fielddisplay(self, 'equatorial_moi', 'mean equatorial moment of inertia [kg m^2]'))
1410+ s = "%s\n%s" % (s, fielddisplay(self, 'polar_moi', 'polar moment of inertia [kg m^2]'))
1411+ s = "%s\n%s" % (s, fielddisplay(self, 'angular_velocity', 'mean rotational velocity of earth [per second]'))
1412+ s = "%s\n%s" % (s, fielddisplay(self, 'ocean_area_scaling', 'correction for model representation of ocean area [default: No correction]'))
1413+ s = "%s\n%s" % (s, fielddisplay(self, 'hydro_rate', 'rate of hydrological expansion [mm / yr]'))
1414+ s = "%s\n%s" % (s, fielddisplay(self, 'geodetic', 'compute geodetic SLR? (in addition to steric?) default 0'))
1415+ s = "%s\n%s" % (s, fielddisplay(self, 'geodetic_run_frequency', 'how many time steps we skip before we run SLR solver during transient (default: 1)'))
1416+ s = "%s\n%s" % (s, fielddisplay(self, 'rigid', 'rigid earth graviational potential perturbation'))
1417+ s = "%s\n%s" % (s, fielddisplay(self, 'elastic', 'elastic earth graviational potential perturbation'))
1418+ s = "%s\n%s" % (s, fielddisplay(self, 'rotation', 'earth rotational potential perturbation'))
1419+ s = "%s\n%s" % (s, fielddisplay(self, 'degacc', 'accuracy (default .01 deg) for numerical discretization of the Green''s functions'))
1420+ s = "%s\n%s" % (s, fielddisplay(self, 'transitions', 'indices into parts of the mesh that will be icecaps'))
1421+ s = "%s\n%s" % (s, fielddisplay(self, 'requested_outputs', 'additional outputs requested'))
1422
1423- return string
1424+ return s
1425 # }}}
1426
1427 def setdefaultparameters(self): # {{{
1428@@ -142,7 +142,7 @@
1429 md = checkfield(md, 'fieldname', 'slr.geodetic_run_frequency', 'size', [1, 1], '>=', 1)
1430 md = checkfield(md, 'fieldname', 'slr.hydro_rate', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
1431 md = checkfield(md, 'fieldname', 'slr.degacc', 'size', [1, 1], '>=', 1e-10)
1432- md = checkfield(md, 'fieldname', 'slr.requested_outputs', 'stringrow', 1)
1433+ md = checkfield(md, 'fieldname', 'slr.requested_outputs', 'srow', 1)
1434 md = checkfield(md, 'fieldname', 'slr.horiz', 'NaN', 1, 'Inf', 1, 'values', [0, 1])
1435
1436 #check that love numbers are provided at the same level of accuracy:
1437@@ -201,5 +201,5 @@
1438 if len(indices) > 0:
1439 outputscopy = outputs[0:max(0, indices[0] - 1)] + self.defaultoutputs(md) + outputs[indices[0] + 1:]
1440 outputs = outputscopy
1441- WriteData(fid, prefix, 'data', outputs, 'name', 'md.slr.requested_outputs', 'format', 'StringArray')
1442+ WriteData(fid, prefix, 'data', outputs, 'name', 'md.slr.requested_outputs', 'format', 'sArray')
1443 # }}}
Note: See TracBrowser for help on using the repository browser.