GridField.h

Go to the documentation of this file.

/* 
   Packing of hard spheres via molecular dynamics
   Developed by Monica Skoge, 2006, Princeton University
   Contact: Aleksandar Donev (adonev@math.princeton.edu) with questions
   This code may be used, modified and distributed freely.
   Please cite:
   
   "Packing Hyperspheres in High-Dimensional Euclidean Spaces"
        M. Skoge, A. Donev, F. H. Stillinger and S. Torquato, 2006
        
   if you use these codes.      
*/


#ifndef GridField_H
#define GridField_H

#include "Vector.h"
namespace algorithms{
namespace ls{
template<int D, class T>
class GridField {
 
 public:
  GridField();
  GridField(const Vector<D, int>&);
  ~GridField();
  T& get(const Vector<D, int>&);
  Vector<D, int> get_size() const;
  void set_size(const Vector<D, int>&);
  void set_size(const int);
  void initialize(const int i);



 private:
  int elements;
  T* f;
  Vector<D, int> size;           // number of grid points for each dimension
  Vector<D, int> offset;
};


// GridField
// ~~~~~~~~~~~~
template<int D, class T>
GridField<D, T>::GridField()
  :  elements(0),f(0)
{
}


// GridField
// ~~~~~~~~~~~~
template<int D, class T>
GridField<D, T>::GridField(const Vector<D, int>& s)
  : f(0)
{
  set_size(s);
}


// ~GridField
// ~~~~~~~~~~~~~
template <int D, class T>
GridField<D, T>::~GridField()
{
  if(f != 0)
    delete[] f;
}


// get_size
// ~~~~~~~~
template<int D, class T>
inline Vector<D, int> GridField<D, T>::get_size() const
{
  return size;
}


// set_size
// ~~~~~~~~
template<int D, class T>
void GridField<D, T>::set_size(const Vector<D, int>& s)
{
  if(f != 0)
    delete[] f;

  size = s;

  elements = 1;
  for(int i=0; i<D; i++) {
    offset[i] = elements;
    elements *= size.x[i];
  }

  f = new T[elements];
}


// set_size
// ~~~~~~~~
template<int D, class T>
void GridField<D, T>::set_size(const int s)
{
  Vector<D, int> square;

  for(int k=0; k<D; k++)
    square[k] = s;

  set_size(square);
}


// get
// ~~~
template<int D, class T>
inline T& GridField<D, T>::get(const Vector<D, int>& pos)
{
  int p=0;
  for(int i=0; i<D; i++)
    p += pos.x[i]*offset[i];

  return f[p];
}


// initialize
// ~~~
template<int D, class T>
void GridField<D, T>::initialize(const int value)
{
  for(int i=0; i<elements; i++)
    f[i] = value;
}

}/* close namespace ls */
}/* close namespace algorithms */

#endif