From 635a87750538fb02ed6fce17a86899a5c84b2f1e Mon Sep 17 00:00:00 2001 From: Daniele Date: Fri, 1 Dec 2023 00:42:19 +0100 Subject: [PATCH] Neighbourlist integer size (#977) * set up test for the NL size problem * updated the test * added a test for 1M atoms * trying to change to size_t the neighbour list * forgot the override on the calculate() function * Squashed commit of the following: commit 1b70cb0928d8fd788cf49639457254058a155173 Author: Daniele Rapetti Date: Tue Oct 10 10:56:35 2023 +0200 Added an exception to catch the memory error commit 19f37da8e3bddfe4b86a7ae678955ba2dc4df1a4 Author: Daniele Rapetti Date: Tue Oct 10 10:26:57 2023 +0200 impreoved readability for NL commit 6c5034127ba8498156760c841bb375628e038019 Author: Daniele Rapetti Date: Tue Oct 10 09:41:14 2023 +0200 test now check that the results for initialization are the same with the "simpler" algortithm commit c469d229a4f09220347c00f89918bd79f2c1b472 Author: Daniele Rapetti Date: Mon Oct 9 16:58:55 2023 +0200 refactored the NL constructors commit 3e8e66ce455a1463838d6c402c9121208f5f83bd Author: Daniele Rapetti Date: Mon Oct 9 16:07:53 2023 +0200 reformatted the test commit 03e5a41844af95082994978f7db99017c55e1fe0 Author: Daniele Rapetti Date: Mon Oct 9 15:12:21 2023 +0200 integer type change and small optimization of the initialization of the NL commit 2d5eb2e8e4f008c5db56b3ccf049a57f3d36bad9 Author: Daniele Rapetti Date: Mon Oct 9 14:09:45 2023 +0200 added a test to control the NL initialization(atomsID) commit 04f762c0add13b17d9e6cea1312a847964396f37 Author: Daniele Rapetti Date: Mon Oct 9 11:45:17 2023 +0200 added a test to control the NL initialization(size) * removed old regtest * changed to one line per thing also int NL.h * Now the test will ignore the stackdump if it is plotted * changed how the "too much memory" for NL on mac works * updated the changelog * added a comment to explain the mac exception * some more explanation on the arbitrary stop * changed the changelog hoping it fast forward to the one in the master branch * changed the changelog hoping it fast forward to the one in the master branch * corrected plumed_merror usage in NL --- .../rt-NeigbourlistInitialization/Makefile | 1 + .../rt-NeigbourlistInitialization/config | 1 + .../rt-NeigbourlistInitialization/main.cpp | 132 +++++++++++++ .../unitTest.reference | 144 ++++++++++++++ .../Makefile | 1 + .../rt-NeigbourlistInitializationError/config | 10 + .../main.cpp | 77 ++++++++ .../unitTest.proc.reference | 6 + src/tools/NeighborList.cpp | 179 ++++++++++++------ src/tools/NeighborList.h | 52 +++-- 10 files changed, 531 insertions(+), 72 deletions(-) create mode 100644 regtest/basic/rt-NeigbourlistInitialization/Makefile create mode 100644 regtest/basic/rt-NeigbourlistInitialization/config create mode 100644 regtest/basic/rt-NeigbourlistInitialization/main.cpp create mode 100644 regtest/basic/rt-NeigbourlistInitialization/unitTest.reference create mode 100644 regtest/basic/rt-NeigbourlistInitializationError/Makefile create mode 100644 regtest/basic/rt-NeigbourlistInitializationError/config create mode 100644 regtest/basic/rt-NeigbourlistInitializationError/main.cpp create mode 100644 regtest/basic/rt-NeigbourlistInitializationError/unitTest.proc.reference diff --git a/regtest/basic/rt-NeigbourlistInitialization/Makefile b/regtest/basic/rt-NeigbourlistInitialization/Makefile new file mode 100644 index 0000000000..3703b27cea --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitialization/Makefile @@ -0,0 +1 @@ +include ../../scripts/test.make diff --git a/regtest/basic/rt-NeigbourlistInitialization/config b/regtest/basic/rt-NeigbourlistInitialization/config new file mode 100644 index 0000000000..df1f95bf3e --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitialization/config @@ -0,0 +1 @@ +type=make diff --git a/regtest/basic/rt-NeigbourlistInitialization/main.cpp b/regtest/basic/rt-NeigbourlistInitialization/main.cpp new file mode 100644 index 0000000000..9d2c253c8d --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitialization/main.cpp @@ -0,0 +1,132 @@ +#include "plumed/tools/AtomNumber.h" +#include "plumed/tools/Communicator.h" +#include "plumed/tools/NeighborList.h" +#include "plumed/tools/Pbc.h" +#include +#include + +using PLMD::AtomNumber; +using PLMD::Communicator; +using PLMD::NeighborList; +using PLMD::Pbc; + +// Testing that the Neigbour list will be intialized with the desired number of +// couples +// We are initializing with distance and stride not set to check the default +// parameters + +#define check(arg) (((arg)) ? "pass\n" : "not pass\n") + +int main(int, char **) { + std::ofstream report("unitTest"); + Pbc pbc{}; + pbc.setBox(PLMD::Tensor({1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0})); + Communicator cm{}; + bool serial = true; + bool do_pbc = false; + for (const size_t nat0 : {100, 500, 1000, 10000}) { + std::vector list0(nat0); + size_t i = 0; + for (auto &an : list0) { + an.setIndex(i); + ++i; + } + { + report << "Single list:\n"; + std::string prepend="["+std::to_string(nat0)+"]"; + size_t expected = ((nat0 - 1) * nat0) / 2; + auto nl = NeighborList(list0, serial, do_pbc, pbc, cm); + + bool expectedcouples = true; + { + size_t cID = 0; + for (size_t i0 = 0; i0 < nat0 && expectedcouples; ++i0) { + for (size_t i1 = i0+1; i1 < nat0 && expectedcouples; ++i1) { + auto couple = nl.getClosePair(cID); + expectedcouples &= couple.first == i0; + expectedcouples &= couple.second == i1; + ++cID; + } + } + } + report << prepend << "Initial number: " + << check(nl.size() == expected); + report << prepend << "getIndexPair(): " + << check(expectedcouples); + report << prepend << "Lastupdate is 0: " + << check(nl.getLastUpdate() == 0); + report << prepend << "Default stride is 0: " + << check(nl.getStride() == 0); + report << "\n"; + } + for (const size_t nat1 : {100, 500, 1000, 10000}) { + + std::vector list1(nat1); + + i = 0; + for (auto &an : list1) { + an.setIndex(i); + ++i; + } + + { + report << "Double list, no pairs:\n"; + std::string prepend="["+std::to_string(nat0) + + ", " + std::to_string(nat1) +"]"; + bool do_pair = false; + size_t expected = nat1 * nat0; + auto nl = NeighborList(list0, list1, serial, do_pair, do_pbc, pbc, cm); + + bool expectedcouples = true; + { + size_t cID = 0; + for (size_t i0 = 0; i0 < nat0 && expectedcouples; ++i0) { + for (size_t i1 = 0; i1 < nat1 && expectedcouples; ++i1) { + auto couple = nl.getClosePair(cID); + //The getIndexPair for non couple input must return this be this + //(cID / nat1); + expectedcouples &= couple.first == i0; + //(cID % nat1 + nat0); + expectedcouples &= couple.second == nat0+i1; + ++cID; + } + } + } + report << prepend << "Initial number: " + << check(nl.size() == expected); + report << prepend << "getIndexPair(): " + << check(expectedcouples); + report << prepend << "Lastupdate is 0: " + << check(nl.getLastUpdate() == 0); + report << prepend << "Default stride is 0: " + << check(nl.getStride() == 0); + report << "\n"; + } + + if (nat1 == nat0) { + report << "Double list, with pairs:\n"; + std::string prepend="["+std::to_string(nat0) + + ", " + std::to_string(nat1) +"]"; + bool do_pair = true; + size_t expected = nat0; + auto nl = NeighborList(list0, list1, serial, do_pair, do_pbc, pbc, cm); + + bool expectedcouples = true; + for (size_t cID = 0; cID < nat0 && expectedcouples; ++cID) { + auto couple = nl.getClosePair(cID); + expectedcouples &= couple.first == cID; + expectedcouples &= couple.second == cID + nat0; + } + report << prepend << "Initial number: " + << check(nl.size() == expected); + report << prepend << "getIndexPair(): " + << check(expectedcouples); + report << prepend << "Lastupdate is 0: " + << check(nl.getLastUpdate() == 0); + report << prepend << "Default stride is 0: " + << check(nl.getStride() == 0); + report << "\n"; + } + } + } +} \ No newline at end of file diff --git a/regtest/basic/rt-NeigbourlistInitialization/unitTest.reference b/regtest/basic/rt-NeigbourlistInitialization/unitTest.reference new file mode 100644 index 0000000000..c0660f998b --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitialization/unitTest.reference @@ -0,0 +1,144 @@ +Single list: +[100]Initial number: pass +[100]getIndexPair(): pass +[100]Lastupdate is 0: pass +[100]Default stride is 0: pass + +Double list, no pairs: +[100, 100]Initial number: pass +[100, 100]getIndexPair(): pass +[100, 100]Lastupdate is 0: pass +[100, 100]Default stride is 0: pass + +Double list, with pairs: +[100, 100]Initial number: pass +[100, 100]getIndexPair(): pass +[100, 100]Lastupdate is 0: pass +[100, 100]Default stride is 0: pass + +Double list, no pairs: +[100, 500]Initial number: pass +[100, 500]getIndexPair(): pass +[100, 500]Lastupdate is 0: pass +[100, 500]Default stride is 0: pass + +Double list, no pairs: +[100, 1000]Initial number: pass +[100, 1000]getIndexPair(): pass +[100, 1000]Lastupdate is 0: pass +[100, 1000]Default stride is 0: pass + +Double list, no pairs: +[100, 10000]Initial number: pass +[100, 10000]getIndexPair(): pass +[100, 10000]Lastupdate is 0: pass +[100, 10000]Default stride is 0: pass + +Single list: +[500]Initial number: pass +[500]getIndexPair(): pass +[500]Lastupdate is 0: pass +[500]Default stride is 0: pass + +Double list, no pairs: +[500, 100]Initial number: pass +[500, 100]getIndexPair(): pass +[500, 100]Lastupdate is 0: pass +[500, 100]Default stride is 0: pass + +Double list, no pairs: +[500, 500]Initial number: pass +[500, 500]getIndexPair(): pass +[500, 500]Lastupdate is 0: pass +[500, 500]Default stride is 0: pass + +Double list, with pairs: +[500, 500]Initial number: pass +[500, 500]getIndexPair(): pass +[500, 500]Lastupdate is 0: pass +[500, 500]Default stride is 0: pass + +Double list, no pairs: +[500, 1000]Initial number: pass +[500, 1000]getIndexPair(): pass +[500, 1000]Lastupdate is 0: pass +[500, 1000]Default stride is 0: pass + +Double list, no pairs: +[500, 10000]Initial number: pass +[500, 10000]getIndexPair(): pass +[500, 10000]Lastupdate is 0: pass +[500, 10000]Default stride is 0: pass + +Single list: +[1000]Initial number: pass +[1000]getIndexPair(): pass +[1000]Lastupdate is 0: pass +[1000]Default stride is 0: pass + +Double list, no pairs: +[1000, 100]Initial number: pass +[1000, 100]getIndexPair(): pass +[1000, 100]Lastupdate is 0: pass +[1000, 100]Default stride is 0: pass + +Double list, no pairs: +[1000, 500]Initial number: pass +[1000, 500]getIndexPair(): pass +[1000, 500]Lastupdate is 0: pass +[1000, 500]Default stride is 0: pass + +Double list, no pairs: +[1000, 1000]Initial number: pass +[1000, 1000]getIndexPair(): pass +[1000, 1000]Lastupdate is 0: pass +[1000, 1000]Default stride is 0: pass + +Double list, with pairs: +[1000, 1000]Initial number: pass +[1000, 1000]getIndexPair(): pass +[1000, 1000]Lastupdate is 0: pass +[1000, 1000]Default stride is 0: pass + +Double list, no pairs: +[1000, 10000]Initial number: pass +[1000, 10000]getIndexPair(): pass +[1000, 10000]Lastupdate is 0: pass +[1000, 10000]Default stride is 0: pass + +Single list: +[10000]Initial number: pass +[10000]getIndexPair(): pass +[10000]Lastupdate is 0: pass +[10000]Default stride is 0: pass + +Double list, no pairs: +[10000, 100]Initial number: pass +[10000, 100]getIndexPair(): pass +[10000, 100]Lastupdate is 0: pass +[10000, 100]Default stride is 0: pass + +Double list, no pairs: +[10000, 500]Initial number: pass +[10000, 500]getIndexPair(): pass +[10000, 500]Lastupdate is 0: pass +[10000, 500]Default stride is 0: pass + +Double list, no pairs: +[10000, 1000]Initial number: pass +[10000, 1000]getIndexPair(): pass +[10000, 1000]Lastupdate is 0: pass +[10000, 1000]Default stride is 0: pass + +Double list, no pairs: +[10000, 10000]Initial number: pass +[10000, 10000]getIndexPair(): pass +[10000, 10000]Lastupdate is 0: pass +[10000, 10000]Default stride is 0: pass + +Double list, with pairs: +[10000, 10000]Initial number: pass +[10000, 10000]getIndexPair(): pass +[10000, 10000]Lastupdate is 0: pass +[10000, 10000]Default stride is 0: pass + diff --git a/regtest/basic/rt-NeigbourlistInitializationError/Makefile b/regtest/basic/rt-NeigbourlistInitializationError/Makefile new file mode 100644 index 0000000000..3703b27cea --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitializationError/Makefile @@ -0,0 +1 @@ +include ../../scripts/test.make diff --git a/regtest/basic/rt-NeigbourlistInitializationError/config b/regtest/basic/rt-NeigbourlistInitializationError/config new file mode 100644 index 0000000000..c51990c1ad --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitializationError/config @@ -0,0 +1,10 @@ +type=make + +plumed_regtest_after(){ + #this discards the lines like + #"(tools/NeighborList.cpp:98) void PLMD::NeighborList::initialize()" + # in this way if NeighborList.cpp is moved or modified this test won't + #trigger a (false) error + awk '/(Single|Double|neighbor) list/{print} + /Exception text/{print}' ./unitTest > unitTest.proc +} diff --git a/regtest/basic/rt-NeigbourlistInitializationError/main.cpp b/regtest/basic/rt-NeigbourlistInitializationError/main.cpp new file mode 100644 index 0000000000..4c74f8066f --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitializationError/main.cpp @@ -0,0 +1,77 @@ +#include "plumed/tools/AtomNumber.h" +#include "plumed/tools/Communicator.h" +#include "plumed/tools/NeighborList.h" +#include "plumed/tools/Pbc.h" +#include +#include + +using PLMD::AtomNumber; +using PLMD::Communicator; +using PLMD::NeighborList; +using PLMD::Pbc; + +// Testing that the Neigbour list must throw an explanatory PLMD::Exception +// when asked to allocate too many pairs + + +#define check(arg) (((arg)) ? "pass\n" : "not pass\n") + +int main(int, char **) { + std::ofstream report("unitTest"); + Pbc pbc{}; + pbc.setBox(PLMD::Tensor({1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0})); + Communicator cm{}; + bool serial = true; + bool do_pbc = false; + // nat0 times nat1 shold be ludicrous big (nat0*nat1*8B~=30~40GB) + const size_t nat0=100000; + const size_t nat1 = 90000; + std::vector list0(nat0); + size_t i = 0; + for (auto &an : list0) { + an.setIndex(i); + ++i; + } + { + report << "Single list:\n"; + std::string prepend="["+std::to_string(nat0)+"]"; + try{ + size_t expected = ((nat0 - 1) * nat0) / 2; + auto nl = NeighborList(list0, serial, do_pbc, pbc, cm); + //I need this line to ensure that nl is not optimized away + report << prepend << "Initial number: " + << check(nl.size() == expected); + } catch ( PLMD::Exception & error ){ + report << prepend <<"Exception text: " + << error.what(); + } + report << "\n"; + } + //doing the same thing with two lists + + std::vector list1(nat1); + i = 0; + for (auto &an : list1) { + an.setIndex(i); + ++i; + } + + { + report << "Double list, no pairs:\n"; + std::string prepend="["+std::to_string(nat0) + + ", " + std::to_string(nat1) +"]"; + bool do_pair = false; + size_t expected = nat1 * nat0; + try{ + auto nl = NeighborList(list0, list1, serial, do_pair, do_pbc, pbc, cm); + report << prepend << "Initial number: " + << check(nl.size() == expected); + + } catch( PLMD::Exception & error) { + report << prepend <<"Exception text: " + << error.what(); + } + report << "\n"; + } +} + diff --git a/regtest/basic/rt-NeigbourlistInitializationError/unitTest.proc.reference b/regtest/basic/rt-NeigbourlistInitializationError/unitTest.proc.reference new file mode 100644 index 0000000000..2883750bba --- /dev/null +++ b/regtest/basic/rt-NeigbourlistInitializationError/unitTest.proc.reference @@ -0,0 +1,6 @@ +Single list: +[100000]Exception text: +An error happened while allocating the neighbor list, please decrease the number of atoms used +Double list, no pairs: +[100000, 90000]Exception text: +An error happened while allocating the neighbor list, please decrease the number of atoms used diff --git a/src/tools/NeighborList.cpp b/src/tools/NeighborList.cpp index 20daf47cff..58fba82a6c 100644 --- a/src/tools/NeighborList.cpp +++ b/src/tools/NeighborList.cpp @@ -30,66 +30,114 @@ #include #include +#ifdef __APPLE__ +//we are using getenv to give the user the opporunity of suppressing +//the too many memory killswitch while compiling on mac +#include +#endif //__APPLE__ namespace PLMD { -NeighborList::NeighborList(const std::vector& list0, const std::vector& list1, - const bool& serial, const bool& do_pair, const bool& do_pbc, const Pbc& pbc, Communicator& cm, - const double& distance, const unsigned& stride): reduced(false), - serial_(serial), do_pair_(do_pair), do_pbc_(do_pbc), pbc_(&pbc), comm(cm), - distance_(distance), stride_(stride) -{ -// store full list of atoms needed - fullatomlist_=list0; +NeighborList::NeighborList(const std::vector& list0, + const std::vector& list1, + const bool& serial, + const bool& do_pair, + const bool& do_pbc, + const Pbc& pbc, + Communicator& cm, + const double& distance, + const unsigned& stride) + : reduced(false), + serial_(serial), + do_pair_(do_pair), + do_pbc_(do_pbc), + twolists_(true), + pbc_(&pbc), + comm(cm), + //copy-initialize fullatomlist_ + fullatomlist_(list0), + distance_(distance), + nlist0_(list0.size()), + nlist1_(list1.size()), + stride_(stride) { + // store the rest of the atoms into fullatomlist_ fullatomlist_.insert(fullatomlist_.end(),list1.begin(),list1.end()); - nlist0_=list0.size(); - nlist1_=list1.size(); - twolists_=true; if(!do_pair) { nallpairs_=nlist0_*nlist1_; } else { - plumed_assert(nlist0_==nlist1_) << "when using PAIR option, the two groups should have the same number of elements\n" - << "the groups you specified have size "<& list0, const bool& serial, const bool& do_pbc, - const Pbc& pbc, Communicator& cm, const double& distance, - const unsigned& stride): reduced(false), - serial_(serial), do_pbc_(do_pbc), pbc_(&pbc), comm(cm), - distance_(distance), stride_(stride) { - fullatomlist_=list0; - nlist0_=list0.size(); - twolists_=false; - nallpairs_=nlist0_*(nlist0_-1)/2; +NeighborList::NeighborList(const std::vector& list0, + const bool& serial, const bool& do_pbc, + const Pbc& pbc, + Communicator& cm, + const double& distance, + const unsigned& stride) + : reduced(false), + serial_(serial), + do_pbc_(do_pbc), + twolists_(false), + pbc_(&pbc), + comm(cm), + //copy-initialize fullatomlist_ + fullatomlist_(list0), + distance_(distance), + nlist0_(list0.size()), + nallpairs_(nlist0_*(nlist0_-1)/2), + stride_(stride) { initialize(); - lastupdate_=0; } +NeighborList::~NeighborList()=default; + void NeighborList::initialize() { - neighbors_.clear(); - for(unsigned int i=0; i 1296000000 ) + plumed_merror("An error happened while allocating the neighbor " + "list, please decrease the number of atoms used"); + } +#endif // __APPLE__ + try { + neighbors_.resize(nallpairs_); + } catch (...) { + plumed_error_nested() << "An error happened while allocating the neighbor " + "list, please decrease the number of atoms used"; } + //TODO: test if this is feasible for accelerating the loop + //#pragma omp parallel for default(shared) + for(unsigned int i=0; i& NeighborList::getFullAtomList() { return fullatomlist_; } -std::pair NeighborList::getIndexPair(unsigned ipair) { - std::pair index; - if(twolists_ && do_pair_) { - index=std::pair(ipair,ipair+nlist0_); - } else if (twolists_ && !do_pair_) { - index=std::pair(ipair/nlist1_,ipair%nlist1_+nlist0_); - } else if (!twolists_) { +NeighborList::pairIDs NeighborList::getIndexPair(const unsigned ipair) { + pairIDs index; + if(twolists_ && do_pair_) + index=pairIDs(ipair,ipair+nlist0_); + else if (twolists_ && !do_pair_) + index=pairIDs(ipair/nlist1_,ipair%nlist1_+nlist0_); + else if (!twolists_) { unsigned ii = nallpairs_-1-ipair; unsigned K = unsigned(std::floor((std::sqrt(double(8*ii+1))+1)/2)); unsigned jj = ii-K*(K-1)/2; - index=std::pair(nlist0_-1-K,nlist0_-1-jj); + index=pairIDs(nlist0_-1-K,nlist0_-1-jj); } return index; } @@ -115,7 +163,7 @@ void NeighborList::update(const std::vector& positions) { std::vector private_flat_nl; #pragma omp for nowait for(unsigned int i=rank; i index=getIndexPair(i); + pairIDs index=getIndexPair(i); unsigned index0=index.first; unsigned index1=index.second; Vector distance; @@ -131,15 +179,21 @@ void NeighborList::update(const std::vector& positions) { } } #pragma omp critical - local_flat_nl.insert(local_flat_nl.end(), private_flat_nl.begin(), private_flat_nl.end()); + local_flat_nl.insert(local_flat_nl.end(), + private_flat_nl.begin(), + private_flat_nl.end()); } // find total dimension of neighborlist std::vector local_nl_size(stride, 0); local_nl_size[rank] = local_flat_nl.size(); - if(!serial_) comm.Sum(&local_nl_size[0], stride); + if(!serial_) + comm.Sum(&local_nl_size[0], stride); int tot_size = std::accumulate(local_nl_size.begin(), local_nl_size.end(), 0); - if(tot_size==0) {setRequestList(); return;} + if(tot_size==0) { + setRequestList(); + return; + } // merge std::vector merge_nl(tot_size, 0); // calculate vector of displacement @@ -151,8 +205,14 @@ void NeighborList::update(const std::vector& positions) { disp[i+1] = rank_size; } // Allgather neighbor list - if(comm.initialized()&&!serial_) comm.Allgatherv((!local_flat_nl.empty()?&local_flat_nl[0]:NULL), local_nl_size[rank], &merge_nl[0], &local_nl_size[0], &disp[0]); - else merge_nl = local_flat_nl; + if(comm.initialized()&&!serial_) { + comm.Allgatherv((!local_flat_nl.empty()?&local_flat_nl[0]:NULL), + local_nl_size[rank], + &merge_nl[0], + &local_nl_size[0], + &disp[0]); + } else + merge_nl = local_flat_nl; // resize neighbor stuff neighbors_.resize(tot_size/2); for(unsigned i=0; i& NeighborList::getReducedAtomList() { - if(!reduced)for(unsigned int i=0; i(newindex0,newindex1); + auto p = std::find(requestlist_.begin(), requestlist_.end(), index0); + plumed_dbg_assert(p!=requestlist_.end()); + unsigned newindex0=p-requestlist_.begin(); + p = std::find(requestlist_.begin(), requestlist_.end(), index1); + plumed_dbg_assert(p!=requestlist_.end()); + unsigned newindex1=p-requestlist_.begin(); + neighbors_[i]=pairIDs(newindex0,newindex1); } + } reduced=true; return requestlist_; } @@ -195,7 +260,7 @@ unsigned NeighborList::getLastUpdate() const { return lastupdate_; } -void NeighborList::setLastUpdate(unsigned step) { +void NeighborList::setLastUpdate(const unsigned step) { lastupdate_=step; } @@ -203,23 +268,27 @@ unsigned NeighborList::size() const { return neighbors_.size(); } -std::pair NeighborList::getClosePair(unsigned i) const { +NeighborList::pairIDs NeighborList::getClosePair(const unsigned i) const { return neighbors_[i]; } -std::pair NeighborList::getClosePairAtomNumber(unsigned i) const { - std::pair Aneigh; - Aneigh=std::pair(fullatomlist_[neighbors_[i].first],fullatomlist_[neighbors_[i].second]); +NeighborList::pairAtomNumbers +NeighborList::getClosePairAtomNumber(const unsigned i) const { + pairAtomNumbers Aneigh=pairAtomNumbers( + fullatomlist_[neighbors_[i].first], + fullatomlist_[neighbors_[i].second]); return Aneigh; } -std::vector NeighborList::getNeighbors(unsigned index) { +std::vector NeighborList::getNeighbors(const unsigned index) { std::vector neighbors; for(unsigned int i=0; i. +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */ -#ifndef __PLUMED_tools_NeighborList_h +#ifndef __PLUMED_tools_NeighborList_h //{ #define __PLUMED_tools_NeighborList_h #include "Vector.h" @@ -34,35 +34,52 @@ class Communicator; /// \ingroup TOOLBOX /// A class that implements neighbor lists from two lists or a single list of atoms -class NeighborList -{ - bool reduced; +class NeighborList { +public: + using pairIDs=std::pair; + using pairAtomNumbers=std::pair; +private: + bool reduced=false; bool serial_; - bool do_pair_,do_pbc_,twolists_; + bool do_pair_; + bool do_pbc_; + bool twolists_; const PLMD::Pbc* pbc_; Communicator& comm; - std::vector fullatomlist_,requestlist_; - std::vector > neighbors_; + std::vector fullatomlist_{}; + std::vector requestlist_{}; + std::vector neighbors_{}; double distance_; - unsigned stride_,nlist0_,nlist1_,nallpairs_,lastupdate_; + size_t nlist0_=0; + size_t nlist1_=0; + size_t nallpairs_; + unsigned stride_=0; + unsigned lastupdate_=0; /// Initialize the neighbor list with all possible pairs void initialize(); /// Return the pair of indexes in the positions array /// of the two atoms forming the i-th pair among all possible pairs - std::pair getIndexPair(unsigned i); + pairIDs getIndexPair(unsigned i); /// Extract the list of atoms from the current list of close pairs void setRequestList(); public: NeighborList(const std::vector& list0, const std::vector& list1, const bool& serial, - const bool& do_pair, const bool& do_pbc, const PLMD::Pbc& pbc, Communicator &cm, - const double& distance=1.0e+30, const unsigned& stride=0); + const bool& do_pair, + const bool& do_pbc, + const PLMD::Pbc& pbc, + Communicator &cm, + const double& distance=1.0e+30, + const unsigned& stride=0); NeighborList(const std::vector& list0, const bool& serial, const bool& do_pbc, - const PLMD::Pbc& pbc, Communicator &cm, const double& distance=1.0e+30, + const PLMD::Pbc& pbc, + Communicator &cm, + const double& distance=1.0e+30, const unsigned& stride=0); + ~NeighborList(); /// Return the list of all atoms. These are needed to rebuild the neighbor list. std::vector& getFullAtomList(); /// Update the indexes in the neighbor list to match the @@ -80,15 +97,16 @@ class NeighborList void setLastUpdate(unsigned step); /// Get the size of the neighbor list unsigned size() const; +/// Get the distance used to create the neighbor list + double distance() const; /// Get the i-th pair of the neighbor list - std::pair getClosePair(unsigned i) const; + pairIDs getClosePair(unsigned i) const; /// Get the list of neighbors of the i-th atom std::vector getNeighbors(unsigned i); - ~NeighborList() {} /// Get the i-th pair of AtomNumbers from the neighbor list - std::pair getClosePairAtomNumber(unsigned i) const; + pairAtomNumbers getClosePairAtomNumber(unsigned i) const; }; -} +} // namespace PLMD -#endif +#endif //__PLUMED_tools_NeighborList_h }