PISM, A Parallel Ice Sheet Model
iceModelVec.hh
Go to the documentation of this file.
1 // Copyright (C) 2008--2021 Ed Bueler, Constantine Khroulev, and David Maxwell
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 #ifndef __IceModelVec_hh
20 #define __IceModelVec_hh
21 
22 #include <initializer_list>
23 #include <memory> // shared_ptr, dynamic_pointer_cast
24 #include <cstdint> // uint64_t
25 #include <set>
26 #include <map>
27 #include <array>
28 
29 #include "Vector2.hh"
30 #include "stencils.hh"
31 #include "pism/util/IceGrid.hh"
32 #include "pism/util/io/IO_Flags.hh"
33 #include "pism/pism_config.hh" // Pism_DEBUG
34 #include "pism/util/error_handling.hh" // RuntimeError
35 
36 namespace pism {
37 
38 class IceGrid;
39 class File;
40 class SpatialVariableMetadata;
41 
42 namespace petsc {
43 class DM;
44 class Vec;
45 class Viewer;
46 } // end of namespace petsc
47 
48 //! What "kind" of a vector to create: with or without ghosts.
50 
52 public:
53  virtual ~PetscAccessible() = default;
54  virtual void begin_access() const = 0;
55  virtual void end_access() const = 0;
56 };
57 
58 //! Makes sure that we call begin_access() and end_access() for all accessed IceModelVecs.
59 class AccessList {
60 public:
61  AccessList();
62  AccessList(std::initializer_list<const PetscAccessible *> vecs);
63  AccessList(const PetscAccessible &v);
64  ~AccessList();
65  void add(const PetscAccessible &v);
66  void add(const std::vector<const PetscAccessible*> &vecs);
67 private:
68  std::vector<const PetscAccessible*> m_vecs;
69 };
70 
71 /*!
72  * Interpolation helper. Does not check if points needed for interpolation are within the current
73  * processor's sub-domain.
74  */
75 template<class F, typename T>
76 T interpolate(const F &field, double x, double y) {
77  auto grid = field.grid();
78 
79  int i_left = 0, i_right = 0, j_bottom = 0, j_top = 0;
80  grid->compute_point_neighbors(x, y, i_left, i_right, j_bottom, j_top);
81 
82  auto w = grid->compute_interp_weights(x, y);
83 
84  return (w[0] * field(i_left, j_bottom) +
85  w[1] * field(i_right, j_bottom) +
86  w[2] * field(i_right, j_top) +
87  w[3] * field(i_left, j_top));
88 }
89 
90 //! \brief Abstract class for reading, writing, allocating, and accessing a
91 //! DA-based PETSc Vec (2D and 3D fields) from within IceModel.
92 /*!
93  @anchor icemodelvec_use
94 
95  This class represents 2D and 3D fields in PISM. Its methods common to all
96  the derived classes can be split (roughly) into six kinds:
97 
98  - memory allocation (create)
99  - point-wise access (begin_access(), end_access())
100  - arithmetic (range(), norm(), add(), shift(), scale(), set(), ...)
101  - setting or reading metadata (set_attrs(), metadata())
102  - file input/output (read, write, regrid)
103  - tracking whether a field was updated (get_state_counter(), inc_state_counter())
104 
105  ## Memory allocation
106 
107  \code
108  IceModelVec2S var(grid, "var_name", WITH_GHOSTS);
109  // var is ready to use
110  \endcode
111 
112  ("WITH_GHOSTS" means "can be used in computations using map-plane neighbors
113  of grid points.)
114 
115  It is usually a good idea to set variable metadata right after creating it.
116  The method set_attrs() is used throughout PISM to set commonly used
117  attributes.
118 
119  ## Point-wise access
120 
121  PETSc performs some pointer arithmetic magic to allow convenient indexing of
122  grid point values. Because of this one needs to surround the code using row,
123  column or level indexes with begin_access() and end_access() calls:
124 
125  \code
126  double foo;
127  int i = 0, j = 0;
128  IceModelVec2S var;
129  // assume that var was allocated
130  ierr = var.begin_access(); CHKERRQ(ierr);
131  foo = var(i,j) * 2;
132  ierr = var.end_access(); CHKERRQ(ierr);
133  \endcode
134 
135  Please see [this page](@ref computational_grid) for a discussion of the
136  organization of PISM's computational grid and examples of for-loops you will
137  probably put between begin_access() and end_access().
138 
139  To ensure that ghost values are up to date add the following call
140  before the code using ghosts:
141 
142  \code
143  ierr = var.update_ghosts(); CHKERRQ(ierr);
144  \endcode
145 
146  ## Reading and writing variables
147 
148  PISM can read variables either from files with data on a grid matching the
149  current grid (read()) or, using bilinear interpolation, from files
150  containing data on a different (but compatible) grid (regrid()).
151 
152  To write a field to a "prepared" NetCDF file, use write(). (A file is prepared
153  if it contains all the necessary dimensions, coordinate variables and global
154  metadata.)
155 
156  If you need to "prepare" a file, do:
157  \code
158  File file(grid.com, PISM_NETCDF3);
159  io::prepare_for_output(file, *grid.ctx());
160  \endcode
161 
162  A note about NetCDF write performance: due to limitations of the NetCDF
163  (classic, version 3) format, it is significantly faster to
164  \code
165  for (all variables)
166  var.define(...);
167 
168  for (all variables)
169  var.write(...);
170  \endcode
171 
172  as opposed to
173 
174  \code
175  for (all variables) {
176  var.define(...);
177  var.write(...);
178  }
179  \endcode
180 
181  IceModelVec::define() is here so that we can use the first approach.
182 
183  ## Tracking if a field changed
184 
185  It is possible to track if a certain field changed with the help of
186  get_state_counter() and inc_state_counter() methods.
187 
188  For example, PISM's SIA code re-computes the smoothed bed only if the bed
189  deformation code updated it:
190 
191  \code
192  if (bed->get_state_counter() > bed_state_counter) {
193  ierr = bed_smoother->preprocess_bed(...); CHKERRQ(ierr);
194  bed_state_counter = bed->get_state_counter();
195  }
196  \endcode
197 
198  The state counter is **not** updated automatically. For the code snippet above
199  to work, a bed deformation model has to call inc_state_counter() after an
200  update.
201 */
202 class IceModelVec : public PetscAccessible {
203 public:
204  virtual ~IceModelVec();
205 
206  typedef std::shared_ptr<IceModelVec> Ptr;
207  typedef std::shared_ptr<const IceModelVec> ConstPtr;
208 
209  //! `dynamic_pointer_cast` wrapper that checks if the cast succeeded.
210  template<class T>
211  static typename T::Ptr cast(IceModelVec::Ptr input) {
212  auto result = std::dynamic_pointer_cast<T,IceModelVec>(input);
213 
214  if (not result) {
215  throw RuntimeError(PISM_ERROR_LOCATION, "dynamic cast failure");
216  }
217 
218  return result;
219  }
220 
221  IceGrid::ConstPtr grid() const;
222  unsigned int ndims() const;
223  std::vector<int> shape() const;
224  //! @brief Returns the number of degrees of freedom per grid point.
225  unsigned int ndof() const;
226  unsigned int stencil_width() const;
227  std::vector<double> levels() const;
228 
229  std::array<double,2> range() const;
230  std::vector<double> norm(int n) const;
231 
232  void add(double alpha, const IceModelVec &x);
233  void shift(double alpha);
234  void scale(double alpha);
235 
236  petsc::Vec& vec() const;
237  std::shared_ptr<petsc::DM> dm() const;
238 
239  void set_name(const std::string &name);
240  const std::string& get_name() const;
241 
242  void set_attrs(const std::string &pism_intent,
243  const std::string &long_name,
244  const std::string &units,
245  const std::string &glaciological_units,
246  const std::string &standard_name,
247  unsigned int component);
248 
249  void define(const File &file, IO_Type default_type = PISM_DOUBLE) const;
250 
251  void read(const std::string &filename, unsigned int time);
252  void read(const File &file, unsigned int time);
253 
254  void write(const std::string &filename) const;
255  void write(const File &file) const;
256 
257  void regrid(const std::string &filename, RegriddingFlag flag,
258  double default_value = 0.0);
259  void regrid(const File &file, RegriddingFlag flag,
260  double default_value = 0.0);
261 
262  virtual void begin_access() const;
263  virtual void end_access() const;
264  void update_ghosts();
265 
266  std::shared_ptr<petsc::Vec> allocate_proc0_copy() const;
267  void put_on_proc0(petsc::Vec &onp0) const;
268  void get_from_proc0(petsc::Vec &onp0);
269 
270  void set(double c);
271 
272  SpatialVariableMetadata& metadata(unsigned int N = 0);
273 
274  const SpatialVariableMetadata& metadata(unsigned int N = 0) const;
275 
276  int state_counter() const;
277  void inc_state_counter();
278  void set_time_independent(bool flag);
279 
280  void view(std::vector<std::shared_ptr<petsc::Viewer> > viewers) const;
281 
282 protected:
284  const std::string &name,
285  IceModelVecKind ghostedp,
286  size_t dof,
287  size_t stencil_width,
288  const std::vector<double> &zlevels);
289  struct Impl;
291 
292  // will be cast to double**, double***, or Vector2** in derived classes
293  // This is not hidden in m_impl to make it possible to inline operator()
294  mutable void *m_array;
295 
296  void set_begin_access_use_dof(bool flag);
297 
298  void read_impl(const File &file, unsigned int time);
299  void regrid_impl(const File &file, RegriddingFlag flag,
300  double default_value = 0.0);
301  void write_impl(const File &file) const;
302 
303  void checkCompatibility(const char *function, const IceModelVec &other) const;
304 
305  //! @brief Check array indices and warn if they are out of range.
306  void check_array_indices(int i, int j, unsigned int k) const;
307 
308  void copy_to_vec(std::shared_ptr<petsc::DM> destination_da, petsc::Vec &destination) const;
309  void get_dof(std::shared_ptr<petsc::DM> da_result, petsc::Vec &result, unsigned int start,
310  unsigned int count=1) const;
311  void set_dof(std::shared_ptr<petsc::DM> da_source, petsc::Vec &source, unsigned int start,
312  unsigned int count=1);
313 private:
314  size_t size() const;
315  // disable copy constructor and the assignment operator:
316  IceModelVec(const IceModelVec &other);
318 public:
319  //! Dump an IceModelVec to a file. *This is for debugging only.*
320  //! Uses const char[] to make it easier to call it from gdb.
321  void dump(const char filename[]) const;
322 
323  uint64_t fletcher64_serial() const;
324  uint64_t fletcher64() const;
325  std::string checksum(bool serial) const;
326  void print_checksum(const char *prefix = "", bool serial=false) const;
327 
329 protected:
330  void put_on_proc0(petsc::Vec &parallel, petsc::Vec &onp0) const;
331  void get_from_proc0(petsc::Vec &onp0, petsc::Vec &parallel) const;
332 };
333 
334 /** A class for storing and accessing scalar 2D fields.
335  IceModelVec2S is just IceModelVec2 with "dof == 1" */
336 class IceModelVec2S : public IceModelVec {
337 public:
338  IceModelVec2S(IceGrid::ConstPtr grid, const std::string &name,
339  IceModelVecKind ghostedp, int width = 1);
340 
341  typedef std::shared_ptr<IceModelVec2S> Ptr;
342  typedef std::shared_ptr<const IceModelVec2S> ConstPtr;
343 
344  // does not need a copy constructor, because it does not add any new data members
345  void copy_from(const IceModelVec2S &source);
346  double** array();
347  double const* const* array() const;
348  void add(double alpha, const IceModelVec2S &x);
349  void add(double alpha, const IceModelVec2S &x, IceModelVec2S &result) const;
350 
351  //! Provides access (both read and write) to the internal double array.
352  /*!
353  Note that i corresponds to the x direction and j to the y.
354  */
355  inline double& operator() (int i, int j);
356  inline const double& operator()(int i, int j) const;
357  inline stencils::Star<double> star(int i, int j) const;
358  inline stencils::Box<double> box(int i, int j) const;
359 };
360 
361 // Finite-difference shortcuts. They may be slower than hard-coding FD approximations of x
362 // and y derivatives. Use with care.
363 double diff_x(const IceModelVec2S &array, int i, int j);
364 double diff_y(const IceModelVec2S &array, int i, int j);
365 
366 // These take grid periodicity into account and use one-sided differences at domain edges.
367 double diff_x_p(const IceModelVec2S &array, int i, int j);
368 double diff_y_p(const IceModelVec2S &array, int i, int j);
369 
370 double sum(const IceModelVec2S &input);
371 double min(const IceModelVec2S &input);
372 double max(const IceModelVec2S &input);
373 double absmax(const IceModelVec2S &input);
374 
375 void apply_mask(const IceModelVec2S &M, double fill, IceModelVec2S &result);
376 
377 void compute_magnitude(const IceModelVec2S &v_x,
378  const IceModelVec2S &v_y,
379  IceModelVec2S &result);
380 
381 class IceModelVec2V;
382 
383 void compute_magnitude(const IceModelVec2V &input, IceModelVec2S &result);
384 
385 std::shared_ptr<IceModelVec2S> duplicate(const IceModelVec2S &source);
386 
387 //! \brief A simple class "hiding" the fact that the mask is stored as
388 //! floating-point scalars (instead of integers).
390 public:
391  IceModelVec2Int(IceGrid::ConstPtr grid, const std::string &name,
392  IceModelVecKind ghostedp, int width = 1);
393 
394  typedef std::shared_ptr<IceModelVec2Int> Ptr;
395  typedef std::shared_ptr<const IceModelVec2Int> ConstPtr;
396 
397  inline int as_int(int i, int j) const;
398  inline stencils::Star<int> star(int i, int j) const;
399  inline stencils::Box<int> box(int i, int j) const;
400 };
401 
402 //! \brief A virtual class collecting methods common to ice and bedrock 3D
403 //! fields.
404 class IceModelVec3 : public IceModelVec {
405 public:
406 
407  // Three-dimensional array with a number of vertical levels
409  const std::string &name,
410  IceModelVecKind ghostedp,
411  const std::vector<double> &levels,
412  unsigned int stencil_width = 1);
413 
414  // A collection of two-dimensional arrays using three-dimensional indexing
416  const std::string &name,
417  IceModelVecKind ghostedp,
418  unsigned int dof,
419  unsigned int stencil_width = 1);
420 
421  virtual ~IceModelVec3() = default;
422 
423  typedef std::shared_ptr<IceModelVec3> Ptr;
424  typedef std::shared_ptr<const IceModelVec3> ConstPtr;
425 
426  void set_column(int i, int j, double c);
427  void set_column(int i, int j, const double *input);
428  double* get_column(int i, int j);
429  const double* get_column(int i, int j) const;
430 
431  double interpolate(int i, int j, double z) const;
432 
433  inline double& operator() (int i, int j, int k);
434  inline const double& operator() (int i, int j, int k) const;
435 
436  void copy_from(const IceModelVec3 &input);
437 };
438 
439 void extract_surface(const IceModelVec3 &data, double z, IceModelVec2S &output);
440 void extract_surface(const IceModelVec3 &data, const IceModelVec2S &z, IceModelVec2S &output);
441 
442 void sum_columns(const IceModelVec3 &data, double A, double B, IceModelVec2S &output);
443 
444 std::shared_ptr<IceModelVec3> duplicate(const IceModelVec3 &source);
445 
446 //! \brief A class for storing and accessing internal staggered-grid 2D fields.
447 //! Uses dof=2 storage. This class is identical to IceModelVec2V, except that
448 //! components are not called `u` and `v` (to avoid confusion).
450 public:
451  IceModelVec2Stag(IceGrid::ConstPtr grid, const std::string &name,
452  IceModelVecKind ghostedp, unsigned int stencil_width = 1);
453 
454  typedef std::shared_ptr<IceModelVec2Stag> Ptr;
455  typedef std::shared_ptr<const IceModelVec2Stag> ConstPtr;
456 
457  //! Returns the values at interfaces of the cell i,j using the staggered grid.
458  /*! The ij member of the return value is set to 0, since it has no meaning in
459  this context.
460  */
461  inline stencils::Star<double> star(int i, int j) const;
462 };
463 
464 std::array<double,2> absmax(const IceModelVec2Stag &input);
465 
466 /**
467  * Convert a PETSc Vec from the units in `from` into units in `to` (in place).
468  *
469  * @param v data to convert
470  * @param system unit system
471  * @param spec1 source unit specification string
472  * @param spec2 destination unit specification string
473  */
474 void convert_vec(petsc::Vec &v, std::shared_ptr<units::System> system,
475  const std::string &spec1, const std::string &spec2);
476 
478 
479 /*!
480  * Average a scalar field from the staggered grid onto the regular grid by considering
481  * only ice-covered grid.
482  *
483  * If `include_floating_ice` is true, include floating ice, otherwise consider grounded
484  * icy cells only.
485  */
486 void staggered_to_regular(const IceModelVec2CellType &cell_type,
487  const IceModelVec2Stag &input,
488  bool include_floating_ice,
489  IceModelVec2S &result);
490 
491 /*!
492  * Average a vector field from the staggered grid onto the regular grid by considering
493  * only ice-covered grid.
494  *
495  * If `include_floating_ice` is true, include floating ice, otherwise consider grounded
496  * icy cells only.
497  */
498 void staggered_to_regular(const IceModelVec2CellType &cell_type,
499  const IceModelVec2Stag &input,
500  bool include_floating_ice,
501  IceModelVec2V &result);
502 
503 } // end of namespace pism
504 
505 // include inline methods; contents are wrapped in namespace pism {...}
506 #include "IceModelVec_inline.hh"
507 
508 #endif /* __IceModelVec_hh */
std::vector< const PetscAccessible * > m_vecs
Definition: iceModelVec.hh:68
void add(const PetscAccessible &v)
Definition: iceModelVec.cc:963
Makes sure that we call begin_access() and end_access() for all accessed IceModelVecs.
Definition: iceModelVec.hh:59
High-level PISM I/O class.
Definition: File.hh:51
std::shared_ptr< const IceGrid > ConstPtr
Definition: IceGrid.hh:233
"Cell type" mask. Adds convenience methods to IceModelVec2Int.
std::shared_ptr< IceModelVec2Int > Ptr
Definition: iceModelVec.hh:394
stencils::Box< int > box(int i, int j) const
int as_int(int i, int j) const
std::shared_ptr< const IceModelVec2Int > ConstPtr
Definition: iceModelVec.hh:395
stencils::Star< int > star(int i, int j) const
IceModelVec2Int(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, int width=1)
A simple class "hiding" the fact that the mask is stored as floating-point scalars (instead of intege...
Definition: iceModelVec.hh:389
void add(double alpha, const IceModelVec2S &x)
IceModelVec2S(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, int width=1)
Definition: iceModelVec2.cc:49
void copy_from(const IceModelVec2S &source)
double & operator()(int i, int j)
Provides access (both read and write) to the internal double array.
stencils::Box< double > box(int i, int j) const
std::shared_ptr< IceModelVec2S > Ptr
Definition: iceModelVec.hh:341
std::shared_ptr< const IceModelVec2S > ConstPtr
Definition: iceModelVec.hh:342
stencils::Star< double > star(int i, int j) const
stencils::Star< double > star(int i, int j) const
Returns the values at interfaces of the cell i,j using the staggered grid.
std::shared_ptr< const IceModelVec2Stag > ConstPtr
Definition: iceModelVec.hh:455
IceModelVec2Stag(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, unsigned int stencil_width=1)
std::shared_ptr< IceModelVec2Stag > Ptr
Definition: iceModelVec.hh:454
A class for storing and accessing internal staggered-grid 2D fields. Uses dof=2 storage....
Definition: iceModelVec.hh:449
double & operator()(int i, int j, int k)
std::shared_ptr< IceModelVec3 > Ptr
Definition: iceModelVec.hh:423
void copy_from(const IceModelVec3 &input)
std::shared_ptr< const IceModelVec3 > ConstPtr
Definition: iceModelVec.hh:424
double * get_column(int i, int j)
IceModelVec3(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, const std::vector< double > &levels, unsigned int stencil_width=1)
Definition: iceModelVec3.cc:43
double interpolate(int i, int j, double z) const
Return value of scalar quantity at level z (m) above base of ice (by linear interpolation).
void set_column(int i, int j, double c)
Set all values of scalar quantity to given a single value in a particular column.
Definition: iceModelVec3.cc:63
virtual ~IceModelVec3()=default
A virtual class collecting methods common to ice and bedrock 3D fields.
Definition: iceModelVec.hh:404
std::vector< int > shape() const
Definition: iceModelVec.cc:156
size_t size() const
Return the total number of elements in the owned part of an array.
Definition: iceModelVec.cc:976
void regrid_impl(const File &file, RegriddingFlag flag, double default_value=0.0)
Gets an IceModelVec from a file file, interpolating onto the current grid.
Definition: iceModelVec.cc:422
void update_ghosts()
Updates ghost points.
Definition: iceModelVec.cc:669
void get_dof(std::shared_ptr< petsc::DM > da_result, petsc::Vec &result, unsigned int start, unsigned int count=1) const
Definition: iceModelVec.cc:277
void read_impl(const File &file, unsigned int time)
Reads appropriate NetCDF variable(s) into an IceModelVec.
Definition: iceModelVec.cc:476
std::string checksum(bool serial) const
void set_name(const std::string &name)
Sets the variable name to name.
Definition: iceModelVec.cc:374
unsigned int stencil_width() const
Get the stencil width of the current IceModelVec. Returns 0 if ghosts are not available.
Definition: iceModelVec.cc:334
SpatialVariableMetadata & metadata(unsigned int N=0)
Returns a reference to the SpatialVariableMetadata object containing metadata for the compoment N.
Definition: iceModelVec.cc:533
void view(std::vector< std::shared_ptr< petsc::Viewer > > viewers) const
View a 2D vector field using existing PETSc viewers.
void dump(const char filename[]) const
Dumps a variable to a file, overwriting this file's contents (for debugging).
Definition: iceModelVec.cc:583
IceModelVec & operator=(const IceModelVec &)
uint64_t fletcher64() const
uint64_t fletcher64_serial() const
void set_attrs(const std::string &pism_intent, const std::string &long_name, const std::string &units, const std::string &glaciological_units, const std::string &standard_name, unsigned int component)
Sets NetCDF attributes of an IceModelVec object.
Definition: iceModelVec.cc:399
pism::AccessList AccessList
Definition: iceModelVec.hh:328
void regrid(const std::string &filename, RegriddingFlag flag, double default_value=0.0)
Definition: iceModelVec.cc:838
void print_checksum(const char *prefix="", bool serial=false) const
std::shared_ptr< IceModelVec > Ptr
Definition: iceModelVec.hh:206
virtual void begin_access() const
Checks if an IceModelVec is allocated and calls DAVecGetArray.
Definition: iceModelVec.cc:622
petsc::Vec & vec() const
Definition: iceModelVec.cc:342
std::shared_ptr< petsc::DM > dm() const
Definition: iceModelVec.cc:356
IceModelVec(const IceModelVec &other)
virtual ~IceModelVec()
Definition: iceModelVec.cc:101
void set(double c)
Result: v[j] <- c for all j.
Definition: iceModelVec.cc:683
void get_from_proc0(petsc::Vec &onp0)
Gets a local IceModelVec2 from processor 0.
void scale(double alpha)
Result: v <- v * alpha. Calls VecScale.
Definition: iceModelVec.cc:252
void read(const std::string &filename, unsigned int time)
Definition: iceModelVec.cc:833
IceGrid::ConstPtr grid() const
Definition: iceModelVec.cc:128
void put_on_proc0(petsc::Vec &onp0) const
Puts a local IceModelVec2S on processor 0.
void check_array_indices(int i, int j, unsigned int k) const
Check array indices and warn if they are out of range.
Definition: iceModelVec.cc:690
void set_dof(std::shared_ptr< petsc::DM > da_source, petsc::Vec &source, unsigned int start, unsigned int count=1)
Definition: iceModelVec.cc:304
static T::Ptr cast(IceModelVec::Ptr input)
dynamic_pointer_cast wrapper that checks if the cast succeeded.
Definition: iceModelVec.hh:211
std::array< double, 2 > range() const
Result: min <- min(v[j]), max <- max(v[j]).
Definition: iceModelVec.cc:193
std::vector< double > levels() const
Definition: iceModelVec.cc:136
void checkCompatibility(const char *function, const IceModelVec &other) const
Checks if two IceModelVecs have compatible sizes, dimensions and numbers of degrees of freedom.
Definition: iceModelVec.cc:599
void define(const File &file, IO_Type default_type=PISM_DOUBLE) const
Define variables corresponding to an IceModelVec in a file opened using file.
Definition: iceModelVec.cc:523
void write(const std::string &filename) const
Definition: iceModelVec.cc:822
std::shared_ptr< const IceModelVec > ConstPtr
Definition: iceModelVec.hh:207
virtual void end_access() const
Checks if an IceModelVec is allocated and calls DAVecRestoreArray.
Definition: iceModelVec.cc:643
void copy_to_vec(std::shared_ptr< petsc::DM > destination_da, petsc::Vec &destination) const
Copies v to a global vector 'destination'. Ghost points are discarded.
Definition: iceModelVec.cc:266
void write_impl(const File &file) const
Writes an IceModelVec to a NetCDF file.
Definition: iceModelVec.cc:544
std::shared_ptr< petsc::Vec > allocate_proc0_copy() const
Definition: iceModelVec.cc:993
void set_begin_access_use_dof(bool flag)
Definition: iceModelVec.cc:181
std::vector< double > norm(int n) const
Computes the norm of all the components of an IceModelVec.
Definition: iceModelVec.cc:769
unsigned int ndims() const
Returns the number of spatial dimensions.
Definition: iceModelVec.cc:152
const std::string & get_name() const
Get the name of an IceModelVec object.
Definition: iceModelVec.cc:386
void inc_state_counter()
Increment the object state counter.
Definition: iceModelVec.cc:147
int state_counter() const
Get the object state counter.
Definition: iceModelVec.cc:124
void shift(double alpha)
Result: v[j] <- v[j] + alpha for all j. Calls VecShift.
Definition: iceModelVec.cc:244
unsigned int ndof() const
Returns the number of degrees of freedom per grid point.
Definition: iceModelVec.cc:132
IceModelVec(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, size_t dof, size_t stencil_width, const std::vector< double > &zlevels)
Definition: iceModelVec.cc:56
void add(double alpha, const IceModelVec &x)
Result: v <- v + alpha * x. Calls VecAXPY.
Definition: iceModelVec.cc:234
void set_time_independent(bool flag)
Set the time independent flag for all variables corresponding to this IceModelVec instance.
Definition: iceModelVec.cc:175
Abstract class for reading, writing, allocating, and accessing a DA-based PETSc Vec (2D and 3D fields...
Definition: iceModelVec.hh:202
virtual ~PetscAccessible()=default
virtual void begin_access() const =0
virtual void end_access() const =0
Spatial NetCDF variable (corresponding to a 2D or 3D scalar field).
#define PISM_ERROR_LOCATION
#define n
Definition: exactTestM.c:37
void apply_mask(const IceModelVec2S &M, double fill, IceModelVec2S &result)
Masks out all the areas where by setting them to fill.
double diff_x(const IceModelVec2S &array, int i, int j)
Returns the x-derivative at i,j approximated using centered finite differences.
double max(const IceModelVec2S &input)
Finds maximum over all the values in an IceModelVec2S object. Ignores ghosts.
static double F(double SL, double B, double H, double alpha)
void compute_magnitude(const IceModelVec2S &v_x, const IceModelVec2S &v_y, IceModelVec2S &result)
Sets an IceModelVec2 to the magnitude of a 2D vector field with components v_x and v_y.
Definition: iceModelVec2.cc:80
IO_Type
Definition: IO_Flags.hh:31
@ PISM_DOUBLE
Definition: IO_Flags.hh:38
double diff_x_p(const IceModelVec2S &array, int i, int j)
Returns the x-derivative at i,j approximated using centered finite differences. Respects grid periodi...
void extract_surface(const IceModelVec3 &data, double z, IceModelVec2S &output)
Copies a horizontal slice at level z of an IceModelVec3 into an IceModelVec2S gslice.
double absmax(const IceModelVec2S &input)
Finds maximum over all the absolute values in an IceModelVec2S object. Ignores ghosts.
static const double k
Definition: exactTestP.cc:45
RegriddingFlag
Definition: IO_Flags.hh:70
T interpolate(const F &field, double x, double y)
Definition: iceModelVec.hh:76
std::shared_ptr< IceModelVec2S > duplicate(const IceModelVec2S &source)
Definition: iceModelVec2.cc:55
void sum_columns(const IceModelVec3 &data, double A, double B, IceModelVec2S &output)
double min(const IceModelVec2S &input)
Finds minimum over all the values in an IceModelVec2S object. Ignores ghosts.
double diff_y_p(const IceModelVec2S &array, int i, int j)
Returns the y-derivative at i,j approximated using centered finite differences. Respects grid periodi...
double diff_y(const IceModelVec2S &array, int i, int j)
Returns the y-derivative at i,j approximated using centered finite differences.
void convert_vec(petsc::Vec &v, units::System::Ptr system, const std::string &spec1, const std::string &spec2)
void staggered_to_regular(const IceModelVec2CellType &cell_type, const IceModelVec2Stag &input, bool include_floating_ice, IceModelVec2S &result)
IceModelVecKind
What "kind" of a vector to create: with or without ghosts.
Definition: iceModelVec.hh:49
@ WITHOUT_GHOSTS
Definition: iceModelVec.hh:49
@ WITH_GHOSTS
Definition: iceModelVec.hh:49
double sum(const IceModelVec2S &input)
Sums up all the values in an IceModelVec2S object. Ignores ghosts.
Star stencil points (in the map-plane).
Definition: stencils.hh:30
int count
Definition: test_cube.c:16