SSAGES  0.1
A MetaDynamics Package
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Groups Pages
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
SSAGES::GridBase< T > Class Template Referenceabstract

Base class for Grids. More...

#include <GridBase.h>

Inheritance diagram for SSAGES::GridBase< T >:
Inheritance graph
[legend]

Public Member Functions

size_t GetDimension () const
 Get the dimension. More...
 
const std::vector< int > & GetNumPoints () const
 Get the number of points for all dimensions. More...
 
int GetNumPoints (size_t dim) const
 Get the number of points for a specific dimension. More...
 
const std::vector< double > & GetLower () const
 Return the lower edges of the Grid. More...
 
double GetLower (size_t dim) const
 Get the lower edge for a specific dimension. More...
 
const std::vector< double > & GetUpper () const
 Return the upper edges of the Grid. More...
 
double GetUpper (size_t dim) const
 Get the upper edge for a specific dimension. More...
 
const std::vector< bool > & GetPeriodic () const
 Return the periodicity of the Grid. More...
 
bool GetPeriodic (size_t dim) const
 Get the periodicity in a specific dimension. More...
 
size_t size () const
 Get the size of the internal storage vector. More...
 
T * data ()
 Get pointer to the internal data storage vector. More...
 
T const * data () const
 Get pointer to const of the internal data storage vector. More...
 
std::vector< int > GetIndices (const std::vector< double > &x) const
 Return the Grid indices for a given point. More...
 
double GetInterpolated (const std::vector< double > &x)
 Return linear interpolation on a coordinate. More...
 
int GetIndex (double x) const
 Return the Grid index for a one-dimensional grid. More...
 
std::vector< double > GetCoordinates (const std::vector< int > &indices)
 Return coordinates of the grid center points. More...
 
double GetCoordinate (int index)
 Return center point of 1d-grid. More...
 
const T & at (const std::vector< int > &indices) const
 Access Grid element read-only. More...
 
T & at (const std::vector< int > &indices)
 Access Grid element read/write. More...
 
template<typename R >
const T & at (std::initializer_list< R > &&x) const
 Const access of Grid element via initializer list. More...
 
template<typename R >
T & at (std::initializer_list< R > &&x)
 Access Grid element via initializer list. More...
 
const T & at (int index) const
 Access 1d Grid by index, read-only. More...
 
T & at (int index)
 Access 1d Grid by index, read-write. More...
 
const T & at (const std::vector< double > &x) const
 Access Grid element pertaining to a specific point – read-only. More...
 
T & at (const std::vector< double > &x)
 Access Grid element pertaining to a specific point – read/write. More...
 
const T & at (double x) const
 Access 1d-Grid by point - read-only. More...
 
T & at (double x)
 Access 1d-Grid by point - read-write. More...
 
const T & operator[] (const std::vector< int > &indices) const
 Access Grid element per [] read-only. More...
 
T & operator[] (const std::vector< int > &indices)
 Access Grid element per [] read-write. More...
 
template<typename R >
const T & operator[] (std::initializer_list< R > &&x) const
 Const access of Grid element via initializer list. More...
 
template<typename R >
T & operator[] (std::initializer_list< R > &&x)
 Access Grid element via initializer list. More...
 
const T & operator[] (int index) const
 Access 1d-Grid per [] operator, read-only. More...
 
T & operator[] (int index)
 Access 1d-Grid per [] operator, read-write. More...
 
const T & operator[] (const std::vector< double > &x) const
 Access Grid element pertaining to a specific point per [] read-only. More...
 
T & operator[] (const std::vector< double > &x)
 Access Grid element pertaining to a specific point per [] read-write. More...
 
const T & operator[] (double x) const
 Access 1d-Grid via specific point, read-only. More...
 
T & operator[] (double x)
 Access 1d-Grid via specific point, read-write. More...
 

Protected Member Functions

std::vector< int > wrapIndices (const std::vector< int > &indices) const
 Wrap the index around periodic boundaries.
 
virtual size_t mapTo1d (const std::vector< int > &indices) const =0
 This function needs to be implemented by child classes. More...
 
 GridBase (std::vector< int > numPoints, std::vector< double > lower, std::vector< double > upper, std::vector< bool > isPeriodic)
 Constructor. More...
 

Protected Attributes

std::vector< T > data_
 Internal storage of the data.
 
size_t dimension_
 Dimension of the grid.
 
std::vector< int > numPoints_
 Number of points in each dimension.
 
std::pair< std::vector< double >
, std::vector< double > > 
edges_
 Edges of the Grid in each dimension.
 
std::vector< bool > isPeriodic_
 Periodicity of the Grid.
 

Detailed Description

template<typename T>
class SSAGES::GridBase< T >

Base class for Grids.

Template Parameters
Typeof data stored on the grid

Base class for all grids. Currently, these are 'Grid' and 'Histogram'.

Definition at line 38 of file GridBase.h.

Constructor & Destructor Documentation

template<typename T>
SSAGES::GridBase< T >::GridBase ( std::vector< int >  numPoints,
std::vector< double >  lower,
std::vector< double >  upper,
std::vector< bool >  isPeriodic 
)
inlineprotected

Constructor.

Parameters
numPointsNumber of grid points in each dimension.
lowerLower edges of the grid.
upperUpper edges of the grid.
isPeriodicBools specifying the periodicity in the respective dimension.

The constructor is protected by design. This makes sure that only child classes of GridBase are constructed.

Definition at line 97 of file GridBase.h.

101  : dimension_(numPoints.size()),
102  numPoints_(numPoints),
103  edges_(std::pair< std::vector<double>, std::vector<double> >(lower, upper)),
104  isPeriodic_(isPeriodic)
105  {
106  // Check that vector sizes are correct
107  if (edges_.first.size() != dimension_ ||
108  edges_.second.size() != dimension_) {
109  throw std::invalid_argument("Size of vector containing upper or "
110  "lower edges, does not match size of vector containing "
111  "number of grid points.");
112  }
113  if (isPeriodic_.size() == 0) {
114  // Default: Non-periodic in all dimensions
115  isPeriodic.resize(dimension_, false);
116  } else if (isPeriodic_.size() != dimension_) {
117  throw std::invalid_argument("Size of vector isPeriodic does not "
118  "match size of vector containing number of grid points.");
119  }
120 
121  }
std::pair< std::vector< double >, std::vector< double > > edges_
Edges of the Grid in each dimension.
Definition: GridBase.h:51
std::vector< int > numPoints_
Number of points in each dimension.
Definition: GridBase.h:48
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45
std::vector< bool > isPeriodic_
Periodicity of the Grid.
Definition: GridBase.h:54

Member Function Documentation

template<typename T>
const T& SSAGES::GridBase< T >::at ( const std::vector< int > &  indices) const
inline

Access Grid element read-only.

Parameters
indicesVector of integers specifying the grid point.
Returns
const reference of the value stored at the given grid point.

In non-periodic dimensions, the index needs to be in the interval [-1, numPoints]. Grid::at(-1) accessed the underflow bin, Grid::at(numPoints) accesses the overflow bin.

In periodic dimensions, the index may take any integer value and will be mapped back to the interval [0, numPoints-1]. Thus, Grid::at(-1) will access the same value as Grid::at(numPoints-1).

Definition at line 452 of file GridBase.h.

Referenced by SSAGES::GridBase< int >::at(), SSAGES::GridBase< int >::GetCoordinate(), SSAGES::GridBase< int >::GetIndex(), SSAGES::GridBase< int >::GetInterpolated(), and SSAGES::GridBase< int >::operator[]().

453  {
454  // Check that indices are in bound.
455  if (indices.size() != GetDimension()) {
456  throw std::invalid_argument("Dimension of indices does not match "
457  "dimension of the grid.");
458  }
459 
460  return data_.at(mapTo1d(indices));
461  }
std::vector< T > data_
Internal storage of the data.
Definition: GridBase.h:42
virtual size_t mapTo1d(const std::vector< int > &indices) const =0
This function needs to be implemented by child classes.
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:128

Here is the caller graph for this function:

template<typename T>
T& SSAGES::GridBase< T >::at ( const std::vector< int > &  indices)
inline

Access Grid element read/write.

Parameters
indicesVector of integers specifying the grid point.
Returns
Reference to value at the specified grid point.

Definition at line 468 of file GridBase.h.

469  {
470  return const_cast<T&>(static_cast<const GridBase<T>* >(this)->at(indices));
471  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
template<typename R >
const T& SSAGES::GridBase< T >::at ( std::initializer_list< R > &&  x) const
inline

Const access of Grid element via initializer list.

Template Parameters
RDatatype in the initializer list
Parameters
xinitializer list
Returns
Const reference to value at the specified point.

This function avoids abiguity if at() is called with a brace-enclosed initializer list. The template parameter makes sure that this function can be called with either ints, specifying a grid point, or doubles, specifying coordinates in space, inside the initializer list.

Definition at line 485 of file GridBase.h.

486  {
487  return at(static_cast<std::vector<R> >(x));
488  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
template<typename R >
T& SSAGES::GridBase< T >::at ( std::initializer_list< R > &&  x)
inline

Access Grid element via initializer list.

Template Parameters
RDatatype in the initializer list
Parameters
xinitializer list
Returns
Reference to value at the specified point.

This function avoids abiguity if at() is called with a brace-enclosed initializer list. The template parameter makes sure that this function can be called with either ints, specifying a grid point, or doubles, specifying coordinates in space, inside the initializer list.

Definition at line 502 of file GridBase.h.

503  {
504  return at(static_cast<std::vector<R> >(x));
505  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
const T& SSAGES::GridBase< T >::at ( int  index) const
inline

Access 1d Grid by index, read-only.

Parameters
indexIndex specifying the grid point.
Returns
Const reference of value at the given index.
Note
This function can only be used for 1d-Grids.

Definition at line 514 of file GridBase.h.

515  {
516  if (dimension_ != 1) {
517  throw std::invalid_argument("Only 1d-Grids can be accessed with a "
518  "single integer as the index.");
519  }
520  return at({index});
521  }
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T& SSAGES::GridBase< T >::at ( int  index)
inline

Access 1d Grid by index, read-write.

Parameters
indexIndex specifying the grid point.
Returns
Reference of value at the given grid point.
Note
This function can only be used for 1d-Grids.

Definition at line 530 of file GridBase.h.

531  {
532  return const_cast<T&>(static_cast<const GridBase<T>* >(this)->at(index));
533  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
const T& SSAGES::GridBase< T >::at ( const std::vector< double > &  x) const
inline

Access Grid element pertaining to a specific point – read-only.

Parameters
xVector of doubles specifying a point.
Returns
Const reference of the value at the given coordinates.

This function is provided for convenience. It is identical to GridBase::at(GridBase::GetIndices(x)).

Definition at line 543 of file GridBase.h.

544  {
545  return at(GetIndices(x));
546  }
std::vector< int > GetIndices(const std::vector< double > &x) const
Return the Grid indices for a given point.
Definition: GridBase.h:296
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T& SSAGES::GridBase< T >::at ( const std::vector< double > &  x)
inline

Access Grid element pertaining to a specific point – read/write.

Parameters
xVector of doubles specifying a point.
Returns
Reference to the value at the given coordinates.

This function is provided for convenience. It is identical to GridBase::at(GridBase::GetIndices(x)).

Definition at line 556 of file GridBase.h.

557  {
558  return at(GetIndices(x));
559  }
std::vector< int > GetIndices(const std::vector< double > &x) const
Return the Grid indices for a given point.
Definition: GridBase.h:296
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
const T& SSAGES::GridBase< T >::at ( double  x) const
inline

Access 1d-Grid by point - read-only.

Parameters
xAccess grid point pertaining to this value.
Returns
Const reference to the value pertaining to the specified coordinate.
Note
This function can only be used for 1d-Grids.

Definition at line 569 of file GridBase.h.

570  {
571  if (dimension_ != 1) {
572  throw std::invalid_argument("Only 1d-Grids can be accessed with a "
573  "single float as the specified point.");
574  }
575  return at({x});
576  }
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T& SSAGES::GridBase< T >::at ( double  x)
inline

Access 1d-Grid by point - read-write.

Parameters
xAccess grid point pertaining to this value.
Returns
Reference to the value pertaining to the specified coordinate.
Note
This function can only be used for 1d Grids.

Definition at line 585 of file GridBase.h.

586  {
587  return const_cast<T&>(static_cast<const GridBase<T>* >(this)->at(x));
588  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T* SSAGES::GridBase< T >::data ( )
inline

Get pointer to the internal data storage vector.

Returns
Pointer to data in the internal storage vector.

It is discouraged to directly access the internal data storage. It might, however be necessary. For example when communicating the data over MPI.

Definition at line 262 of file GridBase.h.

Referenced by SSAGES::Basis::UpdateBias().

263  {
264  return data_.data();
265  }
std::vector< T > data_
Internal storage of the data.
Definition: GridBase.h:42

Here is the caller graph for this function:

template<typename T>
T const* SSAGES::GridBase< T >::data ( ) const
inline

Get pointer to const of the internal data storage vector.

Returns
Const pointer to data in the internal storage vector.

It is discouraged to directly access the internal data storage. It might, however be necessary. For example when communicating data over MPI.

Definition at line 274 of file GridBase.h.

275  {
276  return data_.data();
277  }
std::vector< T > data_
Internal storage of the data.
Definition: GridBase.h:42
template<typename T>
double SSAGES::GridBase< T >::GetCoordinate ( int  index)
inline

Return center point of 1d-grid.

Parameters
indexIndex of the 1d grid.
Returns
Coordinate in real/CV space.
Note
This function is only available for 1d grids.

Definition at line 434 of file GridBase.h.

435  {
436  return GetCoordinates({index}).at(0);
437  }
std::vector< double > GetCoordinates(const std::vector< int > &indices)
Return coordinates of the grid center points.
Definition: GridBase.h:409
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
std::vector<double> SSAGES::GridBase< T >::GetCoordinates ( const std::vector< int > &  indices)
inline

Return coordinates of the grid center points.

Parameters
indicesGrid indices specifying a grid point.
Returns
Vector of double specifying the position of the grid point.

The grid is a discretization of real or cv space. Thus, each grid point is associated with an interval of the underlying space. This function returns the center point of this interval.

Definition at line 409 of file GridBase.h.

Referenced by SSAGES::GridBase< int >::GetCoordinate(), and SSAGES::GridBase< int >::GetInterpolated().

410  {
411  if (indices.size() != dimension_) {
412  throw std::invalid_argument(
413  "Grid indices specified for GetCoordinates() do not have "
414  "the same dimensionality as the grid.");
415  }
416 
417  std::vector<double> v(dimension_);
418 
419  for (size_t i = 0; i < dimension_; ++i) {
420  double spacing = (GetUpper(i) - GetLower(i)) / GetNumPoints(i);
421  v.at(i) = GetLower(i) + (indices[i] + 0.5)*spacing;
422  }
423 
424  return v;
425  }
const std::vector< double > & GetLower() const
Return the lower edges of the Grid.
Definition: GridBase.h:165
const std::vector< double > & GetUpper() const
Return the upper edges of the Grid.
Definition: GridBase.h:191
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45
const std::vector< int > & GetNumPoints() const
Get the number of points for all dimensions.
Definition: GridBase.h:138

Here is the caller graph for this function:

template<typename T>
size_t SSAGES::GridBase< T >::GetDimension ( ) const
inline

Get the dimension.

Returns
Dimensionality of the grid.

Definition at line 128 of file GridBase.h.

Referenced by SSAGES::GridBase< int >::at(), SSAGES::GridBase< int >::GetLower(), SSAGES::GridBase< int >::GetNumPoints(), SSAGES::GridBase< int >::GetPeriodic(), and SSAGES::GridBase< int >::GetUpper().

129  {
130  return dimension_;
131  }
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45

Here is the caller graph for this function:

template<typename T>
int SSAGES::GridBase< T >::GetIndex ( double  x) const
inline

Return the Grid index for a one-dimensional grid.

Parameters
xPoint in space.
Returns
Grid index to which the point pertains.

Return the Grid index pertaining to the given point in space. This function is for convenience when accessing 1d-Grids. For higher-dimensional grids, x needs to be a vector of doubles.

Definition at line 391 of file GridBase.h.

392  {
393  if (dimension_ != 1) {
394  throw std::invalid_argument("1d Grid index can only be accessed for "
395  "1d-Grids can be accessed with a.");
396  }
397  return GetIndices({x}).at(0);
398  }
std::vector< int > GetIndices(const std::vector< double > &x) const
Return the Grid indices for a given point.
Definition: GridBase.h:296
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
std::vector<int> SSAGES::GridBase< T >::GetIndices ( const std::vector< double > &  x) const
inline

Return the Grid indices for a given point.

Parameters
xPoint in space.
Returns
Indices of the grid point to which the point in space pertains.

The grid discretizes the continuous space. For a given point in this continuous space, this function will return the indices of the grid point covering the point in space.

If the grid is non-periodic in a given dimension and x is lower than the lower edge in this dimension, the function will return -1, the index of the underflow bin. Similarly, it will return numPoints, the index of the overflow bin, if x is larger than the upper edge.

In periodic dimensions, the index can take any integer value and will be wrapped to the interval [0, numPoints).

Definition at line 296 of file GridBase.h.

Referenced by SSAGES::GridBase< int >::at(), SSAGES::GridBase< int >::GetIndex(), and SSAGES::GridBase< int >::GetInterpolated().

297  {
298  // Check that input vector has the correct dimensionality
299  if (x.size() != dimension_) {
300  throw std::invalid_argument("Specified point has a larger "
301  "dimensionality than the grid.");
302  }
303 
304  std::vector<int> indices(dimension_);
305  for (size_t i = 0; i < dimension_; ++i) {
306  double xpos = x.at(i);
307  if (!GetPeriodic(i)) {
308  if (xpos < GetLower(i)) {
309  indices.at(i) = -1;
310  continue;
311  } else if (xpos > GetUpper(i)) {
312  indices.at(i) = GetNumPoints(i);
313  continue;
314  }
315  }
316 
317  // To make sure, the value is rounded in the correct direction.
318  double spacing = (GetUpper(i) - GetLower(i)) / GetNumPoints(i);
319  indices.at(i) = std::floor( (xpos - GetLower(i)) / spacing);
320  }
321 
322  return wrapIndices(indices);
323  }
const std::vector< bool > & GetPeriodic() const
Return the periodicity of the Grid.
Definition: GridBase.h:219
std::vector< int > wrapIndices(const std::vector< int > &indices) const
Wrap the index around periodic boundaries.
Definition: GridBase.h:57
const std::vector< double > & GetLower() const
Return the lower edges of the Grid.
Definition: GridBase.h:165
const std::vector< double > & GetUpper() const
Return the upper edges of the Grid.
Definition: GridBase.h:191
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45
const std::vector< int > & GetNumPoints() const
Get the number of points for all dimensions.
Definition: GridBase.h:138

Here is the caller graph for this function:

template<typename T>
double SSAGES::GridBase< T >::GetInterpolated ( const std::vector< double > &  x)
inline

Return linear interpolation on a coordinate.

Parameters
xPoint in space.
Returns
Linearly interpolated value of grid at that point. This function performs a n-linear interpolation on the grid. The formula is bilinear/trilinear interpolation generalized to N-dimensional space.

Definition at line 333 of file GridBase.h.

334  {
335  double interpvalue = 0;
336  int tempindex;
337  for (size_t i = 0 ; i < pow(2,dimension_) ; ++i)
338  {
339 
340  std::vector<double> shiftedvector = x;
341  tempindex = i;
342 
343  double accumulatedweight = 1;
344  for(size_t j = 0; j < dimension_ ; ++j)
345  {
346  double spacing = (GetUpper(j) - GetLower(j)) / GetNumPoints(j);
347  shiftedvector[j] += ((tempindex%2)-0.5)*spacing;
348  tempindex = tempindex/2;
349  }
350 
351  std::vector<int> shiftedindices = GetIndices(shiftedvector);
352  std::vector<double> shiftedcenters = GetCoordinates(shiftedindices);
353 
354  for(size_t j = 0; j < dimension_ ; ++j)
355  {
356  double spacing = (GetUpper(j) - GetLower(j)) / GetNumPoints(j);
357  //Handle Edges
358  if(shiftedindices[j] == -1)
359  {
360  accumulatedweight *= ((std::abs(x[j]-GetCoordinates(GetIndices(x))[j])/spacing));
361  shiftedvector[j] += spacing/2;
362  }
363  else if (shiftedindices[j] == GetNumPoints(j))
364  {
365  accumulatedweight *= ((std::abs(x[j]-GetCoordinates(GetIndices(x))[j])/spacing));
366  shiftedvector[j] -= spacing/2;
367  }
368  else
369  {
370  // Handle Periodicity
371  if(std::abs(x[j]-shiftedcenters[j]) > spacing)
372  accumulatedweight *= (1-(std::abs(std::abs(x[j]-shiftedcenters[j]) - (GetUpper(j) - GetLower(j)))/spacing));
373  else
374  accumulatedweight *= (1-(std::abs(x[j]-shiftedcenters[j])/spacing));
375  }
376  }
377  interpvalue += accumulatedweight*at(shiftedvector);
378  }
379  return interpvalue;
380  }
std::vector< double > GetCoordinates(const std::vector< int > &indices)
Return coordinates of the grid center points.
Definition: GridBase.h:409
std::vector< int > GetIndices(const std::vector< double > &x) const
Return the Grid indices for a given point.
Definition: GridBase.h:296
const std::vector< double > & GetLower() const
Return the lower edges of the Grid.
Definition: GridBase.h:165
const std::vector< double > & GetUpper() const
Return the upper edges of the Grid.
Definition: GridBase.h:191
size_t dimension_
Dimension of the grid.
Definition: GridBase.h:45
const std::vector< int > & GetNumPoints() const
Get the number of points for all dimensions.
Definition: GridBase.h:138
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
const std::vector<double>& SSAGES::GridBase< T >::GetLower ( ) const
inline

Return the lower edges of the Grid.

Returns
Vector containing the lower edges of the grid.

Definition at line 165 of file GridBase.h.

Referenced by SSAGES::GridBase< int >::GetCoordinates(), SSAGES::GridBase< int >::GetIndices(), SSAGES::GridBase< int >::GetInterpolated(), and SSAGES::GridBase< int >::GetLower().

166  {
167  return edges_.first;
168  }
std::pair< std::vector< double >, std::vector< double > > edges_
Edges of the Grid in each dimension.
Definition: GridBase.h:51

Here is the caller graph for this function:

template<typename T>
double SSAGES::GridBase< T >::GetLower ( size_t  dim) const
inline

Get the lower edge for a specific dimension.

Parameters
dimIndex of the dimension.
Returns
Value of the lower edge in the requested dimension.
Note
The first dimension has the index 0.

Definition at line 177 of file GridBase.h.

178  {
179  if (dim >= GetDimension()) {
180  std::cerr << "Warning! Lower edge requested for a dimension larger "
181  "than the grid dimensionality!\n";
182  return 0.0;
183  }
184  return GetLower().at(dim);
185  }
const std::vector< double > & GetLower() const
Return the lower edges of the Grid.
Definition: GridBase.h:165
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:128
template<typename T>
const std::vector<int>& SSAGES::GridBase< T >::GetNumPoints ( ) const
inline

Get the number of points for all dimensions.

Returns
Vector of ints containing the number of grid points for each dimension.

Definition at line 138 of file GridBase.h.

Referenced by SSAGES::Grid< int >::end(), SSAGES::Histogram< int >::end(), SSAGES::GridBase< int >::GetCoordinates(), SSAGES::GridBase< int >::GetIndices(), SSAGES::GridBase< int >::GetInterpolated(), SSAGES::Grid< int >::Grid(), SSAGES::Histogram< int >::Histogram(), SSAGES::Grid< int >::mapTo1d(), SSAGES::Histogram< int >::mapTo1d(), and SSAGES::GridBase< int >::wrapIndices().

139  {
140  return numPoints_;
141  }
std::vector< int > numPoints_
Number of points in each dimension.
Definition: GridBase.h:48

Here is the caller graph for this function:

template<typename T>
int SSAGES::GridBase< T >::GetNumPoints ( size_t  dim) const
inline

Get the number of points for a specific dimension.

Parameters
dimIndex of the dimension.
Returns
Number of grid points in the requested dimension.
Note
The first dimension uses the index 0.

Definition at line 150 of file GridBase.h.

151  {
152  if (dim >= GetDimension()) {
153  std::cerr << "Warning! Grid size requested for a dimension larger "
154  "than the grid dimensionality!\n";
155  return 0;
156  }
157 
158  return numPoints_.at(dim);
159  }
std::vector< int > numPoints_
Number of points in each dimension.
Definition: GridBase.h:48
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:128
template<typename T>
const std::vector<bool>& SSAGES::GridBase< T >::GetPeriodic ( ) const
inline

Return the periodicity of the Grid.

Returns
Vector of bools. The values are True (False ) if the grid is periodic (non-periodic) in the given dimension.

Definition at line 219 of file GridBase.h.

Referenced by SSAGES::GridBase< int >::GetIndices(), SSAGES::GridBase< int >::GetPeriodic(), SSAGES::Histogram< int >::mapTo1d(), and SSAGES::GridBase< int >::wrapIndices().

220  {
221  return isPeriodic_;
222  }
std::vector< bool > isPeriodic_
Periodicity of the Grid.
Definition: GridBase.h:54

Here is the caller graph for this function:

template<typename T>
bool SSAGES::GridBase< T >::GetPeriodic ( size_t  dim) const
inline

Get the periodicity in a specific dimension.

Parameters
dimIndex of the dimension.
Returns
True (False ) if the grid is periodic (non-periodic) in the specified dimension.
Note
The dimensions are indexed starting with 0.

Definition at line 232 of file GridBase.h.

233  {
234  if (dim >= GetDimension()) {
235  std::cerr << "Warning! Periodicity requested for a dimension larger "
236  "than the grid dimensionality!\n";
237  return false;
238  }
239  return GetPeriodic().at(dim);
240  }
const std::vector< bool > & GetPeriodic() const
Return the periodicity of the Grid.
Definition: GridBase.h:219
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:128
template<typename T>
const std::vector<double>& SSAGES::GridBase< T >::GetUpper ( ) const
inline

Return the upper edges of the Grid.

Returns
Vector containing the upper edges of the grid.

Definition at line 191 of file GridBase.h.

Referenced by SSAGES::GridBase< int >::GetCoordinates(), SSAGES::GridBase< int >::GetIndices(), SSAGES::GridBase< int >::GetInterpolated(), and SSAGES::GridBase< int >::GetUpper().

192  {
193  return edges_.second;
194  }
std::pair< std::vector< double >, std::vector< double > > edges_
Edges of the Grid in each dimension.
Definition: GridBase.h:51

Here is the caller graph for this function:

template<typename T>
double SSAGES::GridBase< T >::GetUpper ( size_t  dim) const
inline

Get the upper edge for a specific dimension.

Parameters
dimIndex of the dimension.
Returns
Value of the upper edge in the given dimension.
Note
The dimensions are indexed starting with 0.

Definition at line 204 of file GridBase.h.

205  {
206  if (dim >= GetDimension()) {
207  std::cerr << "Warning! Upper edge requested for a dimension larger "
208  "than the grid dimensionality!\n";
209  return 0.0;
210  }
211  return GetUpper().at(dim);
212  }
const std::vector< double > & GetUpper() const
Return the upper edges of the Grid.
Definition: GridBase.h:191
size_t GetDimension() const
Get the dimension.
Definition: GridBase.h:128
template<typename T>
virtual size_t SSAGES::GridBase< T >::mapTo1d ( const std::vector< int > &  indices) const
protectedpure virtual

This function needs to be implemented by child classes.

Parameters
indicesThe indices specifying the grid point.
Returns
Index of the grid point in the 1d storage vector.

mapTo1d maps the indices onto a single index to access the data stored in the data_ vector. This mapping will be different for different grid implementations.

Implemented in SSAGES::Histogram< T >, SSAGES::Histogram< int >, SSAGES::Grid< T >, SSAGES::Grid< double >, SSAGES::Grid< Vector >, and SSAGES::Grid< int >.

Referenced by SSAGES::GridBase< int >::at().

Here is the caller graph for this function:

template<typename T>
const T& SSAGES::GridBase< T >::operator[] ( const std::vector< int > &  indices) const
inline

Access Grid element per [] read-only.

Parameters
indicesVector of integers specifying the grid point.
Returns
Const reference to value to the given grid point.

Definition at line 595 of file GridBase.h.

596  {
597  return at(indices);
598  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T& SSAGES::GridBase< T >::operator[] ( const std::vector< int > &  indices)
inline

Access Grid element per [] read-write.

Parameters
indicesVector of integers specifying the grid point.
Returns
Reference to value at the specified grid point.

Definition at line 605 of file GridBase.h.

606  {
607  return at(indices);
608  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
template<typename R >
const T& SSAGES::GridBase< T >::operator[] ( std::initializer_list< R > &&  x) const
inline

Const access of Grid element via initializer list.

Template Parameters
RDatatype in the initializer list
Parameters
xinitializer list
Returns
Const reference to the value at the specified point.

This function avoids abiguity if operator[] is called with a brace-enclosed initializer list.

Example: grid[{0,1}] or grid[{-1.23, 4.2, 0.0}]

Definition at line 622 of file GridBase.h.

623  {
624  return at(static_cast<std::vector<R> >(x));
625  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
template<typename R >
T& SSAGES::GridBase< T >::operator[] ( std::initializer_list< R > &&  x)
inline

Access Grid element via initializer list.

Template Parameters
RDatatype in the initializer list
Parameters
xinitializer list
Returns
Reference to the value at the specified point.

This function avoids abiguity if operator[] is called with a brace-enclosed initializer list.

Example: grid[{0,1}]

Definition at line 639 of file GridBase.h.

640  {
641  return at(static_cast<std::vector<R> >(x));
642  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
const T& SSAGES::GridBase< T >::operator[] ( int  index) const
inline

Access 1d-Grid per [] operator, read-only.

Parameters
indexIndex of the grid point.
Returns
Const reference to the value at the given grid point.
Note
This function can only be used for 1d grids.

Definition at line 651 of file GridBase.h.

652  {
653  return at(index);
654  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T& SSAGES::GridBase< T >::operator[] ( int  index)
inline

Access 1d-Grid per [] operator, read-write.

Parameters
indexIndex of the grid point.
Returns
Reference to the value at the given grid point.
Note
This function can only be used for 1d grids.

Definition at line 663 of file GridBase.h.

664  {
665  return at(index);
666  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
const T& SSAGES::GridBase< T >::operator[] ( const std::vector< double > &  x) const
inline

Access Grid element pertaining to a specific point per [] read-only.

Parameters
xVector of doubles specifying the point in space.
Returns
Const reference to the value pertaining to the given coordinates.

Definition at line 673 of file GridBase.h.

674  {
675  return at(x);
676  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T& SSAGES::GridBase< T >::operator[] ( const std::vector< double > &  x)
inline

Access Grid element pertaining to a specific point per [] read-write.

Parameters
xVector of doubles specifying the point in space.
Returns
Reference to the value pertaining to the given coordinates.

Definition at line 683 of file GridBase.h.

684  {
685  return at(x);
686  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
const T& SSAGES::GridBase< T >::operator[] ( double  x) const
inline

Access 1d-Grid via specific point, read-only.

Parameters
xPoint specifying the desired Grid point.
Returns
Const reference to the value at the specified coordinate.
Note
This function can only be used for 1d grids.

Definition at line 695 of file GridBase.h.

696  {
697  return at(x);
698  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
T& SSAGES::GridBase< T >::operator[] ( double  x)
inline

Access 1d-Grid via specific point, read-write.

Parameters
xPoint specifying the desired Grid point.
Returns
Reference to value at the specified coordinate.
Note
This function can only be used for 1d grids.

Definition at line 707 of file GridBase.h.

708  {
709  return at(x);
710  }
const T & at(const std::vector< int > &indices) const
Access Grid element read-only.
Definition: GridBase.h:452
template<typename T>
size_t SSAGES::GridBase< T >::size ( ) const
inline

Get the size of the internal storage vector.

Returns
Size of the internal storage vector.

This function returns the size of the internal storage vector. This is also the total number of grid points including the over/underflow bins in case of a histogram.

Definition at line 250 of file GridBase.h.

251  {
252  return data_.size();
253  }
std::vector< T > data_
Internal storage of the data.
Definition: GridBase.h:42

The documentation for this class was generated from the following file: