ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Finish stencil classes

Browse Source
This commit is contained in:
Joel Clemmer
2020-11-10 20:03:00 -07:00
parent af57879416
commit dd1cce1da5
15 changed files with 779 additions and 707 deletions
Download Patch File Download Diff File Expand all files Collapse all files

167
src/nstencil.cpp
View File

@ -31,33 +31,33 @@ using namespace LAMMPS_NS;
sx,sy,sz = bin bounds = furthest the stencil could possibly extend sx,sy,sz = bin bounds = furthest the stencil could possibly extend
calculated below in create_setup() calculated below in create_setup()
3d creates xyz stencil, 2d creates xy stencil 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 stencil is all surrounding bins including self
regardless of triclinic 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 is bins to the "upper right" of central bin
stencil does not include self stencil does not include self
no versions that allow ghost on (no callers need it?) no versions that allow ghost on (no callers need it?)
for half list with newton on and triclinic: 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 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 in 2d is all bins in y-plane of self and above, but not below
stencil includes self stencil includes self
no versions that allow ghost on (no callers need it?) 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: for multi:
create one stencil for each atom type create one stencil for each atom type
stencil follows same rules for half/full, newton on/off, triclinic stencil follows same rules for half/full, newton on/off, triclinic
cutoff is not cutneighmaxsq, but max cutoff for that atom type cutoff is not cutneighmaxsq, but max cutoff for that atom type
no versions that allow ghost on (any need for it?) no versions that allow ghost on (any need for it?)
for multi/tiered: for multi2:
create one stencil for each itype-jtype pairing create one stencil for each itype-jtype pairing
stencils do not generally follow the same rules for half/full or newton on/off 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 full stencils including all surrounding bins are always used except
for same-type stencils with newton on which uses a split stencil for same-type stencils with newton on which uses a half stencil
for orthogonal boxes, a split stencil includes bins to the "upper right" of central bin for orthogonal boxes, a half 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 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 cutoff is not cutneighmaxsq, but max cutoff for that atom type
no versions that allow ghost on (any need for it?) no versions that allow ghost on (any need for it?)
------------------------------------------------------------------------- */ ------------------------------------------------------------------------- */
@ -74,7 +74,7 @@ NStencil::NStencil(LAMMPS *lmp) : Pointers(lmp)
stencil_multi = nullptr; stencil_multi = nullptr;
distsq_multi = nullptr; distsq_multi = nullptr;
stencil_split = nullptr; stencil_half = nullptr;
stencil_skip = nullptr; stencil_skip = nullptr;
stencil_bin_type = nullptr; stencil_bin_type = nullptr;
stencil_cut = nullptr; stencil_cut = nullptr;
@ -101,31 +101,31 @@ NStencil::~NStencil()
delete [] distsq_multi; delete [] distsq_multi;
} }
if (stencil_multi_tiered) { if (stencil_multi2) {
int n = atom->ntypes; int n = atom->ntypes;
memory->destroy(nstencil_multi_tiered); memory->destroy(nstencil_multi2);
for (int i = 1; i <= n; i++) { for (int i = 1; i <= n; i++) {
for (int j = 0; j <= n; j++) { for (int j = 0; j <= n; j++) {
if (! stencil_skip[i][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; delete [] stencil_multi2;
memory->destroy(maxstencil_multi_tiered); memory->destroy(maxstencil_multi2);
memory->destroy(stencil_split); memory->destroy(stencil_half);
memory->destroy(stencil_skip); memory->destroy(stencil_skip);
memory->destroy(stencil_bin_type); memory->destroy(stencil_bin_type);
memory->destroy(stencil_cut); memory->destroy(stencil_cut);
memory->destroy(sx_multi_tiered); memory->destroy(sx_multi2);
memory->destroy(sy_multi_tiered); memory->destroy(sy_multi2);
memory->destroy(sz_multi_tiered); memory->destroy(sz_multi2);
memory->destroy(binsizex_multi_tiered); memory->destroy(binsizex_multi2);
memory->destroy(binsizey_multi_tiered); memory->destroy(binsizey_multi2);
memory->destroy(binsizez_multi_tiered); 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 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]; mbinx = nb->mbinx2[type];
mbiny = nb->mbiny_tiered[type]; mbiny = nb->mbiny2[type];
mbinz = nb->mbinz_tiered[type]; mbinz = nb->mbinz2[type];
binsizex = nb->binsizex_tiered[type]; binsizex = nb->binsizex2[type];
binsizey = nb->binsizey_tiered[type]; binsizey = nb->binsizey2[type];
binsizez = nb->binsizez_tiered[type]; binsizez = nb->binsizez2[type];
bininvx = nb->bininvx_tiered[type]; bininvx = nb->bininvx2[type];
bininvy = nb->bininvy_tiered[type]; bininvy = nb->bininvy2[type];
bininvz = nb->bininvz_tiered[type]; bininvz = nb->bininvz2[type];
} }
/* ---------------------------------------------------------------------- /* ----------------------------------------------------------------------
@ -200,7 +200,7 @@ void NStencil::copy_bin_info_multi_tiered(int type)
void NStencil::create_setup() void NStencil::create_setup()
{ {
if (neighstyle != Neighbor::MULTI_TIERED){ if (neighstyle != Neighbor::MULTI2){
if (nb) copy_bin_info(); if (nb) copy_bin_info();
last_stencil = update->ntimestep; last_stencil = update->ntimestep;
@ -263,34 +263,53 @@ void NStencil::create_setup()
int n = atom->ntypes; int n = atom->ntypes;
// Allocate arrays to store stencil information // Allocate arrays to store stencil information
memory->create(stencil_split, n, n, memory->create(stencil_half, n+1, n+1,
"neighstencil:stencil_split");" "neighstencil:stencil_half");
memory->create(stencil_skip, n, n, memory->create(stencil_skip, n+1, n+1,
"neighstencil:stencil_skip");" "neighstencil:stencil_skip");
memory->create(stencil_bin_type, n, n, memory->create(stencil_bin_type, n+1, n+1,
"neighstencil:stencil_bin_type");" "neighstencil:stencil_bin_type");
memory->create(stencil_cut, n, n, memory->create(stencil_cut, n+1, n+1,
"neighstencil:stencil_cut");" "neighstencil:stencil_cut");
memory->create(sx_multi_tiered, n, n, memory->create(sx_multi2, n+1, n+1, "neighstencil:sx_multi2");
"neighstencil:sx_multi_tiered");" memory->create(sy_multi2, n+1, n+1, "neighstencil:sy_multi2");
memory->create(sy_multi_tiered, n, n, memory->create(sz_multi2, n+1, n+1, "neighstencil:sz_multi2");
"neighstencil:sy_multi_tiered");"
memory->create(sz_multi_tiered, n, n,
"neighstencil:sz_multi_tiered");"
memory->create(binsizex_multi_tiered, n, n, memory->create(binsizex_multi2, n+1, n+1,
"neighstencil:binsizex_multi_tiered");" "neighstencil:binsizex_multi2");
memory->create(binsizey_multi_tiered, n, n, memory->create(binsizey_multi2, n+1, n+1,
"neighstencil:binsizey_multi_tiered");" "neighstencil:binsizey_multi2");
memory->create(binsizez_multi_tiered, n, n, memory->create(binsizez_multi2, n+1, n+1,
"neighstencil:binsizez_multi_tiered");" "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 // Determine which stencils need to be built
set_stencil_properties(); set_stencil_properties();
// 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) { for (i = 1; i <= n; ++i) {
stencil_multi2[i] = new int*[n+1]();
for (j = 1; j <= n; ++j) { 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 // Skip creation of unused stencils
if (stencil_skip[i][j]) continue; if (stencil_skip[i][j]) continue;
@ -299,9 +318,9 @@ void NStencil::create_setup()
bin_type = stencil_bin_type[i][j]; bin_type = stencil_bin_type[i][j];
copy_bin_info_bytype(bin_type); copy_bin_info_bytype(bin_type);
binsizex_multi_tiered[i][j] = binsizex; binsizex_multi2[i][j] = binsizex;
binsizey_multi_tiered[i][j] = binsizey; binsizey_multi2[i][j] = binsizey;
binsizez_multi_tiered[i][j] = binsizez; binsizez_multi2[i][j] = binsizez;
stencil_range = stencil_cut[i][j]; stencil_range = stencil_cut[i][j];
@ -312,17 +331,17 @@ void NStencil::create_setup()
sz = static_cast<int> (stencil_range*bininvz); sz = static_cast<int> (stencil_range*bininvz);
if (sz*binsizez < stencil_range) sz++; if (sz*binsizez < stencil_range) sz++;
sx_multi_tiered[i][j] = sx; sx_multi2[i][j] = sx;
sy_multi_tiered[i][j] = sy; sy_multi2[i][j] = sy;
sz_multi_tiered[i][j] = sz; sz_multi2[i][j] = sz;
smax = ((2*sx+1) * (2*sy+1) * (2*sz+1)); smax = ((2*sx+1) * (2*sy+1) * (2*sz+1));
if (smax > maxstencil_multi_tiered[i][j]) { if (smax > maxstencil_multi2[i][j]) {
maxstencil_multi_tiered[i][j] = smax; maxstencil_multi2[i][j] = smax;
memory->destroy(stencil_multi_tiered[i][j]); memory->destroy(stencil_multi2[i][j]);
memory->create(stencil_multi_tiered[i][j], smax, memory->create(stencil_multi2[i][j], smax,
"neighstencil::stencil_multi_tiered"); "neighstencil::stencil_multi2");
} }
} }
} }
@ -363,10 +382,18 @@ double NStencil::memory_usage()
} else if (neighstyle == Neighbor::MULTI) { } else if (neighstyle == Neighbor::MULTI) {
bytes += atom->ntypes*maxstencil_multi * sizeof(int); bytes += atom->ntypes*maxstencil_multi * sizeof(int);
bytes += atom->ntypes*maxstencil_multi * sizeof(double); bytes += atom->ntypes*maxstencil_multi * sizeof(double);
} else if (neighstyle == Neighbor::MULTI_TIERED) { } else if (neighstyle == Neighbor::MULTI2) {
bytes += atom->ntypes*maxstencil_multi * sizeof(int); int n = atom->ntypes;
bytes += atom->ntypes*maxstencil_multi * sizeof(int); for (i = 1; i <= n; ++i) {
bytes += atom->ntypes*maxstencil_multi * sizeof(double); 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; return bytes;
} }

23
src/nstencil.h
View File

@ -29,14 +29,17 @@ class NStencil : protected Pointers {
int **stencilxyz; // bin offsets in xyz dims int **stencilxyz; // bin offsets in xyz dims
int *nstencil_multi; // # bins in each type-based multi stencil int *nstencil_multi; // # bins in each type-based multi stencil
int **stencil_multi; // list of bin offsets in each 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 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,sy,sz; // extent of stencil in each dim
int **sx_multi_tiered; // analogs for multi tiered int **sx_multi2; // analogs for multi tiered
int **sy_multi_tiered; int **sy_multi2;
int **sz_multi_tiered; int **sz_multi2;
double cutoff_custom; // cutoff set by requestor double cutoff_custom; // cutoff set by requestor
@ -75,9 +78,9 @@ class NStencil : protected Pointers {
// analogs for multi-tiered // analogs for multi-tiered
double **binsizex_multi_tiered; double **binsizex_multi2;
double **binsizey_multi_tiered; double **binsizey_multi2;
double **binsizez_multi_tiered; double **binsizez_multi2;
// data common to all NStencil variants // data common to all NStencil variants
@ -94,7 +97,7 @@ class NStencil : protected Pointers {
// methods for multi/tiered NStencil // 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 void set_stencil_properties(); // determine which stencils to build and how
}; };

97
src/nstencil_full_multi2_2d.cpp Normal file
View File

@ -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 <math.h>
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;
}
}
}

23
src/nstencil_half_multi2_3d_delete_noff.h → src/nstencil_full_multi2_2d.h
View File

@ -13,32 +13,27 @@
#ifdef NSTENCIL_CLASS #ifdef NSTENCIL_CLASS
NStencilStyle(half/bytype/3d/newtoff, NStencilStyle(full/multi2/2d,
NStencilHalfBytype3dNewtoff, NStencilFullMulti22d, NS_FULL | NS_Multi2 | NS_2D | NS_ORTHO | NS_TRI)
NS_HALF | NS_BYTYPE | NS_3D | NS_NEWTOFF | NS_ORTHO | NS_TRI)
#else #else
#ifndef LMP_NSTENCIL_HALF_BYTYPE_3D_NEWTOFF_H #ifndef LMP_NSTENCIL_FULL_MULTI2_2D_H
#define LMP_NSTENCIL_HALF_BYTYPE_3D_NEWTOFF_H #define LMP_NSTENCIL_FULL_MULTI2_2D_H
#include "nstencil.h" #include "nstencil.h"
namespace LAMMPS_NS { namespace LAMMPS_NS {
class NStencilHalfBytype3dNewtoff : public NStencil { class NStencilFullMulti22d : public NStencil {
public: public:
NStencilHalfBytype3dNewtoff(class LAMMPS *); NStencilFullMulti22d(class LAMMPS *);
~NStencilHalfBytype3dNewtoff(); ~NStencilFullMulti22d();
void create_setup();
void create(); void create();
private: protected:
int ** maxstencil_type; void set_stencil_properties();
void copy_bin_info_bytype(int);
int copy_neigh_info_bytype(int);
void create_newtoff(int, int, double);
}; };
} }

187
src/nstencil_full_multi2_3d.cpp
View File

@ -11,7 +11,7 @@
See the README file in the top-level LAMMPS directory. 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 "neighbor.h"
#include "neigh_list.h" #include "neigh_list.h"
#include "nbin.h" #include "nbin.h"
@ -23,164 +23,79 @@ using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
NStencilFullBytype3d::NStencilFullBytype3d(LAMMPS *lmp) : NStencilFullMulti23d::NStencilFullMulti23d(LAMMPS *lmp) : NStencil(lmp) {}
NStencil(lmp)
/* ---------------------------------------------------------------------- */
void NStencilFullMulti23d::set_stencil_properties()
{ {
maxstencil_type = NULL; int n = atom->ntypes;
} int i, j;
/* ---------------------------------------------------------------------- */ // like -> like => use half stencil
for (i = 1; i <= n; i++) {
NStencilFullBytype3d::~NStencilFullBytype3d() { stencil_half[i][i] = 0;
stencil_skip[i][i] = 0;
if (maxstencil_type) { stencil_bin_type[i][i] = i;
memory->destroy(nstencil_type); stencil_cut[i][i] = sqrt(cutneighsq[i][i]);
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<int> (cutneighmax*bininvx);
if (sx*binsizex < cutneighmax) sx++;
sy = static_cast<int> (cutneighmax*bininvy);
if (sy*binsizey < cutneighmax) sy++;
sz = static_cast<int> (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);
} }
// smaller -> larger => use existing newtoff stencil in larger bin // smaller -> larger => use existing newtoff stencil in larger bin
// larger -> smaller => use multi-like stencil for small-large in smaller bin // larger -> smaller => use multi-like stencil for small-large in smaller bin
// If types are same cutoff, use existing like-like stencil. // If types are same cutoff, use existing like-like stencil.
for (itype = 1; itype <= maxtypes; ++itype) { for (i = 1; i <= n; i++) {
for (jtype = 1; jtype <= maxtypes; ++jtype) { for (j = 1; j <= n; j++) {
if (itype == jtype) continue; if(i == j) 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);
cutneighmaxsq = cuttypesq[jtype]; Does it need to be this big? Can't I use cutneighsq[i][j]? stencil_half[i][j] = 0;
cutneighmax = sqrt(cutneighmaxsq); stencil_skip[i][j] = 0;
sx = static_cast<int> (cutneighmax*bininvx); stencil_cut[i][j] = sqrt(cutneighsq[i][j]);
if (sx*binsizex < cutneighmax) sx++;
sy = static_cast<int> (cutneighmax*bininvy);
if (sy*binsizey < cutneighmax) sy++;
sz = static_cast<int> (cutneighmax*bininvz);
if (sz*binsizez < cutneighmax) sz++;
smax = (2*sx+1) * (2*sy+1) * (2*sz+1); if(cuttypesq[i] <= cuttypesq[j]){
if (smax > maxstencil_type[itype][jtype]) { stencil_bin_type[i][j] = i;
maxstencil_type[itype][jtype] = smax; } else {
memory->destroy(stencil_type[itype][jtype]); stencil_bin_type[i][j] = j;
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]);
// }
// }
} }
/* ---------------------------------------------------------------------- */ /* ----------------------------------------------------------------------
create stencils based on bin geometry and cutoff
------------------------------------------------------------------------- */
void NStencilFullBytype3d::create_newtoff(int itype, int jtype, double cutsq) { void NStencilFullMulti23d::create()
{
int itype, jtype, i, j, k, ns;
int n = atom->ntypes;
double cutsq;
int i, j, k, ns;
for (itype = 1; itype <= n; itype++) {
for (jtype = 1; jtype <= n; jtype++) {
if (stencil_skip[itype][jtype]) continue;
ns = 0; 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 (k = -sz; k <= sz; k++)
for (j = -sy; j <= sy; j++) for (j = -sy; j <= sy; j++)
for (i = -sx; i <= sx; i++) for (i = -sx; i <= sx; i++)
if (bin_distance(i,j,k) < cutsq) if (bin_distance(i,j,k) < cutsq)
stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; stencil_type[itype][jtype][ns++] =
k*mbiny*mbinx + j*mbinx + i;
nstencil_type[itype][jtype] = ns; nstencil_multi2[itype][jtype] = ns;
} }
}
/* ----------------------------------------------------------------------
create stencil based on bin geometry and cutoff
------------------------------------------------------------------------- */
void NStencilFullBytype3d::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;
} }

22
src/nstencil_full_multi2_3d.h
View File

@ -13,30 +13,26 @@
#ifdef NSTENCIL_CLASS #ifdef NSTENCIL_CLASS
NStencilStyle(full/multi/tiered/3d, NStencilStyle(full/multi2/3d,
NStencilFullMultiTiered3d, NS_FULL | NS_Multi_Tiered | NS_3D | NS_ORTHO | NS_TRI) NStencilFullMulti23d, NS_FULL | NS_Multi2 | NS_3D | NS_ORTHO | NS_TRI)
#else #else
#ifndef LMP_NSTENCIL_FULL_MULTI_TIERED_3D_H #ifndef LMP_NSTENCIL_FULL_MULTI2_3D_H
#define LMP_NSTENCIL_FULL_MULTI_TIERED_3D_H #define LMP_NSTENCIL_FULL_MULTI2_3D_H
#include "nstencil.h" #include "nstencil.h"
namespace LAMMPS_NS { namespace LAMMPS_NS {
class NStencilFullBytype3d : public NStencil { class NStencilFullMulti23d : public NStencil {
public: public:
NStencilFullBytype3d(class LAMMPS *); NStencilFullMulti23d(class LAMMPS *);
~NStencilFullBytype3d(); ~NStencilFullMulti23d();
void create(); void create();
void create_setup();
private: protected:
int ** maxstencil_type; void set_stencil_properties();
int copy_neigh_info_bytype(int);
void create_newtoff(int, int, double);
}; };

111
src/nstencil_half_multi2_2d.cpp Normal file
View File

@ -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 <math.h>
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;
}
}
}

46
src/nstencil_half_multi2_2d.h Normal file
View File

@ -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:
*/

108
src/nstencil_half_multi2_2d_tri.cpp Executable file
View File

@ -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 <math.h>
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;
}
}
}

45
src/nstencil_half_multi2_2d_tri.h Executable file
View File

@ -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:
*/

210
src/nstencil_half_multi2_3d.cpp
View File

@ -11,7 +11,7 @@
See the README file in the top-level LAMMPS directory. 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 "neighbor.h"
#include "neigh_list.h" #include "neigh_list.h"
#include "nbin.h" #include "nbin.h"
@ -23,183 +23,95 @@ using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
NStencilHalfBytype3dNewton::NStencilHalfBytype3dNewton(LAMMPS *lmp) : NStencilHalfMulti23d::NStencilHalfMulti23d(LAMMPS *lmp) :
NStencil(lmp) NStencil(lmp) {}
/* ---------------------------------------------------------------------- */
void NStencilHalfMulti23d::set_stencil_properties()
{ {
maxstencil_type = NULL; int n = atom->ntypes;
} int i, j;
/* ---------------------------------------------------------------------- */ // like -> like => use half stencil
for (i = 1; i <= n; i++) {
NStencilHalfBytype3dNewton::~NStencilHalfBytype3dNewton() { stencil_half[i][i] = 1;
stencil_skip[i][i] = 0;
memory->destroy(nstencil_type); stencil_bin_type[i][i] = i;
for (int i = 1; i <= atom->ntypes; i++) { stencil_cut[i][i] = sqrt(cutneighsq[i][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<int> (cutneighmax*bininvx);
if (sx*binsizex < cutneighmax) sx++;
sy = static_cast<int> (cutneighmax*bininvy);
if (sy*binsizey < cutneighmax) sy++;
sz = static_cast<int> (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;
}
}
} }
// like -> like => use standard Newton stencil at bin // Cross types: use full stencil, looking one way through hierarchy
// smaller -> larger => use full stencil in larger 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
// larger -> smaller => no nstecil required // larger -> smaller => no nstecil required
// If cut offs are same, use existing type-type stencil // If cut offs are same, use existing type-type stencil
for (itype = 1; itype <= maxtypes; ++itype) { for (i = 1; i <= n; i++) {
for (jtype = 1; jtype <= maxtypes; ++jtype) { for (j = 1; j <= n; j++) {
if (itype == jtype) continue; if(i == j) continue;
if (cuttypesq[itype] == cuttypesq[jtype]) { if(cuttypesq[i] > cuttypesq[j]) continue;
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);
cutneighmaxsq = cuttypesq[jtype]; stencil_skip[i][j] = 0;
cutneighmax = sqrt(cutneighmaxsq); stencil_cut[i][j] = sqrt(cutneighsq[i][j]);
sx = static_cast<int> (cutneighmax*bininvx);
if (sx*binsizex < cutneighmax) sx++;
sy = static_cast<int> (cutneighmax*bininvy);
if (sy*binsizey < cutneighmax) sy++;
sz = static_cast<int> (cutneighmax*bininvz);
if (sz*binsizez < cutneighmax) sz++;
smax = (2*sx+1) * (2*sy+1) * (2*sz+1); if(cuttypesq[i] == cuttypesq[j]){
if (smax > maxstencil_type[itype][jtype]) { stencil_half[i][j] = 1;
maxstencil_type[itype][jtype] = smax; stencil_bin_type[i][j] = i;
memory->destroy(stencil_type[itype][jtype]); } else {
memory->create(stencil_type[itype][jtype], smax, "stencil_type[]"); stencil_half[i][j] = 0;
} stencil_bin_type[i][j] = j;
create_newtoff(itype, jtype, cuttypesq[jtype]);
} }
} }
} }
} }
/* ---------------------------------------------------------------------- */ /* ----------------------------------------------------------------------
create stencils based on bin geometry and cutoff
------------------------------------------------------------------------- */
void NStencilHalfBytype3dNewton::create_newton(int itype, int jtype, double cutsq) { void NStencilHalfMulti23d::create()
{
int itype, jtype, i, j, k, ns;
int n = atom->ntypes;
double cutsq;
int i, j, k, ns;
for (itype = 1; itype <= n; itype++) {
for (jtype = 1; jtype <= n; jtype++) {
if (stencil_skip[itype][jtype]) continue;
ns = 0; 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 (k = 0; k <= sz; k++)
for (j = -sy; j <= sy; j++) for (j = -sy; j <= sy; j++)
for (i = -sx; i <= sx; i++) for (i = -sx; i <= sx; i++)
if (k > 0 || j > 0 || (j == 0 && i > 0)) { if (k > 0 || j > 0 || (j == 0 && i > 0)) {
if (bin_distance(i,j,k) < cutsq) if (bin_distance(i,j,k) < cutsq)
stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; stencil_type[itype][jtype][ns++] =
k*mbiny*mbinx + j*mbinx + i;
} }
nstencil_type[itype][jtype] = ns; } else {
}
/* ---------------------------------------------------------------------- */
void NStencilHalfBytype3dNewton::create_newtoff(int itype, int jtype, double cutsq) {
int i, j, k, ns;
ns = 0;
for (k = -sz; k <= sz; k++) for (k = -sz; k <= sz; k++)
for (j = -sy; j <= sy; j++) for (j = -sy; j <= sy; j++)
for (i = -sx; i <= sx; i++) for (i = -sx; i <= sx; i++)
if (bin_distance(i,j,k) < cutsq) if (bin_distance(i,j,k) < cutsq)
stencil_type[itype][jtype][ns++] = k*mbiny*mbinx + j*mbinx + i; stencil_type[itype][jtype][ns++] =
k*mbiny*mbinx + j*mbinx + i;
}
nstencil_type[itype][jtype] = ns; nstencil_multi2[itype][jtype] = ns;
}
}
} }
/* ----------------------------------------------------------------------
create stencil based on bin geometry and cutoff
------------------------------------------------------------------------- */
void NStencilHalfBytype3dNewton::create()
{
// KS. Move "creation" here.
}

23
src/nstencil_half_multi2_3d.h
View File

@ -13,32 +13,27 @@
#ifdef NSTENCIL_CLASS #ifdef NSTENCIL_CLASS
NStencilStyle(half/bytype/3d, NStencilStyle(half/multi2/3d,
NStencilHalfBytype3dNewton, NS_HALF | NS_BYTYPE | NS_3D | NS_ORTHO) NStencilHalfMulti23d, NS_HALF | NS_MULTI2 | NS_3D | NS_ORTHO)
#else #else
#ifndef LMP_NSTENCIL_HALF_BYTYPE_3D_H #ifndef LMP_NSTENCIL_HALF_MULTI2_3D_H
#define LMP_NSTENCIL_HALF_BYTYPE_3D_H #define LMP_NSTENCIL_HALF_MULTI2_3D_H
#include "nstencil.h" #include "nstencil.h"
namespace LAMMPS_NS { namespace LAMMPS_NS {
class NStencilHalfBytype3dNewton : public NStencil { class NStencilHalfMulti23d : public NStencil {
public: public:
NStencilHalfBytype3dNewton(class LAMMPS *); NStencilHalfMulti23d(class LAMMPS *);
~NStencilHalfBytype3dNewton(); ~NStencilHalfMulti23d();
void create_setup();
void create(); void create();
private: protected:
int ** maxstencil_type; 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);
}; };
} }

198
src/nstencil_half_multi2_3d_delete_noff.cpp
View File

@ -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 <math.h>
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<int> (cutneighmax*bininvx);
if (sx*binsizex < cutneighmax) sx++;
sy = static_cast<int> (cutneighmax*bininvy);
if (sy*binsizey < cutneighmax) sy++;
sz = static_cast<int> (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<int> (cutneighmax*bininvx);
if (sx*binsizex < cutneighmax) sx++;
sy = static_cast<int> (cutneighmax*bininvy);
if (sy*binsizey < cutneighmax) sy++;
sz = static_cast<int> (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;
}

146
src/nstencil_half_multi2_3d_tri.cpp
View File

@ -11,7 +11,7 @@
See the README file in the top-level LAMMPS directory. 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 "neighbor.h"
#include "neigh_list.h" #include "neigh_list.h"
#include "nbin.h" #include "nbin.h"
@ -23,31 +23,100 @@ using namespace LAMMPS_NS;
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
NStencilHalfBytype3dNewton::NStencilHalfBytype3dNewton(LAMMPS *lmp) : NStencilHalfMulti23dTri::NStencilHalfMulti23d(LAMMPS *lmp) :
NStencil(lmp) NStencil(lmp) {}
{
maxstencil_type = NULL;
}
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
NStencilHalfBytype3dNewton::~NStencilHalfBytype3dNewton() { void NStencilHalfMulti23d::set_stencil_properties()
{
int n = atom->ntypes;
int i, j;
memory->destroy(nstencil_type); // like -> like => use half stencil
for (int i = 1; i <= atom->ntypes; i++) { for (i = 1; i <= n; i++) {
for (int j = 0; j <= atom->ntypes; j++) { stencil_half[i][i] = 1;
if (maxstencil_type[i][j] > 0) memory->destroy(stencil_type[i][j]); 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[i];
} }
delete [] stencil_type;
memory->destroy(maxstencil_type);
} }
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
// KS To superclass // KS To superclass
void NStencilHalfBytype3dNewton::copy_bin_info_bytype(int itype) { void NStencilHalfMulti23d::copy_bin_info_bytype(int itype) {
mbinx = nb->mbinx_type[itype]; mbinx = nb->mbinx_type[itype];
mbiny = nb->mbiny_type[itype]; mbiny = nb->mbiny_type[itype];
@ -63,7 +132,7 @@ void NStencilHalfBytype3dNewton::copy_bin_info_bytype(int itype) {
/* ---------------------------------------------------------------------- */ /* ---------------------------------------------------------------------- */
// KS To superclass? // KS To superclass?
int NStencilHalfBytype3dNewton::copy_neigh_info_bytype(int itype) { int NStencilHalfMulti23d::copy_neigh_info_bytype(int itype) {
cutneighmaxsq = neighbor->cutneighsq[itype][itype]; cutneighmaxsq = neighbor->cutneighsq[itype][itype];
cutneighmax = sqrt(cutneighmaxsq); 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; 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.
} }

24
src/nstencil_half_multi2_3d_tri.h
View File

@ -13,32 +13,26 @@
#ifdef NSTENCIL_CLASS #ifdef NSTENCIL_CLASS
NStencilStyle(half/bytype/3d, NStencilStyle(half/multi2/3d/tri,
NStencilHalfBytype3d, NS_HALF | NS_BYTYPE | NS_3D | NS_TRI) NStencilHalfMulti23dTri, NS_HALF | NS_MULTI2 | NS_3D | NS_TRI)
#else #else
#ifndef LMP_NSTENCIL_HALF_BYTYPE_3D_H #ifndef LMP_NSTENCIL_HALF_MULTI2_3D_TRI_H
#define LMP_NSTENCIL_HALF_BYTYPE_3D_H #define LMP_NSTENCIL_HALF_MULTI2_3D_TRI_H
#include "nstencil.h" #include "nstencil.h"
namespace LAMMPS_NS { namespace LAMMPS_NS {
class NStencilHalfBytype3dNewton : public NStencil { class NStencilHalfMulti23dTri : public NStencil {
public: public:
NStencilHalfBytype3dNewton(class LAMMPS *); NStencilHalfMulti23dTri(class LAMMPS *);
~NStencilHalfBytype3dNewton(); ~NStencilHalfMulti23dTri();
void create_setup();
void create(); void create();
private: protected:
int ** maxstencil_type; 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);
}; };
} }