473a3ebeef
an alternate route to address this issue would be to allow an "ANY" setting for neighbor list requests and then query the neighbor list for newton setting instead of the force class.
2076 lines
64 KiB
C++
2076 lines
64 KiB
C++
/* ----------------------------------------------------------------------
|
|
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
|
http://lammps.sandia.gov, Sandia National Laboratories
|
|
Steve Plimpton, sjplimp@sandia.gov
|
|
|
|
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
|
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
|
certain rights in this software. This software is distributed under
|
|
the GNU General Public License.
|
|
|
|
See the README file in the top-level LAMMPS directory.
|
|
------------------------------------------------------------------------- */
|
|
|
|
/* ----------------------------------------------------------------------
|
|
Contributing author (triclinic and multi-neigh) : Pieter in 't Veld (SNL)
|
|
------------------------------------------------------------------------- */
|
|
|
|
#include <mpi.h>
|
|
#include <math.h>
|
|
#include <stdlib.h>
|
|
#include <string.h>
|
|
#include "neighbor.h"
|
|
#include "neigh_list.h"
|
|
#include "neigh_request.h"
|
|
#include "style_nbin.h"
|
|
#include "style_nstencil.h"
|
|
#include "style_npair.h"
|
|
#include "style_ntopo.h"
|
|
#include "atom.h"
|
|
#include "atom_vec.h"
|
|
#include "comm.h"
|
|
#include "force.h"
|
|
#include "pair.h"
|
|
#include "domain.h"
|
|
#include "group.h"
|
|
#include "modify.h"
|
|
#include "fix.h"
|
|
#include "compute.h"
|
|
#include "update.h"
|
|
#include "respa.h"
|
|
#include "output.h"
|
|
#include "citeme.h"
|
|
#include "memory.h"
|
|
#include "error.h"
|
|
|
|
#include <map>
|
|
|
|
using namespace LAMMPS_NS;
|
|
using namespace NeighConst;
|
|
|
|
#define RQDELTA 1
|
|
#define EXDELTA 1
|
|
|
|
#define BIG 1.0e20
|
|
|
|
enum{NSQ,BIN,MULTI}; // also in NBin, NeighList, NStencil
|
|
enum{NONE,ALL,PARTIAL,TEMPLATE};
|
|
|
|
static const char cite_neigh_multi[] =
|
|
"neighbor multi command:\n\n"
|
|
"@Article{Intveld08,\n"
|
|
" author = {P.{\\,}J.~in{\\,}'t~Veld and S.{\\,}J.~Plimpton"
|
|
" and G.{\\,}S.~Grest},\n"
|
|
" title = {Accurate and Efficient Methods for Modeling Colloidal\n"
|
|
" Mixtures in an Explicit Solvent using Molecular Dynamics},\n"
|
|
" journal = {Comp.~Phys.~Comm.},\n"
|
|
" year = 2008,\n"
|
|
" volume = 179,\n"
|
|
" pages = {320--329}\n"
|
|
"}\n\n";
|
|
|
|
//#define NEIGH_LIST_DEBUG 1
|
|
|
|
/* ---------------------------------------------------------------------- */
|
|
|
|
Neighbor::Neighbor(LAMMPS *lmp) : Pointers(lmp),
|
|
pairclass(NULL), pairnames(NULL), pairmasks(NULL)
|
|
{
|
|
MPI_Comm_rank(world,&me);
|
|
MPI_Comm_size(world,&nprocs);
|
|
|
|
firsttime = 1;
|
|
|
|
style = BIN;
|
|
every = 1;
|
|
delay = 10;
|
|
dist_check = 1;
|
|
pgsize = 100000;
|
|
oneatom = 2000;
|
|
binsizeflag = 0;
|
|
build_once = 0;
|
|
cluster_check = 0;
|
|
ago = -1;
|
|
|
|
cutneighmax = 0.0;
|
|
cutneighsq = NULL;
|
|
cutneighghostsq = NULL;
|
|
cuttype = NULL;
|
|
cuttypesq = NULL;
|
|
fixchecklist = NULL;
|
|
|
|
// pairwise neighbor lists and associated data structs
|
|
|
|
nlist = 0;
|
|
lists = NULL;
|
|
|
|
nbin = 0;
|
|
neigh_bin = NULL;
|
|
|
|
nstencil = 0;
|
|
neigh_stencil = NULL;
|
|
|
|
neigh_pair = NULL;
|
|
|
|
nstencil_perpetual = 0;
|
|
slist = NULL;
|
|
|
|
npair_perpetual = 0;
|
|
plist = NULL;
|
|
|
|
nrequest = maxrequest = 0;
|
|
requests = NULL;
|
|
|
|
old_nrequest = 0;
|
|
old_requests = NULL;
|
|
|
|
old_style = style;
|
|
old_triclinic = 0;
|
|
old_pgsize = pgsize;
|
|
old_oneatom = oneatom;
|
|
|
|
zeroes = NULL;
|
|
|
|
binclass = NULL;
|
|
binnames = NULL;
|
|
binmasks = NULL;
|
|
stencilclass = NULL;
|
|
stencilnames = NULL;
|
|
stencilmasks = NULL;
|
|
|
|
// topology lists
|
|
|
|
bondwhich = anglewhich = dihedralwhich = improperwhich = NONE;
|
|
|
|
neigh_bond = NULL;
|
|
neigh_angle = NULL;
|
|
neigh_dihedral = NULL;
|
|
neigh_improper = NULL;
|
|
|
|
// coords at last neighboring
|
|
|
|
maxhold = 0;
|
|
xhold = NULL;
|
|
lastcall = -1;
|
|
last_setup_bins = -1;
|
|
|
|
// pair exclusion list info
|
|
|
|
includegroup = 0;
|
|
|
|
nex_type = maxex_type = 0;
|
|
ex1_type = ex2_type = NULL;
|
|
ex_type = NULL;
|
|
|
|
nex_group = maxex_group = 0;
|
|
ex1_group = ex2_group = ex1_bit = ex2_bit = NULL;
|
|
|
|
nex_mol = maxex_mol = 0;
|
|
ex_mol_group = ex_mol_bit = NULL;
|
|
|
|
// Kokkos setting
|
|
|
|
copymode = 0;
|
|
}
|
|
|
|
/* ---------------------------------------------------------------------- */
|
|
|
|
Neighbor::~Neighbor()
|
|
{
|
|
if (copymode) return;
|
|
|
|
memory->destroy(cutneighsq);
|
|
memory->destroy(cutneighghostsq);
|
|
delete [] cuttype;
|
|
delete [] cuttypesq;
|
|
delete [] fixchecklist;
|
|
|
|
for (int i = 0; i < nlist; i++) delete lists[i];
|
|
for (int i = 0; i < nbin; i++) delete neigh_bin[i];
|
|
for (int i = 0; i < nstencil; i++) delete neigh_stencil[i];
|
|
for (int i = 0; i < nlist; i++) delete neigh_pair[i];
|
|
delete [] lists;
|
|
delete [] neigh_bin;
|
|
delete [] neigh_stencil;
|
|
delete [] neigh_pair;
|
|
|
|
delete [] slist;
|
|
delete [] plist;
|
|
|
|
for (int i = 0; i < nrequest; i++) delete requests[i];
|
|
memory->sfree(requests);
|
|
for (int i = 0; i < old_nrequest; i++) delete old_requests[i];
|
|
memory->sfree(old_requests);
|
|
|
|
delete [] zeroes;
|
|
|
|
delete [] binclass;
|
|
delete [] binnames;
|
|
delete [] binmasks;
|
|
delete [] stencilclass;
|
|
delete [] stencilnames;
|
|
delete [] stencilmasks;
|
|
delete [] pairclass;
|
|
delete [] pairnames;
|
|
delete [] pairmasks;
|
|
|
|
delete neigh_bond;
|
|
delete neigh_angle;
|
|
delete neigh_dihedral;
|
|
delete neigh_improper;
|
|
|
|
memory->destroy(xhold);
|
|
|
|
memory->destroy(ex1_type);
|
|
memory->destroy(ex2_type);
|
|
memory->destroy(ex_type);
|
|
|
|
memory->destroy(ex1_group);
|
|
memory->destroy(ex2_group);
|
|
delete [] ex1_bit;
|
|
delete [] ex2_bit;
|
|
|
|
memory->destroy(ex_mol_group);
|
|
delete [] ex_mol_bit;
|
|
}
|
|
|
|
/* ---------------------------------------------------------------------- */
|
|
|
|
void Neighbor::init()
|
|
{
|
|
int i,j,n;
|
|
|
|
ncalls = ndanger = 0;
|
|
dimension = domain->dimension;
|
|
triclinic = domain->triclinic;
|
|
newton_pair = force->newton_pair;
|
|
|
|
// error check
|
|
|
|
if (delay > 0 && (delay % every) != 0)
|
|
error->all(FLERR,"Neighbor delay must be 0 or multiple of every setting");
|
|
|
|
if (pgsize < 10*oneatom)
|
|
error->all(FLERR,"Neighbor page size must be >= 10x the one atom setting");
|
|
|
|
// ------------------------------------------------------------------
|
|
// settings
|
|
|
|
// bbox lo/hi ptrs = bounding box of entire domain, stored by Domain
|
|
|
|
if (triclinic == 0) {
|
|
bboxlo = domain->boxlo;
|
|
bboxhi = domain->boxhi;
|
|
} else {
|
|
bboxlo = domain->boxlo_bound;
|
|
bboxhi = domain->boxhi_bound;
|
|
}
|
|
|
|
// set neighbor cutoffs (force cutoff + skin)
|
|
// trigger determines when atoms migrate and neighbor lists are rebuilt
|
|
// needs to be non-zero for migration distance check
|
|
// even if pair = NULL and no neighbor lists are used
|
|
// cutneigh = force cutoff + skin if cutforce > 0, else cutneigh = 0
|
|
// cutneighghost = pair cutghost if it requests it, else same as cutneigh
|
|
|
|
triggersq = 0.25*skin*skin;
|
|
boxcheck = 0;
|
|
if (domain->box_change && (domain->xperiodic || domain->yperiodic ||
|
|
(dimension == 3 && domain->zperiodic)))
|
|
boxcheck = 1;
|
|
|
|
n = atom->ntypes;
|
|
if (cutneighsq == NULL) {
|
|
if (lmp->kokkos) init_cutneighsq_kokkos(n);
|
|
else memory->create(cutneighsq,n+1,n+1,"neigh:cutneighsq");
|
|
memory->create(cutneighghostsq,n+1,n+1,"neigh:cutneighghostsq");
|
|
cuttype = new double[n+1];
|
|
cuttypesq = new double[n+1];
|
|
}
|
|
|
|
double cutoff,delta,cut;
|
|
cutneighmin = BIG;
|
|
cutneighmax = 0.0;
|
|
|
|
for (i = 1; i <= n; i++) {
|
|
cuttype[i] = cuttypesq[i] = 0.0;
|
|
for (j = 1; j <= n; j++) {
|
|
if (force->pair) cutoff = sqrt(force->pair->cutsq[i][j]);
|
|
else cutoff = 0.0;
|
|
if (cutoff > 0.0) delta = skin;
|
|
else delta = 0.0;
|
|
cut = cutoff + delta;
|
|
|
|
cutneighsq[i][j] = cut*cut;
|
|
cuttype[i] = MAX(cuttype[i],cut);
|
|
cuttypesq[i] = MAX(cuttypesq[i],cut*cut);
|
|
cutneighmin = MIN(cutneighmin,cut);
|
|
cutneighmax = MAX(cutneighmax,cut);
|
|
|
|
if (force->pair && force->pair->ghostneigh) {
|
|
cut = force->pair->cutghost[i][j] + skin;
|
|
cutneighghostsq[i][j] = cut*cut;
|
|
} else cutneighghostsq[i][j] = cut*cut;
|
|
}
|
|
}
|
|
cutneighmaxsq = cutneighmax * cutneighmax;
|
|
|
|
// rRESPA cutoffs
|
|
|
|
int respa = 0;
|
|
if (update->whichflag == 1 && strstr(update->integrate_style,"respa")) {
|
|
if (((Respa *) update->integrate)->level_inner >= 0) respa = 1;
|
|
if (((Respa *) update->integrate)->level_middle >= 0) respa = 2;
|
|
}
|
|
|
|
if (respa) {
|
|
double *cut_respa = ((Respa *) update->integrate)->cutoff;
|
|
cut_inner_sq = (cut_respa[1] + skin) * (cut_respa[1] + skin);
|
|
cut_middle_sq = (cut_respa[3] + skin) * (cut_respa[3] + skin);
|
|
cut_middle_inside_sq = (cut_respa[0] - skin) * (cut_respa[0] - skin);
|
|
if (cut_respa[0]-skin < 0) cut_middle_inside_sq = 0.0;
|
|
}
|
|
|
|
// fixchecklist = other classes that can induce reneighboring in decide()
|
|
|
|
restart_check = 0;
|
|
if (output->restart_flag) restart_check = 1;
|
|
|
|
delete [] fixchecklist;
|
|
fixchecklist = NULL;
|
|
fixchecklist = new int[modify->nfix];
|
|
|
|
fix_check = 0;
|
|
for (i = 0; i < modify->nfix; i++)
|
|
if (modify->fix[i]->force_reneighbor)
|
|
fixchecklist[fix_check++] = i;
|
|
|
|
must_check = 0;
|
|
if (restart_check || fix_check) must_check = 1;
|
|
|
|
// set special_flag for 1-2, 1-3, 1-4 neighbors
|
|
// flag[0] is not used, flag[1] = 1-2, flag[2] = 1-3, flag[3] = 1-4
|
|
// flag = 0 if both LJ/Coulomb special values are 0.0
|
|
// flag = 1 if both LJ/Coulomb special values are 1.0
|
|
// flag = 2 otherwise or if KSpace solver is enabled
|
|
// pairwise portion of KSpace solver uses all 1-2,1-3,1-4 neighbors
|
|
// or selected Coulomb-approixmation pair styles require it
|
|
|
|
if (force->special_lj[1] == 0.0 && force->special_coul[1] == 0.0)
|
|
special_flag[1] = 0;
|
|
else if (force->special_lj[1] == 1.0 && force->special_coul[1] == 1.0)
|
|
special_flag[1] = 1;
|
|
else special_flag[1] = 2;
|
|
|
|
if (force->special_lj[2] == 0.0 && force->special_coul[2] == 0.0)
|
|
special_flag[2] = 0;
|
|
else if (force->special_lj[2] == 1.0 && force->special_coul[2] == 1.0)
|
|
special_flag[2] = 1;
|
|
else special_flag[2] = 2;
|
|
|
|
if (force->special_lj[3] == 0.0 && force->special_coul[3] == 0.0)
|
|
special_flag[3] = 0;
|
|
else if (force->special_lj[3] == 1.0 && force->special_coul[3] == 1.0)
|
|
special_flag[3] = 1;
|
|
else special_flag[3] = 2;
|
|
|
|
if (force->kspace || force->pair_match("coul/wolf",0) ||
|
|
force->pair_match("coul/dsf",0) || force->pair_match("thole",0))
|
|
special_flag[1] = special_flag[2] = special_flag[3] = 2;
|
|
|
|
// maxwt = max multiplicative factor on atom indices stored in neigh list
|
|
|
|
maxwt = 0;
|
|
if (special_flag[1] == 2) maxwt = 2;
|
|
if (special_flag[2] == 2) maxwt = 3;
|
|
if (special_flag[3] == 2) maxwt = 4;
|
|
|
|
// ------------------------------------------------------------------
|
|
// xhold array
|
|
|
|
// free if not needed for this run
|
|
|
|
if (dist_check == 0) {
|
|
memory->destroy(xhold);
|
|
maxhold = 0;
|
|
xhold = NULL;
|
|
}
|
|
|
|
// first time allocation
|
|
|
|
if (dist_check) {
|
|
if (maxhold == 0) {
|
|
maxhold = atom->nmax;
|
|
memory->create(xhold,maxhold,3,"neigh:xhold");
|
|
}
|
|
}
|
|
|
|
// ------------------------------------------------------------------
|
|
// exclusion lists
|
|
|
|
// depend on type, group, molecule settings from neigh_modify
|
|
// warn if exclusions used with KSpace solver
|
|
|
|
n = atom->ntypes;
|
|
|
|
if (nex_type == 0 && nex_group == 0 && nex_mol == 0) exclude = 0;
|
|
else exclude = 1;
|
|
|
|
if (nex_type) {
|
|
if (lmp->kokkos)
|
|
init_ex_type_kokkos(n);
|
|
else {
|
|
memory->destroy(ex_type);
|
|
memory->create(ex_type,n+1,n+1,"neigh:ex_type");
|
|
}
|
|
|
|
for (i = 1; i <= n; i++)
|
|
for (j = 1; j <= n; j++)
|
|
ex_type[i][j] = 0;
|
|
|
|
for (i = 0; i < nex_type; i++) {
|
|
if (ex1_type[i] <= 0 || ex1_type[i] > n ||
|
|
ex2_type[i] <= 0 || ex2_type[i] > n)
|
|
error->all(FLERR,"Invalid atom type in neighbor exclusion list");
|
|
ex_type[ex1_type[i]][ex2_type[i]] = 1;
|
|
ex_type[ex2_type[i]][ex1_type[i]] = 1;
|
|
}
|
|
}
|
|
|
|
if (nex_group) {
|
|
if (lmp->kokkos)
|
|
init_ex_bit_kokkos();
|
|
else {
|
|
delete [] ex1_bit;
|
|
delete [] ex2_bit;
|
|
ex1_bit = new int[nex_group];
|
|
ex2_bit = new int[nex_group];
|
|
}
|
|
|
|
for (i = 0; i < nex_group; i++) {
|
|
ex1_bit[i] = group->bitmask[ex1_group[i]];
|
|
ex2_bit[i] = group->bitmask[ex2_group[i]];
|
|
}
|
|
}
|
|
|
|
if (nex_mol) {
|
|
if (lmp->kokkos)
|
|
init_ex_mol_bit_kokkos();
|
|
else {
|
|
delete [] ex_mol_bit;
|
|
ex_mol_bit = new int[nex_mol];
|
|
}
|
|
|
|
for (i = 0; i < nex_mol; i++)
|
|
ex_mol_bit[i] = group->bitmask[ex_mol_group[i]];
|
|
}
|
|
|
|
if (exclude && force->kspace && me == 0)
|
|
error->warning(FLERR,"Neighbor exclusions used with KSpace solver "
|
|
"may give inconsistent Coulombic energies");
|
|
|
|
// ------------------------------------------------------------------
|
|
// create pairwise lists
|
|
// one-time call to init_styles() to scan style files and setup
|
|
// init_pair() creates auxiliary classes: NBin, NStencil, NPair
|
|
|
|
if (firsttime) init_styles();
|
|
firsttime = 0;
|
|
|
|
init_pair();
|
|
|
|
// invoke copy_neighbor_info() in Bin,Stencil,Pair classes
|
|
// copied once per run in case any cutoff, exclusion, special info changed
|
|
|
|
for (i = 0; i < nbin; i++) neigh_bin[i]->copy_neighbor_info();
|
|
for (i = 0; i < nstencil; i++) neigh_stencil[i]->copy_neighbor_info();
|
|
for (i = 0; i < nlist; i++)
|
|
if (neigh_pair[i]) neigh_pair[i]->copy_neighbor_info();
|
|
|
|
if (!same && comm->me == 0) print_pairwise_info();
|
|
requests_new2old();
|
|
|
|
// ------------------------------------------------------------------
|
|
// create topology lists
|
|
// instantiated topo styles can change from run to run
|
|
|
|
init_topology();
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
create and initialize lists of Nbin, Nstencil, NPair classes
|
|
lists have info on all classes in 3 style*.h files
|
|
cannot do this in constructor, b/c too early to instantiate classes
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::init_styles()
|
|
{
|
|
// extract info from NBin classes listed in style_nbin.h
|
|
|
|
nbclass = 0;
|
|
|
|
#define NBIN_CLASS
|
|
#define NBinStyle(key,Class,bitmasks) nbclass++;
|
|
#include "style_nbin.h"
|
|
#undef NBinStyle
|
|
#undef NBIN_CLASS
|
|
|
|
binclass = new BinCreator[nbclass];
|
|
binnames = new char*[nbclass];
|
|
binmasks = new int[nbclass];
|
|
nbclass = 0;
|
|
|
|
#define NBIN_CLASS
|
|
#define NBinStyle(key,Class,bitmasks) \
|
|
binnames[nbclass] = (char *) #key; \
|
|
binclass[nbclass] = &bin_creator<Class>; \
|
|
binmasks[nbclass++] = bitmasks;
|
|
#include "style_nbin.h"
|
|
#undef NBinStyle
|
|
#undef NBIN_CLASS
|
|
|
|
// extract info from NStencil classes listed in style_nstencil.h
|
|
|
|
nsclass = 0;
|
|
|
|
#define NSTENCIL_CLASS
|
|
#define NStencilStyle(key,Class,bitmasks) nsclass++;
|
|
#include "style_nstencil.h"
|
|
#undef NStencilStyle
|
|
#undef NSTENCIL_CLASS
|
|
|
|
stencilclass = new StencilCreator[nsclass];
|
|
stencilnames = new char*[nsclass];
|
|
stencilmasks = new int[nsclass];
|
|
nsclass = 0;
|
|
|
|
#define NSTENCIL_CLASS
|
|
#define NStencilStyle(key,Class,bitmasks) \
|
|
stencilnames[nsclass] = (char *) #key; \
|
|
stencilclass[nsclass] = &stencil_creator<Class>; \
|
|
stencilmasks[nsclass++] = bitmasks;
|
|
#include "style_nstencil.h"
|
|
#undef NStencilStyle
|
|
#undef NSTENCIL_CLASS
|
|
|
|
// extract info from NPair classes listed in style_npair.h
|
|
|
|
npclass = 0;
|
|
|
|
#define NPAIR_CLASS
|
|
#define NPairStyle(key,Class,bitmasks) npclass++;
|
|
#include "style_npair.h"
|
|
#undef NPairStyle
|
|
#undef NPAIR_CLASS
|
|
|
|
pairclass = new PairCreator[npclass];
|
|
pairnames = new char*[npclass];
|
|
pairmasks = new int[npclass];
|
|
npclass = 0;
|
|
|
|
#define NPAIR_CLASS
|
|
#define NPairStyle(key,Class,bitmasks) \
|
|
pairnames[npclass] = (char *) #key; \
|
|
pairclass[npclass] = &pair_creator<Class>; \
|
|
pairmasks[npclass++] = bitmasks;
|
|
#include "style_npair.h"
|
|
#undef NPairStyle
|
|
#undef NPAIR_CLASS
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
create and initialize NPair classes
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::init_pair()
|
|
{
|
|
int i,j,k,m;
|
|
|
|
// test if pairwise lists need to be re-created
|
|
// no need to re-create if:
|
|
// neigh style, triclinic, pgsize, oneatom have not changed
|
|
// current requests = old requests
|
|
// first archive request params for current requests
|
|
// before possibly changing them below
|
|
|
|
for (i = 0; i < nrequest; i++) requests[i]->archive();
|
|
|
|
same = 1;
|
|
if (style != old_style) same = 0;
|
|
if (triclinic != old_triclinic) same = 0;
|
|
if (pgsize != old_pgsize) same = 0;
|
|
if (oneatom != old_oneatom) same = 0;
|
|
if (nrequest != old_nrequest) same = 0;
|
|
else
|
|
for (i = 0; i < nrequest; i++)
|
|
if (requests[i]->identical(old_requests[i]) == 0) same = 0;
|
|
|
|
#ifdef NEIGH_LIST_DEBUG
|
|
if (comm->me == 0) printf("SAME flag %d\n",same);
|
|
#endif
|
|
|
|
if (same) return;
|
|
|
|
// delete old lists and create new ones
|
|
|
|
for (i = 0; i < nlist; i++) delete lists[i];
|
|
for (i = 0; i < nbin; i++) delete neigh_bin[i];
|
|
for (i = 0; i < nstencil; i++) delete neigh_stencil[i];
|
|
for (i = 0; i < nlist; i++) delete neigh_pair[i];
|
|
delete [] lists;
|
|
delete [] neigh_bin;
|
|
delete [] neigh_stencil;
|
|
delete [] neigh_pair;
|
|
|
|
nlist = nrequest;
|
|
|
|
lists = new NeighList*[nrequest];
|
|
neigh_bin = new NBin*[nrequest];
|
|
neigh_stencil = new NStencil*[nrequest];
|
|
neigh_pair = new NPair*[nrequest];
|
|
|
|
// create individual lists, one per request
|
|
// pass list ptr back to requestor (except for Command class)
|
|
// wait to allocate initial pages until copy lists are detected
|
|
|
|
for (i = 0; i < nrequest; i++) {
|
|
if (requests[i]->kokkos_host || requests[i]->kokkos_device)
|
|
create_kokkos_list(i);
|
|
else
|
|
lists[i] = new NeighList(lmp);
|
|
lists[i]->index = i;
|
|
|
|
if (requests[i]->pair) {
|
|
Pair *pair = (Pair *) requests[i]->requestor;
|
|
pair->init_list(requests[i]->id,lists[i]);
|
|
} else if (requests[i]->fix) {
|
|
Fix *fix = (Fix *) requests[i]->requestor;
|
|
fix->init_list(requests[i]->id,lists[i]);
|
|
} else if (requests[i]->compute) {
|
|
Compute *compute = (Compute *) requests[i]->requestor;
|
|
compute->init_list(requests[i]->id,lists[i]);
|
|
}
|
|
}
|
|
|
|
// morph requests via A,B,C rules
|
|
// this is to avoid duplicate or inefficient builds
|
|
// update both request and list when morph
|
|
|
|
// (A) rule:
|
|
// invoke post_constructor() for all lists
|
|
// processes copy,skip,half_from_full,granhistory,respaouter lists
|
|
// error checks and resets internal ptrs to other lists that now exist
|
|
|
|
for (i = 0; i < nrequest; i++)
|
|
lists[i]->post_constructor(requests[i]);
|
|
|
|
// (B) rule:
|
|
// if request = pair, half, newton != 2
|
|
// and full perpetual non-skip/copy list exists,
|
|
// then morph to half_from_full of matching parent list
|
|
// NOTE: should be OK if parent is skip list?
|
|
// see build method comments
|
|
// parent can be pair or fix, so long as perpetual fix
|
|
// NOTE: could remove newton != 2 restriction if added
|
|
// half_from_full_newtoff_ghost NPair class
|
|
// this would require full list having ghost info
|
|
// would be useful when reax/c used in hybrid mode, e.g. with airebo
|
|
|
|
for (i = 0; i < nrequest; i++) {
|
|
if (requests[i]->pair && requests[i]->half && requests[i]->newton != 2) {
|
|
for (j = 0; j < nrequest; j++) {
|
|
// Kokkos doesn't yet support half from full
|
|
if (requests[i]->kokkos_device || requests[j]->kokkos_device) continue;
|
|
if (requests[i]->kokkos_host || requests[j]->kokkos_host) continue;
|
|
|
|
if (requests[j]->full && requests[j]->occasional == 0 &&
|
|
!requests[j]->skip && !requests[j]->copy) break;
|
|
}
|
|
if (j < nrequest) {
|
|
requests[i]->half = 0;
|
|
requests[i]->half_from_full = 1;
|
|
lists[i]->listfull = lists[j];
|
|
}
|
|
}
|
|
}
|
|
|
|
// (C) rule:
|
|
// for fix/compute requests, occasional or not does not matter
|
|
// 1st check:
|
|
// if request = half and non-skip/copy pair half/respaouter request exists,
|
|
// or if request = full and non-skip/copy pair full request exists,
|
|
// or if request = gran and non-skip/copy pair gran request exists,
|
|
// then morph to copy of the matching parent list
|
|
// 2nd check: only if no match to 1st check
|
|
// if request = half and non-skip/copy pair full request exists,
|
|
// then morph to half_from_full of the matching parent list
|
|
// for 1st or 2nd check, parent can be copy list or pair or fix
|
|
|
|
for (i = 0; i < nrequest; i++) {
|
|
if (!requests[i]->fix && !requests[i]->compute) continue;
|
|
for (j = 0; j < nrequest; j++) {
|
|
// Kokkos flags must match
|
|
if (requests[i]->kokkos_device != requests[j]->kokkos_device) continue;
|
|
if (requests[i]->kokkos_host != requests[j]->kokkos_host) continue;
|
|
|
|
if (requests[i]->ssa != requests[j]->ssa) continue;
|
|
// newton 2 and newton 0 both are newton off
|
|
if ((requests[i]->newton & 2) != (requests[j]->newton & 2)) continue;
|
|
|
|
if (requests[i]->half && requests[j]->pair &&
|
|
!requests[j]->skip && requests[j]->half && !requests[j]->copy)
|
|
break;
|
|
if (requests[i]->half && requests[j]->pair &&
|
|
!requests[j]->skip && requests[j]->respaouter && !requests[j]->copy)
|
|
break;
|
|
if (requests[i]->full && requests[j]->pair &&
|
|
!requests[j]->skip && requests[j]->full && !requests[j]->copy)
|
|
break;
|
|
if (requests[i]->gran && requests[j]->pair &&
|
|
!requests[j]->skip && requests[j]->gran && !requests[j]->copy)
|
|
break;
|
|
}
|
|
if (j < nrequest) {
|
|
requests[i]->copy = 1;
|
|
requests[i]->otherlist = j;
|
|
lists[i]->copy = 1;
|
|
lists[i]->listcopy = lists[j];
|
|
continue;
|
|
}
|
|
for (j = 0; j < nrequest; j++) {
|
|
// Kokkos doesn't yet support half from full
|
|
if (requests[i]->kokkos_device || requests[j]->kokkos_device) continue;
|
|
if (requests[i]->kokkos_host || requests[j]->kokkos_host) continue;
|
|
|
|
if (requests[i]->half && requests[j]->pair &&
|
|
!requests[j]->skip && requests[j]->full && !requests[j]->copy)
|
|
break;
|
|
}
|
|
if (j < nrequest) {
|
|
requests[i]->half = 0;
|
|
requests[i]->half_from_full = 1;
|
|
lists[i]->listfull = lists[j];
|
|
}
|
|
}
|
|
|
|
// assign Bin,Stencil,Pair style to each list
|
|
|
|
int flag;
|
|
for (i = 0; i < nrequest; i++) {
|
|
flag = choose_bin(requests[i]);
|
|
lists[i]->bin_method = flag;
|
|
if (flag < 0)
|
|
error->all(FLERR,"Requested neighbor bin option does not exist");
|
|
|
|
flag = choose_stencil(requests[i]);
|
|
lists[i]->stencil_method = flag;
|
|
if (flag < 0)
|
|
error->all(FLERR,"Requested neighbor stencil method does not exist");
|
|
|
|
flag = choose_pair(requests[i]);
|
|
lists[i]->pair_method = flag;
|
|
if (flag < 0)
|
|
error->all(FLERR,"Requested neighbor pair method does not exist");
|
|
}
|
|
|
|
// instantiate unique Bin,Stencil classes in neigh_bin & neigh_stencil vecs
|
|
// instantiate one Pair class per list in neigh_pair vec
|
|
|
|
nbin = 0;
|
|
for (i = 0; i < nrequest; i++) {
|
|
flag = lists[i]->bin_method;
|
|
if (flag == 0) continue;
|
|
for (j = 0; j < nbin; j++)
|
|
if (neigh_bin[j]->istyle == flag) break;
|
|
if (j < nbin) continue;
|
|
BinCreator bin_creator = binclass[flag-1];
|
|
neigh_bin[nbin] = bin_creator(lmp);
|
|
neigh_bin[nbin]->istyle = flag;
|
|
nbin++;
|
|
}
|
|
|
|
nstencil = 0;
|
|
for (i = 0; i < nrequest; i++) {
|
|
flag = lists[i]->stencil_method;
|
|
if (flag == 0) continue;
|
|
for (j = 0; j < nstencil; j++)
|
|
if (neigh_stencil[j]->istyle == flag) break;
|
|
if (j < nstencil) continue;
|
|
StencilCreator stencil_creator = stencilclass[flag-1];
|
|
neigh_stencil[nstencil] = stencil_creator(lmp);
|
|
neigh_stencil[nstencil]->istyle = flag;
|
|
int bin_method = lists[i]->bin_method;
|
|
for (k = 0; k < nbin; k++) {
|
|
if (neigh_bin[k]->istyle == bin_method) {
|
|
neigh_stencil[nstencil]->nb = neigh_bin[k];
|
|
break;
|
|
}
|
|
}
|
|
if (k == nbin)
|
|
error->all(FLERR,"Could not assign bin method to neighbor stencil");
|
|
nstencil++;
|
|
}
|
|
|
|
for (i = 0; i < nrequest; i++) {
|
|
flag = lists[i]->pair_method;
|
|
if (flag == 0) {
|
|
neigh_pair[i] = NULL;
|
|
continue;
|
|
}
|
|
PairCreator pair_creator = pairclass[flag-1];
|
|
neigh_pair[i] = pair_creator(lmp);
|
|
neigh_pair[i]->istyle = flag;
|
|
|
|
int bin_method = lists[i]->bin_method;
|
|
if (bin_method == 0) neigh_pair[i]->nb = NULL;
|
|
else {
|
|
for (k = 0; k < nbin; k++) {
|
|
if (neigh_bin[k]->istyle == bin_method) {
|
|
neigh_pair[i]->nb = neigh_bin[k];
|
|
break;
|
|
}
|
|
}
|
|
if (k == nbin)
|
|
error->all(FLERR,"Could not assign bin method to neighbor pair");
|
|
}
|
|
int stencil_method = lists[i]->stencil_method;
|
|
if (stencil_method == 0) neigh_pair[i]->ns = NULL;
|
|
else {
|
|
for (k = 0; k < nstencil; k++) {
|
|
if (neigh_stencil[k]->istyle == stencil_method) {
|
|
neigh_pair[i]->ns = neigh_stencil[k];
|
|
break;
|
|
}
|
|
}
|
|
if (k == nstencil)
|
|
error->all(FLERR,"Could not assign stencil method to neighbor pair");
|
|
}
|
|
}
|
|
|
|
// allocate initial pages for each list, except if copy flag set
|
|
// allocate dnum vector of zeroes if set
|
|
|
|
int dnummax = 0;
|
|
for (i = 0; i < nlist; i++) {
|
|
if (lists[i]->copy) continue;
|
|
lists[i]->setup_pages(pgsize,oneatom);
|
|
dnummax = MAX(dnummax,lists[i]->dnum);
|
|
}
|
|
|
|
if (dnummax) {
|
|
delete [] zeroes;
|
|
zeroes = new double[dnummax];
|
|
for (i = 0; i < dnummax; i++) zeroes[i] = 0.0;
|
|
}
|
|
|
|
// first-time allocation of per-atom data for lists that are built and store
|
|
// lists that are not built: granhistory, respa inner/middle (no neigh_pair)
|
|
// lists that do not store: copy
|
|
// use atom->nmax for both grow() args
|
|
// i.e. grow first time to expanded size to avoid future reallocs
|
|
// also Kokkos list initialization
|
|
|
|
int maxatom = atom->nmax;
|
|
for (i = 0; i < nlist; i++)
|
|
if (neigh_pair[i] && !lists[i]->copy) lists[i]->grow(maxatom,maxatom);
|
|
|
|
// plist = indices of perpetual NPair classes
|
|
// perpetual = non-occasional, re-built at every reneighboring
|
|
// slist = indices of perpetual NStencil classes
|
|
// perpetual = used by any perpetual NPair class
|
|
|
|
delete [] slist;
|
|
delete [] plist;
|
|
nstencil_perpetual = npair_perpetual = 0;
|
|
slist = new int[nstencil];
|
|
plist = new int[nlist];
|
|
|
|
for (i = 0; i < nlist; i++) {
|
|
if (lists[i]->occasional == 0 && lists[i]->pair_method)
|
|
plist[npair_perpetual++] = i;
|
|
}
|
|
|
|
for (i = 0; i < nstencil; i++) {
|
|
flag = 0;
|
|
for (j = 0; j < npair_perpetual; j++)
|
|
if (lists[plist[j]]->stencil_method == neigh_stencil[i]->istyle)
|
|
flag = 1;
|
|
if (flag) slist[nstencil_perpetual++] = i;
|
|
}
|
|
|
|
// reorder plist vector if necessary
|
|
// relevant for lists that copy/skip/half-full from parent
|
|
// the child index must appear in plist after the parent index
|
|
// swap two indices within plist when dependency is mis-ordered
|
|
// done when entire pass thru plist results in no swaps
|
|
|
|
NeighList *ptr;
|
|
|
|
int done = 0;
|
|
while (!done) {
|
|
done = 1;
|
|
for (i = 0; i < npair_perpetual; i++) {
|
|
ptr = NULL;
|
|
if (lists[plist[i]]->listcopy) ptr = lists[plist[i]]->listcopy;
|
|
if (lists[plist[i]]->listskip) ptr = lists[plist[i]]->listskip;
|
|
if (lists[plist[i]]->listfull) ptr = lists[plist[i]]->listfull;
|
|
if (ptr == NULL) continue;
|
|
for (m = 0; m < nrequest; m++)
|
|
if (ptr == lists[m]) break;
|
|
for (j = 0; j < npair_perpetual; j++)
|
|
if (m == plist[j]) break;
|
|
if (j < i) continue;
|
|
int tmp = plist[i]; // swap I,J indices
|
|
plist[i] = plist[j];
|
|
plist[j] = tmp;
|
|
done = 0;
|
|
break;
|
|
}
|
|
}
|
|
|
|
// debug output
|
|
|
|
#ifdef NEIGH_LIST_DEBUG
|
|
for (i = 0; i < nrequest; i++) lists[i]->print_attributes();
|
|
#endif
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
create and initialize NTopo classes
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::init_topology()
|
|
{
|
|
int i,m;
|
|
|
|
if (!atom->molecular) return;
|
|
|
|
// set flags that determine which topology neighbor classes to use
|
|
// these settings could change from run to run, depending on fixes defined
|
|
// bonds,etc can only be broken for atom->molecular = 1, not 2
|
|
// SHAKE sets bonds and angles negative
|
|
// gcmc sets all bonds, angles, etc negative
|
|
// bond_quartic sets bonds to 0
|
|
// delete_bonds sets all interactions negative
|
|
|
|
int bond_off = 0;
|
|
int angle_off = 0;
|
|
for (i = 0; i < modify->nfix; i++)
|
|
if ((strcmp(modify->fix[i]->style,"shake") == 0)
|
|
|| (strcmp(modify->fix[i]->style,"rattle") == 0))
|
|
bond_off = angle_off = 1;
|
|
if (force->bond && force->bond_match("quartic")) bond_off = 1;
|
|
|
|
if (atom->avec->bonds_allow && atom->molecular == 1) {
|
|
for (i = 0; i < atom->nlocal; i++) {
|
|
if (bond_off) break;
|
|
for (m = 0; m < atom->num_bond[i]; m++)
|
|
if (atom->bond_type[i][m] <= 0) bond_off = 1;
|
|
}
|
|
}
|
|
|
|
if (atom->avec->angles_allow && atom->molecular == 1) {
|
|
for (i = 0; i < atom->nlocal; i++) {
|
|
if (angle_off) break;
|
|
for (m = 0; m < atom->num_angle[i]; m++)
|
|
if (atom->angle_type[i][m] <= 0) angle_off = 1;
|
|
}
|
|
}
|
|
|
|
int dihedral_off = 0;
|
|
if (atom->avec->dihedrals_allow && atom->molecular == 1) {
|
|
for (i = 0; i < atom->nlocal; i++) {
|
|
if (dihedral_off) break;
|
|
for (m = 0; m < atom->num_dihedral[i]; m++)
|
|
if (atom->dihedral_type[i][m] <= 0) dihedral_off = 1;
|
|
}
|
|
}
|
|
|
|
int improper_off = 0;
|
|
if (atom->avec->impropers_allow && atom->molecular == 1) {
|
|
for (i = 0; i < atom->nlocal; i++) {
|
|
if (improper_off) break;
|
|
for (m = 0; m < atom->num_improper[i]; m++)
|
|
if (atom->improper_type[i][m] <= 0) improper_off = 1;
|
|
}
|
|
}
|
|
|
|
for (i = 0; i < modify->nfix; i++)
|
|
if ((strcmp(modify->fix[i]->style,"gcmc") == 0))
|
|
bond_off = angle_off = dihedral_off = improper_off = 1;
|
|
|
|
// sync on/off settings across all procs
|
|
|
|
int onoff = bond_off;
|
|
MPI_Allreduce(&onoff,&bond_off,1,MPI_INT,MPI_MAX,world);
|
|
onoff = angle_off;
|
|
MPI_Allreduce(&onoff,&angle_off,1,MPI_INT,MPI_MAX,world);
|
|
onoff = dihedral_off;
|
|
MPI_Allreduce(&onoff,&dihedral_off,1,MPI_INT,MPI_MAX,world);
|
|
onoff = improper_off;
|
|
MPI_Allreduce(&onoff,&improper_off,1,MPI_INT,MPI_MAX,world);
|
|
|
|
// instantiate NTopo classes
|
|
|
|
if (atom->avec->bonds_allow) {
|
|
int old_bondwhich = bondwhich;
|
|
if (atom->molecular == 2) bondwhich = TEMPLATE;
|
|
else if (bond_off) bondwhich = PARTIAL;
|
|
else bondwhich = ALL;
|
|
if (!neigh_bond || bondwhich != old_bondwhich) {
|
|
delete neigh_bond;
|
|
if (bondwhich == ALL)
|
|
neigh_bond = new NTopoBondAll(lmp);
|
|
else if (bondwhich == PARTIAL)
|
|
neigh_bond = new NTopoBondPartial(lmp);
|
|
else if (bondwhich == TEMPLATE)
|
|
neigh_bond = new NTopoBondTemplate(lmp);
|
|
}
|
|
}
|
|
|
|
if (atom->avec->angles_allow) {
|
|
int old_anglewhich = anglewhich;
|
|
if (atom->molecular == 2) anglewhich = TEMPLATE;
|
|
else if (angle_off) anglewhich = PARTIAL;
|
|
else anglewhich = ALL;
|
|
if (!neigh_angle || anglewhich != old_anglewhich) {
|
|
delete neigh_angle;
|
|
if (anglewhich == ALL)
|
|
neigh_angle = new NTopoAngleAll(lmp);
|
|
else if (anglewhich == PARTIAL)
|
|
neigh_angle = new NTopoAnglePartial(lmp);
|
|
else if (anglewhich == TEMPLATE)
|
|
neigh_angle = new NTopoAngleTemplate(lmp);
|
|
}
|
|
}
|
|
|
|
if (atom->avec->dihedrals_allow) {
|
|
int old_dihedralwhich = dihedralwhich;
|
|
if (atom->molecular == 2) dihedralwhich = TEMPLATE;
|
|
else if (dihedral_off) dihedralwhich = PARTIAL;
|
|
else dihedralwhich = ALL;
|
|
if (!neigh_dihedral || dihedralwhich != old_dihedralwhich) {
|
|
delete neigh_dihedral;
|
|
if (dihedralwhich == ALL)
|
|
neigh_dihedral = new NTopoDihedralAll(lmp);
|
|
else if (dihedralwhich == PARTIAL)
|
|
neigh_dihedral = new NTopoDihedralPartial(lmp);
|
|
else if (dihedralwhich == TEMPLATE)
|
|
neigh_dihedral = new NTopoDihedralTemplate(lmp);
|
|
}
|
|
}
|
|
|
|
if (atom->avec->impropers_allow) {
|
|
int old_improperwhich = improperwhich;
|
|
if (atom->molecular == 2) improperwhich = TEMPLATE;
|
|
else if (improper_off) improperwhich = PARTIAL;
|
|
else improperwhich = ALL;
|
|
if (!neigh_improper || improperwhich != old_improperwhich) {
|
|
delete neigh_improper;
|
|
if (improperwhich == ALL)
|
|
neigh_improper = new NTopoImproperAll(lmp);
|
|
else if (improperwhich == PARTIAL)
|
|
neigh_improper = new NTopoImproperPartial(lmp);
|
|
else if (improperwhich == TEMPLATE)
|
|
neigh_improper = new NTopoImproperTemplate(lmp);
|
|
}
|
|
}
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
output summary of pairwise neighbor list info
|
|
only called by proc 0
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::print_pairwise_info()
|
|
{
|
|
int i,j,m;
|
|
char str[128];
|
|
const char *kind;
|
|
FILE *out;
|
|
|
|
const double cutghost = MAX(cutneighmax,comm->cutghostuser);
|
|
|
|
double binsize, bbox[3];
|
|
bbox[0] = bboxhi[0]-bboxlo[0];
|
|
bbox[1] = bboxhi[1]-bboxlo[1];
|
|
bbox[2] = bboxhi[2]-bboxlo[2];
|
|
if (binsizeflag) binsize = binsize_user;
|
|
else if (style == BIN) binsize = 0.5*cutneighmax;
|
|
else binsize = 0.5*cutneighmin;
|
|
if (binsize == 0.0) binsize = bbox[0];
|
|
|
|
int nperpetual = 0;
|
|
int noccasional = 0;
|
|
int nextra = 0;
|
|
for (i = 0; i < nlist; i++) {
|
|
if (lists[i]->pair_method == 0) nextra++;
|
|
else if (lists[i]->occasional) noccasional++;
|
|
else nperpetual++;
|
|
}
|
|
|
|
for (m = 0; m < 2; m++) {
|
|
if (m == 0) out = screen;
|
|
else out = logfile;
|
|
|
|
if (out) {
|
|
fprintf(out,"Neighbor list info ...\n");
|
|
fprintf(out," update every %d steps, delay %d steps, check %s\n",
|
|
every,delay,dist_check ? "yes" : "no");
|
|
fprintf(out," max neighbors/atom: %d, page size: %d\n",
|
|
oneatom, pgsize);
|
|
fprintf(out," master list distance cutoff = %g\n",cutneighmax);
|
|
fprintf(out," ghost atom cutoff = %g\n",cutghost);
|
|
if (style != NSQ)
|
|
fprintf(out," binsize = %g, bins = %g %g %g\n",binsize,
|
|
ceil(bbox[0]/binsize), ceil(bbox[1]/binsize),
|
|
ceil(bbox[2]/binsize));
|
|
|
|
fprintf(out," %d neighbor lists, "
|
|
"perpetual/occasional/extra = %d %d %d\n",
|
|
nlist,nperpetual,noccasional,nextra);
|
|
|
|
for (i = 0; i < nlist; i++) {
|
|
if (requests[i]->pair) {
|
|
char *pname = force->pair_match_ptr((Pair *) requests[i]->requestor);
|
|
sprintf(str," (%d) pair %s",i+1,pname);
|
|
} else if (requests[i]->fix) {
|
|
sprintf(str," (%d) fix %s",i+1,
|
|
((Fix *) requests[i]->requestor)->style);
|
|
} else if (requests[i]->compute) {
|
|
sprintf(str," (%d) compute %s",i+1,
|
|
((Compute *) requests[i]->requestor)->style);
|
|
} else {
|
|
sprintf(str," (%d) command %s",i+1,requests[i]->command_style);
|
|
}
|
|
fprintf(out,"%s",str);
|
|
|
|
if (requests[i]->half) kind = "half";
|
|
else if (requests[i]->full) kind = "full";
|
|
else if (requests[i]->gran) kind = "size";
|
|
else if (requests[i]->granhistory) kind = "size/history";
|
|
else if (requests[i]->respainner) kind = "respa/inner";
|
|
else if (requests[i]->respamiddle) kind = "respa/middle";
|
|
else if (requests[i]->respaouter) kind = "respa/outer";
|
|
else if (requests[i]->half_from_full) kind = "half/from/full";
|
|
if (requests[i]->occasional) fprintf(out,", occasional");
|
|
else fprintf(out,", perpetual");
|
|
if (requests[i]->ghost) fprintf(out,", ghost");
|
|
if (requests[i]->ssa) fprintf(out,", ssa");
|
|
if (requests[i]->omp) fprintf(out,", omp");
|
|
if (requests[i]->intel) fprintf(out,", intel");
|
|
if (requests[i]->kokkos_device) fprintf(out,", kokkos_device");
|
|
if (requests[i]->kokkos_host) fprintf(out,", kokkos_host");
|
|
if (requests[i]->copy)
|
|
fprintf(out,", copy from (%d)",requests[i]->otherlist+1);
|
|
if (requests[i]->skip)
|
|
fprintf(out,", skip from (%d)",requests[i]->otherlist+1);
|
|
if (requests[i]->off2on) fprintf(out,", off2on");
|
|
fprintf(out,"\n");
|
|
|
|
if (lists[i]->pair_method == 0) fprintf(out," pair build: none\n");
|
|
else fprintf(out," pair build: %s\n",
|
|
pairnames[lists[i]->pair_method-1]);
|
|
|
|
if (lists[i]->stencil_method == 0) fprintf(out," stencil: none\n");
|
|
else fprintf(out," stencil: %s\n",
|
|
stencilnames[lists[i]->stencil_method-1]);
|
|
|
|
if (lists[i]->bin_method == 0) fprintf(out," bin: none\n");
|
|
else fprintf(out," bin: %s\n",binnames[lists[i]->bin_method-1]);
|
|
}
|
|
|
|
/*
|
|
fprintf(out," %d stencil methods\n",nstencil);
|
|
for (i = 0; i < nstencil; i++)
|
|
fprintf(out," (%d) %s\n",
|
|
i+1,stencilnames[neigh_stencil[i]->istyle-1]);
|
|
|
|
fprintf(out," %d bin methods\n",nbin);
|
|
for (i = 0; i < nbin; i++)
|
|
fprintf(out," (%d) %s\n",i+1,binnames[neigh_bin[i]->istyle-1]);
|
|
*/
|
|
}
|
|
}
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
delete old NeighRequests
|
|
copy current requests and params to old for next run
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::requests_new2old()
|
|
{
|
|
for (int i = 0; i < old_nrequest; i++) delete old_requests[i];
|
|
memory->sfree(old_requests);
|
|
|
|
old_nrequest = nrequest;
|
|
old_requests = requests;
|
|
nrequest = maxrequest = 0;
|
|
requests = NULL;
|
|
old_style = style;
|
|
old_triclinic = triclinic;
|
|
old_pgsize = pgsize;
|
|
old_oneatom = oneatom;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
assign NBin class to a NeighList
|
|
use neigh request settings to build mask
|
|
match mask to list of masks of known Nbin classes
|
|
return index+1 of match in list of masks
|
|
return 0 for no binning
|
|
return -1 if no match
|
|
------------------------------------------------------------------------- */
|
|
|
|
int Neighbor::choose_bin(NeighRequest *rq)
|
|
{
|
|
// no binning needed
|
|
|
|
if (style == NSQ) return 0;
|
|
if (rq->skip || rq->copy || rq->half_from_full) return 0;
|
|
if (rq->granhistory) return 0;
|
|
if (rq->respainner || rq->respamiddle) return 0;
|
|
|
|
// flags for settings the request + system requires of NBin class
|
|
// ssaflag = no/yes ssa request
|
|
// intelflag = no/yes intel request
|
|
// kokkos_device_flag = no/yes kokkos device request
|
|
// kokkos_host_flag = no/yes kokkos host request
|
|
|
|
int ssaflag,intelflag,kokkos_device_flag,kokkos_host_flag;
|
|
|
|
ssaflag = intelflag = kokkos_device_flag = kokkos_host_flag = 0;
|
|
|
|
if (rq->ssa) ssaflag = NB_SSA;
|
|
if (rq->intel) intelflag = NB_INTEL;
|
|
if (rq->kokkos_device) kokkos_device_flag = NB_KOKKOS_DEVICE;
|
|
if (rq->kokkos_host) kokkos_host_flag = NB_KOKKOS_HOST;
|
|
|
|
// use flags to match exactly one of NBin class masks, bit by bit
|
|
|
|
int mask;
|
|
|
|
for (int i = 0; i < nbclass; i++) {
|
|
mask = binmasks[i];
|
|
|
|
if (ssaflag != (mask & NB_SSA)) continue;
|
|
if (intelflag != (mask & NB_INTEL)) continue;
|
|
if (kokkos_device_flag != (mask & NB_KOKKOS_DEVICE)) continue;
|
|
if (kokkos_host_flag != (mask & NB_KOKKOS_HOST)) continue;
|
|
|
|
return i+1;
|
|
}
|
|
|
|
// error return if matched none
|
|
|
|
return -1;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
assign NStencil class to a NeighList
|
|
use neigh request settings to build mask
|
|
match mask to list of masks of known NStencil classes
|
|
return index+1 of match in list of masks
|
|
return 0 for no binning
|
|
return -1 if no match
|
|
------------------------------------------------------------------------- */
|
|
|
|
int Neighbor::choose_stencil(NeighRequest *rq)
|
|
{
|
|
// no stencil creation needed
|
|
|
|
if (style == NSQ) return 0;
|
|
if (rq->skip || rq->copy || rq->half_from_full) return 0;
|
|
if (rq->granhistory) return 0;
|
|
if (rq->respainner || rq->respamiddle) return 0;
|
|
|
|
// flags for settings the request + system requires of NStencil class
|
|
// halfflag = half request (gran and respa are also half lists)
|
|
// fullflag = full request
|
|
// ghostflag = no/yes ghost request
|
|
// ssaflag = no/yes ssa request
|
|
// dimension = 2d/3d
|
|
// newtflag = newton off/on request
|
|
// triclinic = orthgonal/triclinic box
|
|
|
|
int halfflag,fullflag,ghostflag,ssaflag;
|
|
|
|
halfflag = fullflag = ghostflag = ssaflag = 0;
|
|
|
|
if (rq->half) halfflag = 1;
|
|
if (rq->full) fullflag = 1;
|
|
if (rq->gran) halfflag = 1;
|
|
if (rq->respaouter) halfflag = 1;
|
|
|
|
if (rq->ghost) ghostflag = NS_GHOST;
|
|
if (rq->ssa) ssaflag = NS_SSA;
|
|
|
|
int newtflag;
|
|
if (rq->newton == 0 && newton_pair) newtflag = 1;
|
|
else if (rq->newton == 0 && !newton_pair) newtflag = 0;
|
|
else if (rq->newton == 1) newtflag = 1;
|
|
else if (rq->newton == 2) newtflag = 0;
|
|
|
|
// use flags to match exactly one of NStencil class masks, bit by bit
|
|
// exactly one of halfflag,fullflag is set and thus must match
|
|
|
|
int mask;
|
|
|
|
for (int i = 0; i < nsclass; i++) {
|
|
mask = stencilmasks[i];
|
|
|
|
if (halfflag) {
|
|
if (!(mask & NS_HALF)) continue;
|
|
} else if (fullflag) {
|
|
if (!(mask & NS_FULL)) continue;
|
|
}
|
|
|
|
if (ghostflag != (mask & NS_GHOST)) continue;
|
|
if (ssaflag != (mask & NS_SSA)) continue;
|
|
|
|
if (style == BIN && !(mask & NS_BIN)) continue;
|
|
if (style == MULTI && !(mask & NS_MULTI)) continue;
|
|
|
|
if (dimension == 2 && !(mask & NS_2D)) continue;
|
|
if (dimension == 3 && !(mask & NS_3D)) continue;
|
|
|
|
if (newtflag && !(mask & NS_NEWTON)) continue;
|
|
if (!newtflag && !(mask & NS_NEWTOFF)) continue;
|
|
|
|
if (!triclinic && !(mask & NS_ORTHO)) continue;
|
|
if (triclinic && !(mask & NS_TRI)) continue;
|
|
|
|
return i+1;
|
|
}
|
|
|
|
// error return if matched none
|
|
|
|
return -1;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
assign NPair class to a NeighList
|
|
use neigh request settings to build mask
|
|
match mask to list of masks of known NPair classes
|
|
return index+1 of match in list of masks
|
|
return 0 for no binning
|
|
return -1 if no match
|
|
------------------------------------------------------------------------- */
|
|
|
|
int Neighbor::choose_pair(NeighRequest *rq)
|
|
{
|
|
// no NPair build performed
|
|
|
|
if (rq->granhistory) return 0;
|
|
if (rq->respainner || rq->respamiddle) return 0;
|
|
|
|
// error check for includegroup with ghost neighbor request
|
|
|
|
if (includegroup && rq->ghost)
|
|
error->all(FLERR,"Neighbor include group not allowed "
|
|
"with ghost neighbors");
|
|
|
|
// flags for settings the request + system requires of NPair class
|
|
// some are set to 0/1, others are set to mask bit
|
|
// comparisons below in loop over classes reflect that
|
|
// copyflag = no/yes copy request
|
|
// skipflag = no/yes skip request
|
|
// halfflag = half request (gran and respa are also half lists)
|
|
// fullflag = full request
|
|
// halffullflag = half_from_full request
|
|
// sizeflag = no/yes gran request for finite-size particles
|
|
// ghostflag = no/yes ghost request
|
|
// respaflag = no/yes respa request
|
|
// off2onflag = no/yes off2on request
|
|
// onesideflag = no/yes granonesided request
|
|
// ssaflag = no/yes request
|
|
// ompflag = no/yes omp request
|
|
// intelflag = no/yes intel request
|
|
// kokkos_device_flag = no/yes Kokkos device request
|
|
// kokkos_host_flag = no/yes Kokkos host request
|
|
// newtflag = newton off/on request
|
|
// style = NSQ/BIN/MULTI neighbor style
|
|
// triclinic = orthgonal/triclinic box
|
|
|
|
int copyflag,skipflag,halfflag,fullflag,halffullflag,sizeflag,respaflag,
|
|
ghostflag,off2onflag,onesideflag,ssaflag,ompflag,intelflag,
|
|
kokkos_device_flag,kokkos_host_flag;
|
|
|
|
copyflag = skipflag = halfflag = fullflag = halffullflag = sizeflag =
|
|
ghostflag = respaflag = off2onflag = onesideflag = ssaflag =
|
|
ompflag = intelflag = kokkos_device_flag = kokkos_host_flag = 0;
|
|
|
|
if (rq->copy) copyflag = NP_COPY;
|
|
if (rq->skip) skipflag = NP_SKIP;
|
|
|
|
// NOTE: exactly one of these request flags is set (see neigh_request.h)
|
|
// this requires gran/respaouter also set halfflag
|
|
// can simplify this logic, if follow NOTE in neigh_request.h
|
|
// why do size/off2on and size/off2on/oneside set NP_HALF
|
|
// either should set both half & full, or half should be in file name
|
|
// to be consistent with how other NP classes use "half"
|
|
|
|
if (rq->half) halfflag = 1;
|
|
if (rq->full) fullflag = 1;
|
|
if (rq->half_from_full) halffullflag = 1;
|
|
if (rq->gran) {
|
|
sizeflag = NP_SIZE;
|
|
halfflag = 1;
|
|
}
|
|
if (rq->respaouter) {
|
|
respaflag = NP_RESPA;
|
|
halfflag = 1;
|
|
}
|
|
|
|
if (rq->ghost) ghostflag = NP_GHOST;
|
|
if (rq->off2on) off2onflag = NP_OFF2ON;
|
|
if (rq->granonesided) onesideflag = NP_ONESIDE;
|
|
if (rq->ssa) ssaflag = NP_SSA;
|
|
if (rq->omp) ompflag = NP_OMP;
|
|
if (rq->intel) intelflag = NP_INTEL;
|
|
if (rq->kokkos_device) kokkos_device_flag = NP_KOKKOS_DEVICE;
|
|
if (rq->kokkos_host) kokkos_host_flag = NP_KOKKOS_HOST;
|
|
|
|
int newtflag;
|
|
if (rq->newton == 0 && newton_pair) newtflag = 1;
|
|
else if (rq->newton == 0 && !newton_pair) newtflag = 0;
|
|
else if (rq->newton == 1) newtflag = 1;
|
|
else if (rq->newton == 2) newtflag = 0;
|
|
|
|
// use flags to match exactly one of NPair class masks
|
|
// sequence of checks is bit by bit in NeighConst
|
|
|
|
int mask;
|
|
|
|
//printf("FLAGS: %d %d %d %d %d %d %d %d %d %d %d %d %d %d\n",
|
|
// copyflag,skipflag,halfflag,fullflag,halffullflag,
|
|
// sizeflag,respaflag,ghostflag,off2onflag,onesideflag,ssaflag,
|
|
// ompflag,intelflag,newtflag);
|
|
|
|
for (int i = 0; i < npclass; i++) {
|
|
mask = pairmasks[i];
|
|
|
|
// if copyflag set, return or continue with no further checks
|
|
|
|
if (copyflag) {
|
|
if (!(mask & NP_COPY)) continue;
|
|
if (kokkos_device_flag != (mask & NP_KOKKOS_DEVICE)) continue;
|
|
if (kokkos_host_flag != (mask & NP_KOKKOS_HOST)) continue;
|
|
return i+1;
|
|
}
|
|
|
|
// skipflag must match along with other flags, so do not return
|
|
|
|
if (skipflag != (mask & NP_SKIP)) continue;
|
|
|
|
// exactly one of halfflag,fullflag,halffullflag is set and must match
|
|
|
|
if (halfflag) {
|
|
if (!(mask & NP_HALF)) continue;
|
|
} else if (fullflag) {
|
|
if (!(mask & NP_FULL)) continue;
|
|
} else if (halffullflag) {
|
|
if (!(mask & NP_HALFFULL)) continue;
|
|
}
|
|
|
|
if (sizeflag != (mask & NP_SIZE)) continue;
|
|
if (respaflag != (mask & NP_RESPA)) continue;
|
|
if (ghostflag != (mask & NP_GHOST)) continue;
|
|
if (off2onflag != (mask & NP_OFF2ON)) continue;
|
|
if (onesideflag != (mask & NP_ONESIDE)) continue;
|
|
if (ssaflag != (mask & NP_SSA)) continue;
|
|
if (ompflag != (mask & NP_OMP)) continue;
|
|
if (intelflag != (mask & NP_INTEL)) continue;
|
|
|
|
// style is one of NSQ,BIN,MULTI and must match
|
|
|
|
if (style == NSQ) {
|
|
if (!(mask & NP_NSQ)) continue;
|
|
} else if (style == BIN) {
|
|
if (!(mask & NP_BIN)) continue;
|
|
} else if (style == MULTI) {
|
|
if (!(mask & NP_MULTI)) continue;
|
|
}
|
|
|
|
// newtflag is on or off and must match
|
|
|
|
if (newtflag) {
|
|
if (!(mask & NP_NEWTON)) continue;
|
|
} else if (!newtflag) {
|
|
if (!(mask & NP_NEWTOFF)) continue;
|
|
}
|
|
|
|
// triclinic flag is on or off and must match
|
|
|
|
if (triclinic) {
|
|
if (!(mask & NP_TRI)) continue;
|
|
} else if (!triclinic) {
|
|
if (!(mask & NP_ORTHO)) continue;
|
|
}
|
|
|
|
// Kokkos flags
|
|
|
|
if (kokkos_device_flag != (mask & NP_KOKKOS_DEVICE)) continue;
|
|
if (kokkos_host_flag != (mask & NP_KOKKOS_HOST)) continue;
|
|
|
|
return i+1;
|
|
}
|
|
|
|
//printf("NO MATCH\n");
|
|
|
|
// error return if matched none
|
|
|
|
return -1;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
called by other classes to request a pairwise neighbor list
|
|
------------------------------------------------------------------------- */
|
|
|
|
int Neighbor::request(void *requestor, int instance)
|
|
{
|
|
if (nrequest == maxrequest) {
|
|
maxrequest += RQDELTA;
|
|
requests = (NeighRequest **)
|
|
memory->srealloc(requests,maxrequest*sizeof(NeighRequest *),
|
|
"neighbor:requests");
|
|
}
|
|
|
|
requests[nrequest] = new NeighRequest(lmp);
|
|
requests[nrequest]->index = nrequest;
|
|
requests[nrequest]->requestor = requestor;
|
|
requests[nrequest]->requestor_instance = instance;
|
|
nrequest++;
|
|
return nrequest-1;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
one instance per entry in style_neigh_bin.h
|
|
------------------------------------------------------------------------- */
|
|
|
|
template <typename T>
|
|
NBin *Neighbor::bin_creator(LAMMPS *lmp)
|
|
{
|
|
return new T(lmp);
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
one instance per entry in style_neigh_stencil.h
|
|
------------------------------------------------------------------------- */
|
|
|
|
template <typename T>
|
|
NStencil *Neighbor::stencil_creator(LAMMPS *lmp)
|
|
{
|
|
return new T(lmp);
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
one instance per entry in style_neigh_pair.h
|
|
------------------------------------------------------------------------- */
|
|
|
|
template <typename T>
|
|
NPair *Neighbor::pair_creator(LAMMPS *lmp)
|
|
{
|
|
return new T(lmp);
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
setup neighbor binning and neighbor stencils
|
|
called before run and every reneighbor if box size/shape changes
|
|
only operates on perpetual lists
|
|
build_one() operates on occasional lists
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::setup_bins()
|
|
{
|
|
// invoke setup_bins() for all NBin
|
|
// actual binning is performed in build()
|
|
|
|
for (int i = 0; i < nbin; i++)
|
|
neigh_bin[i]->setup_bins(style);
|
|
|
|
// invoke create_setup() and create() for all perpetual NStencil
|
|
// same ops performed for occasional lists in build_one()
|
|
|
|
for (int i = 0; i < nstencil_perpetual; i++) {
|
|
neigh_stencil[slist[i]]->create_setup();
|
|
neigh_stencil[slist[i]]->create();
|
|
}
|
|
|
|
last_setup_bins = update->ntimestep;
|
|
}
|
|
|
|
/* ---------------------------------------------------------------------- */
|
|
|
|
int Neighbor::decide()
|
|
{
|
|
if (must_check) {
|
|
bigint n = update->ntimestep;
|
|
if (restart_check && n == output->next_restart) return 1;
|
|
for (int i = 0; i < fix_check; i++)
|
|
if (n == modify->fix[fixchecklist[i]]->next_reneighbor) return 1;
|
|
}
|
|
|
|
ago++;
|
|
if (ago >= delay && ago % every == 0) {
|
|
if (build_once) return 0;
|
|
if (dist_check == 0) return 1;
|
|
return check_distance();
|
|
} else return 0;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
if any atom moved trigger distance (half of neighbor skin) return 1
|
|
shrink trigger distance if box size has changed
|
|
conservative shrink procedure:
|
|
compute distance each of 8 corners of box has moved since last reneighbor
|
|
reduce skin distance by sum of 2 largest of the 8 values
|
|
new trigger = 1/2 of reduced skin distance
|
|
for orthogonal box, only need 2 lo/hi corners
|
|
for triclinic, need all 8 corners since deformations can displace all 8
|
|
------------------------------------------------------------------------- */
|
|
|
|
int Neighbor::check_distance()
|
|
{
|
|
double delx,dely,delz,rsq;
|
|
double delta,deltasq,delta1,delta2;
|
|
|
|
if (boxcheck) {
|
|
if (triclinic == 0) {
|
|
delx = bboxlo[0] - boxlo_hold[0];
|
|
dely = bboxlo[1] - boxlo_hold[1];
|
|
delz = bboxlo[2] - boxlo_hold[2];
|
|
delta1 = sqrt(delx*delx + dely*dely + delz*delz);
|
|
delx = bboxhi[0] - boxhi_hold[0];
|
|
dely = bboxhi[1] - boxhi_hold[1];
|
|
delz = bboxhi[2] - boxhi_hold[2];
|
|
delta2 = sqrt(delx*delx + dely*dely + delz*delz);
|
|
delta = 0.5 * (skin - (delta1+delta2));
|
|
deltasq = delta*delta;
|
|
} else {
|
|
domain->box_corners();
|
|
delta1 = delta2 = 0.0;
|
|
for (int i = 0; i < 8; i++) {
|
|
delx = corners[i][0] - corners_hold[i][0];
|
|
dely = corners[i][1] - corners_hold[i][1];
|
|
delz = corners[i][2] - corners_hold[i][2];
|
|
delta = sqrt(delx*delx + dely*dely + delz*delz);
|
|
if (delta > delta1) delta1 = delta;
|
|
else if (delta > delta2) delta2 = delta;
|
|
}
|
|
delta = 0.5 * (skin - (delta1+delta2));
|
|
deltasq = delta*delta;
|
|
}
|
|
} else deltasq = triggersq;
|
|
|
|
double **x = atom->x;
|
|
int nlocal = atom->nlocal;
|
|
if (includegroup) nlocal = atom->nfirst;
|
|
|
|
int flag = 0;
|
|
for (int i = 0; i < nlocal; i++) {
|
|
delx = x[i][0] - xhold[i][0];
|
|
dely = x[i][1] - xhold[i][1];
|
|
delz = x[i][2] - xhold[i][2];
|
|
rsq = delx*delx + dely*dely + delz*delz;
|
|
if (rsq > deltasq) flag = 1;
|
|
}
|
|
|
|
int flagall;
|
|
MPI_Allreduce(&flag,&flagall,1,MPI_INT,MPI_MAX,world);
|
|
if (flagall && ago == MAX(every,delay)) ndanger++;
|
|
return flagall;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
build perpetual neighbor lists
|
|
called at setup and every few timesteps during run or minimization
|
|
topology lists also built if topoflag = 1 (Kokkos calls with topoflag=0)
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::build(int topoflag)
|
|
{
|
|
int i,m;
|
|
|
|
ago = 0;
|
|
ncalls++;
|
|
lastcall = update->ntimestep;
|
|
|
|
int nlocal = atom->nlocal;
|
|
int nall = nlocal + atom->nghost;
|
|
|
|
// check that using special bond flags will not overflow neigh lists
|
|
|
|
if (nall > NEIGHMASK)
|
|
error->one(FLERR,"Too many local+ghost atoms for neighbor list");
|
|
|
|
// store current atom positions and box size if needed
|
|
|
|
if (dist_check) {
|
|
double **x = atom->x;
|
|
if (includegroup) nlocal = atom->nfirst;
|
|
if (atom->nmax > maxhold) {
|
|
maxhold = atom->nmax;
|
|
memory->destroy(xhold);
|
|
memory->create(xhold,maxhold,3,"neigh:xhold");
|
|
}
|
|
for (i = 0; i < nlocal; i++) {
|
|
xhold[i][0] = x[i][0];
|
|
xhold[i][1] = x[i][1];
|
|
xhold[i][2] = x[i][2];
|
|
}
|
|
if (boxcheck) {
|
|
if (triclinic == 0) {
|
|
boxlo_hold[0] = bboxlo[0];
|
|
boxlo_hold[1] = bboxlo[1];
|
|
boxlo_hold[2] = bboxlo[2];
|
|
boxhi_hold[0] = bboxhi[0];
|
|
boxhi_hold[1] = bboxhi[1];
|
|
boxhi_hold[2] = bboxhi[2];
|
|
} else {
|
|
domain->box_corners();
|
|
corners = domain->corners;
|
|
for (i = 0; i < 8; i++) {
|
|
corners_hold[i][0] = corners[i][0];
|
|
corners_hold[i][1] = corners[i][1];
|
|
corners_hold[i][2] = corners[i][2];
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
// bin atoms for all NBin instances
|
|
// not just NBin associated with perpetual lists
|
|
// b/c cannot wait to bin occasional lists in build_one() call
|
|
// if bin then, atoms may have moved outside of proc domain & bin extent,
|
|
// leading to errors or even a crash
|
|
|
|
if (style != NSQ) {
|
|
for (int i = 0; i < nbin; i++) {
|
|
neigh_bin[i]->bin_atoms_setup(nall);
|
|
neigh_bin[i]->bin_atoms();
|
|
}
|
|
}
|
|
|
|
// build pairwise lists for all perpetual NPair/NeighList
|
|
// grow() with nlocal/nall args so that only realloc if have to
|
|
|
|
for (i = 0; i < npair_perpetual; i++) {
|
|
m = plist[i];
|
|
if (!lists[m]->copy) lists[m]->grow(nlocal,nall);
|
|
neigh_pair[m]->build_setup();
|
|
neigh_pair[m]->build(lists[m]);
|
|
}
|
|
|
|
// build topology lists for bonds/angles/etc
|
|
|
|
if (atom->molecular && topoflag) build_topology();
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
build topology neighbor lists: bond, angle, dihedral, improper
|
|
copy their list info back to Neighbor for access by bond/angle/etc classes
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::build_topology()
|
|
{
|
|
if (force->bond) {
|
|
neigh_bond->build();
|
|
nbondlist = neigh_bond->nbondlist;
|
|
bondlist = neigh_bond->bondlist;
|
|
}
|
|
if (force->angle) {
|
|
neigh_angle->build();
|
|
nanglelist = neigh_angle->nanglelist;
|
|
anglelist = neigh_angle->anglelist;
|
|
}
|
|
if (force->dihedral) {
|
|
neigh_dihedral->build();
|
|
ndihedrallist = neigh_dihedral->ndihedrallist;
|
|
dihedrallist = neigh_dihedral->dihedrallist;
|
|
}
|
|
if (force->improper) {
|
|
neigh_improper->build();
|
|
nimproperlist = neigh_improper->nimproperlist;
|
|
improperlist = neigh_improper->improperlist;
|
|
}
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
build a single occasional pairwise neighbor list indexed by I
|
|
called by other classes
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::build_one(class NeighList *mylist, int preflag)
|
|
{
|
|
// check if list structure is initialized
|
|
|
|
if (mylist == NULL)
|
|
error->all(FLERR,"Trying to build an occasional neighbor list "
|
|
"before initialization completed");
|
|
|
|
// build_one() should never be invoked on a perpetual list
|
|
|
|
if (!mylist->occasional)
|
|
error->all(FLERR,"Neighbor build one invoked on perpetual list");
|
|
|
|
// no need to build if already built since last re-neighbor
|
|
// preflag is set by fix bond/create and fix bond/swap
|
|
// b/c they invoke build_one() on same step neigh list is re-built,
|
|
// but before re-build, so need to use ">" instead of ">="
|
|
|
|
NPair *np = neigh_pair[mylist->index];
|
|
|
|
if (preflag) {
|
|
if (np->last_build > lastcall) return;
|
|
} else {
|
|
if (np->last_build >= lastcall) return;
|
|
}
|
|
|
|
// if this is copy list and parent is occasional list,
|
|
// or this is half_from_full and parent is occasional list,
|
|
// insure parent is current
|
|
|
|
if (mylist->listcopy && mylist->listcopy->occasional)
|
|
build_one(mylist->listcopy,preflag);
|
|
if (mylist->listfull && mylist->listfull->occasional)
|
|
build_one(mylist->listfull,preflag);
|
|
|
|
// create stencil if hasn't been created since last setup_bins() call
|
|
|
|
NStencil *ns = np->ns;
|
|
if (ns && ns->last_create < last_setup_bins) {
|
|
ns->create_setup();
|
|
ns->create();
|
|
}
|
|
|
|
// build the list
|
|
|
|
np->build_setup();
|
|
np->build(mylist);
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
set neighbor style and skin distance
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::set(int narg, char **arg)
|
|
{
|
|
if (narg != 2) error->all(FLERR,"Illegal neighbor command");
|
|
|
|
skin = force->numeric(FLERR,arg[0]);
|
|
if (skin < 0.0) error->all(FLERR,"Illegal neighbor command");
|
|
|
|
if (strcmp(arg[1],"nsq") == 0) style = NSQ;
|
|
else if (strcmp(arg[1],"bin") == 0) style = BIN;
|
|
else if (strcmp(arg[1],"multi") == 0) style = MULTI;
|
|
else error->all(FLERR,"Illegal neighbor command");
|
|
|
|
if (style == MULTI && lmp->citeme) lmp->citeme->add(cite_neigh_multi);
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
reset timestamps in all NeignBin, NStencil, NPair classes
|
|
so that neighbor lists will rebuild properly with timestep change
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::reset_timestep(bigint ntimestep)
|
|
{
|
|
for (int i = 0; i < nbin; i++) {
|
|
neigh_bin[i]->last_setup = -1;
|
|
neigh_bin[i]->last_bin = -1;
|
|
neigh_bin[i]->last_bin_memory = -1;
|
|
}
|
|
|
|
for (int i = 0; i < nstencil; i++) {
|
|
neigh_stencil[i]->last_create = -1;
|
|
neigh_stencil[i]->last_stencil_memory = -1;
|
|
neigh_stencil[i]->last_copy_bin = -1;
|
|
}
|
|
|
|
for (int i = 0; i < nlist; i++) {
|
|
if (!neigh_pair[i]) continue;
|
|
neigh_pair[i]->last_build = -1;
|
|
neigh_pair[i]->last_copy_bin_setup = -1;
|
|
neigh_pair[i]->last_copy_bin = -1;
|
|
neigh_pair[i]->last_copy_stencil = -1;
|
|
}
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
modify parameters of the pair-wise neighbor build
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::modify_params(int narg, char **arg)
|
|
{
|
|
int iarg = 0;
|
|
while (iarg < narg) {
|
|
if (strcmp(arg[iarg],"every") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
every = force->inumeric(FLERR,arg[iarg+1]);
|
|
if (every <= 0) error->all(FLERR,"Illegal neigh_modify command");
|
|
iarg += 2;
|
|
} else if (strcmp(arg[iarg],"delay") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
delay = force->inumeric(FLERR,arg[iarg+1]);
|
|
if (delay < 0) error->all(FLERR,"Illegal neigh_modify command");
|
|
iarg += 2;
|
|
} else if (strcmp(arg[iarg],"check") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
if (strcmp(arg[iarg+1],"yes") == 0) dist_check = 1;
|
|
else if (strcmp(arg[iarg+1],"no") == 0) dist_check = 0;
|
|
else error->all(FLERR,"Illegal neigh_modify command");
|
|
iarg += 2;
|
|
} else if (strcmp(arg[iarg],"once") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
if (strcmp(arg[iarg+1],"yes") == 0) build_once = 1;
|
|
else if (strcmp(arg[iarg+1],"no") == 0) build_once = 0;
|
|
else error->all(FLERR,"Illegal neigh_modify command");
|
|
iarg += 2;
|
|
} else if (strcmp(arg[iarg],"page") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
old_pgsize = pgsize;
|
|
pgsize = force->inumeric(FLERR,arg[iarg+1]);
|
|
iarg += 2;
|
|
} else if (strcmp(arg[iarg],"one") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
old_oneatom = oneatom;
|
|
oneatom = force->inumeric(FLERR,arg[iarg+1]);
|
|
iarg += 2;
|
|
} else if (strcmp(arg[iarg],"binsize") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
binsize_user = force->numeric(FLERR,arg[iarg+1]);
|
|
if (binsize_user <= 0.0) binsizeflag = 0;
|
|
else binsizeflag = 1;
|
|
iarg += 2;
|
|
} else if (strcmp(arg[iarg],"cluster") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
if (strcmp(arg[iarg+1],"yes") == 0) cluster_check = 1;
|
|
else if (strcmp(arg[iarg+1],"no") == 0) cluster_check = 0;
|
|
else error->all(FLERR,"Illegal neigh_modify command");
|
|
iarg += 2;
|
|
|
|
} else if (strcmp(arg[iarg],"include") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
includegroup = group->find(arg[iarg+1]);
|
|
if (includegroup < 0)
|
|
error->all(FLERR,"Invalid group ID in neigh_modify command");
|
|
if (includegroup && (atom->firstgroupname == NULL ||
|
|
strcmp(arg[iarg+1],atom->firstgroupname) != 0))
|
|
error->all(FLERR,
|
|
"Neigh_modify include group != atom_modify first group");
|
|
iarg += 2;
|
|
|
|
} else if (strcmp(arg[iarg],"exclude") == 0) {
|
|
if (iarg+2 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
|
|
if (strcmp(arg[iarg+1],"type") == 0) {
|
|
if (iarg+4 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
if (nex_type == maxex_type) {
|
|
maxex_type += EXDELTA;
|
|
memory->grow(ex1_type,maxex_type,"neigh:ex1_type");
|
|
memory->grow(ex2_type,maxex_type,"neigh:ex2_type");
|
|
}
|
|
ex1_type[nex_type] = force->inumeric(FLERR,arg[iarg+2]);
|
|
ex2_type[nex_type] = force->inumeric(FLERR,arg[iarg+3]);
|
|
nex_type++;
|
|
iarg += 4;
|
|
|
|
} else if (strcmp(arg[iarg+1],"group") == 0) {
|
|
if (iarg+4 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
if (nex_group == maxex_group) {
|
|
maxex_group += EXDELTA;
|
|
memory->grow(ex1_group,maxex_group,"neigh:ex1_group");
|
|
memory->grow(ex2_group,maxex_group,"neigh:ex2_group");
|
|
}
|
|
ex1_group[nex_group] = group->find(arg[iarg+2]);
|
|
ex2_group[nex_group] = group->find(arg[iarg+3]);
|
|
if (ex1_group[nex_group] == -1 || ex2_group[nex_group] == -1)
|
|
error->all(FLERR,"Invalid group ID in neigh_modify command");
|
|
nex_group++;
|
|
iarg += 4;
|
|
|
|
} else if (strcmp(arg[iarg+1],"molecule") == 0) {
|
|
if (iarg+3 > narg) error->all(FLERR,"Illegal neigh_modify command");
|
|
if (atom->molecule_flag == 0)
|
|
error->all(FLERR,"Neigh_modify exclude molecule "
|
|
"requires atom attribute molecule");
|
|
if (nex_mol == maxex_mol) {
|
|
maxex_mol += EXDELTA;
|
|
memory->grow(ex_mol_group,maxex_mol,"neigh:ex_mol_group");
|
|
}
|
|
ex_mol_group[nex_mol] = group->find(arg[iarg+2]);
|
|
if (ex_mol_group[nex_mol] == -1)
|
|
error->all(FLERR,"Invalid group ID in neigh_modify command");
|
|
nex_mol++;
|
|
iarg += 3;
|
|
|
|
} else if (strcmp(arg[iarg+1],"none") == 0) {
|
|
nex_type = nex_group = nex_mol = 0;
|
|
iarg += 2;
|
|
|
|
} else error->all(FLERR,"Illegal neigh_modify command");
|
|
|
|
} else error->all(FLERR,"Illegal neigh_modify command");
|
|
}
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
remove the first group-group exclusion matching group1, group2
|
|
------------------------------------------------------------------------- */
|
|
|
|
void Neighbor::exclusion_group_group_delete(int group1, int group2)
|
|
{
|
|
int m, mlast;
|
|
for (m = 0; m < nex_group; m++)
|
|
if (ex1_group[m] == group1 && ex2_group[m] == group2 )
|
|
break;
|
|
|
|
mlast = m;
|
|
if (mlast == nex_group)
|
|
error->all(FLERR,"Unable to find group-group exclusion");
|
|
|
|
for (m = mlast+1; m < nex_group; m++) {
|
|
ex1_group[m-1] = ex1_group[m];
|
|
ex2_group[m-1] = ex2_group[m];
|
|
ex1_bit[m-1] = ex1_bit[m];
|
|
ex2_bit[m-1] = ex2_bit[m];
|
|
}
|
|
nex_group--;
|
|
}
|
|
|
|
|
|
/* ----------------------------------------------------------------------
|
|
return the value of exclude - used to check compatibility with GPU
|
|
------------------------------------------------------------------------- */
|
|
|
|
int Neighbor::exclude_setting()
|
|
{
|
|
return exclude;
|
|
}
|
|
|
|
/* ----------------------------------------------------------------------
|
|
return # of bytes of allocated memory
|
|
------------------------------------------------------------------------- */
|
|
|
|
bigint Neighbor::memory_usage()
|
|
{
|
|
bigint bytes = 0;
|
|
bytes += memory->usage(xhold,maxhold,3);
|
|
|
|
for (int i = 0; i < nlist; i++)
|
|
if (lists[i]) bytes += lists[i]->memory_usage();
|
|
for (int i = 0; i < nstencil; i++)
|
|
bytes += neigh_stencil[i]->memory_usage();
|
|
for (int i = 0; i < nbin; i++)
|
|
bytes += neigh_bin[i]->memory_usage();
|
|
|
|
if (neigh_bond) bytes += neigh_bond->memory_usage();
|
|
if (neigh_angle) bytes += neigh_angle->memory_usage();
|
|
if (neigh_dihedral) bytes += neigh_dihedral->memory_usage();
|
|
if (neigh_improper) bytes += neigh_improper->memory_usage();
|
|
|
|
return bytes;
|
|
}
|