tmodel.py - pism - [fork] customized build of PISM, the parallel ice sheet model (tillflux branch)
(HTM) git clone git://src.adamsgaard.dk/pism
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) LICENSE
---
tmodel.py (40185B)
---
1 # Copyright (C) 2011, 2014, 2015, 2016, 2017, 2018, 2019 David Maxwell and Constantine Khroulev
2 #
3 # This file is part of PISM.
4 #
5 # PISM is free software; you can redistribute it and/or modify it under the
6 # terms of the GNU General Public License as published by the Free Software
7 # Foundation; either version 3 of the License, or (at your option) any later
8 # version.
9 #
10 # PISM is distributed in the hope that it will be useful, but WITHOUT ANY
11 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 # FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
13 # details.
14 #
15 # You should have received a copy of the GNU General Public License
16 # along with PISM; if not, write to the Free Software
17 # Foundation, Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
18
19 """Module containing classes and functions for managing computational
20 grids, model physics, and creating standard vectors."""
21
22 import PISM
23 from PISM.util import convert
24
25 class ModelVecs(object):
26
27 """A dictionary of :cpp:class:`IceModelVec`'s containing model data.
28 This plays a similar role as the C++ :cpp:class:`Vars` with the key
29 difference that it keeps track of fields that may need to be written
30 to a file.::
31
32 # Construct a tau_c vector with intrinsic name 'tauc'
33 yieldstress = PISM.IceModelVec2S(grid, 'tauc', PISM.WITHOUT_GHOSTS)
34 yieldstress.set_attrs("diagnostic", desc, "Pa", "Pa", "", 0)
35
36 # Construct a thickness vector with intrinsic name 'thk'
37 thickness = PISM.IceModelVec2S(grid, 'thk', PISM.WITHOUT_GHOSTS)
38 thickness.set_attrs("model_state", "land ice thickness", "m", "m", "land_ice_thickness", 0)
39 thickness.metadata().set_number("valid_min", 0.0)
40
41 vecs = PISM.model.ModelVecs()
42
43 # Add 'yieldstress' to be accessed via the given name 'basal_yield_stress'
44 vecs.add(yieldstress, 'basal_yield_stress')
45
46 # Add 'thickness' to be accessed by its intrinsic name, (i.e. 'thk')
47 vecs.add(thickness)
48
49 # Access the vectors
50 v1 = vecs.basal_yield_stress
51 v2 = vecs.thk
52
53 See also functions :func:`createYieldStressVec` and friends for one
54 step creation of standard PISM vecs.
55
56 """
57
58 def __init__(self, variables):
59 self.needs_writing = set()
60 self._vecs = variables
61
62 def __getattr__(self, name):
63 """Allows access to contents via dot notation, e.g. `vecs.thickness`."""
64 try:
65 return self.get(name)
66 except RuntimeError as e:
67 print("#### ModelVecs.__getattr__ got a PISM-level exception:", e)
68 print("# The contents of this ModelVecs instance:")
69 print(self)
70 print("#### End of the contents")
71 raise AttributeError(name)
72
73 def __repr__(self):
74 """:returns: a string representation of the :class:`ModelVecs` listing its contents."""
75 s = 'ModelVecs:\n'
76 items = [(key, self._vecs.get(key)) for key in list(self._vecs.keys())]
77 for (key, value) in items:
78 short_name = value.metadata().get_string("short_name")
79 long_name = value.metadata().get_string("long_name")
80 standard_name = value.metadata().get_string("standard_name")
81
82 if key != short_name or len(standard_name) == 0:
83 # alternative name was given, so we can't look up by standard name
84 s += " %s: %s\n" % (key, long_name)
85 else:
86 s += " %s or %s: %s\n" % (short_name, standard_name, long_name)
87 return s
88
89 def __str__(self):
90 return self.__repr__()
91
92 def get(self, name):
93 """:returns: the vector with the given *name*, raising an :exc:AttributeError if one is not found."""
94 return self._vecs.get(name)
95
96 def add(self, var, name=None, writing=False):
97 """Adds a vector to the collection.
98
99 :param var: The :cpp:class:`IceModelVec` to add.
100 :param name: The name to associate with `var`. If `None`, then
101 the `var`'s `name` attribute is used.
102 :param writing: A boolean. `True` indicates the vector is among
103 those that should be written to a file if
104 :meth:`write` is called.
105
106 """
107
108 if name is None:
109 self._vecs.add(var)
110 else:
111 self._vecs.add(var, name)
112
113 if writing:
114 self.needs_writing.add(var)
115
116 def has(self, name):
117 """:returns: `True` if the dictionary contains a vector with the given name."""
118 return self._vecs.is_available(name)
119
120 def write(self, output_filename):
121 """Writes any member vectors that have been flagged as need writing to the given :file:`.nc` filename."""
122 vlist = [v for v in self.needs_writing]
123 vlist.sort(key=lambda v: v.get_name())
124 for v in vlist:
125 v.write(output_filename)
126
127 def writeall(self, output_filename):
128 """Writes all member vectors to the given :file:`.nc` filename."""
129 vlist = [self._vecs.get(var) for var in list(self._vecs.keys())]
130 vlist.sort(key=lambda v: v.get_name())
131
132 try:
133 com = vlist[0].get_grid().ctx().com()
134 except:
135 com = PISM.PETSc.COMM_WORLD
136
137 f = PISM.File(com, output_filename, PISM.PISM_NETCDF3, PISM.PISM_READWRITE)
138 for v in vlist:
139 v.define(f, PISM.PISM_DOUBLE)
140 for v in vlist:
141 v.write(f)
142 f.close()
143
144 def markForWriting(self, var):
145 """Marks a given vector as needing writing.
146
147 :param `var`: Either the name of a member vector, or the vector itself."""
148
149 if isinstance(var, str):
150 var = self.get(var)
151 self.needs_writing.add(var)
152
153
154 class ModelData(object):
155
156 """A collection of data and physics comprising a PISM model. The
157 collection consists of
158
159 * ``grid`` the computation grid
160 * ``config`` the associated :cpp:class:`Config`
161 * ``basal`` a :cpp:class:`YieldStress`.
162 * ``enthalpy`` an :cpp:class:`EnthalpyConverter`.
163 * ``vecs`` a :class:`ModelVecs` containing vector data.
164
165 These member variables are intended to be accessed as public attributes.
166
167 """
168
169 def __init__(self, grid, config=None):
170 self.grid = grid
171 if config is None:
172 config = grid.ctx().config()
173 self.config = config
174
175 # Components for the physics
176 self.basal = None #: Attribute for a :cpp:class:`YieldStress`
177 self.enthalpyconverter = None #: Attribute for an :cpp:class:`EnthalpyConverter`
178
179 self.vecs = ModelVecs(grid.variables()) #: The list
180
181 def setPhysics(self, enthalpyconverter):
182 """Sets the physics components of the model.
183
184 :param enthalpyconverter: An instance of :cpp:class:`EnthalpyConverter`
185 """
186
187 self.enthalpyconverter = enthalpyconverter
188
189
190 def initGrid(ctx, Lx, Ly, Lz, Mx, My, Mz, r):
191 """Initialize a :cpp:class:`IceGrid` intended for 3-d computations.
192
193 :param ctx: The execution context.
194 :param Lx: Half width of the ice model grid in x-direction (m)
195 :param Ly: Half width of the ice model grid in y-direction (m)
196 :param Lz: Extent of the grid in the vertical direction (m)
197 :param Mx: number of grid points in the x-direction
198 :param My: number of grid points in the y-direction
199 :param Mz: number of grid points in the z-direction
200 :param r: grid registration (one of ``PISM.CELL_CENTER``, ``PISM.CELL_CORNER``)
201 """
202 P = PISM.GridParameters(ctx.config)
203
204 P.Lx = Lx
205 P.Ly = Ly
206 P.Mx = Mx
207 P.My = My
208 P.registration = r
209 z = PISM.IceGrid.compute_vertical_levels(Lz, Mz, PISM.EQUAL)
210 P.z = PISM.DoubleVector(z)
211 P.horizontal_extent_from_options()
212 P.ownership_ranges_from_options(ctx.size)
213
214 return PISM.IceGrid(ctx.ctx, P)
215
216 def _stencil_width(grid, ghost_type, width):
217 "Compute stencil width for a field."
218 if ghost_type == PISM.WITH_GHOSTS:
219 if width is None:
220 return int(grid.ctx().config().get_number("grid.max_stencil_width"))
221 else:
222 return width
223 else:
224 return 0
225
226
227 def createIceSurfaceVec(grid, name='usurf', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
228 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
229 surface elevation data.
230
231 :param grid: The grid to associate with the vector.
232 :param name: The intrinsic name to give the vector.
233
234 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
235 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
236 is ghosted.
237
238 :param stencil_width: The size of the ghost stencil. Ignored if
239 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
240 ``stencil_width`` is ``None`` and the vector
241 is ghosted, the grid's maximum stencil width
242 is used.
243
244 """
245 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
246
247 surface = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
248 surface.set_attrs("diagnostic", "ice upper surface elevation", "m", "m", "surface_altitude", 0)
249 return surface
250
251 def createSeaLevelVec(grid, name='sea_level', ghost_type=PISM.WITH_GHOSTS, stencil_width=1):
252 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
253 sea level data.
254
255 :param grid: The grid to associate with the vector.
256 :param name: The intrinsic name to give the vector.
257
258 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
259 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
260 is ghosted.
261
262 :param stencil_width: The size of the ghost stencil. Ignored if
263 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
264 ``stencil_width`` is ``None`` and the vector
265 is ghosted, the grid's maximum stencil width
266 is used.
267
268 """
269 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
270
271 sea_level = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
272 sea_level.set_attrs("diagnostic",
273 "sea level elevation above reference ellipsoid",
274 "m", "m",
275 "sea_surface_height_above_reference_ellipsoid", 0)
276 return sea_level
277
278
279 def createIceThicknessVec(grid, name='thk', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
280 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for ice thickness data.
281
282 :param grid: The grid to associate with the vector.
283 :param name: The intrinsic name to give the vector.
284
285 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
286 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
287 is ghosted.
288
289 :param stencil_width: The size of the ghost stencil. Ignored if
290 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
291 ``stencil_width`` is ``None`` and the vector
292 is ghosted, the grid's maximum stencil width
293 is used.
294
295 """
296 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
297
298 thickness = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
299 thickness.set_attrs("model_state", "land ice thickness", "m", "m", "land_ice_thickness", 0)
300 thickness.metadata().set_number("valid_min", 0.0)
301 return thickness
302
303
304 def createIceSurfaceStoreVec(grid, name='usurfstore', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
305 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
306 saved surface elevation data in regional models.
307
308 :param grid: The grid to associate with the vector.
309 :param name: The intrinsic name to give the vector.
310
311 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
312 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
313 is ghosted.
314
315 :param stencil_width: The size of the ghost stencil. Ignored if
316 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
317 ``stencil_width`` is ``None`` and the vector
318 is ghosted, the grid's maximum stencil width
319 is used.
320
321 """
322 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
323
324 usurfstore = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
325 usurfstore.set_attrs("model_state",
326 "saved surface elevation for use to keep surface gradient constant in no_model strip",
327 "m", "m", "", 0)
328 return usurfstore
329
330
331 def createIceThicknessStoreVec(grid, name='thkstore', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
332 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
333 saved ice thickness data in regional models.
334
335 :param grid: The grid to associate with the vector.
336 :param name: The intrinsic name to give the vector.
337
338 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
339 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
340 is ghosted.
341
342 :param stencil_width: The size of the ghost stencil. Ignored if
343 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
344 ``stencil_width`` is ``None`` and the vector
345 is ghosted, the grid's maximum stencil width
346 is used.
347
348 """
349 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
350
351 thkstore = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
352 thkstore.set_attrs("model_state",
353 "saved ice thickness for use to keep driving stress constant in no_model strip",
354 "m", "m", "", 0)
355 return thkstore
356
357
358 def createBedrockElevationVec(grid, name='topg', desc="bedrock surface elevation",
359 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
360 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for bedrock elevation data.
361
362 :param grid: The grid to associate with the vector.
363 :param name: The intrinsic name to give the vector.
364
365 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
366 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
367 is ghosted.
368
369 :param stencil_width: The size of the ghost stencil. Ignored if
370 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
371 ``stencil_width`` is ``None`` and the vector
372 is ghosted, the grid's maximum stencil width
373 is used.
374
375 """
376
377 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
378
379 bed = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
380 bed.set_attrs("model_state", desc, "m", "m", "bedrock_altitude", 0)
381 return bed
382
383
384 def createYieldStressVec(grid, name='tauc', desc="yield stress for basal till (plastic or pseudo-plastic model)",
385 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
386 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for basal yield stress data.
387
388 :param grid: The grid to associate with the vector.
389 :param name: The intrinsic name to give the vector.
390
391 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
392 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
393 is ghosted.
394
395 :param stencil_width: The size of the ghost stencil. Ignored if
396 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
397 ``stencil_width`` is ``None`` and the vector
398 is ghosted, the grid's maximum stencil width
399 is used.
400
401 """
402 # yield stress for basal till (plastic or pseudo-plastic model)
403 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
404
405 tauc = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
406 tauc.set_attrs("diagnostic", desc, "Pa", "Pa", "", 0)
407 return tauc
408
409
410 def createAveragedHardnessVec(grid, name='hardav', desc="vertical average of ice hardness",
411 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
412 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
413 vertically-averaged ice-hardness data.
414
415 :param grid: The grid to associate with the vector.
416 :param name: The intrinsic name to give the vector.
417
418 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
419 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
420 is ghosted.
421
422 :param stencil_width: The size of the ghost stencil. Ignored if
423 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
424 ``stencil_width`` is ``None`` and the vector
425 is ghosted, the grid's maximum stencil width
426 is used.
427
428 """
429 # yield stress for basal till (plastic or pseudo-plastic model)
430 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
431
432 hardav = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
433
434 power = 1.0 / grid.ctx().config().get_number("stress_balance.ssa.Glen_exponent")
435 units = "Pa s%f" % power
436
437 hardav.set_attrs("diagnostic", desc, units, units, "", 0)
438
439 hardav.metadata().set_number("valid_min", 0.0)
440
441 return hardav
442
443
444 def createEnthalpyVec(grid, name='enthalpy', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
445 """Returns an :cpp:class:`IceModelVec3` with attributes setup for enthalpy data.
446
447 :param grid: The grid to associate with the vector.
448 :param name: The intrinsic name to give the vector.
449
450 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
451 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
452 is ghosted.
453
454 :param stencil_width: The size of the ghost stencil. Ignored if
455 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
456 ``stencil_width`` is ``None`` and the vector
457 is ghosted, the grid's maximum stencil width
458 is used.
459
460 """
461
462 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
463 enthalpy = PISM.IceModelVec3(grid, name, ghost_type, stencil_width)
464
465 enthalpy.set_attrs("model_state",
466 "ice enthalpy (includes sensible heat, latent heat, pressure)",
467 "J kg-1", "J kg-1", "", 0)
468 return enthalpy
469
470
471 def createStrainHeatingVec(grid, name='strain_heating',
472 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
473 """Returns an :cpp:class:`IceModelVec3` with attributes setup for
474 strain heating.
475
476 """
477
478 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
479 result = PISM.IceModelVec3(grid, name, ghost_type, stencil_width)
480
481 result.set_attrs("internal",
482 "rate of strain heating in ice (dissipation heating)",
483 "W m-3", "W m-3", "", 0)
484 return result
485
486
487 def createAgeVec(grid, name='age', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
488 """Returns an :cpp:class:`IceModelVec3` with attributes setup for age data.
489
490 :param grid: The grid to associate with the vector.
491 :param name: The intrinsic name to give the vector.
492
493 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
494 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
495 is ghosted.
496
497 :param stencil_width: The size of the ghost stencil. Ignored if
498 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
499 ``stencil_width`` is ``None`` and the vector
500 is ghosted, the grid's maximum stencil width
501 is used.
502
503 """
504
505 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
506 age = PISM.IceModelVec3(grid, name, ghost_type, stencil_width)
507
508 age.set_attrs("model_state", "age of the ice", "s", "s", "", 0)
509 return age
510
511
512 def create3DVelocityVecs(grid, ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
513 """Returns a size 3 tuple with instances of :cpp:class:`IceModelVec3`
514 with attributes setup for the x, y, and z-components of the 3D
515 velocity.
516
517 """
518 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
519
520 u = PISM.IceModelVec3(grid, "uvel", ghost_type, stencil_width)
521 u.set_attrs("diagnostic", "x-component of the ice velocity", "m s-1", "m year-1", "", 0)
522
523 v = PISM.IceModelVec3(grid, "vvel", ghost_type, stencil_width)
524 v.set_attrs("diagnostic", "y-component of the ice velocity", "m s-1", "m year-1", "", 0)
525
526 w = PISM.IceModelVec3(grid, "wvel_rel", ghost_type, stencil_width)
527 w.set_attrs("diagnostic",
528 "z-component of the ice velocity relative to the base directly below",
529 "m s-1", "m year-1", "", 0)
530
531 return u, v, w
532
533
534 def createBasalMeltRateVec(grid, name='bmelt', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
535 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for basal melt rate data.
536
537 :param grid: The grid to associate with the vector.
538 :param name: The intrinsic name to give the vector.
539
540 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
541 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
542 is ghosted.
543
544 :param stencil_width: The size of the ghost stencil. Ignored if
545 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
546 ``stencil_width`` is ``None`` and the vector
547 is ghosted, the grid's maximum stencil width
548 is used.
549
550 """
551 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
552
553 bmr = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
554 bmr.set_attrs("model_state",
555 "ice basal melt rate in ice thickness per time",
556 "m s-1", "m year-1", "land_ice_basal_melt_rate", 0)
557
558 bmr.metadata().set_string("comment", "positive basal melt rate corresponds to ice loss")
559 return bmr
560
561
562 def createTillPhiVec(grid, name='tillphi', desc="friction angle for till under grounded ice sheet",
563 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
564 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
565 till friction angle data.
566
567 :param grid: The grid to associate with the vector.
568 :param name: The intrinsic name to give the vector.
569
570 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
571 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
572 is ghosted.
573
574 :param stencil_width: The size of the ghost stencil. Ignored if
575 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
576 ``stencil_width`` is ``None`` and the vector
577 is ghosted, the grid's maximum stencil width
578 is used.
579
580 """
581 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
582
583 tillphi = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
584 # // ghosted to allow the "redundant" computation of tauc
585 # // PROPOSED standard_name = land_ice_basal_material_friction_angle
586 tillphi.set_attrs("climate_steady", desc,
587 "degrees", "degrees", "", 0)
588 tillphi.time_independent = True
589 return tillphi
590
591
592 def createBasalWaterVec(grid, name='tillwat', desc="effective thickness of subglacial melt water",
593 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
594 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for basal melt water thickness data.
595
596 :param grid: The grid to associate with the vector.
597 :param name: The intrinsic name to give the vector.
598
599 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
600 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
601 is ghosted.
602
603 :param stencil_width: The size of the ghost stencil. Ignored if
604 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
605 ``stencil_width`` is ``None`` and the vector
606 is ghosted, the grid's maximum stencil width
607 is used.
608
609 """
610
611 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
612
613 tillwat = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
614 tillwat.set_attrs("model_state", desc, "m", "m", "", 0)
615
616 valid_max = PISM.Context().config.get_number("hydrology.tillwat_max")
617 tillwat.metadata().set_number("valid_min", 0.0)
618 tillwat.metadata().set_number("valid_max", valid_max)
619 return tillwat
620
621
622 def createGroundingLineMask(grid, name='gl_mask', desc="fractional floatation mask",
623 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
624 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for the fractional floatation mask.
625
626 :param grid: The grid to associate with the vector.
627 :param name: The intrinsic name to give the vector.
628
629 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
630 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
631 is ghosted.
632
633 :param stencil_width: The size of the ghost stencil. Ignored if
634 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
635 ``stencil_width`` is ``None`` and the vector
636 is ghosted, the grid's maximum stencil width
637 is used.
638
639 """
640
641 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
642
643 gl_mask = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
644 gl_mask.set_attrs("model_state", desc, "1", "1", "", 0)
645 gl_mask.metadata().set_number("valid_min", 0.0)
646 gl_mask.metadata().set_number("valid_max", 1.0)
647 return gl_mask
648
649
650 def create2dVelocityVec(grid, name="", desc="", intent="", ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
651 """Returns an :cpp:class:`IceModelVec2V` with attributes setup for horizontal velocity data.
652
653 :param grid: The grid to associate with the vector.
654 :param name: The intrinsic name to give the vector.
655
656 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
657 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
658 is ghosted.
659
660 :param stencil_width: The size of the ghost stencil. Ignored if
661 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
662 ``stencil_width`` is ``None`` and the vector
663 is ghosted, the grid's maximum stencil width
664 is used.
665
666 """
667 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
668
669 vel = PISM.IceModelVec2V(grid, name, ghost_type, stencil_width)
670 vel.set_attrs(intent, "%s%s" % ("X-component of the ", desc), "m s-1", "m year-1", "", 0)
671 vel.set_attrs(intent, "%s%s" % ("Y-component of the ", desc), "m s-1", "m year-1", "", 1)
672
673 huge_vel = convert(1e10, "m/year", "m/second")
674 attrs = [("valid_min", -huge_vel), ("valid_max", huge_vel), ("_FillValue", 2 * huge_vel)]
675 for a in attrs:
676 for component in [0, 1]:
677 vel.metadata(component).set_number(a[0], a[1])
678 vel.set(2 * huge_vel)
679 return vel
680
681
682 def createDrivingStressXVec(grid, name="ssa_driving_stress_x", desc="driving stress", intent="model_state",
683 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
684 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
685 driving stresses int the X-direction.
686
687 :param grid: The grid to associate with the vector.
688 :param name: The intrinsic name to give the vector.
689
690 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
691 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
692 is ghosted.
693
694 :param stencil_width: The size of the ghost stencil. Ignored if
695 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
696 ``stencil_width`` is ``None`` and the vector
697 is ghosted, the grid's maximum stencil width
698 is used.
699
700 """
701
702 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
703
704 stress = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
705 stress.set_attrs(intent, "%s%s" % ("X-component of the ", desc), "Pa", "Pa", "", 0)
706 return stress
707
708
709 def createDrivingStressYVec(grid, name="ssa_driving_stress_y", desc="driving stress", intent="model_state",
710 ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
711 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
712 driving stresses int the Y-direction.
713
714 :param grid: The grid to associate with the vector.
715 :param name: The intrinsic name to give the vector.
716
717 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
718 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
719 is ghosted.
720
721 :param stencil_width: The size of the ghost stencil. Ignored if
722 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
723 ``stencil_width`` is ``None`` and the vector
724 is ghosted, the grid's maximum stencil width
725 is used.
726
727 """
728
729 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
730
731 stress = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
732 stress.set_attrs(intent, "%s%s" % ("X-component of the ", desc), "Pa", "Pa", "", 0)
733 return stress
734
735
736 def createVelocityMisfitWeightVec(grid, name="vel_misfit_weight", ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
737 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
738 weights for inversion velocity misfit functionals.
739
740 :param grid: The grid to associate with the vector.
741 :param name: The intrinsic name to give the vector.
742
743 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
744 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
745 is ghosted.
746
747 :param stencil_width: The size of the ghost stencil. Ignored if
748 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
749 ``stencil_width`` is ``None`` and the vector
750 is ghosted, the grid's maximum stencil width
751 is used.
752
753 """
754 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
755
756 vel_misfit_weight = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
757 vel_misfit_weight.set_attrs("diagnostic",
758 "weight for surface velocity misfit functional", "", "", "", 0)
759 return vel_misfit_weight
760
761
762 def createCBarVec(grid, name="velbar_mag", ghost_type=PISM.WITHOUT_GHOSTS, stencil_width=None):
763 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for horizontal velocity magnitudes.
764
765 :param grid: The grid to associate with the vector.
766 :param name: The intrinsic name to give the vector.
767
768 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
769 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
770 is ghosted.
771
772 :param stencil_width: The size of the ghost stencil. Ignored if
773 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
774 ``stencil_width`` is ``None`` and the vector
775 is ghosted, the grid's maximum stencil width
776 is used.
777
778 """
779
780 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
781
782 velbar_mag = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
783
784 velbar_mag.set_attrs("diagnostic",
785 "magnitude of vertically-integrated horizontal velocity of ice",
786 "m s-1", "m year-1", "", 0)
787
788 velbar_mag.metadata().set_number("valid_min", 0.0)
789
790 return velbar_mag
791
792
793 def createIceMaskVec(grid, name='mask', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
794 """Returns an :cpp:class:`IceModelVec2CellType` with attributes setup for PISM's mask determining ice
795 mode at grid points (grounded vs. floating and icy vs. ice-free).
796
797 :param grid: The grid to associate with the vector.
798 :param name: The intrinsic name to give the vector.
799
800 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
801 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
802 is ghosted.
803
804 :param stencil_width: The size of the ghost stencil. Ignored if
805 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
806 ``stencil_width`` is ``None`` and the vector
807 is ghosted, the grid's maximum stencil width
808 is used.
809
810 """
811 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
812
813 ice_mask = PISM.IceModelVec2CellType()
814 ice_mask.create(grid, name, ghost_type, stencil_width)
815
816 ice_mask.set_attrs("model_state", "grounded_dragging_floating integer mask", "", "", "", 0)
817 mask_values = [PISM.MASK_ICE_FREE_BEDROCK, PISM.MASK_GROUNDED, PISM.MASK_FLOATING,
818 PISM.MASK_ICE_FREE_OCEAN]
819 ice_mask.metadata().set_numbers("flag_values", mask_values)
820 ice_mask.metadata().set_string("flag_meanings",
821 "ice_free_bedrock dragging_sheet floating ice_free_ocean")
822 return ice_mask
823
824
825 def createBCMaskVec(grid, name='bc_mask', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
826 """Returns an :cpp:class:`IceModelVec2Int` with attributes setup for
827 masks indicating where Dirichlet data should be applied.
828
829 :param grid: The grid to associate with the vector.
830 :param name: The intrinsic name to give the vector.
831
832 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
833 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
834 is ghosted.
835
836 :param stencil_width: The size of the ghost stencil. Ignored if
837 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
838 ``stencil_width`` is ``None`` and the vector
839 is ghosted, the grid's maximum stencil width
840 is used.
841
842 """
843 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
844
845 bc_mask = PISM.IceModelVec2Int(grid, name, ghost_type, stencil_width)
846 bc_mask.set_attrs("model_state", "grounded_dragging_floating integer mask", "", "", "", 0)
847 mask_values = [0, 1]
848 bc_mask.metadata().set_numbers("flag_values", mask_values)
849 bc_mask.metadata().set_string("flag_meanings", "no_data ssa.dirichlet_bc_location")
850
851 return bc_mask
852
853
854 def createNoModelMaskVec(grid, name='no_model_mask', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
855 """Returns an :cpp:class:`IceModelVec2Int` with attributes setup for
856 regional model flags indicating grid points outside the primary
857 domain, (i.e. where ice does does not evolve in time).
858
859 :param grid: The grid to associate with the vector.
860 :param name: The intrinsic name to give the vector.
861
862 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
863 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
864 is ghosted.
865
866 :param stencil_width: The size of the ghost stencil. Ignored if
867 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
868 ``stencil_width`` is ``None`` and the vector
869 is ghosted, the grid's maximum stencil width
870 is used.
871
872 """
873
874 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
875
876 no_model_mask = PISM.IceModelVec2Int(grid, name, ghost_type, stencil_width)
877
878 no_model_mask.set_attrs("model_state",
879 "mask: zeros (modeling domain) and ones (no-model buffer near grid edges)",
880 "", "", "", 0)
881
882 mask_values = [0, 1]
883 no_model_mask.metadata().set_numbers("flag_values", mask_values)
884 no_model_mask.metadata().set_string("flag_meanings", "normal special_treatment")
885 no_model_mask.time_independent = True
886 no_model_mask.set(0)
887 return no_model_mask
888
889
890 def createZetaFixedMaskVec(grid, name='zeta_fixed_mask', ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
891 """Returns an :cpp:class:`IceModelVec2s` with attributes for the mask
892 determining where parameterized design variables have fixed, known
893 values.
894
895 :param grid: The grid to associate with the vector.
896 :param name: The intrinsic name to give the vector.
897
898 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
899 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
900 is ghosted.
901
902 :param stencil_width: The size of the ghost stencil. Ignored if
903 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
904 ``stencil_width`` is ``None`` and the vector
905 is ghosted, the grid's maximum stencil width
906 is used.
907
908 """
909
910 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
911
912 zeta_fixed_mask = PISM.IceModelVec2Int(grid, name, ghost_type, stencil_width)
913
914 zeta_fixed_mask.set_attrs("model_state", "tauc_unchanging integer mask", "", "", "", 0)
915 mask_values = [0, 1]
916 zeta_fixed_mask.metadata().set_numbers("flag_values", mask_values)
917 zeta_fixed_mask.metadata().set_string("flag_meanings", "tauc_changable tauc_unchangeable")
918
919 return zeta_fixed_mask
920
921
922 def createLongitudeVec(grid, name="lon", ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
923 """Returns an :cpp:class:`IceModelVec2S` with attributes setup for
924 longitude data.
925
926 :param grid: The grid to associate with the vector.
927 :param name: The intrinsic name to give the vector.
928
929 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
930 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
931 is ghosted.
932
933 :param stencil_width: The size of the ghost stencil. Ignored if
934 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
935 ``stencil_width`` is ``None`` and the vector
936 is ghosted, the grid's maximum stencil width
937 is used.
938
939 """
940 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
941 longitude = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
942 longitude.set_attrs("mapping", "longitude", "degree_east", "degree_east", "longitude", 0)
943 longitude.time_independent = True
944 longitude.metadata().set_string("coordinates", "")
945 longitude.metadata().set_string("grid_mapping", "")
946 return longitude
947
948
949 def createLatitudeVec(grid, name="lat", ghost_type=PISM.WITH_GHOSTS, stencil_width=None):
950 """Returns an :cpp:class:`IceModelVec2s` with attributes setup for
951 latitude data.
952
953 :param grid: The grid to associate with the vector.
954 :param name: The intrinsic name to give the vector.
955
956 :param ghost_type: One of ``PISM.WITH_GHOSTS`` or
957 ``PISM.WITHOUT_GHOSTS`` indicating if the vector
958 is ghosted.
959
960 :param stencil_width: The size of the ghost stencil. Ignored if
961 ``ghost_type`` is ``PISM.WITHOUT_GHOSTS``. If
962 ``stencil_width`` is ``None`` and the vector
963 is ghosted, the grid's maximum stencil width
964 is used.
965
966 """
967 stencil_width = _stencil_width(grid, ghost_type, stencil_width)
968
969 latitude = PISM.IceModelVec2S(grid, name, ghost_type, stencil_width)
970
971 latitude.set_attrs("mapping", "latitude", "degree_north", "degree_north", "latitude", 0)
972 latitude.time_independent = True
973 latitude.metadata().set_string("coordinates", "")
974 latitude.metadata().set_string("grid_mapping", "")
975 return latitude