19 #ifndef __IceModelVec_hh
20 #define __IceModelVec_hh
22 #include <initializer_list>
31 #include "pism/util/IceGrid.hh"
32 #include "pism/util/io/IO_Flags.hh"
33 #include "pism/pism_config.hh"
34 #include "pism/util/error_handling.hh"
40 class SpatialVariableMetadata;
62 AccessList(std::initializer_list<const PetscAccessible *> vecs);
66 void add(
const std::vector<const PetscAccessible*> &vecs);
68 std::vector<const PetscAccessible*>
m_vecs;
75 template<
class F,
typename T>
77 auto grid = field.grid();
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);
82 auto w = grid->compute_interp_weights(x, y);
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));
206 typedef std::shared_ptr<IceModelVec>
Ptr;
207 typedef std::shared_ptr<const IceModelVec>
ConstPtr;
212 auto result = std::dynamic_pointer_cast<T,IceModelVec>(input);
222 unsigned int ndims()
const;
223 std::vector<int>
shape()
const;
225 unsigned int ndof()
const;
227 std::vector<double>
levels()
const;
229 std::array<double,2>
range()
const;
230 std::vector<double>
norm(
int n)
const;
233 void shift(
double alpha);
234 void scale(
double alpha);
237 std::shared_ptr<petsc::DM>
dm()
const;
239 void set_name(
const std::string &name);
240 const std::string&
get_name()
const;
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);
251 void read(
const std::string &filename,
unsigned int time);
252 void read(
const File &file,
unsigned int time);
254 void write(
const std::string &filename)
const;
258 double default_value = 0.0);
260 double default_value = 0.0);
280 void view(std::vector<std::shared_ptr<petsc::Viewer> > viewers)
const;
284 const std::string &name,
288 const std::vector<double> &zlevels);
300 double default_value = 0.0);
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);
321 void dump(
const char filename[])
const;
325 std::string
checksum(
bool serial)
const;
326 void print_checksum(
const char *prefix =
"",
bool serial=
false)
const;
341 typedef std::shared_ptr<IceModelVec2S>
Ptr;
342 typedef std::shared_ptr<const IceModelVec2S>
ConstPtr;
347 double const*
const*
array()
const;
356 inline const double&
operator()(
int i,
int j)
const;
394 typedef std::shared_ptr<IceModelVec2Int>
Ptr;
395 typedef std::shared_ptr<const IceModelVec2Int>
ConstPtr;
397 inline int as_int(
int i,
int j)
const;
409 const std::string &name,
411 const std::vector<double> &
levels,
416 const std::string &name,
423 typedef std::shared_ptr<IceModelVec3>
Ptr;
424 typedef std::shared_ptr<const IceModelVec3>
ConstPtr;
427 void set_column(
int i,
int j,
const double *input);
434 inline const double&
operator() (
int i,
int j,
int k)
const;
454 typedef std::shared_ptr<IceModelVec2Stag>
Ptr;
455 typedef std::shared_ptr<const IceModelVec2Stag>
ConstPtr;
475 const std::string &spec1,
const std::string &spec2);
488 bool include_floating_ice,
500 bool include_floating_ice,
std::vector< const PetscAccessible * > m_vecs
void add(const PetscAccessible &v)
Makes sure that we call begin_access() and end_access() for all accessed IceModelVecs.
High-level PISM I/O class.
std::shared_ptr< const IceGrid > ConstPtr
"Cell type" mask. Adds convenience methods to IceModelVec2Int.
std::shared_ptr< IceModelVec2Int > Ptr
stencils::Box< int > box(int i, int j) const
int as_int(int i, int j) const
std::shared_ptr< const IceModelVec2Int > ConstPtr
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...
void add(double alpha, const IceModelVec2S &x)
IceModelVec2S(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, int width=1)
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
std::shared_ptr< const IceModelVec2S > ConstPtr
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
IceModelVec2Stag(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, unsigned int stencil_width=1)
std::shared_ptr< IceModelVec2Stag > Ptr
A class for storing and accessing internal staggered-grid 2D fields. Uses dof=2 storage....
double & operator()(int i, int j, int k)
std::shared_ptr< IceModelVec3 > Ptr
void copy_from(const IceModelVec3 &input)
std::shared_ptr< const IceModelVec3 > ConstPtr
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)
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.
virtual ~IceModelVec3()=default
A virtual class collecting methods common to ice and bedrock 3D fields.
std::vector< int > shape() const
size_t size() const
Return the total number of elements in the owned part of an array.
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.
void update_ghosts()
Updates ghost points.
void get_dof(std::shared_ptr< petsc::DM > da_result, petsc::Vec &result, unsigned int start, unsigned int count=1) const
void read_impl(const File &file, unsigned int time)
Reads appropriate NetCDF variable(s) into an IceModelVec.
std::string checksum(bool serial) const
void set_name(const std::string &name)
Sets the variable name to name.
unsigned int stencil_width() const
Get the stencil width of the current IceModelVec. Returns 0 if ghosts are not available.
SpatialVariableMetadata & metadata(unsigned int N=0)
Returns a reference to the SpatialVariableMetadata object containing metadata for the compoment N.
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).
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.
pism::AccessList AccessList
void regrid(const std::string &filename, RegriddingFlag flag, double default_value=0.0)
void print_checksum(const char *prefix="", bool serial=false) const
std::shared_ptr< IceModelVec > Ptr
virtual void begin_access() const
Checks if an IceModelVec is allocated and calls DAVecGetArray.
std::shared_ptr< petsc::DM > dm() const
IceModelVec(const IceModelVec &other)
void set(double c)
Result: v[j] <- c for all j.
void get_from_proc0(petsc::Vec &onp0)
Gets a local IceModelVec2 from processor 0.
void scale(double alpha)
Result: v <- v * alpha. Calls VecScale.
void read(const std::string &filename, unsigned int time)
IceGrid::ConstPtr grid() const
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.
void set_dof(std::shared_ptr< petsc::DM > da_source, petsc::Vec &source, unsigned int start, unsigned int count=1)
static T::Ptr cast(IceModelVec::Ptr input)
dynamic_pointer_cast wrapper that checks if the cast succeeded.
std::array< double, 2 > range() const
Result: min <- min(v[j]), max <- max(v[j]).
std::vector< double > levels() const
void checkCompatibility(const char *function, const IceModelVec &other) const
Checks if two IceModelVecs have compatible sizes, dimensions and numbers of degrees of freedom.
void define(const File &file, IO_Type default_type=PISM_DOUBLE) const
Define variables corresponding to an IceModelVec in a file opened using file.
void write(const std::string &filename) const
std::shared_ptr< const IceModelVec > ConstPtr
virtual void end_access() const
Checks if an IceModelVec is allocated and calls DAVecRestoreArray.
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.
void write_impl(const File &file) const
Writes an IceModelVec to a NetCDF file.
std::shared_ptr< petsc::Vec > allocate_proc0_copy() const
void set_begin_access_use_dof(bool flag)
std::vector< double > norm(int n) const
Computes the norm of all the components of an IceModelVec.
unsigned int ndims() const
Returns the number of spatial dimensions.
const std::string & get_name() const
Get the name of an IceModelVec object.
void inc_state_counter()
Increment the object state counter.
int state_counter() const
Get the object state counter.
void shift(double alpha)
Result: v[j] <- v[j] + alpha for all j. Calls VecShift.
unsigned int ndof() const
Returns the number of degrees of freedom per grid point.
IceModelVec(IceGrid::ConstPtr grid, const std::string &name, IceModelVecKind ghostedp, size_t dof, size_t stencil_width, const std::vector< double > &zlevels)
void add(double alpha, const IceModelVec &x)
Result: v <- v + alpha * x. Calls VecAXPY.
void set_time_independent(bool flag)
Set the time independent flag for all variables corresponding to this IceModelVec instance.
Abstract class for reading, writing, allocating, and accessing a DA-based PETSc Vec (2D and 3D fields...
virtual ~PetscAccessible()=default
virtual void begin_access() const =0
virtual void end_access() const =0
#define PISM_ERROR_LOCATION
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.
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.
T interpolate(const F &field, double x, double y)
std::shared_ptr< IceModelVec2S > duplicate(const IceModelVec2S &source)
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.
double sum(const IceModelVec2S &input)
Sums up all the values in an IceModelVec2S object. Ignores ghosts.
Star stencil points (in the map-plane).