source: src/LinkedCell/LinkedCell_View.cpp@ 2390e6

/*
 * Project: MoleCuilder
 * Description: creates and alters molecular systems
 * Copyright (C)  2012 University of Bonn. All rights reserved.
 * 
 *
 *   This file is part of MoleCuilder.
 *
 *    MoleCuilder 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.
 *
 *    MoleCuilder 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 MoleCuilder.  If not, see <http://www.gnu.org/licenses/>.
 */

/*
 * LinkedCell_View.cpp
 *
 *  Created on: Nov 15, 2011
 *      Author: heber
 */

// include config.h
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif

//#include "CodePatterns/MemDebug.hpp"

#include "LinkedCell_View.hpp"

#include <algorithm>

#include "CodePatterns/Assert.hpp"
#include "CodePatterns/Log.hpp"

#include "Atom/TesselPoint.hpp"
#include "Box.hpp"
#include "LinearAlgebra/Vector.hpp"
#include "LinkedCell.hpp"
#include "LinkedCell_Model.hpp"
#include "LinkedCell_View_ModelWrapper.hpp"
#include "tripleIndex.hpp"

using namespace LinkedCell;

// static multimap instance
LinkedCell_View::ModelInstanceMap LinkedCell_View::RAIIMap;

/** Constructor of class LinkedCell_View.
 *
 * We add ourselves to RAIIMap here.
 */
LinkedCell_View::LinkedCell_View(const LinkedCell_Model &_LC) :
    LC( new LinkedCell_View_ModelWrapper(&_LC) )
{
  // add us to RAII counting map
#ifndef NDEBUG
  std::pair< ModelInstanceMap::iterator, bool> inserter =
#endif
      RAIIMap.insert( this );
  ASSERT( inserter.second,
      "LinkedCell_View::LinkedCell_View() - we "+toString(this)+" are already present in RAIIMap.");
  LOG(3, "INFO: Placing instance "+toString(this)+" into RAIIMap.");
}

/** Copy Constructor of class LinkedCell_View.
 *
 * We add ourselves to RAIIMap here.
 */
LinkedCell_View::LinkedCell_View(const LinkedCell_View &_view) :
    LC( new LinkedCell_View_ModelWrapper(_view.LC->getModel()) )
{
  if (this != &_view) {
    // add us to RAII counting map
  #ifndef NDEBUG
    std::pair< ModelInstanceMap::iterator, bool> inserter =
  #endif
        RAIIMap.insert( this );
    ASSERT( inserter.second,
        "LinkedCell_View::LinkedCell_View(&view) - we "+toString(this)+" are already present in RAIIMap.");
    LOG(3, "INFO: Placing instance "+toString(this)+" copied from "+toString(&_view)+" into RAIIMap.");
  }
}

/** Destructor of class LinkedCell_View.
 *
 * We remove ourselves from RAIIMap here.
 */
LinkedCell_View::~LinkedCell_View()
{
  ModelInstanceMap::iterator iter = RAIIMap.find(this);
  ASSERT( iter != RAIIMap.end(),
      "LinkedCell_View::~LinkedCell_View() - there is no instance "
      +toString(this)+" in RAIIMap.");
  if (iter != RAIIMap.end()) {
    RAIIMap.erase(iter);
    LOG(3, "INFO: Removing instance "+toString(this)+" from RAIIMap.");
  } else {
    ELOG(1, "Failed to remove instance "+toString(this)+" from RAIIMap.");
  }
  delete LC;
}

/** Return at least as many points as are inside a sphere of \a radius around \a center.
 *
 * \sa LinkedCell_View::getPointsInsideSphere()
 *
 * @param radius radius of sphere
 * @param center center of sphere
 * @return a list containing at least all points inside described sphere
 */
LinkedList LinkedCell_View::getAllNeighbors(const double radius, const Vector &center) const
{
  LinkedList TesselList;  // we do not need a set, as nodes are uniquely associated to a cell.

  const LinkedCell_Model * const LCmodel = LC->getModel();  // get quick ref to model
  // get relative bounds
  const tripleIndex step = LCmodel->getStep(radius);
  const tripleIndex index = LCmodel->getIndexToVector(center);
  LinkedCell_Model::LinkedCellNeighborhoodBounds neighbors =
      LCmodel->getNeighborhoodBounds(index, step);

  tripleIndex n;
  for (n[0] = 0; n[0] < neighbors.second[0]; n[0]++)
    for (n[1] = 0; n[1] < neighbors.second[1]; n[1]++)
      for (n[2] = 0; n[2] < neighbors.second[2]; n[2]++) {
        tripleIndex absolute_n = neighbors.first + n;
        if (!LCmodel->checkArrayBounds(absolute_n))
          LCmodel->applyBoundaryConditions(absolute_n);
        const LinkedCell &List = LCmodel->getCell(absolute_n);
        LOG(4, "INFO: Current cell is " << neighbors.first << " plus " << n << ", yielding " << absolute_n << ".");
        for (LinkedCell::const_iterator Runner = List.begin(); Runner != List.end(); Runner++) {
          TesselList.push_back(*Runner);
          LOG(5, "INFO: Inserting " << **Runner);
        }
      }
  return TesselList;
}

LinkedList LinkedCell_View::getPointsInsideSphere(const double radius, const Vector &center) const
{
  // get overly much points
  const LinkedList TesselList = getAllNeighbors(radius, center);
  LinkedList ReturnList;

  // remove all unnecessary ones
  const Box& domain = LC->getModel()->getDomain();
  for (LinkedList::const_iterator iter = TesselList.begin(); iter != TesselList.end(); ++iter) {
    if (domain.periodicDistanceSquared(center, (*iter)->getPosition()) <= radius*radius)
      ReturnList.push_back(*iter);
  }

  return ReturnList;
}