ddMd/neighbor/CellList.h

#ifndef DDMD_CELL_LIST_H
#define DDMD_CELL_LIST_H

/*
* Simpatico - Simulation Package for Polymeric and Molecular Liquids
*
* Copyright 2010 - 2017, The Regents of the University of Minnesota
* Distributed under the terms of the GNU General Public License.
*/

#include "Cell.h"
#include <ddMd/chemistry/Atom.h>
#include <simp/boundary/Boundary.h>
#include <util/space/Grid.h>
#include <util/containers/DArray.h>
#include <util/containers/GArray.h>
#include <util/global.h>

namespace DdMd
{

   using namespace Util;
   using namespace Simp;

   class CellList
   {

   public:

      // Public methods

      CellList();

      virtual ~CellList();

      void allocate(int atomCapacity, const Vector& lower, const Vector& upper, 
                    const Vector& cutoffs, int nCellCut = 1);

      void allocate(int atomCapacity, const Vector& lower, const Vector& upper, 
                    double cutoff, int nCellCut = 1);

      void 
      makeGrid(const Vector& lower, const Vector& upper, const Vector& cutoffs, 
               int nCellCut = 1);

      void placeAtom(Atom &atom);

      void build();

      void update();

      void clear();

      const Grid& grid() const;

      double cellLength(int i) const;

      int cellIndexFromPosition(const Vector& position) const;

      const Cell* begin() const;

      const Cell& cell(int i) const;

      int nAtom() const;

      int nReject() const;

      #ifdef UTIL_DEBUG

      int maxNAtomCell() const;
      #endif

      int atomCapacity() const;

      int cellCapacity() const;

      bool isAllocated() const;

      bool isBuilt() const;

      bool isValid() const;

   private:

      /*
      * Temporary storage for atom pointers, before copying to cells. 
      */
      struct Tag {
         Atom* ptr;
         int cellRank;
      };

      Cell::OffsetArray offsets_;

      Grid grid_;

      DArray<Tag> tags_;

      DArray<CellAtom> atoms_;

      GArray<Cell> cells_;

      Vector lower_; 

      Vector upper_; 

      Vector cellLengths_; 

      Vector lowerOuter_; 

      Vector upperOuter_; 

      Cell* begin_;

      int nAtom_;

      int nReject_;

      #ifdef UTIL_DEBUG
      int maxNAtomCell_;
      #endif

      bool isBuilt_;

      void setGridDimensions(const Vector& lower, const Vector& upper, 
                             const Vector& cutoffs, int nCellCut);

      bool isValidAtomId(int atomId);

   }; 

   // Public inline method definitions:

   /*
   * Identify the cell for an Atom, based on its position.
   */
   inline int CellList::cellIndexFromPosition(const Vector& position) const
   {
      IntVector r;
      for (int i = 0; i < Dimension; ++i) {
         if (position[i] <= lowerOuter_[i]) {
            return -1;
         }
         if (position[i] >= upperOuter_[i]) {
            return -1;
         }
         r[i] = int( (position[i] - lowerOuter_[i])/ cellLengths_[i] );
         assert(r[i] < grid_.dimension(i));
      }
      return grid_.rank(r);
   }

   /*
   * Add an Atom to the appropriate cell, based on its position.
   */
   inline void CellList::placeAtom(Atom &atom)
   {
      // Preconditon
      assert(nAtom_ < tags_.capacity());

      int rank = cellIndexFromPosition(atom.position());
      if (rank >= 0) {
         tags_[nAtom_].cellRank = rank;
         tags_[nAtom_].ptr = &atom;
         cells_[rank].incrementCapacity();
         ++nAtom_;
      } else {
         ++nReject_;
      }
   }

   /*
   * Return true iff atomId is valid, i.e., if 0 <= 0 < atomCapacity.
   */
   inline bool CellList::isValidAtomId(int atomId)
   { return ( (0 <= atomId) && (atomId < tags_.capacity()) ); }

   /*
   * Return associated Grid object.
   */
   inline const Grid& CellList::grid() const
   {  return grid_; }

   /*
   * Return length of each cell in direction i.
   */
   inline double CellList::cellLength(int i) const
   {  return cellLengths_[i]; }

   /*
   * Return reference to cell number i.
   */
   inline const Cell& CellList::cell(int i) const
   {
      assert(i < cells_.size());  
      return cells_[i]; 
   }

   /*
   * Return pointer to first Cell.
   */
   inline const Cell* CellList::begin() const
   {  return begin_; }

   /*
   * Is this CellList allocated?
   */
   inline bool CellList::isAllocated() const
   {  return (cells_.capacity() > 0); }

}
#endif