/home/kitware/dashboards/avogadro/libavogadro/src/neighborlist.h

/*********************************************************************
  NeighborList - NeighborList class

  Copyright (C) 2009 by Tim Vandermeersch

  This file is part of the Avogadro molecular editor project.
  For more information, see <http://avogadro.openmolecules.net/>

  Avogadro is free software; you can redistribute it and/or modify
  it under the terms of the GNU General Public License as published by
  the Free Software Foundation; either version 2 of the License, or
  (at your option) any later version.

  Avogadro is distributed in the hope that it will be useful,
  but WITHOUT ANY WARRANTY; without even the implied warranty of
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  GNU General Public License for more details.

  You should have received a copy of the GNU General Public License
  along with this program; if not, write to the Free Software
  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
  02110-1301, USA.
***********************************************************************/

#ifndef NEIGHBORLIST_H
#define NEIGHBORLIST_H

#include <avogadro/molecule.h>

#include <Eigen/Core>

namespace Avogadro
{
  class Atom;

  class A_EXPORT NeighborList
  {
    private:
      typedef std::vector<Atom*>::iterator atom_iter;

    public:
      NeighborList(Molecule *mol, double rcut, bool periodic = false, int boxSize = 1);
      NeighborList(const QList<Atom*> &atoms, double rcut, bool periodic = false, int boxSize = 1);

      void update();
      QList<Atom*> nbrs(Atom *atom, bool uniqueOnly = true);
      QList<Atom*> nbrs(const Eigen::Vector3f *pos);
      inline double r2(unsigned int index) const
      {
        return m_r2.at(index);
      }

    private:
      inline unsigned int ghostIndex(int i, int j, int k) const
      {
        i += m_boxSize + 1;
        j += m_boxSize + 1;
        k += m_boxSize + 1;
        return i + j * m_ghostX + k * m_ghostXY;
      }

      inline unsigned int ghostIndex(const Eigen::Vector3i &index) const
      {
        return ghostIndex(index.x(), index.y(), index.z());
      }

      inline unsigned int cellIndex(int i, int j, int k) const
      {
        return i + j * m_dim.x() + k * m_xyDim;
      }

      inline unsigned int cellIndex(const Eigen::Vector3i &index) const
      {
        return index.x() + index.y() * m_dim.x() + index.z() * m_xyDim;
      }


      inline unsigned int cellIndex(const Eigen::Vector3d &pos) const
      {
        return cellIndex( int(floor( (pos.x() - m_min.x()) / m_edgeLength )),
                          int(floor( (pos.y() - m_min.y()) / m_edgeLength )),
                          int(floor( (pos.z() - m_min.z()) / m_edgeLength )) );
      }

      inline Eigen::Vector3i cellIndexes(const Eigen::Vector3d *pos) const
      {
        Eigen::Vector3i index;
        index.x() = int(floor( (pos->x() - m_min.x()) / m_edgeLength ));
        index.y() = int(floor( (pos->y() - m_min.y()) / m_edgeLength ));
        index.z() = int(floor( (pos->z() - m_min.z()) / m_edgeLength ));
        return index;
      }

      inline bool IsOneTwo(unsigned int i, unsigned int j) const
      {
        std::vector<unsigned int>::const_iterator iter;
        for (iter = m_oneTwo.at(i).begin(); iter != m_oneTwo.at(i).end(); ++iter)
          if (*iter == j)
            return true;

        return false;
      }

      inline bool IsOneThree(unsigned int i, unsigned int j) const
      {
        std::vector<unsigned int>::const_iterator iter;
        for (iter = m_oneThree.at(i).begin(); iter != m_oneThree.at(i).end(); ++iter)
          if (*iter == j)
            return true;

        return false;
      }

      void initOneTwo();
      void initCells();
      void updateCells();
      void initOffsetMap();
      void initGhostMap(bool periodic = false);
      bool insideShpere(const Eigen::Vector3i &index);

      QList<Atom*>                        m_atoms;
      double                              m_rcut, m_rcut2;
      double                              m_edgeLength;
      int                                 m_boxSize;
      int                                 m_updateCounter;


      Eigen::Vector3d                     m_min, m_max;
      Eigen::Vector3i                     m_dim;
      int                                 m_xyDim;
      std::vector<std::vector<Atom*> >    m_cells;

      std::vector<Eigen::Vector3i>        m_offsetMap;
      std::vector<Eigen::Vector3i>        m_ghostMap;
      int                                 m_ghostX;
      int                                 m_ghostXY;

      std::vector<double>                 m_r2;

      std::vector<std::vector<unsigned int> > m_oneTwo;
      std::vector<std::vector<unsigned int> > m_oneThree;
  };

} // end namespace OpenBabel


#endif