22#include <initializer_list>
26#include "pism/util/error_handling.hh"
67 AccessScope(std::initializer_list<const PetscAccessible *> vecs);
71 void add(
const std::vector<const PetscAccessible*> &vecs);
73 std::vector<const PetscAccessible*>
m_vecs;
80template<
class F,
typename T>
82 auto grid = field.grid();
84 int i_left = 0, i_right = 0, j_bottom = 0, j_top = 0;
85 grid->compute_point_neighbors(x, y, i_left, i_right, j_bottom, j_top);
87 auto w = grid->compute_interp_weights(x, y);
89 return (w[0] * field(i_left, j_bottom) +
90 w[1] * field(i_right, j_bottom) +
91 w[2] * field(i_right, j_top) +
92 w[3] * field(i_left, j_top));
211 std::shared_ptr<const Grid>
grid()
const;
212 unsigned int ndims()
const;
213 std::vector<int>
shape()
const;
215 unsigned int ndof()
const;
218 const std::vector<double>&
levels()
const;
220 std::vector<double>
norm(
int n)
const;
222 void add(
double alpha,
const Array &x);
223 void shift(
double alpha);
224 void scale(
double alpha);
227 std::shared_ptr<petsc::DM>
dm()
const;
229 void set_name(
const std::string &name);
230 const std::string&
get_name()
const;
234 void read(
const std::string &filename,
unsigned int time);
235 void read(
const File &file,
unsigned int time);
237 void write(
const std::string &filename)
const;
262 void view(std::vector<std::shared_ptr<petsc::Viewer> > viewers)
const;
266 const std::string &name,
270 const std::vector<double> &zlevels);
290 void get_dof(std::shared_ptr<petsc::DM> da_result,
petsc::Vec &result,
unsigned int start,
291 unsigned int count=1)
const;
292 void set_dof(std::shared_ptr<petsc::DM> da_source,
petsc::Vec &source,
unsigned int start,
293 unsigned int count=1);
302 void dump(
const char filename[])
const;
306 std::string
checksum(
bool serial)
const;
307 void print_checksum(
const char *prefix =
"",
bool serial=
false)
const;
316static typename std::shared_ptr<T>
cast(std::shared_ptr<Array> input) {
317 auto result = std::dynamic_pointer_cast<T, Array>(input);
316static typename std::shared_ptr<T>
cast(std::shared_ptr<Array> input) {
…}
337 const std::string &spec1,
const std::string &spec2);
High-level PISM I/O class.
Describes the PISM grid and the distribution of data across processors.
virtual ~PetscAccessible()=default
virtual void begin_access() const =0
virtual void end_access() const =0
void add(const PetscAccessible &v)
std::vector< const PetscAccessible * > m_vecs
Makes sure that we call begin_access() and end_access() for all accessed array::Arrays.
void read(const std::string &filename, unsigned int time)
void write_impl(const File &file) const
Writes an Array to a NetCDF file.
virtual void end_access() const
Checks if an Array is allocated and calls DAVecRestoreArray.
virtual void regrid_impl(const File &file, io::Default default_value)
Gets an Array from a file file, interpolating onto the current grid.
unsigned int ndof() const
Returns the number of degrees of freedom per grid point.
void set_interpolation_type(InterpolationType type)
void set_name(const std::string &name)
Sets the variable name to name.
void dump(const char filename[]) const
Dumps a variable to a file, overwriting this file's contents (for debugging).
Array(const Array &other)
void define(const File &file, io::Type default_type) const
Define variables corresponding to an Array in a file opened using file.
virtual void begin_access() const
Checks if an Array is allocated and calls DAVecGetArray.
void shift(double alpha)
Result: v[j] <- v[j] + alpha for all j. Calls VecShift.
const std::vector< double > & levels() const
size_t size() const
Return the total number of elements in the owned part of an array.
void get_from_proc0(petsc::Vec &onp0)
Gets a local Array2 from processor 0.
void write(const std::string &filename) const
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.
const std::string & get_name() const
Get the name of an Array object.
uint64_t fletcher64_serial() const
void get_dof(std::shared_ptr< petsc::DM > da_result, petsc::Vec &result, unsigned int start, unsigned int count=1) const
void scale(double alpha)
Result: v <- v * alpha. Calls VecScale.
Array & operator=(const Array &)
void print_checksum(const char *prefix="", bool serial=false) const
int state_counter() const
Get the object state counter.
std::vector< int > shape() const
std::string checksum(bool serial) const
void view(std::vector< std::shared_ptr< petsc::Viewer > > viewers) const
View a 2D vector field using existing PETSc viewers.
std::shared_ptr< const Grid > grid() const
void set(double c)
Result: v[j] <- c for all j.
void set_dof(std::shared_ptr< petsc::DM > da_source, petsc::Vec &source, unsigned int start, unsigned int count=1)
void inc_state_counter()
Increment the object state counter.
std::shared_ptr< petsc::DM > dm() const
unsigned int ndims() const
Returns the number of spatial dimensions.
void regrid(const std::string &filename, io::Default default_value)
std::shared_ptr< petsc::Vec > allocate_proc0_copy() const
uint64_t fletcher64() const
void put_on_proc0(petsc::Vec &onp0) const
Puts a local array::Scalar on processor 0.
void set_begin_access_use_dof(bool flag)
void add(double alpha, const Array &x)
Result: v <- v + alpha * x. Calls VecAXPY.
std::vector< double > norm(int n) const
Computes the norm of all the components of an Array.
void update_ghosts()
Updates ghost points.
unsigned int stencil_width() const
Get the stencil width of the current Array. Returns 0 if ghosts are not available.
void checkCompatibility(const char *function, const Array &other) const
Checks if two Arrays have compatible sizes, dimensions and numbers of degrees of freedom.
void read_impl(const File &file, unsigned int time)
Reads appropriate NetCDF variable(s) into an Array.
void check_array_indices(int i, int j, unsigned int k) const
Check array indices and warn if they are out of range.
SpatialVariableMetadata & metadata(unsigned int N=0)
Returns a reference to the SpatialVariableMetadata object containing metadata for the compoment N.
Abstract class for reading, writing, allocating, and accessing a DA-based PETSc Vec (2D and 3D fields...
#define PISM_ERROR_LOCATION
static std::shared_ptr< T > cast(std::shared_ptr< Array > input)
std::dynamic_pointer_cast wrapper that checks if the cast succeeded.
Kind
What "kind" of a vector to create: with or without ghosts.
T interpolate(const F &field, double x, double y)
static double F(double SL, double B, double H, double alpha)
void convert_vec(petsc::Vec &v, units::System::Ptr system, const std::string &spec1, const std::string &spec2)