diff --git a/src/nstencil.cpp b/src/nstencil.cpp index 3518b7db55..525369d573 100644 --- a/src/nstencil.cpp +++ b/src/nstencil.cpp @@ -31,33 +31,33 @@ using namespace LAMMPS_NS; sx,sy,sz = bin bounds = furthest the stencil could possibly extend calculated below in create_setup() 3d creates xyz stencil, 2d creates xy stencil - for half list with newton off: + for full list or half list with newton off + use a full stencil stencil is all surrounding bins including self regardless of triclinic - for half list with newton on: + for half list with newton on + use a half stencil stencil is bins to the "upper right" of central bin stencil does not include self no versions that allow ghost on (no callers need it?) for half list with newton on and triclinic: + use a half stencil stencil is all bins in z-plane of self and above, but not below in 2d is all bins in y-plane of self and above, but not below stencil includes self no versions that allow ghost on (no callers need it?) - for full list: - stencil is all surrounding bins including self - regardless of newton on/off or triclinic for multi: create one stencil for each atom type stencil follows same rules for half/full, newton on/off, triclinic cutoff is not cutneighmaxsq, but max cutoff for that atom type no versions that allow ghost on (any need for it?) - for multi/tiered: + for multi2: create one stencil for each itype-jtype pairing stencils do not generally follow the same rules for half/full or newton on/off - whole stencils including all surrounding bins are always used except - for same-type stencils with newton on which uses a split stencil - for orthogonal boxes, a split stencil includes bins to the "upper right" of central bin - for triclinic, a split stencil includes bins in the z (3D) or y (2D) plane of self and above + full stencils including all surrounding bins are always used except + for same-type stencils with newton on which uses a half stencil + for orthogonal boxes, a half stencil includes bins to the "upper right" of central bin + for triclinic, a half stencil includes bins in the z (3D) or y (2D) plane of self and above cutoff is not cutneighmaxsq, but max cutoff for that atom type no versions that allow ghost on (any need for it?) ------------------------------------------------------------------------- */ @@ -74,7 +74,7 @@ NStencil::NStencil(LAMMPS *lmp) : Pointers(lmp) stencil_multi = nullptr; distsq_multi = nullptr; - stencil_split = nullptr; + stencil_half = nullptr; stencil_skip = nullptr; stencil_bin_type = nullptr; stencil_cut = nullptr; @@ -101,31 +101,31 @@ NStencil::~NStencil() delete [] distsq_multi; } - if (stencil_multi_tiered) { + if (stencil_multi2) { int n = atom->ntypes; - memory->destroy(nstencil_multi_tiered); + memory->destroy(nstencil_multi2); for (int i = 1; i <= n; i++) { for (int j = 0; j <= n; j++) { if (! stencil_skip[i][j]) - memory->destroy(stencil_multi_tiered[i][j]); + memory->destroy(stencil_multi2[i][j]); } - delete [] stencil_multi_tiered[i]; + delete [] stencil_multi2[i]; } - delete [] stencil_multi_tiered; - memory->destroy(maxstencil_multi_tiered); - memory->destroy(stencil_split); + delete [] stencil_multi2; + memory->destroy(maxstencil_multi2); + memory->destroy(stencil_half); memory->destroy(stencil_skip); memory->destroy(stencil_bin_type); memory->destroy(stencil_cut); - memory->destroy(sx_multi_tiered); - memory->destroy(sy_multi_tiered); - memory->destroy(sz_multi_tiered); + memory->destroy(sx_multi2); + memory->destroy(sy_multi2); + memory->destroy(sz_multi2); - memory->destroy(binsizex_multi_tiered); - memory->destroy(binsizey_multi_tiered); - memory->destroy(binsizez_multi_tiered); + memory->destroy(binsizex_multi2); + memory->destroy(binsizey_multi2); + memory->destroy(binsizez_multi2); } } @@ -179,17 +179,17 @@ void NStencil::copy_bin_info() copy needed info for a given type from NBin class to this stencil class ------------------------------------------------------------------------- */ -void NStencil::copy_bin_info_multi_tiered(int type) +void NStencil::copy_bin_info_multi2(int type) { - mbinx = nb->mbinx_tiered[type]; - mbiny = nb->mbiny_tiered[type]; - mbinz = nb->mbinz_tiered[type]; - binsizex = nb->binsizex_tiered[type]; - binsizey = nb->binsizey_tiered[type]; - binsizez = nb->binsizez_tiered[type]; - bininvx = nb->bininvx_tiered[type]; - bininvy = nb->bininvy_tiered[type]; - bininvz = nb->bininvz_tiered[type]; + mbinx = nb->mbinx2[type]; + mbiny = nb->mbiny2[type]; + mbinz = nb->mbinz2[type]; + binsizex = nb->binsizex2[type]; + binsizey = nb->binsizey2[type]; + binsizez = nb->binsizez2[type]; + bininvx = nb->bininvx2[type]; + bininvy = nb->bininvy2[type]; + bininvz = nb->bininvz2[type]; } /* ---------------------------------------------------------------------- @@ -200,7 +200,7 @@ void NStencil::copy_bin_info_multi_tiered(int type) void NStencil::create_setup() { - if (neighstyle != Neighbor::MULTI_TIERED){ + if (neighstyle != Neighbor::MULTI2){ if (nb) copy_bin_info(); last_stencil = update->ntimestep; @@ -263,34 +263,53 @@ void NStencil::create_setup() int n = atom->ntypes; // Allocate arrays to store stencil information - memory->create(stencil_split, n, n, - "neighstencil:stencil_split");" - memory->create(stencil_skip, n, n, - "neighstencil:stencil_skip");" - memory->create(stencil_bin_type, n, n, - "neighstencil:stencil_bin_type");" - memory->create(stencil_cut, n, n, - "neighstencil:stencil_cut");" + memory->create(stencil_half, n+1, n+1, + "neighstencil:stencil_half"); + memory->create(stencil_skip, n+1, n+1, + "neighstencil:stencil_skip"); + memory->create(stencil_bin_type, n+1, n+1, + "neighstencil:stencil_bin_type"); + memory->create(stencil_cut, n+1, n+1, + "neighstencil:stencil_cut"); - memory->create(sx_multi_tiered, n, n, - "neighstencil:sx_multi_tiered");" - memory->create(sy_multi_tiered, n, n, - "neighstencil:sy_multi_tiered");" - memory->create(sz_multi_tiered, n, n, - "neighstencil:sz_multi_tiered");" + memory->create(sx_multi2, n+1, n+1, "neighstencil:sx_multi2"); + memory->create(sy_multi2, n+1, n+1, "neighstencil:sy_multi2"); + memory->create(sz_multi2, n+1, n+1, "neighstencil:sz_multi2"); - memory->create(binsizex_multi_tiered, n, n, - "neighstencil:binsizex_multi_tiered");" - memory->create(binsizey_multi_tiered, n, n, - "neighstencil:binsizey_multi_tiered");" - memory->create(binsizez_multi_tiered, n, n, - "neighstencil:binsizez_multi_tiered");" - + memory->create(binsizex_multi2, n+1, n+1, + "neighstencil:binsizex_multi2"); + memory->create(binsizey_multi2, n+1, n+1, + "neighstencil:binsizey_multi2"); + memory->create(binsizez_multi2, n+1, n+1, + "neighstencil:binsizez_multi2"); + + // Skip all stencils by default, initialize smax + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { + stencil_skip[i][j] = 1; + maxstencil_multi2[i][j] = 0; + } + } + // Determine which stencils need to be built set_stencil_properties(); - for (i = 1; i <= n; ++i) { - for (j = 1; j <= n; ++j) { + // Allocate arrays to store stencils + + if (!maxstencil_multi2) { + memory->create(maxstencil_multi2, n+1, n+1, "neighstencil::stencil_multi2"); + memory->create(nstencil_multi2, n+1, n+1, "neighstencil::nstencil_multi2"); + stencil_multi2 = new int**[n+1](); + for (i = 1; i <= n; ++i) { + stencil_multi2[i] = new int*[n+1](); + for (j = 1; j <= n; ++j) { + maxstencil_multi2[i][j] = 0; + } + } + } + + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { // Skip creation of unused stencils if (stencil_skip[i][j]) continue; @@ -299,9 +318,9 @@ void NStencil::create_setup() bin_type = stencil_bin_type[i][j]; copy_bin_info_bytype(bin_type); - binsizex_multi_tiered[i][j] = binsizex; - binsizey_multi_tiered[i][j] = binsizey; - binsizez_multi_tiered[i][j] = binsizez; + binsizex_multi2[i][j] = binsizex; + binsizey_multi2[i][j] = binsizey; + binsizez_multi2[i][j] = binsizez; stencil_range = stencil_cut[i][j]; @@ -312,17 +331,17 @@ void NStencil::create_setup() sz = static_cast (stencil_range*bininvz); if (sz*binsizez < stencil_range) sz++; - sx_multi_tiered[i][j] = sx; - sy_multi_tiered[i][j] = sy; - sz_multi_tiered[i][j] = sz; + sx_multi2[i][j] = sx; + sy_multi2[i][j] = sy; + sz_multi2[i][j] = sz; smax = ((2*sx+1) * (2*sy+1) * (2*sz+1)); - if (smax > maxstencil_multi_tiered[i][j]) { - maxstencil_multi_tiered[i][j] = smax; - memory->destroy(stencil_multi_tiered[i][j]); - memory->create(stencil_multi_tiered[i][j], smax, - "neighstencil::stencil_multi_tiered"); + if (smax > maxstencil_multi2[i][j]) { + maxstencil_multi2[i][j] = smax; + memory->destroy(stencil_multi2[i][j]); + memory->create(stencil_multi2[i][j], smax, + "neighstencil::stencil_multi2"); } } } @@ -363,10 +382,18 @@ double NStencil::memory_usage() } else if (neighstyle == Neighbor::MULTI) { bytes += atom->ntypes*maxstencil_multi * sizeof(int); bytes += atom->ntypes*maxstencil_multi * sizeof(double); - } else if (neighstyle == Neighbor::MULTI_TIERED) { - bytes += atom->ntypes*maxstencil_multi * sizeof(int); - bytes += atom->ntypes*maxstencil_multi * sizeof(int); - bytes += atom->ntypes*maxstencil_multi * sizeof(double); + } else if (neighstyle == Neighbor::MULTI2) { + int n = atom->ntypes; + for (i = 1; i <= n; ++i) { + for (j = 1; j <= n; ++j) { + bytes += maxstencil_multi2[i][j] * sizeof(int); + bytes += maxstencil_multi2[i][j] * sizeof(int); + bytes += maxstencil_multi2[i][j] * sizeof(double); + } + } + bytes += 2 * n * n * sizeof(bool); + bytes += 6 * n * n * sizeof(int); + bytes += 4 * n * n * sizeof(double); } return bytes; } diff --git a/src/nstencil.h b/src/nstencil.h index 58c39ee13c..7c6a04a56a 100644 --- a/src/nstencil.h +++ b/src/nstencil.h @@ -29,14 +29,17 @@ class NStencil : protected Pointers { int **stencilxyz; // bin offsets in xyz dims int *nstencil_multi; // # bins in each type-based multi stencil int **stencil_multi; // list of bin offsets in each stencil - int ** nstencil_multi_tiered; // # bins bins in each itype-jtype tiered multi stencil - int *** stencil_multi_tiered; // list of bin offsets in each tiered multi stencil - int ** maxstencil_type; // max double **distsq_multi; // sq distances to bins in each stencil + int ** nstencil_multi2; // # bins bins in each itype-jtype multi2 stencil + int *** stencil_multi2; // list of bin offsets in each multi2 stencil + int ** maxstencil_multi2; // max stencil size for each multi2 stencil + + ^do i need a multi 2 analog to distsq_multi? + int sx,sy,sz; // extent of stencil in each dim - int **sx_multi_tiered; // analogs for multi tiered - int **sy_multi_tiered; - int **sz_multi_tiered; + int **sx_multi2; // analogs for multi tiered + int **sy_multi2; + int **sz_multi2; double cutoff_custom; // cutoff set by requestor @@ -75,9 +78,9 @@ class NStencil : protected Pointers { // analogs for multi-tiered - double **binsizex_multi_tiered; - double **binsizey_multi_tiered; - double **binsizez_multi_tiered; + double **binsizex_multi2; + double **binsizey_multi2; + double **binsizez_multi2; // data common to all NStencil variants @@ -94,7 +97,7 @@ class NStencil : protected Pointers { // methods for multi/tiered NStencil - void copy_bin_info_multi_tiered(int); // copy mult/tiered info from NBin class + void copy_bin_info_multi2(int); // copy mult/tiered info from NBin class void set_stencil_properties(); // determine which stencils to build and how }; diff --git a/src/nstencil_full_multi2_2d.cpp b/src/nstencil_full_multi2_2d.cpp new file mode 100644 index 0000000000..f62cddb0b5 --- /dev/null +++ b/src/nstencil_full_multi2_2d.cpp @@ -0,0 +1,97 @@ +/* ---------------------------------------------------------------------- + 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. +------------------------------------------------------------------------- */ + +#include "nstencil_full_multi2_2d.h" +#include "neighbor.h" +#include "neigh_list.h" +#include "nbin.h" +#include "memory.h" +#include "atom.h" +#include + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NStencilFullMulti22d::NStencilFullMulti22d(LAMMPS *lmp) : NStencil(lmp) {} + +/* ---------------------------------------------------------------------- */ + +void NStencilFullMulti22d::set_stencil_properties() +{ + int n = atom->ntypes; + int i, j; + + // like -> like => use half stencil + for (i = 1; i <= n; i++) { + stencil_half[i][i] = 0; + stencil_skip[i][i] = 0; + stencil_bin_type[i][i] = i; + stencil_cut[i][i] = sqrt(cutneighsq[i][i]); + } + + // smaller -> larger => use existing newtoff stencil in larger bin + // larger -> smaller => use multi-like stencil for small-large in smaller bin + // If types are same cutoff, use existing like-like stencil. + + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { + if(i == j) continue; + + stencil_half[i][j] = 0; + stencil_skip[i][j] = 0; + stencil_cut[i][j] = sqrt(cutneighsq[i][j]); + + if(cuttypesq[i] <= cuttypesq[j]){ + stencil_bin_type[i][j] = i; + } else { + stencil_bin_type[i][j] = j; + } + } + } +} + +/* ---------------------------------------------------------------------- + create stencils based on bin geometry and cutoff +------------------------------------------------------------------------- */ + +void NStencilFullMulti22d::create() +{ + int itype, jtype, i, j, k, ns; + int n = atom->ntypes; + double cutsq; + + + for (itype = 1; itype <= n; itype++) { + for (jtype = 1; jtype <= n; jtype++) { + if (stencil_skip[itype][jtype]) continue; + + ns = 0; + + sx = sx_multi2[itype][jtype]; + sy = sy_multi2[itype][jtype]; + + mbinx = mbinx_multi2[itype][jtype]; + mbiny = mbiny_multi2[itype][jtype]; + + cutsq = stencil_cut[itype][jtype]; + + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,0) < cutsq) + stencil_type[itype][jtype][ns++] = j*mbinx + i; + + nstencil_multi2[itype][jtype] = ns; + } + } +} diff --git a/src/nstencil_half_multi2_3d_delete_noff.h b/src/nstencil_full_multi2_2d.h similarity index 58% rename from src/nstencil_half_multi2_3d_delete_noff.h rename to src/nstencil_full_multi2_2d.h index bf83a31e9a..1dfc5ae4d6 100644 --- a/src/nstencil_half_multi2_3d_delete_noff.h +++ b/src/nstencil_full_multi2_2d.h @@ -13,32 +13,27 @@ #ifdef NSTENCIL_CLASS -NStencilStyle(half/bytype/3d/newtoff, - NStencilHalfBytype3dNewtoff, - NS_HALF | NS_BYTYPE | NS_3D | NS_NEWTOFF | NS_ORTHO | NS_TRI) +NStencilStyle(full/multi2/2d, + NStencilFullMulti22d, NS_FULL | NS_Multi2 | NS_2D | NS_ORTHO | NS_TRI) #else -#ifndef LMP_NSTENCIL_HALF_BYTYPE_3D_NEWTOFF_H -#define LMP_NSTENCIL_HALF_BYTYPE_3D_NEWTOFF_H +#ifndef LMP_NSTENCIL_FULL_MULTI2_2D_H +#define LMP_NSTENCIL_FULL_MULTI2_2D_H #include "nstencil.h" namespace LAMMPS_NS { -class NStencilHalfBytype3dNewtoff : public NStencil { +class NStencilFullMulti22d : public NStencil { public: - NStencilHalfBytype3dNewtoff(class LAMMPS *); - ~NStencilHalfBytype3dNewtoff(); - void create_setup(); + NStencilFullMulti22d(class LAMMPS *); + ~NStencilFullMulti22d(); void create(); - private: - int ** maxstencil_type; + protected: + void set_stencil_properties(); - void copy_bin_info_bytype(int); - int copy_neigh_info_bytype(int); - void create_newtoff(int, int, double); }; } diff --git a/src/nstencil_full_multi2_3d.cpp b/src/nstencil_full_multi2_3d.cpp index 099e2adeb8..77ae62738c 100644 --- a/src/nstencil_full_multi2_3d.cpp +++ b/src/nstencil_full_multi2_3d.cpp @@ -11,7 +11,7 @@ See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ -#include "nstencil_full_bytype_3d.h" +#include "nstencil_full_multi2_3d.h" #include "neighbor.h" #include "neigh_list.h" #include "nbin.h" @@ -23,164 +23,79 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -NStencilFullBytype3d::NStencilFullBytype3d(LAMMPS *lmp) : - NStencil(lmp) +NStencilFullMulti23d::NStencilFullMulti23d(LAMMPS *lmp) : NStencil(lmp) {} + +/* ---------------------------------------------------------------------- */ + +void NStencilFullMulti23d::set_stencil_properties() { - maxstencil_type = NULL; -} - -/* ---------------------------------------------------------------------- */ - -NStencilFullBytype3d::~NStencilFullBytype3d() { - - if (maxstencil_type) { - memory->destroy(nstencil_type); - for (int i = 1; i <= atom->ntypes; i++) { - for (int j = 0; j <= atom->ntypes; j++) { - if (maxstencil_type[i][j] > 0) memory->destroy(stencil_type[i][j]); - } - delete [] stencil_type[i]; - } - delete [] stencil_type; - memory->destroy(maxstencil_type); - } -} - - -/* ---------------------------------------------------------------------- */ -INCORPORTE INTO CREATE THEN DELETE, NOTE NAME OF CUTNEIGHMAX ETC -int NStencilFullBytype3d::copy_neigh_info_bytype(int itype) { - - cutneighmaxsq = neighbor->cutneighsq[itype][itype]; - cutneighmax = sqrt(cutneighmaxsq); - cuttypesq = neighbor->cuttypesq; - - // sx,sy,sz = max range of stencil in each dim - // smax = max possible size of entire 3d stencil - // stencil will be empty if cutneighmax = 0.0 - - sx = static_cast (cutneighmax*bininvx); - if (sx*binsizex < cutneighmax) sx++; - sy = static_cast (cutneighmax*bininvy); - if (sy*binsizey < cutneighmax) sy++; - sz = static_cast (cutneighmax*bininvz); - if (sz*binsizez < cutneighmax) sz++; - - return ((2*sx+1) * (2*sy+1) * (2*sz+1)); -} - -/* ---------------------------------------------------------------------- */ - -void NStencilFullBytype3d::create_setup() -{ - - int itype, jtype; - int maxtypes; - int smax; - - maxtypes = atom->ntypes; - - if (maxstencil_type == NULL) { - memory->create(maxstencil_type, maxtypes+1, maxtypes+1, "BAD A"); - memory->create(nstencil_type, maxtypes+1, maxtypes+1, "BAD B"); - stencil_type = new int**[maxtypes+1](); - for (itype = 1; itype <= maxtypes; ++itype) { - stencil_type[itype] = new int*[maxtypes+1](); - for (jtype = 1; jtype <= maxtypes; ++jtype) { - maxstencil_type[itype][jtype] = 0; - } - } - } MOVE TO PARENT CLASS - - // like -> like => use standard newtoff stencil at bin - - for (itype = 1; itype <= maxtypes; ++itype) { - copy_bin_info_bytype(itype); - smax = copy_neigh_info_bytype(itype); uses cutneighsq[itype][itype] to create s's - if (smax > maxstencil_type[itype][itype]) { - maxstencil_type[itype][itype] = smax; - memory->destroy(stencil_type[itype][itype]); - memory->create(stencil_type[itype][itype], smax, - "NStencilFullBytype::create_steup() stencil"); - } - create_newtoff(itype, itype, cutneighmaxsq); + int n = atom->ntypes; + int i, j; + + // like -> like => use half stencil + for (i = 1; i <= n; i++) { + stencil_half[i][i] = 0; + stencil_skip[i][i] = 0; + stencil_bin_type[i][i] = i; + stencil_cut[i][i] = sqrt(cutneighsq[i][i]); } // smaller -> larger => use existing newtoff stencil in larger bin // larger -> smaller => use multi-like stencil for small-large in smaller bin // If types are same cutoff, use existing like-like stencil. - for (itype = 1; itype <= maxtypes; ++itype) { - for (jtype = 1; jtype <= maxtypes; ++jtype) { - if (itype == jtype) continue; - if (cuttypesq[itype] <= cuttypesq[jtype]) { This does work to say which prticle is smller - // Potential destroy/create problem? - nstencil_type[itype][jtype] = nstencil_type[jtype][jtype]; - stencil_type[itype][jtype] = stencil_type[jtype][jtype]; - } - else { - copy_bin_info_bytype(jtype); - // smax = copy_neigh_info_bytype(jtype); + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { + if(i == j) continue; - cutneighmaxsq = cuttypesq[jtype]; Does it need to be this big? Can't I use cutneighsq[i][j]? - cutneighmax = sqrt(cutneighmaxsq); - sx = static_cast (cutneighmax*bininvx); - if (sx*binsizex < cutneighmax) sx++; - sy = static_cast (cutneighmax*bininvy); - if (sy*binsizey < cutneighmax) sy++; - sz = static_cast (cutneighmax*bininvz); - if (sz*binsizez < cutneighmax) sz++; - - smax = (2*sx+1) * (2*sy+1) * (2*sz+1); - if (smax > maxstencil_type[itype][jtype]) { - maxstencil_type[itype][jtype] = smax; - memory->destroy(stencil_type[itype][jtype]); - memory->create(stencil_type[itype][jtype], smax, "Bad C"); - } - create_newtoff(itype, jtype, cuttypesq[jtype]); + stencil_half[i][j] = 0; + stencil_skip[i][j] = 0; + stencil_cut[i][j] = sqrt(cutneighsq[i][j]); + + if(cuttypesq[i] <= cuttypesq[j]){ + stencil_bin_type[i][j] = i; + } else { + stencil_bin_type[i][j] = j; } } } - - //for (itype = 1; itype <= maxtypes; itype++) { - // for (jtype = 1; jtype <= maxtypes; jtype++) { - // printf("i j n %d %d %d\n", itype, jtype, nstencil_type[itype][jtype]); - // printf("i j n %d %d %d\n", itype, jtype, maxstencil_type[itype][jtype]); - // } - // } - -} - -/* ---------------------------------------------------------------------- */ - -void NStencilFullBytype3d::create_newtoff(int itype, int jtype, double cutsq) { - - int i, j, k, ns; - - ns = 0; - - for (k = -sz; k <= sz; k++) - for (j = -sy; j <= sy; j++) - for (i = -sx; i <= sx; i++) - if (bin_distance(i,j,k) < cutsq) - stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; - - nstencil_type[itype][jtype] = ns; } /* ---------------------------------------------------------------------- - create stencil based on bin geometry and cutoff + create stencils based on bin geometry and cutoff ------------------------------------------------------------------------- */ -void NStencilFullBytype3d::create() +void NStencilFullMulti23d::create() { - //int i,j,k; - - //nstencil = 0; - - //for (k = -sz; k <= sz; k++) - // for (j = -sy; j <= sy; j++) - // for (i = -sx; i <= sx; i++) - // if (bin_distance(i,j,k) < cutneighmaxsq) - // stencil[nstencil++] = k*mbiny*mbinx + j*mbinx + i; + int itype, jtype, i, j, k, ns; + int n = atom->ntypes; + double cutsq; + + + for (itype = 1; itype <= n; itype++) { + for (jtype = 1; jtype <= n; jtype++) { + if (stencil_skip[itype][jtype]) continue; + + ns = 0; + + sx = sx_multi2[itype][jtype]; + sy = sy_multi2[itype][jtype]; + sz = sz_multi2[itype][jtype]; + + mbinx = mbinx_multi2[itype][jtype]; + mbiny = mbiny_multi2[itype][jtype]; + mbinz = mbinz_multi2[itype][jtype]; + + cutsq = stencil_cut[itype][jtype]; + + for (k = -sz; k <= sz; k++) + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,k) < cutsq) + stencil_type[itype][jtype][ns++] = + k*mbiny*mbinx + j*mbinx + i; + + nstencil_multi2[itype][jtype] = ns; + } + } } diff --git a/src/nstencil_full_multi2_3d.h b/src/nstencil_full_multi2_3d.h index 48bdfd2611..f7baefc8e9 100644 --- a/src/nstencil_full_multi2_3d.h +++ b/src/nstencil_full_multi2_3d.h @@ -13,30 +13,26 @@ #ifdef NSTENCIL_CLASS -NStencilStyle(full/multi/tiered/3d, - NStencilFullMultiTiered3d, NS_FULL | NS_Multi_Tiered | NS_3D | NS_ORTHO | NS_TRI) +NStencilStyle(full/multi2/3d, + NStencilFullMulti23d, NS_FULL | NS_Multi2 | NS_3D | NS_ORTHO | NS_TRI) #else -#ifndef LMP_NSTENCIL_FULL_MULTI_TIERED_3D_H -#define LMP_NSTENCIL_FULL_MULTI_TIERED_3D_H +#ifndef LMP_NSTENCIL_FULL_MULTI2_3D_H +#define LMP_NSTENCIL_FULL_MULTI2_3D_H #include "nstencil.h" namespace LAMMPS_NS { -class NStencilFullBytype3d : public NStencil { +class NStencilFullMulti23d : public NStencil { public: - NStencilFullBytype3d(class LAMMPS *); - ~NStencilFullBytype3d(); + NStencilFullMulti23d(class LAMMPS *); + ~NStencilFullMulti23d(); void create(); - void create_setup(); -private: - int ** maxstencil_type; - - int copy_neigh_info_bytype(int); - void create_newtoff(int, int, double); + protected: + void set_stencil_properties(); }; diff --git a/src/nstencil_half_multi2_2d.cpp b/src/nstencil_half_multi2_2d.cpp new file mode 100644 index 0000000000..86da6ea91e --- /dev/null +++ b/src/nstencil_half_multi2_2d.cpp @@ -0,0 +1,111 @@ +/* ---------------------------------------------------------------------- + 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. +------------------------------------------------------------------------- */ + +#include "nstencil_half_multi2_2d.h" +#include "neighbor.h" +#include "neigh_list.h" +#include "nbin.h" +#include "memory.h" +#include "atom.h" +#include + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NStencilHalfMulti22d::NStencilHalfMulti23d(LAMMPS *lmp) : + NStencil(lmp) {} + +/* ---------------------------------------------------------------------- */ + +void NStencilHalfMulti23d::set_stencil_properties() +{ + int n = atom->ntypes; + int i, j; + + // like -> like => use half stencil + for (i = 1; i <= n; i++) { + stencil_half[i][i] = 1; + stencil_skip[i][i] = 0; + stencil_bin_type[i][i] = i; + stencil_cut[i][i] = sqrt(cutneighsq[i][i]); + } + + // Cross types: use full stencil, looking one way through hierarchy + // smaller -> larger => use full stencil in larger bin + // larger -> smaller => no nstecil required + // If cut offs are same, use existing type-type stencil + + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { + if(i == j) continue; + if(cuttypesq[i] > cuttypesq[j]) continue; + + stencil_skip[i][j] = 0; + stencil_cut[i][j] = sqrt(cutneighsq[i][j]); + + if(cuttypesq[i] == cuttypesq[j]){ + stencil_half[i][j] = 1; + stencil_bin_type[i][j] = i; + } else { + stencil_half[i][j] = 0; + stencil_bin_type[i][j] = j; + } + } + } +} + +/* ---------------------------------------------------------------------- + create stencils based on bin geometry and cutoff +------------------------------------------------------------------------- */ + +void NStencilHalfMulti23d::create() +{ + int itype, jtype, i, j, ns; + int n = atom->ntypes; + double cutsq; + + + for (itype = 1; itype <= n; itype++) { + for (jtype = 1; jtype <= n; jtype++) { + if (stencil_skip[itype][jtype]) continue; + + ns = 0; + + sx = sx_multi2[itype][jtype]; + sy = sy_multi2[itype][jtype]; + + mbinx = mbinx_multi2[itype][jtype]; + mbiny = mbiny_multi2[itype][jtype]; + + cutsq = stencil_cut[itype][jtype]; + + if (stencil_half[itype][jtype]) { + for (j = 0; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (j > 0 || (j == 0 && i > 0)) { + if (bin_distance(i,j,0) < cutsq) + stencil_type[itype][jtype][ns++] = j*mbinx + i; + } + } else { + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,0) < cutsq) + stencil_type[itype][jtype][ns++] = j*mbinx + i; + } + + nstencil_multi2[itype][jtype] = ns; + } + } +} + diff --git a/src/nstencil_half_multi2_2d.h b/src/nstencil_half_multi2_2d.h new file mode 100644 index 0000000000..ba3a471833 --- /dev/null +++ b/src/nstencil_half_multi2_2d.h @@ -0,0 +1,46 @@ +/* -*- 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. +------------------------------------------------------------------------- */ + +#ifdef NSTENCIL_CLASS + +NStencilStyle(half/multi2/2d, + NStencilHalfMulti22d, NS_HALF | NS_MULTI2 | NS_2D | NS_ORTHO) + +#else + +#ifndef LMP_NSTENCIL_HALF_MULTI2_2D_H +#define LMP_NSTENCIL_HALF_MULTI2_2D_H + +#include "nstencil.h" + +namespace LAMMPS_NS { + +class NStencilHalfMulti22d : public NStencil { + public: + NStencilHalfMulti22d(class LAMMPS *); + ~NStencilHalfMulti22d(); + void create(); + + protected: + void set_stencil_properties(); + +}; + +} + +#endif +#endif + +/* ERROR/WARNING messages: + +*/ diff --git a/src/nstencil_half_multi2_2d_tri.cpp b/src/nstencil_half_multi2_2d_tri.cpp new file mode 100755 index 0000000000..46f55300ce --- /dev/null +++ b/src/nstencil_half_multi2_2d_tri.cpp @@ -0,0 +1,108 @@ +/* ---------------------------------------------------------------------- + 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. +------------------------------------------------------------------------- */ + +#include "nstencil_half_multi2_2d_tri.h" +#include "neighbor.h" +#include "neigh_list.h" +#include "nbin.h" +#include "memory.h" +#include "atom.h" +#include + +using namespace LAMMPS_NS; + +/* ---------------------------------------------------------------------- */ + +NStencilHalfMulti22dTri::NStencilHalfMulti22d(LAMMPS *lmp) : + NStencil(lmp) {} + +/* ---------------------------------------------------------------------- */ + +void NStencilHalfMulti22d::set_stencil_properties() +{ + int n = atom->ntypes; + int i, j; + + // like -> like => use half stencil + for (i = 1; i <= n; i++) { + stencil_half[i][i] = 1; + stencil_skip[i][i] = 0; + stencil_bin_type[i][i] = i; + stencil_cut[i][i] = sqrt(cutneighsq[i][i]); + } + + // Cross types: use full stencil, looking one way through hierarchy + // smaller -> larger => use full stencil in larger bin + // larger -> smaller => no nstecil required + // If cut offs are same, use existing type-type stencil + + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { + if(i == j) continue; + if(cuttypesq[i] > cuttypesq[j]) continue; + + stencil_skip[i][j] = 0; + stencil_cut[i][j] = sqrt(cutneighsq[i][j]); + + if(cuttypesq[i] == cuttypesq[j]){ + stencil_half[i][j] = 1; + stencil_bin_type[i][j] = i; + } else { + stencil_half[i][j] = 0; + stencil_bin_type[i][j] = j; + } + } + } +} + +/* ---------------------------------------------------------------------- + create stencils based on bin geometry and cutoff +------------------------------------------------------------------------- */ + +void NStencilHalfMulti22d::create() +{ + int itype, jtype, i, j, ns; + int n = atom->ntypes; + double cutsq; + + + for (itype = 1; itype <= n; itype++) { + for (jtype = 1; jtype <= n; jtype++) { + if (stencil_skip[itype][jtype]) continue; + + ns = 0; + + sx = sx_multi2[itype][jtype]; + sy = sy_multi2[itype][jtype]; + + mbinx = mbinx_multi2[itype][jtype]; + mbiny = mbiny_multi2[itype][jtype]; + + cutsq = stencil_cut[itype][jtype]; + + if (stencil_half[itype][jtype]) { + for (j = 0; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,0) < cutsq) + stencil_type[itype][jtype][ns++] = j*mbinx + i; + } else { + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,0) < cutsq) + stencil_type[itype][jtype][ns++] = j*mbinx + i; + } + + nstencil_multi2[itype][jtype] = ns; + } + } +} diff --git a/src/nstencil_half_multi2_2d_tri.h b/src/nstencil_half_multi2_2d_tri.h new file mode 100755 index 0000000000..1d25e56776 --- /dev/null +++ b/src/nstencil_half_multi2_2d_tri.h @@ -0,0 +1,45 @@ +/* -*- 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. +------------------------------------------------------------------------- */ + +#ifdef NSTENCIL_CLASS + +NStencilStyle(half/multi2/2d/tri, + NStencilHalfMulti22dTri, NS_HALF | NS_MULTI2 | NS_2D | NS_TRI) + +#else + +#ifndef LMP_NSTENCIL_HALF_MULTI2_2D_TRI_H +#define LMP_NSTENCIL_HALF_MULTI2_2D_TRI_H + +#include "nstencil.h" + +namespace LAMMPS_NS { + +class NStencilHalfMulti22dTri : public NStencil { + public: + NStencilHalfMulti22dTri(class LAMMPS *); + ~NStencilHalfMulti22dTri(); + void create(); + + protected: + void set_stencil_properties(); +}; + +} + +#endif +#endif + +/* ERROR/WARNING messages: + +*/ diff --git a/src/nstencil_half_multi2_3d.cpp b/src/nstencil_half_multi2_3d.cpp index 7fc53fe762..02a65bd3a1 100644 --- a/src/nstencil_half_multi2_3d.cpp +++ b/src/nstencil_half_multi2_3d.cpp @@ -11,7 +11,7 @@ See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ -#include "nstencil_half_bytype_3d_newton.h" +#include "nstencil_half_multi2_3d.h" #include "neighbor.h" #include "neigh_list.h" #include "nbin.h" @@ -23,183 +23,95 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -NStencilHalfBytype3dNewton::NStencilHalfBytype3dNewton(LAMMPS *lmp) : - NStencil(lmp) +NStencilHalfMulti23d::NStencilHalfMulti23d(LAMMPS *lmp) : + NStencil(lmp) {} + +/* ---------------------------------------------------------------------- */ + +void NStencilHalfMulti23d::set_stencil_properties() { - maxstencil_type = NULL; -} - -/* ---------------------------------------------------------------------- */ - -NStencilHalfBytype3dNewton::~NStencilHalfBytype3dNewton() { - - memory->destroy(nstencil_type); - for (int i = 1; i <= atom->ntypes; i++) { - for (int j = 0; j <= atom->ntypes; j++) { - if (maxstencil_type[i][j] > 0) memory->destroy(stencil_type[i][j]); - } - delete [] stencil_type[i]; - } - delete [] stencil_type; - memory->destroy(maxstencil_type); -} - -/* ---------------------------------------------------------------------- */ - -// KS To superclass -void NStencilHalfBytype3dNewton::copy_bin_info_bytype(int itype) { - - mbinx = nb->mbinx_type[itype]; - mbiny = nb->mbiny_type[itype]; - mbinz = nb->mbinz_type[itype]; - binsizex = nb->binsizex_type[itype]; - binsizey = nb->binsizey_type[itype]; - binsizez = nb->binsizez_type[itype]; - bininvx = nb->bininvx_type[itype]; - bininvy = nb->bininvy_type[itype]; - bininvz = nb->bininvz_type[itype]; -} - -/* ---------------------------------------------------------------------- */ - -// KS To superclass? -int NStencilHalfBytype3dNewton::copy_neigh_info_bytype(int itype) { - - cutneighmaxsq = neighbor->cutneighsq[itype][itype]; - cutneighmax = sqrt(cutneighmaxsq); - cuttypesq = neighbor->cuttypesq; - - // sx,sy,sz = max range of stencil in each dim - // smax = max possible size of entire 3d stencil - // stencil will be empty if cutneighmax = 0.0 - - sx = static_cast (cutneighmax*bininvx); - if (sx*binsizex < cutneighmax) sx++; - sy = static_cast (cutneighmax*bininvy); - if (sy*binsizey < cutneighmax) sy++; - sz = static_cast (cutneighmax*bininvz); - if (sz*binsizez < cutneighmax) sz++; - - return ((2*sx+1) * (2*sy+1) * (2*sz+1)); -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewton::create_setup() -{ - - int itype, jtype; - int maxtypes; - int smax; - - // maxstencil_type to superclass? - maxtypes = atom->ntypes; - - if (maxstencil_type == NULL) { - memory->create(maxstencil_type, maxtypes+1, maxtypes+1, "maxstencil_type"); - memory->create(nstencil_type, maxtypes+1, maxtypes+1, "nstencil_type"); - stencil_type = new int**[maxtypes+1](); - for (itype = 1; itype <= maxtypes; ++itype) { - stencil_type[itype] = new int*[maxtypes+1](); - for (jtype = 1; jtype <= maxtypes; ++jtype) { - maxstencil_type[itype][jtype] = 0; - nstencil_type[itype][jtype] = 0; - } - } + int n = atom->ntypes; + int i, j; + + // like -> like => use half stencil + for (i = 1; i <= n; i++) { + stencil_half[i][i] = 1; + stencil_skip[i][i] = 0; + stencil_bin_type[i][i] = i; + stencil_cut[i][i] = sqrt(cutneighsq[i][i]); } - // like -> like => use standard Newton stencil at bin - - for (itype = 1; itype <= maxtypes; ++itype) { - copy_bin_info_bytype(itype); - smax = copy_neigh_info_bytype(itype); - if (smax > maxstencil_type[itype][itype]) { - maxstencil_type[itype][itype] = smax; - memory->destroy(stencil_type[itype][itype]); - memory->create(stencil_type[itype][itype], smax, - "NStencilHalfBytypeNewton::create_steup() stencil"); - } - create_newton(itype, itype, cutneighmaxsq); - } - - // Cross types: "Newton on" reached by using Newton off stencil and - // looking one way through hierarchy - // smaller -> larger => use Newton off stencil in larger bin + // Cross types: use full stencil, looking one way through hierarchy + // smaller -> larger => use full stencil in larger bin // larger -> smaller => no nstecil required // If cut offs are same, use existing type-type stencil - for (itype = 1; itype <= maxtypes; ++itype) { - for (jtype = 1; jtype <= maxtypes; ++jtype) { - if (itype == jtype) continue; - if (cuttypesq[itype] == cuttypesq[jtype]) { - nstencil_type[itype][jtype] = nstencil_type[jtype][jtype]; - stencil_type[itype][jtype] = stencil_type[jtype][jtype]; - } - else if (cuttypesq[itype] < cuttypesq[jtype]) { - copy_bin_info_bytype(jtype); + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { + if(i == j) continue; + if(cuttypesq[i] > cuttypesq[j]) continue; - cutneighmaxsq = cuttypesq[jtype]; - cutneighmax = sqrt(cutneighmaxsq); - sx = static_cast (cutneighmax*bininvx); - if (sx*binsizex < cutneighmax) sx++; - sy = static_cast (cutneighmax*bininvy); - if (sy*binsizey < cutneighmax) sy++; - sz = static_cast (cutneighmax*bininvz); - if (sz*binsizez < cutneighmax) sz++; - - smax = (2*sx+1) * (2*sy+1) * (2*sz+1); - if (smax > maxstencil_type[itype][jtype]) { - maxstencil_type[itype][jtype] = smax; - memory->destroy(stencil_type[itype][jtype]); - memory->create(stencil_type[itype][jtype], smax, "stencil_type[]"); - } - create_newtoff(itype, jtype, cuttypesq[jtype]); + stencil_skip[i][j] = 0; + stencil_cut[i][j] = sqrt(cutneighsq[i][j]); + + if(cuttypesq[i] == cuttypesq[j]){ + stencil_half[i][j] = 1; + stencil_bin_type[i][j] = i; + } else { + stencil_half[i][j] = 0; + stencil_bin_type[i][j] = j; } } } - -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewton::create_newton(int itype, int jtype, double cutsq) { - - int i, j, k, ns; - - ns = 0; - - for (k = 0; k <= sz; k++) - for (j = -sy; j <= sy; j++) - for (i = -sx; i <= sx; i++) - if (k > 0 || j > 0 || (j == 0 && i > 0)) { - if (bin_distance(i,j,k) < cutsq) - stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; - } - nstencil_type[itype][jtype] = ns; -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewton::create_newtoff(int itype, int jtype, double cutsq) { - - int i, j, k, ns; - - ns = 0; - - for (k = -sz; k <= sz; k++) - for (j = -sy; j <= sy; j++) - for (i = -sx; i <= sx; i++) - if (bin_distance(i,j,k) < cutsq) - stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; - - nstencil_type[itype][jtype] = ns; } /* ---------------------------------------------------------------------- - create stencil based on bin geometry and cutoff + create stencils based on bin geometry and cutoff ------------------------------------------------------------------------- */ -void NStencilHalfBytype3dNewton::create() +void NStencilHalfMulti23d::create() { - // KS. Move "creation" here. + int itype, jtype, i, j, k, ns; + int n = atom->ntypes; + double cutsq; + + + for (itype = 1; itype <= n; itype++) { + for (jtype = 1; jtype <= n; jtype++) { + if (stencil_skip[itype][jtype]) continue; + + ns = 0; + + sx = sx_multi2[itype][jtype]; + sy = sy_multi2[itype][jtype]; + sz = sz_multi2[itype][jtype]; + + mbinx = mbinx_multi2[itype][jtype]; + mbiny = mbiny_multi2[itype][jtype]; + mbinz = mbinz_multi2[itype][jtype]; + + cutsq = stencil_cut[itype][jtype]; + + if (stencil_half[itype][jtype]) { + for (k = 0; k <= sz; k++) + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (k > 0 || j > 0 || (j == 0 && i > 0)) { + if (bin_distance(i,j,k) < cutsq) + stencil_type[itype][jtype][ns++] = + k*mbiny*mbinx + j*mbinx + i; + } + } else { + for (k = -sz; k <= sz; k++) + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,k) < cutsq) + stencil_type[itype][jtype][ns++] = + k*mbiny*mbinx + j*mbinx + i; + } + + nstencil_multi2[itype][jtype] = ns; + } + } } + diff --git a/src/nstencil_half_multi2_3d.h b/src/nstencil_half_multi2_3d.h index e39db1c1e1..9cae8269cb 100644 --- a/src/nstencil_half_multi2_3d.h +++ b/src/nstencil_half_multi2_3d.h @@ -13,32 +13,27 @@ #ifdef NSTENCIL_CLASS -NStencilStyle(half/bytype/3d, - NStencilHalfBytype3dNewton, NS_HALF | NS_BYTYPE | NS_3D | NS_ORTHO) +NStencilStyle(half/multi2/3d, + NStencilHalfMulti23d, NS_HALF | NS_MULTI2 | NS_3D | NS_ORTHO) #else -#ifndef LMP_NSTENCIL_HALF_BYTYPE_3D_H -#define LMP_NSTENCIL_HALF_BYTYPE_3D_H +#ifndef LMP_NSTENCIL_HALF_MULTI2_3D_H +#define LMP_NSTENCIL_HALF_MULTI2_3D_H #include "nstencil.h" namespace LAMMPS_NS { -class NStencilHalfBytype3dNewton : public NStencil { +class NStencilHalfMulti23d : public NStencil { public: - NStencilHalfBytype3dNewton(class LAMMPS *); - ~NStencilHalfBytype3dNewton(); - void create_setup(); + NStencilHalfMulti23d(class LAMMPS *); + ~NStencilHalfMulti23d(); void create(); - private: - int ** maxstencil_type; + protected: + void set_stencil_properties(); - void copy_bin_info_bytype(int); - int copy_neigh_info_bytype(int); - void create_newton(int, int, double); - void create_newtoff(int, int, double); }; } diff --git a/src/nstencil_half_multi2_3d_delete_noff.cpp b/src/nstencil_half_multi2_3d_delete_noff.cpp deleted file mode 100644 index fedc8060ce..0000000000 --- a/src/nstencil_half_multi2_3d_delete_noff.cpp +++ /dev/null @@ -1,198 +0,0 @@ -/* ---------------------------------------------------------------------- - 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. -------------------------------------------------------------------------- */ - -#include "nstencil_half_bytype_3d_newtoff.h" -#include "neighbor.h" -#include "neigh_list.h" -#include "nbin.h" -#include "memory.h" -#include "atom.h" -#include - -using namespace LAMMPS_NS; - -/* ---------------------------------------------------------------------- */ - -NStencilHalfBytype3dNewtoff::NStencilHalfBytype3dNewtoff(LAMMPS *lmp) : - NStencil(lmp) -{ - maxstencil_type = NULL; -} - -/* ---------------------------------------------------------------------- */ - -NStencilHalfBytype3dNewtoff::~NStencilHalfBytype3dNewtoff() { - - memory->destroy(nstencil_type); - for (int i = 1; i <= atom->ntypes; i++) { - for (int j = 0; j <= atom->ntypes; j++) { - if (maxstencil_type[i][j] > 0) memory->destroy(stencil_type[i][j]); - } - delete [] stencil_type[i]; - } - delete [] stencil_type; - memory->destroy(maxstencil_type); -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewtoff::copy_bin_info_bytype(int itype) { - - mbinx = nb->mbinx_type[itype]; - mbiny = nb->mbiny_type[itype]; - mbinz = nb->mbinz_type[itype]; - binsizex = nb->binsizex_type[itype]; - binsizey = nb->binsizey_type[itype]; - binsizez = nb->binsizez_type[itype]; - bininvx = nb->bininvx_type[itype]; - bininvy = nb->bininvy_type[itype]; - bininvz = nb->bininvz_type[itype]; -} - -/* ---------------------------------------------------------------------- */ - -int NStencilHalfBytype3dNewtoff::copy_neigh_info_bytype(int itype) { - - cutneighmaxsq = neighbor->cutneighsq[itype][itype]; - cutneighmax = sqrt(cutneighmaxsq); - cuttypesq = neighbor->cuttypesq; - - // sx,sy,sz = max range of stencil in each dim - // smax = max possible size of entire 3d stencil - // stencil will be empty if cutneighmax = 0.0 - - sx = static_cast (cutneighmax*bininvx); - if (sx*binsizex < cutneighmax) sx++; - sy = static_cast (cutneighmax*bininvy); - if (sy*binsizey < cutneighmax) sy++; - sz = static_cast (cutneighmax*bininvz); - if (sz*binsizez < cutneighmax) sz++; - - return ((2*sx+1) * (2*sy+1) * (2*sz+1)); -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewtoff::create_setup() -{ - - int itype, jtype; - int maxtypes; - int smax; - - maxtypes = atom->ntypes; - - if (maxstencil_type == NULL) { - memory->create(maxstencil_type, maxtypes+1, maxtypes+1, "BAD A"); - memory->create(nstencil_type, maxtypes+1, maxtypes+1, "BAD B"); - stencil_type = new int**[maxtypes+1](); - for (itype = 1; itype <= maxtypes; ++itype) { - stencil_type[itype] = new int*[maxtypes+1](); - for (jtype = 1; jtype <= maxtypes; ++jtype) { - maxstencil_type[itype][jtype] = 0; - } - } - } - - // like -> like => use standard newtoff stencil at bin - - for (itype = 1; itype <= maxtypes; ++itype) { - copy_bin_info_bytype(itype); - smax = copy_neigh_info_bytype(itype); - if (smax > maxstencil_type[itype][itype]) { - maxstencil_type[itype][itype] = smax; - memory->destroy(stencil_type[itype][itype]); - memory->create(stencil_type[itype][itype], smax, - "NStencilHalfBytypeNewtoff::create_steup() stencil"); - } - create_newtoff(itype, itype, cutneighmaxsq); - } - - // smaller -> larger => use existing newtoff stencil in larger bin - // larger -> smaller => use multi-like stencil for small-large in smaller bin - // If types are same cutoff, use existing like-like stencil. - - for (itype = 1; itype <= maxtypes; ++itype) { - for (jtype = 1; jtype <= maxtypes; ++jtype) { - if (itype == jtype) continue; - if (cuttypesq[itype] <= cuttypesq[jtype]) { - // Potential destroy/create problem? - nstencil_type[itype][jtype] = nstencil_type[jtype][jtype]; - stencil_type[itype][jtype] = stencil_type[jtype][jtype]; - } - else { - copy_bin_info_bytype(jtype); - // smax = copy_neigh_info_bytype(jtype); - - cutneighmaxsq = cuttypesq[jtype]; - cutneighmax = sqrt(cutneighmaxsq); - sx = static_cast (cutneighmax*bininvx); - if (sx*binsizex < cutneighmax) sx++; - sy = static_cast (cutneighmax*bininvy); - if (sy*binsizey < cutneighmax) sy++; - sz = static_cast (cutneighmax*bininvz); - if (sz*binsizez < cutneighmax) sz++; - - smax = (2*sx+1) * (2*sy+1) * (2*sz+1); - if (smax > maxstencil_type[itype][jtype]) { - maxstencil_type[itype][jtype] = smax; - memory->destroy(stencil_type[itype][jtype]); - memory->create(stencil_type[itype][jtype], smax, "Bad C"); - } - create_newtoff(itype, jtype, cuttypesq[jtype]); - } - } - } - - //for (itype = 1; itype <= maxtypes; itype++) { - // for (jtype = 1; jtype <= maxtypes; jtype++) { - // printf("i j n %d %d %d\n", itype, jtype, nstencil_type[itype][jtype]); - // printf("i j n %d %d %d\n", itype, jtype, maxstencil_type[itype][jtype]); - // } - // } - -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewtoff::create_newtoff(int itype, int jtype, double cutsq) { - - int i, j, k, ns; - - ns = 0; - - for (k = -sz; k <= sz; k++) - for (j = -sy; j <= sy; j++) - for (i = -sx; i <= sx; i++) - if (bin_distance(i,j,k) < cutsq) - stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; - - nstencil_type[itype][jtype] = ns; -} - -/* ---------------------------------------------------------------------- - create stencil based on bin geometry and cutoff -------------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewtoff::create() -{ - //int i,j,k; - - //nstencil = 0; - - //for (k = -sz; k <= sz; k++) - // for (j = -sy; j <= sy; j++) - // for (i = -sx; i <= sx; i++) - // if (bin_distance(i,j,k) < cutneighmaxsq) - // stencil[nstencil++] = k*mbiny*mbinx + j*mbinx + i; -} diff --git a/src/nstencil_half_multi2_3d_tri.cpp b/src/nstencil_half_multi2_3d_tri.cpp index f1ab713725..bcaf22abbd 100755 --- a/src/nstencil_half_multi2_3d_tri.cpp +++ b/src/nstencil_half_multi2_3d_tri.cpp @@ -11,7 +11,7 @@ See the README file in the top-level LAMMPS directory. ------------------------------------------------------------------------- */ -#include "nstencil_half_bytype_3d_newton.h" +#include "nstencil_half_multi2_3d_tri.h" #include "neighbor.h" #include "neigh_list.h" #include "nbin.h" @@ -23,31 +23,100 @@ using namespace LAMMPS_NS; /* ---------------------------------------------------------------------- */ -NStencilHalfBytype3dNewton::NStencilHalfBytype3dNewton(LAMMPS *lmp) : - NStencil(lmp) -{ - maxstencil_type = NULL; -} +NStencilHalfMulti23dTri::NStencilHalfMulti23d(LAMMPS *lmp) : + NStencil(lmp) {} /* ---------------------------------------------------------------------- */ -NStencilHalfBytype3dNewton::~NStencilHalfBytype3dNewton() { - - memory->destroy(nstencil_type); - for (int i = 1; i <= atom->ntypes; i++) { - for (int j = 0; j <= atom->ntypes; j++) { - if (maxstencil_type[i][j] > 0) memory->destroy(stencil_type[i][j]); - } - delete [] stencil_type[i]; +void NStencilHalfMulti23d::set_stencil_properties() +{ + int n = atom->ntypes; + int i, j; + + // like -> like => use half stencil + for (i = 1; i <= n; i++) { + stencil_half[i][i] = 1; + stencil_skip[i][i] = 0; + stencil_bin_type[i][i] = i; + stencil_cut[i][i] = sqrt(cutneighsq[i][i]); + } + + // Cross types: use full stencil, looking one way through hierarchy + // smaller -> larger => use full stencil in larger bin + // larger -> smaller => no nstecil required + // If cut offs are same, use existing type-type stencil + + for (i = 1; i <= n; i++) { + for (j = 1; j <= n; j++) { + if(i == j) continue; + if(cuttypesq[i] > cuttypesq[j]) continue; + + stencil_skip[i][j] = 0; + stencil_cut[i][j] = sqrt(cutneighsq[i][j]); + + if(cuttypesq[i] == cuttypesq[j]){ + stencil_half[i][j] = 1; + stencil_bin_type[i][j] = i; + } else { + stencil_half[i][j] = 0; + stencil_bin_type[i][j] = j; + } + } + } +} + +/* ---------------------------------------------------------------------- + create stencils based on bin geometry and cutoff +------------------------------------------------------------------------- */ + +void NStencilHalfMulti23d::create() +{ + int itype, jtype, i, j, k, ns; + int n = atom->ntypes; + double cutsq; + + + for (itype = 1; itype <= n; itype++) { + for (jtype = 1; jtype <= n; jtype++) { + if (stencil_skip[itype][jtype]) continue; + + ns = 0; + + sx = sx_multi2[itype][jtype]; + sy = sy_multi2[itype][jtype]; + sz = sz_multi2[itype][jtype]; + + mbinx = mbinx_multi2[itype][jtype]; + mbiny = mbiny_multi2[itype][jtype]; + mbinz = mbinz_multi2[itype][jtype]; + + cutsq = stencil_cut[itype][jtype]; + + if (stencil_half[itype][jtype]) { + for (k = 0; k <= sz; k++) + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,k) < cutsq) + stencil_type[itype][jtype][ns++] = + k*mbiny*mbinx + j*mbinx + i; + } else { + for (k = -sz; k <= sz; k++) + for (j = -sy; j <= sy; j++) + for (i = -sx; i <= sx; i++) + if (bin_distance(i,j,k) < cutsq) + stencil_type[itype][jtype][ns++] = + k*mbiny*mbinx + j*mbinx + i; + } + + nstencil_multi2[itype][jtype] = ns; + } } - delete [] stencil_type; - memory->destroy(maxstencil_type); } /* ---------------------------------------------------------------------- */ // KS To superclass -void NStencilHalfBytype3dNewton::copy_bin_info_bytype(int itype) { +void NStencilHalfMulti23d::copy_bin_info_bytype(int itype) { mbinx = nb->mbinx_type[itype]; mbiny = nb->mbiny_type[itype]; @@ -63,7 +132,7 @@ void NStencilHalfBytype3dNewton::copy_bin_info_bytype(int itype) { /* ---------------------------------------------------------------------- */ // KS To superclass? -int NStencilHalfBytype3dNewton::copy_neigh_info_bytype(int itype) { +int NStencilHalfMulti23d::copy_neigh_info_bytype(int itype) { cutneighmaxsq = neighbor->cutneighsq[itype][itype]; cutneighmax = sqrt(cutneighmaxsq); @@ -85,7 +154,7 @@ int NStencilHalfBytype3dNewton::copy_neigh_info_bytype(int itype) { /* ---------------------------------------------------------------------- */ -void NStencilHalfBytype3dNewton::create_setup() +void NStencilHalfMulti23d::create_setup() { int itype, jtype; @@ -157,47 +226,4 @@ void NStencilHalfBytype3dNewton::create_setup() } } } - -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewton::create_newton(int itype, int jtype, double cutsq) { - - int i, j, k, ns; - - ns = 0; - - for (k = 0; k <= sz; k++) - for (j = -sy; j <= sy; j++) - for (i = -sx; i <= sx; i++) - if (bin_distance(i,j,k) < cutsq) - stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; - nstencil_type[itype][jtype] = ns; -} - -/* ---------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewton::create_newtoff(int itype, int jtype, double cutsq) { - - int i, j, k, ns; - - ns = 0; - - for (k = -sz; k <= sz; k++) - for (j = -sy; j <= sy; j++) - for (i = -sx; i <= sx; i++) - if (bin_distance(i,j,k) < cutsq) - stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; - - nstencil_type[itype][jtype] = ns; -} - -/* ---------------------------------------------------------------------- - create stencil based on bin geometry and cutoff -------------------------------------------------------------------------- */ - -void NStencilHalfBytype3dNewton::create() -{ - // KS. Move "creation" here. } diff --git a/src/nstencil_half_multi2_3d_tri.h b/src/nstencil_half_multi2_3d_tri.h index 5424d73f11..f4b37ebebe 100755 --- a/src/nstencil_half_multi2_3d_tri.h +++ b/src/nstencil_half_multi2_3d_tri.h @@ -13,32 +13,26 @@ #ifdef NSTENCIL_CLASS -NStencilStyle(half/bytype/3d, - NStencilHalfBytype3d, NS_HALF | NS_BYTYPE | NS_3D | NS_TRI) +NStencilStyle(half/multi2/3d/tri, + NStencilHalfMulti23dTri, NS_HALF | NS_MULTI2 | NS_3D | NS_TRI) #else -#ifndef LMP_NSTENCIL_HALF_BYTYPE_3D_H -#define LMP_NSTENCIL_HALF_BYTYPE_3D_H +#ifndef LMP_NSTENCIL_HALF_MULTI2_3D_TRI_H +#define LMP_NSTENCIL_HALF_MULTI2_3D_TRI_H #include "nstencil.h" namespace LAMMPS_NS { -class NStencilHalfBytype3dNewton : public NStencil { +class NStencilHalfMulti23dTri : public NStencil { public: - NStencilHalfBytype3dNewton(class LAMMPS *); - ~NStencilHalfBytype3dNewton(); - void create_setup(); + NStencilHalfMulti23dTri(class LAMMPS *); + ~NStencilHalfMulti23dTri(); void create(); - private: - int ** maxstencil_type; - - void copy_bin_info_bytype(int); - int copy_neigh_info_bytype(int); - void create_newton(int, int, double); - void create_newtoff(int, int, double); + protected: + void set_stencil_properties(); }; }