1 // -*- C++ -*-
2 //
3 // KimNeighborLocator.cpp: Common base class for KIM interface neighbor locators.
4 //
5 // Copyright (C) 2012-2013 Jakob Schiotz and the Department of Physics,
6 // Technical University of Denmark. Email: schiotz@fysik.dtu.dk
7 //
8 // This file is part of Asap version 3.
9 // Asap is released under the GNU Lesser Public License (LGPL) version 3.
10 // However, the parts of Asap distributed within the OpenKIM project
11 // (including this file) are also released under the Common Development
12 // and Distribution License (CDDL) version 1.0.
13 //
14 // This program is free software: you can redistribute it and/or
15 // modify it under the terms of the GNU Lesser General Public License
16 // version 3 as published by the Free Software Foundation. Permission
17 // to use other versions of the GNU Lesser General Public License may
18 // granted by Jakob Schiotz or the head of department of the
19 // Department of Physics, Technical University of Denmark, as
20 // described in section 14 of the GNU General Public License.
21 //
22 // This program is distributed in the hope that it will be useful,
23 // but WITHOUT ANY WARRANTY; without even the implied warranty of
24 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25 // GNU General Public License for more details.
26 //
27 // You should have received a copy of the GNU General Public License
28 // and the GNU Lesser Public License along with this program. If not,
29 // see <http://www.gnu.org/licenses/>.
30
31 #include "KimNeighborLocator.h"
32 #include "KIM_ModelHeaders.hpp"
33 #include "Asap.h"
34 #include "Debug.h"
35
36 namespace ASAPSPACE {
37
PyAsap_NewKimNeighborLocator(KimAtoms * atoms,double rCut)38 PyAsap_NeighborLocatorObject *PyAsap_NewKimNeighborLocator(KimAtoms *atoms,
39 double rCut)
40 {
41 PyAsap_NeighborLocatorObject *self;
42
43 self = (PyAsap_NeighborLocatorObject*) malloc(sizeof(PyAsap_NeighborLocatorObject));
44 if (self == NULL)
45 throw AsapError("malloc failed.");
46
47 self->ob_refcnt = 1;
48 self->weakrefs = NULL;
49 self->fulllist = false;
50 self->cobj = new KimNeighborLocator(atoms, rCut);
51 if (self->cobj == NULL)
52 {
53 CHECKREF(self);
54 Py_DECREF(self);
55 throw AsapError("Failed to create a new NeighborList object.");
56 }
57 return self;
58 }
59
60 } // end namespace
61
62
63
KimNeighborLocator(KimAtoms * atoms,double rCut)64 KimNeighborLocator::KimNeighborLocator(KimAtoms *atoms, double rCut)
65 {
66 CONSTRUCTOR;
67 this->atoms = atoms;
68 AsapAtoms_INCREF(atoms);
69 nAtoms = nGhosts = 0;
70 this->rcut = rCut;
71 rcut2 = rCut*rCut;
72 }
73
~KimNeighborLocator()74 KimNeighborLocator::~KimNeighborLocator()
75 {
76 DESTRUCTOR;
77 AsapAtoms_DECREF(atoms);
78 }
79
CheckNeighborList()80 bool KimNeighborLocator::CheckNeighborList()
81 {
82 bool result = (nAtoms != atoms->GetNumberOfAtoms());
83 UpdateNeighborList();
84 nAtoms = atoms->GetNumberOfAtoms();
85 nGhosts = atoms->GetNumberOfAtoms();
86 return result;
87 }
88
GetFullNeighbors(int n,int * neighbors,Vec * diffs,double * diffs2,int & size,double r) const89 int KimNeighborLocator::GetFullNeighbors(int n, int *neighbors, Vec *diffs, double *diffs2,
90 int& size, double r /* = -1.0 */ ) const
91 {
92 int numberOfNeighbors;
93 const int *rawneighbors;
94 KIM::ModelComputeArguments const * const modelComputeArguments = atoms->GetModelComputeArgumentsObject();
95 assert(modelComputeArguments != NULL);
96 int error = modelComputeArguments->GetNeighborList(0, n, &numberOfNeighbors, &rawneighbors);
97 if (error)
98 throw AsapError("modelComputeArguments->GetNeighborLists failed ") << __FILE__ << ":" << __LINE__;
99 // Now construct the list of distance vectors
100 const Vec *pos = atoms->GetPositions();
101 const Vec *thispos = pos + n;
102 int numnb = 0;
103 double rcut2 = this->rcut2;
104 if (r > 0)
105 rcut2 = r*r;
106 for (int i = 0; i < numberOfNeighbors; i++)
107 {
108 const Vec *otherpos = pos + rawneighbors[i];
109 diffs[numnb] = *otherpos - *thispos;
110 diffs2[numnb] = diffs[numnb] * diffs[numnb];
111 if (diffs2[numnb] <= rcut2)
112 {
113 neighbors[numnb] = rawneighbors[i];
114 numnb++;
115 }
116 }
117 assert(numnb <= size);
118 size -= numnb;
119 return numnb;
120 }
121