source: issm/oecreview/Archive/25834-26739/ISSM-26012-26013.diff

Last change on this file was 26740, checked in by Mathieu Morlighem, 3 years ago

CHG: added 25834-26739
File size: 10.3 KB
RevLine 
[26740]1Index: ../trunk-jpl/src/m/classes/initialization.py
2===================================================================
3--- ../trunk-jpl/src/m/classes/initialization.py (revision 26012)
4+++ ../trunk-jpl/src/m/classes/initialization.py (revision 26013)
5@@ -28,6 +28,7 @@
6 self.epl_thickness = np.nan
7 self.hydraulic_potential = np.nan
8 self.channelarea = np.nan
9+ self.sample = np.nan
10
11 #set defaults
12 self.setdefaultparameters()
13@@ -49,6 +50,7 @@
14 s += '{}\n'.format(fielddisplay(self, 'epl_thickness', 'thickness of the epl [m]'))
15 s += '{}\n'.format(fielddisplay(self, 'hydraulic_potential', 'Hydraulic potential (for GlaDS) [Pa]'))
16 s += '{}\n'.format(fielddisplay(self, 'channelarea', 'subglaciale water channel area (for GlaDS) [m2]'))
17+ s += '{}\n'.format(fielddisplay(self,'sample','Realization of a Gaussian random field'))
18 return s
19 #}}}
20
21@@ -123,6 +125,9 @@
22 md = checkfield(md, 'fieldname', 'initialization.watercolumn', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
23 md = checkfield(md, 'fieldname', 'initialization.hydraulic_potential', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
24 md = checkfield(md, 'fieldname', 'initialization.channelarea', 'NaN', 1, 'Inf', 1, '>=', 0, 'size', [md.mesh.numberofelements])
25+ if 'SamplingAnalysis' in analyses and not solution == 'TransientSolution' and not md.transient.issampling:
26+ if np.any(np.isnan(md.initialization.sample)):
27+ md = checkfield(md,'fieldname','initialization.sample','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1])
28 return md
29 # }}}
30
31@@ -140,6 +145,7 @@
32 WriteData(fid, prefix, 'object', self, 'fieldname', 'watercolumn', 'format', 'DoubleMat', 'mattype', 1)
33 WriteData(fid, prefix, 'object', self, 'fieldname', 'channelarea', 'format', 'DoubleMat', 'mattype', 1)
34 WriteData(fid, prefix, 'object', self, 'fieldname', 'hydraulic_potential', 'format', 'DoubleMat', 'mattype', 1)
35+ WriteData(fid,prefix,'object',self,'fieldname','sample','format','DoubleMat','mattype',1)
36 if md.thermal.isenthalpy:
37 if (np.size(self.enthalpy) <= 1):
38 tpmp = md.materials.meltingpoint - md.materials.beta * md.initialization.pressure
39Index: ../trunk-jpl/src/m/classes/model.py
40===================================================================
41--- ../trunk-jpl/src/m/classes/model.py (revision 26012)
42+++ ../trunk-jpl/src/m/classes/model.py (revision 26013)
43@@ -75,6 +75,7 @@
44 from ElementConnectivity import ElementConnectivity
45 from contourenvelope import contourenvelope
46 from DepthAverage import DepthAverage
47+from sampling import sampling
48 #}}}
49
50
51@@ -115,6 +116,7 @@
52 self.gia = None
53 self.love = None
54 self.esa = None
55+ self.sampling = None
56 self.autodiff = None
57 self.inversion = None
58 self.qmu = None
59@@ -172,6 +174,7 @@
60 'gia',
61 'love',
62 'esa',
63+ 'sampling',
64 'autodiff',
65 'inversion',
66 'qmu',
67@@ -223,6 +226,7 @@
68 s = "%s\n%s" % (s, "%19s: %-22s -- %s" % ("frontalforcings", "[%s %s]" % ("1x1", obj.frontalforcings.__class__.__name__), "parameters for frontalforcings"))
69 s = "%s\n%s" % (s, "%19s: %-22s -- %s" % ("gia", "[%s %s]" % ("1x1", obj.gia.__class__.__name__), "parameters for gia solution"))
70 s = "%s\n%s" % (s, '%19s: %-22s -- %s' % ("esa", "[%s %s]" % ("1x1", obj.esa.__class__.__name__), "parameters for elastic adjustment solution"))
71+ s = "%s\n%s" % (s, '%19s: %-22s -- %s' % ("sampling", "[%s %s]" % ("1x1", obj.sampling.__class__.__name__), "parameters for stochastic sampler"))
72 s = "%s\n%s" % (s, '%19s: %-22s -- %s' % ("love", "[%s %s]" % ("1x1", obj.love.__class__.__name__), "parameters for love solution"))
73 s = "%s\n%s" % (s, "%19s: %-22s -- %s" % ("autodiff", "[%s %s]" % ("1x1", obj.autodiff.__class__.__name__), "automatic differentiation parameters"))
74 s = "%s\n%s" % (s, "%19s: %-22s -- %s" % ("inversion", "[%s %s]" % ("1x1", obj.inversion.__class__.__name__), "parameters for inverse methods"))
75@@ -270,6 +274,7 @@
76 self.gia = giamme()
77 self.love = fourierlove()
78 self.esa = esa()
79+ self.sampling = sampling()
80 self.autodiff = autodiff()
81 self.inversion = inversion()
82 self.qmu = qmu()
83Index: ../trunk-jpl/src/m/consistency/ismodelselfconsistent.py
84===================================================================
85--- ../trunk-jpl/src/m/consistency/ismodelselfconsistent.py (revision 26012)
86+++ ../trunk-jpl/src/m/consistency/ismodelselfconsistent.py (revision 26013)
87@@ -78,6 +78,8 @@
88 analyses = ['L2ProjectionBaseAnalysis', 'HydrologyShreveAnalysis', 'HydrologyDCInefficientAnalysis', 'HydrologyDCEfficientAnalysis']
89 elif 'DamageEvolutionSolution':
90 analyses = ['DamageEvolutionAnalysis']
91+ elseif strcmp(solutiontype,'SamplingSolution')
92+ analyses=['SamplingAnalysis']
93 else:
94 raise TypeError('solution type: {} not supported yet!'.format(solutiontype))
95
96Index: ../trunk-jpl/src/m/sampling.py
97===================================================================
98--- ../trunk-jpl/src/m/sampling.py (revision 26012)
99+++ ../trunk-jpl/src/m/sampling.py (revision 26013)
100@@ -54,9 +54,9 @@
101 # Temporal correlation factor
102 self.phi = 0
103 # Exponent in fraction SPDE (default=2, biLaplacian covariance operator)
104- self.alpha=2;
105+ self.alpha = 2
106 # Seed for pseudorandom number generator (default -1 for random seed)
107- self.alpha=-1;
108+ self.alpha = -1
109 # Default output
110 self.requested_outputs = ['default']
111 return self
112@@ -67,26 +67,26 @@
113 if (SamplingAnalysis' not in analyses):
114 return md
115
116- md = checkfield(md,'fieldname','sampling.kappa','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1],'>',0);
117- md = checkfield(md,'fieldname','sampling.tau','NaN',1,'Inf',1,'numel',1,'>',0);
118- md = checkfield(md,'fieldname','sampling.beta','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1],'>',0);
119- md = checkfield(md,'fieldname','sampling.phi','NaN',1,'Inf',1,'numel',1,'>=',0);
120- md = checkfield(md,'fieldname','sampling.alpha','NaN',1,'Inf',1,'numel',1,'>',0);
121- md = checkfield(md,'fieldname','sampling.robin','numel',1,'values',[0 1]);
122- md = checkfield(md,'fieldname','sampling.seed','NaN',1,'Inf',1,'numel',1);
123- md = checkfield(md,'fieldname','sampling.requested_outputs','stringrow',1);
124+ md = checkfield(md,'fieldname','sampling.kappa','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1],'>',0)
125+ md = checkfield(md,'fieldname','sampling.tau','NaN',1,'Inf',1,'numel',1,'>',0)
126+ md = checkfield(md,'fieldname','sampling.beta','NaN',1,'Inf',1,'size',[md.mesh.numberofvertices 1],'>',0)
127+ md = checkfield(md,'fieldname','sampling.phi','NaN',1,'Inf',1,'numel',1,'>=',0)
128+ md = checkfield(md,'fieldname','sampling.alpha','NaN',1,'Inf',1,'numel',1,'>',0)
129+ md = checkfield(md,'fieldname','sampling.robin','numel',1,'values',[0 1])
130+ md = checkfield(md,'fieldname','sampling.seed','NaN',1,'Inf',1,'numel',1)
131+ md = checkfield(md,'fieldname','sampling.requested_outputs','stringrow',1)
132
133 return md
134 # }}}
135
136 def marshall(self, prefix, md, fid): # {{{
137- WriteData(fid,prefix,'object',self,'fieldname','kappa','format','DoubleMat','mattype',1);
138- WriteData(fid,prefix,'object',self,'fieldname','tau','format','Double');
139- WriteData(fid,prefix,'object',self,'fieldname','beta','format','DoubleMat','mattype',1);
140- WriteData(fid,prefix,'object',self,'fieldname','phi','format','Double');
141- WriteData(fid,prefix,'object',self,'fieldname','alpha','format','Integer');
142- WriteData(fid,prefix,'object',self,'fieldname','robin','format','Boolean');
143- WriteData(fid,prefix,'object',self,'fieldname','seed','format','Integer');
144+ WriteData(fid,prefix,'object',self,'fieldname','kappa','format','DoubleMat','mattype',1)
145+ WriteData(fid,prefix,'object',self,'fieldname','tau','format','Double')
146+ WriteData(fid,prefix,'object',self,'fieldname','beta','format','DoubleMat','mattype',1)
147+ WriteData(fid,prefix,'object',self,'fieldname','phi','format','Double')
148+ WriteData(fid,prefix,'object',self,'fieldname','alpha','format','Integer')
149+ WriteData(fid,prefix,'object',self,'fieldname','robin','format','Boolean')
150+ WriteData(fid,prefix,'object',self,'fieldname','seed','format','Integer')
151
152 # Process requested outputs
153 outputs = self.requested_outputs
154@@ -97,10 +97,10 @@
155 WriteData(fid, prefix, 'data', outputs, 'name', 'md.masstransport.requested_outputs', 'format', 'StringArray')
156
157 # Process requested outputs
158- outputs = self.requested_outputs;
159+ outputs = self.requested_outputs
160 indices = [i for i, x in enumerate(outputs) if x == 'default']
161 if len(indices) > 0:
162 outputscopy = outputs[0:max(0, indices[0] - 1)] + self.defaultoutputs(md) + outputs[indices[0] + 1:]
163 outputs = outputscopy
164- WriteData(fid,prefix,'data',outputs,'name','md.sampling.requested_outputs','format','StringArray');
165+ WriteData(fid,prefix,'data',outputs,'name','md.sampling.requested_outputs','format','StringArray')
166 # }}}
167\ No newline at end of file
168Index: ../trunk-jpl/src/m/solve/solve.py
169===================================================================
170--- ../trunk-jpl/src/m/solve/solve.py (revision 26012)
171+++ ../trunk-jpl/src/m/solve/solve.py (revision 26013)
172@@ -34,6 +34,7 @@
173 - 'Esa' or 'esa'
174 - 'Sealevelchange' or 'slc'
175 - 'Love' or 'lv'
176+ - 'Sampling' or 'smp'
177
178 Extra options:
179 - loadonly : does not solve. only load results
180@@ -79,6 +80,8 @@
181 solutionstring = 'EsaSolution'
182 elif solutionstring.lower() == 'slc' or solutionstring.lower() == 'sealevelchange':
183 solutionstring = 'SealevelchangeSolution'
184+ elif solutionstring.lower() == 'smp' or solutionstring.lower() == 'sampling':
185+ solutionstring = 'SamplingSolution'
186 else:
187 raise ValueError("solutionstring '%s' not supported!" % solutionstring)
188 options = pairoptions('solutionstring', solutionstring, *args)
Note: See TracBrowser for help on using the repository browser.