Skip to content

Commit

Permalink
Merge pull request lammps#4015 from srtee/develop
Browse files Browse the repository at this point in the history 
add border information to Atom::add_custom and Atom::find_custom
  • Loading branch information
@akohlmey
akohlmey authored Jan 30, 2024
2 parents 7d53f8d + 6927934 commit 52de76d
Show file tree
Hide file tree
Showing 6 changed files with 200 additions and 107 deletions.
208 changes: 112 additions & 96 deletions src/AMOEBA/pair_amoeba.cpp
View file
Original file line number Diff line number Diff line change
@@ -1,4 +1,3 @@
// clang-format off
/* ----------------------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
https://www.lammps.org/ Sandia National Laboratories
Expand Down Expand Up @@ -39,17 +38,19 @@ using namespace LAMMPS_NS;

using MathSpecial::powint;

enum{INDUCE,RSD,SETUP_AMOEBA,SETUP_HIPPO,KMPOLE,AMGROUP,PVAL}; // forward comm
enum{FIELD,ZRSD,TORQUE,UFLD}; // reverse comm
enum{ARITHMETIC,GEOMETRIC,CUBIC_MEAN,R_MIN,SIGMA,DIAMETER,HARMONIC,HHG,W_H};
enum{HAL,REPULSE,QFER,DISP,MPOLE,POLAR,USOLV,DISP_LONG,MPOLE_LONG,POLAR_LONG};
enum{MPOLE_GRID,POLAR_GRID,POLAR_GRIDC,DISP_GRID,INDUCE_GRID,INDUCE_GRIDC};
enum{MUTUAL,OPT,TCG,DIRECT};
enum{GEAR,ASPC,LSQR};
enum { INDUCE, RSD, SETUP_AMOEBA, SETUP_HIPPO, KMPOLE, AMGROUP, PVAL }; // forward comm
enum { FIELD, ZRSD, TORQUE, UFLD }; // reverse comm
enum { ARITHMETIC, GEOMETRIC, CUBIC_MEAN, R_MIN, SIGMA, DIAMETER, HARMONIC, HHG, W_H };
enum { HAL, REPULSE, QFER, DISP, MPOLE, POLAR, USOLV, DISP_LONG, MPOLE_LONG, POLAR_LONG };
enum { MPOLE_GRID, POLAR_GRID, POLAR_GRIDC, DISP_GRID, INDUCE_GRID, INDUCE_GRIDC };
enum { MUTUAL, OPT, TCG, DIRECT };
enum { GEAR, ASPC, LSQR };

static constexpr int DELTASTACK = 16;
#define DEBUG_AMOEBA 0

// clang-format off

/* ---------------------------------------------------------------------- */

PairAmoeba::PairAmoeba(LAMMPS *lmp) : Pair(lmp)
Expand Down Expand Up @@ -429,7 +430,7 @@ void PairAmoeba::compute(int eflag, int vflag)
// output FF settings to screen and logfile
// delay until here because RMS force accuracy is computed based on rpole

if (first_flag_compute && (comm->me == 0)) print_settings();
if (first_flag_compute) print_settings();
first_flag_compute = 0;

if (amoeba) pbc_xred();
Expand Down Expand Up @@ -827,28 +828,36 @@ void PairAmoeba::init_style()

// check if all custom atom arrays were set via fix property/atom

int flag,cols;

index_amtype = atom->find_custom("amtype",flag,cols);
if (index_amtype < 0 || flag || cols)
error->all(FLERR,"Pair {} amtype is not defined", mystyle);
index_amgroup = atom->find_custom("amgroup",flag,cols);
if (index_amgroup < 0 || flag || cols)
error->all(FLERR,"Pair {} amgroup is not defined", mystyle);

index_redID = atom->find_custom("redID",flag,cols);
if (index_redID < 0 || !flag || cols)
error->all(FLERR,"Pair {} redID is not defined", mystyle);
index_xyzaxis = atom->find_custom("xyzaxis",flag,cols);
if (index_xyzaxis < 0 || !flag || cols == 0)
error->all(FLERR,"Pair {} xyzaxis is not defined", mystyle);

index_polaxe = atom->find_custom("polaxe",flag,cols);
if (index_polaxe < 0 || flag || cols)
error->all(FLERR,"Pair {} polaxe is not defined", mystyle);
index_pval = atom->find_custom("pval",flag,cols);
if (index_pval < 0 || !flag || cols)
error->all(FLERR,"Pair {} pval is not defined", mystyle);
// clang-format on
const char *names[6] = {"amtype", "amgroup", "redID", "xyzaxis", "polaxe", "pval"};
const int flag_check[6] = {0, 0, 1, 1, 0, 1}; // correct type (0 int, 1 dbl)
const int cols_check[6] = {0, 0, 0, 3, 0, 0}; // xyzaxis 3 cols, all others 0
const int ghost_check[6] = {1, 1, 1, 0, 0, 1}; // which types need ghost; TO-DO: check
int flag, cols, ghost, index[6];

// clang-format off

for (int i = 0; i < 6; i++) {
if (ghost_check[i]) {
index[i] = atom->find_custom_ghost(names[i], flag, cols, ghost);
} else {
index[i] = atom->find_custom(names[i], flag, cols);
}
std::string err = "";
if (index[i] < 0) err = "was not defined";
else if (flag_check[i] != flag) err = "has the wrong type";
else if (cols_check[i] != cols) err = "has the wrong number of columns";
else if (ghost_check[i] && !ghost) err = "must be set by fix property/atom with ghost yes";
if (err != "")
error->all(FLERR,"Pair {} per-atom variable {} {}", mystyle, names[i], err);
}

index_amtype = index[0];
index_amgroup = index[1];
index_redID = index[2];
index_xyzaxis = index[3];
index_polaxe = index[4];
index_pval = index[5];

// -------------------------------------------------------------------
// one-time initializations
Expand Down Expand Up @@ -1069,79 +1078,86 @@ void PairAmoeba::init_style()
void PairAmoeba::print_settings()
{
std::string mesg = utils::uppercase(mystyle) + " force field settings\n";
double estimated_mpole_accuracy = 0.0;

if (amoeba) {
choose(HAL);
mesg += fmt::format(" hal: cut {} taper {} vscale {} {} {} {}\n", sqrt(off2),sqrt(cut2),
special_hal[1],special_hal[2],special_hal[3],special_hal[4]);
} else {
choose(REPULSE);
mesg += fmt::format(" repulsion: cut {} taper {} rscale {} {} {} {}\n", sqrt(off2),sqrt(cut2),
special_repel[1],special_repel[2],special_repel[3],special_repel[4]);
if (use_ewald) {
choose(MPOLE_LONG);
estimated_mpole_accuracy = final_accuracy_mpole();
}

choose(QFER);
mesg += fmt::format(" qxfer: cut {} taper {} mscale {} {} {} {}\n", sqrt(off2),sqrt(cut2),
special_mpole[1],special_mpole[2],special_mpole[3],special_mpole[4]);
if (comm->me == 0) {
if (amoeba) {
choose(HAL);
mesg += fmt::format(" hal: cut {} taper {} vscale {} {} {} {}\n", sqrt(off2),sqrt(cut2),
special_hal[1],special_hal[2],special_hal[3],special_hal[4]);
} else {
choose(REPULSE);
mesg += fmt::format(" repulsion: cut {} taper {} rscale {} {} {} {}\n", sqrt(off2),sqrt(cut2),
special_repel[1],special_repel[2],special_repel[3],special_repel[4]);

choose(QFER);
mesg += fmt::format(" qxfer: cut {} taper {} mscale {} {} {} {}\n", sqrt(off2),sqrt(cut2),
special_mpole[1],special_mpole[2],special_mpole[3],special_mpole[4]);

if (use_dewald) {
choose(DISP_LONG);
mesg += fmt::format(" dispersion: cut {} aewald {} bsorder {} FFT {} {} {} "
"dspscale {} {} {} {}\n", sqrt(off2),aewald,bsdorder,ndfft1,ndfft2,ndfft3,
special_disp[1],special_disp[2],special_disp[3],special_disp[4]);
} else {
choose(DISP);
mesg += fmt::format(" dispersion: cut {} aewald {} dspscale {} {} {} {}\n",
sqrt(off2),aewald,special_disp[1],
special_disp[2],special_disp[3],special_disp[4]);
}
}

if (use_dewald) {
choose(DISP_LONG);
mesg += fmt::format(" dispersion: cut {} aewald {} bsorder {} FFT {} {} {} "
"dspscale {} {} {} {}\n", sqrt(off2),aewald,bsdorder,ndfft1,ndfft2,ndfft3,
special_disp[1],special_disp[2],special_disp[3],special_disp[4]);
if (use_ewald) {
choose(MPOLE_LONG);
mesg += fmt::format(" multipole: cut {} aewald {} bsorder {} FFT {} {} {}\n"
" estimated absolute RMS force accuracy = {:.8g}\n"
" estimated relative RMS force accuracy = {:.8g}\n"
" mscale {} {} {} {}\n",
sqrt(off2),aewald,bseorder,nefft1,nefft2,nefft3,
estimated_mpole_accuracy,estimated_mpole_accuracy/two_charge_force,
special_mpole[1],special_mpole[2],special_mpole[3],special_mpole[4]);
} else {
choose(DISP);
mesg += fmt::format(" dispersion: cut {} aewald {} dspscale {} {} {} {}\n",
sqrt(off2),aewald,special_disp[1],
special_disp[2],special_disp[3],special_disp[4]);
choose(MPOLE);
mesg += fmt::format(" multipole: cut {} aewald {} mscale {} {} {} {}\n", sqrt(off2),aewald,
special_mpole[1],special_mpole[2],special_mpole[3],special_mpole[4]);
}
}

if (use_ewald) {
choose(MPOLE_LONG);
double estimated_accuracy = final_accuracy_mpole();
mesg += fmt::format(" multipole: cut {} aewald {} bsorder {} FFT {} {} {}; "
"estimated absolute RMS force accuracy = {:.8g}; "
"estimated relative RMS force accuracy = {:.8g}; "
"mscale {} {} {} {}\n",
sqrt(off2),aewald,bseorder,nefft1,nefft2,nefft3,
estimated_accuracy,estimated_accuracy/two_charge_force,
special_mpole[1],special_mpole[2],special_mpole[3],special_mpole[4]);
} else {
choose(MPOLE);
mesg += fmt::format(" multipole: cut {} aewald {} mscale {} {} {} {}\n", sqrt(off2),aewald,
special_mpole[1],special_mpole[2],special_mpole[3],special_mpole[4]);
}
if (use_ewald) {
choose(POLAR_LONG);
mesg += fmt::format(" polar: cut {} aewald {} bsorder {} FFT {} {} {}\n",
sqrt(off2),aewald,bsporder,nefft1,nefft2,nefft3);
mesg += fmt::format(" pscale {} {} {} {} piscale {} {} {} {} "
"wscale {} {} {} {} d/u scale {} {}\n",
special_polar_pscale[1],special_polar_pscale[2],
special_polar_pscale[3],special_polar_pscale[4],
special_polar_piscale[1],special_polar_piscale[2],
special_polar_piscale[3],special_polar_piscale[4],
special_polar_wscale[1],special_polar_wscale[2],
special_polar_wscale[3],special_polar_wscale[4],
polar_dscale,polar_uscale);
} else {
choose(POLAR);
mesg += fmt::format(" polar: cut {} aewald {}\n",sqrt(off2),aewald);
mesg += fmt::format(" pscale {} {} {} {} piscale {} {} {} {} "
"wscale {} {} {} {} d/u scale {} {}\n",
special_polar_pscale[1],special_polar_pscale[2],
special_polar_pscale[3],special_polar_pscale[4],
special_polar_piscale[1],special_polar_piscale[2],
special_polar_piscale[3],special_polar_piscale[4],
special_polar_wscale[1],special_polar_wscale[2],
special_polar_wscale[3],special_polar_wscale[4],
polar_dscale,polar_uscale);
}

if (use_ewald) {
choose(POLAR_LONG);
mesg += fmt::format(" polar: cut {} aewald {} bsorder {} FFT {} {} {}\n",
sqrt(off2),aewald,bsporder,nefft1,nefft2,nefft3);
mesg += fmt::format(" pscale {} {} {} {} piscale {} {} {} {} "
"wscale {} {} {} {} d/u scale {} {}\n",
special_polar_pscale[1],special_polar_pscale[2],
special_polar_pscale[3],special_polar_pscale[4],
special_polar_piscale[1],special_polar_piscale[2],
special_polar_piscale[3],special_polar_piscale[4],
special_polar_wscale[1],special_polar_wscale[2],
special_polar_wscale[3],special_polar_wscale[4],
polar_dscale,polar_uscale);
} else {
choose(POLAR);
mesg += fmt::format(" polar: cut {} aewald {}\n",sqrt(off2),aewald);
mesg += fmt::format(" pscale {} {} {} {} piscale {} {} {} {} "
"wscale {} {} {} {} d/u scale {} {}\n",
special_polar_pscale[1],special_polar_pscale[2],
special_polar_pscale[3],special_polar_pscale[4],
special_polar_piscale[1],special_polar_piscale[2],
special_polar_piscale[3],special_polar_piscale[4],
special_polar_wscale[1],special_polar_wscale[2],
special_polar_wscale[3],special_polar_wscale[4],
polar_dscale,polar_uscale);
choose(USOLV);
mesg += fmt::format(" precondition: cut {}\n",sqrt(off2));
utils::logmesg(lmp, mesg);
}

choose(USOLV);
mesg += fmt::format(" precondition: cut {}\n",sqrt(off2));
utils::logmesg(lmp, mesg);
}

/* ----------------------------------------------------------------------
Expand Down
13 changes: 12 additions & 1 deletion src/KOKKOS/atom_kokkos.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -300,7 +300,7 @@ void AtomKokkos::grow(unsigned int mask)
return index in ivector or dvector of its location
------------------------------------------------------------------------- */

int AtomKokkos::add_custom(const char *name, int flag, int cols)
int AtomKokkos::add_custom(const char *name, int flag, int cols, int ghost)
{
int index = -1;

Expand All @@ -309,6 +309,8 @@ int AtomKokkos::add_custom(const char *name, int flag, int cols)
nivector++;
ivname = (char **) memory->srealloc(ivname, nivector * sizeof(char *), "atom:ivname");
ivname[index] = utils::strdup(name);
ivghost = (int *) memory->srealloc(ivghost,nivector * sizeof(int),"atom:ivghost");
ivghost[index] = ghost;
ivector = (int **) memory->srealloc(ivector, nivector * sizeof(int *), "atom:ivector");
memory->create(ivector[index], nmax, "atom:ivector");

Expand All @@ -317,6 +319,8 @@ int AtomKokkos::add_custom(const char *name, int flag, int cols)
ndvector++;
dvname = (char **) memory->srealloc(dvname, ndvector * sizeof(char *), "atom:dvname");
dvname[index] = utils::strdup(name);
dvghost = (int *) memory->srealloc(dvghost, ndvector * sizeof(int), "atom:dvghost");
dvghost[index] = ghost;
dvector = (double **) memory->srealloc(dvector, ndvector * sizeof(double *), "atom:dvector");
this->sync(Device, DVECTOR_MASK);
memoryKK->grow_kokkos(k_dvector, dvector, ndvector, nmax, "atom:dvector");
Expand All @@ -327,6 +331,8 @@ int AtomKokkos::add_custom(const char *name, int flag, int cols)
niarray++;
ianame = (char **) memory->srealloc(ianame, niarray * sizeof(char *), "atom:ianame");
ianame[index] = utils::strdup(name);
iaghost = (int *) memory->srealloc(iaghost, niarray * sizeof(int), "atom:iaghost");
iaghost[index] = ghost;
iarray = (int ***) memory->srealloc(iarray, niarray * sizeof(int **), "atom:iarray");
memory->create(iarray[index], nmax, cols, "atom:iarray");

Expand All @@ -338,13 +344,18 @@ int AtomKokkos::add_custom(const char *name, int flag, int cols)
ndarray++;
daname = (char **) memory->srealloc(daname, ndarray * sizeof(char *), "atom:daname");
daname[index] = utils::strdup(name);
daghost = (int *) memory->srealloc(daghost, ndarray * sizeof(int), "atom:daghost");
daghost[index] = ghost;
darray = (double ***) memory->srealloc(darray, ndarray * sizeof(double **), "atom:darray");
memory->create(darray[index], nmax, cols, "atom:darray");

dcols = (int *) memory->srealloc(dcols, ndarray * sizeof(int), "atom:dcols");
dcols[index] = cols;
}

if (index < 0)
error->all(FLERR,"Invalid call to AtomKokkos::add_custom()");

return index;
}

Expand Down
2 changes: 1 addition & 1 deletion src/KOKKOS/atom_kokkos.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -154,7 +154,7 @@ class AtomKokkos : public Atom {
void sync_overlapping_device(const ExecutionSpace space, unsigned int mask);
void sort() override;
virtual void grow(unsigned int mask);
int add_custom(const char *, int, int) override;
int add_custom(const char *, int, int, int border = 0) override;
void remove_custom(int, int, int) override;
virtual void deallocate_topology();
private:
Expand Down
Loading

0 comments on commit 52de76d

Please sign in to comment.