From b4b6ba91e02bcec48b372ad5f508ecb0a58befb4 Mon Sep 17 00:00:00 2001 From: Aidan Thompson Date: Wed, 20 Apr 2022 18:01:24 -0600 Subject: [PATCH 01/16] Slight adjustment to definition of Sinner and Dinner --- doc/src/compute_sna_atom.rst | 29 ++++++++++++++++------------- doc/src/pair_snap.rst | 10 +++++----- src/ML-SNAP/pair_snap.cpp | 1 - src/ML-SNAP/sna.cpp | 12 ++++++------ 4 files changed, 27 insertions(+), 25 deletions(-) diff --git a/doc/src/compute_sna_atom.rst b/doc/src/compute_sna_atom.rst index e43e61cd00..992c87c816 100644 --- a/doc/src/compute_sna_atom.rst +++ b/doc/src/compute_sna_atom.rst @@ -59,9 +59,9 @@ Syntax *bikflag* value = *0* or *1* (only implemented for compute snap) *0* = per-atom bispectrum descriptors are summed over atoms *1* = per-atom bispectrum descriptors are not summed over atoms - *switchinnerflag* values = *rinnerlist* *drinnerlist* - *rinnerlist* = *ntypes* values of rinner (distance units) - *drinnerlist* = *ntypes* values of drinner (distance units) + *switchinnerflag* values = *sinnerlist* *dinnerlist* + *sinnerlist* = *ntypes* values of *Sinner* (distance units) + *dinnerlist* = *ntypes* values of *Dinner* (distance units) Examples """""""" @@ -316,21 +316,24 @@ The keyword *switchinnerflag* activates an additional radial switching function similar to :math:`f_c(r)` above, but acting to switch off smoothly contributions from neighbor atoms at short separation distances. This is useful when SNAP is used in combination with a simple -repulsive potential. The keyword is followed by the *ntypes* -values for :math:`r_{inner}` and the *ntypes* -values for :math:`\Delta r_{inner}`. For a neighbor atom at +repulsive potential. For a neighbor atom at distance :math:`r`, its contribution is scaled by a multiplicative factor :math:`f_{inner}(r)` defined as follows: .. math:: - = & 0, r \leq r_{inner} \\ - f_{inner}(r) = & \frac{1}{2}(1 - \cos(\pi \frac{r-r_{inner}}{\Delta r_{inner}})), r_{inner} < r \leq r_{inner} + \Delta r_{inner} \\ - = & 1, r > r_{inner} + \Delta r_{inner} + = & 0, r \leq S_{inner} - D_{inner} \\ + f_{inner}(r) = & \frac{1}{2}(1 - \cos(\frac{\pi}{2} (1 + \frac{r-S_{inner}}{D_{inner}})), S_{inner} - D_{inner} < r \leq S_{inner} + D_{inner} \\ + = & 1, r > S_{inner} + D_{inner} -The values of :math:`r_{inner}` and :math:`\Delta r_{inner}` are -the arithmetic means of the values for the central atom of type I -and the neighbor atom of type J. +where the switching region is centered at :math:`S_{inner}` and it extends a distance :math:`D_{inner}` +to the left and to the right of this. +The keyword is followed by the *ntypes* +values for :math:`S_{inner}` and the *ntypes* +values for :math:`D_{inner}`. +When the central atom and the neighbor atom have different types, +the values of :math:`S_{inner}` and :math:`D_{inner}` are +the arithmetic means of the values for both types. .. note:: @@ -450,7 +453,7 @@ Default The optional keyword defaults are *rmin0* = 0, *switchflag* = 1, *bzeroflag* = 1, *quadraticflag* = 0, -*bnormflag* = 0, *wselfallflag* = 0 +*bnormflag* = 0, *wselfallflag* = 0, *switchinnerflag* = off, ---------- diff --git a/doc/src/pair_snap.rst b/doc/src/pair_snap.rst index 350561fe4d..f7cc7ff731 100644 --- a/doc/src/pair_snap.rst +++ b/doc/src/pair_snap.rst @@ -136,7 +136,7 @@ keyword/value pair. The required keywords are *rcutfac* and *twojmax*\ . Optional keywords are *rfac0*, *rmin0*, *switchflag*, *bzeroflag*, *quadraticflag*, *chemflag*, *bnormflag*, *wselfallflag*, *switchinnerflag*, -*rinner*, *drinner*, *chunksize*, and *parallelthresh*\ . +*Sinner*, *Dinner*, *chunksize*, and *parallelthresh*\ . The default values for these keywords are @@ -152,6 +152,9 @@ The default values for these keywords are * *chunksize* = 32768 * *parallelthresh* = 8192 +For detailed definitions of all of these keywords, +see the :doc:`compute sna/atom ` doc page. + If *quadraticflag* is set to 1, then the SNAP energy expression includes additional quadratic terms that have been shown to increase the overall accuracy of the potential without much increase in computational cost @@ -194,7 +197,7 @@ pair_coeff command, to avoid ambiguity in the number of coefficients. The keyword *switchinnerflag* activates an additional switching function that smoothly turns off contributions to the SNAP potential from neighbor atoms at short separations. If *switchinnerflag* is set to 1 then -the additional keywords *rinner* and *drinner* must also be provided. +the additional keywords *Sinner* and *Dinner* must also be provided. Each of these is followed by *nelements* values, where *nelements* is the number of unique elements appearing in appearing in the LAMMPS pair_coeff command. The element order should correspond to the order @@ -217,9 +220,6 @@ already large enough to saturate the GPU threads. Extra parallelism will be performed if the *chunksize* (or total number of atoms per GPU) is smaller than *parallelthresh*. -Detailed definitions for all the other keywords -are given on the :doc:`compute sna/atom ` doc page. - .. note:: The previously used *diagonalstyle* keyword was removed in 2019, diff --git a/src/ML-SNAP/pair_snap.cpp b/src/ML-SNAP/pair_snap.cpp index 6a108ad3dd..3ba14c5872 100644 --- a/src/ML-SNAP/pair_snap.cpp +++ b/src/ML-SNAP/pair_snap.cpp @@ -703,7 +703,6 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) if (keywd == "rinner") { keyval = words[iword]; for (int ielem = 0; ielem < nelements; ielem++) { - printf("rinnerelem = %p ielem = %d nelements = %d iword = %d nwords = %d\n",rinnerelem, ielem, nelements, iword, nwords); rinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); iword++; } diff --git a/src/ML-SNAP/sna.cpp b/src/ML-SNAP/sna.cpp index 9d32f4ab9c..717d0014f5 100644 --- a/src/ML-SNAP/sna.cpp +++ b/src/ML-SNAP/sna.cpp @@ -1582,10 +1582,10 @@ double SNA::compute_sfac_inner(double r, double rinner, double drinner) double sfac; if (switch_inner_flag == 0) sfac = 1.0; else if (r >= rinner + drinner) sfac = 1.0; - else if (r <= rinner) sfac = 0.0; + else if (r <= rinner - drinner) sfac = 0.0; else { - double rcutfac = MY_PI / drinner; - sfac = 0.5 * (1.0 - cos((r - rinner) * rcutfac)); + double rcutfac = MY_PI2 / drinner; + sfac = 0.5 * (1.0 - cos(MY_PI2 + (r - rinner) * rcutfac)); } return sfac; } @@ -1597,10 +1597,10 @@ double SNA::compute_dsfac_inner(double r, double rinner, double drinner) double dsfac; if (switch_inner_flag == 0) dsfac = 0.0; else if (r >= rinner + drinner) dsfac = 0.0; - else if (r <= rinner) dsfac = 0.0; + else if (r <= rinner - drinner) dsfac = 0.0; else { - double rcutfac = MY_PI / drinner; - dsfac = 0.5 * sin((r - rinner) * rcutfac) * rcutfac; + double rcutfac = MY_PI2 / drinner; + dsfac = 0.5 * rcutfac * sin(MY_PI2 + (r - rinner) * rcutfac); } return dsfac; } From e8cacc438001d57cd530aae4642cd305de0e6eb5 Mon Sep 17 00:00:00 2001 From: Aidan Thompson Date: Wed, 20 Apr 2022 18:17:23 -0600 Subject: [PATCH 02/16] Updated example syntax --- doc/src/compute_sna_atom.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/doc/src/compute_sna_atom.rst b/doc/src/compute_sna_atom.rst index 992c87c816..eecb3295ab 100644 --- a/doc/src/compute_sna_atom.rst +++ b/doc/src/compute_sna_atom.rst @@ -73,7 +73,7 @@ Examples compute vb all sna/atom 1.4 0.95 6 2.0 1.0 compute snap all snap 1.4 0.95 6 2.0 1.0 compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 chem 2 0 1 - compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 switchinnerflag 1.1 1.3 0.5 0.6 + compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 switchinnerflag 1.35 1.6 0.25 0.3 Description """"""""""" From 4de9ab85ce0aa348606e1293734b0bcf1d076618 Mon Sep 17 00:00:00 2001 From: Aidan Thompson Date: Thu, 21 Apr 2022 17:30:12 -0600 Subject: [PATCH 03/16] Completed inner cutoff, KOKKOS still in progress --- doc/src/compute_sna_atom.rst | 21 +++--- doc/src/pair_snap.rst | 4 +- src/KOKKOS/pair_snap_kokkos_impl.h | 2 +- src/KOKKOS/sna_kokkos.h | 20 ++++-- src/KOKKOS/sna_kokkos_impl.h | 52 ++++++++++---- src/ML-SNAP/compute_sna_atom.cpp | 43 +++++++++--- src/ML-SNAP/compute_sna_atom.h | 4 +- src/ML-SNAP/compute_snad_atom.cpp | 52 ++++++++++---- src/ML-SNAP/compute_snad_atom.h | 4 +- src/ML-SNAP/compute_snap.cpp | 45 ++++++++---- src/ML-SNAP/compute_snap.h | 4 +- src/ML-SNAP/compute_snav_atom.cpp | 51 ++++++++++---- src/ML-SNAP/compute_snav_atom.h | 4 +- src/ML-SNAP/pair_snap.cpp | 46 ++++++------- src/ML-SNAP/pair_snap.h | 4 +- src/ML-SNAP/sna.cpp | 107 ++++++++++++----------------- src/ML-SNAP/sna.h | 11 ++- 17 files changed, 286 insertions(+), 188 deletions(-) diff --git a/doc/src/compute_sna_atom.rst b/doc/src/compute_sna_atom.rst index eecb3295ab..e34c9a8889 100644 --- a/doc/src/compute_sna_atom.rst +++ b/doc/src/compute_sna_atom.rst @@ -33,7 +33,7 @@ Syntax * R_1, R_2,... = list of cutoff radii, one for each type (distance units) * w_1, w_2,... = list of neighbor weights, one for each type * zero or more keyword/value pairs may be appended -* keyword = *rmin0* or *switchflag* or *bzeroflag* or *quadraticflag* or *chem* or *bnormflag* or *wselfallflag* or *bikflag* or *switchinnerflag* +* keyword = *rmin0* or *switchflag* or *bzeroflag* or *quadraticflag* or *chem* or *bnormflag* or *wselfallflag* or *bikflag* or *switchinnerflag* or *sinner* or *dinner* .. parsed-literal:: @@ -59,8 +59,12 @@ Syntax *bikflag* value = *0* or *1* (only implemented for compute snap) *0* = per-atom bispectrum descriptors are summed over atoms *1* = per-atom bispectrum descriptors are not summed over atoms - *switchinnerflag* values = *sinnerlist* *dinnerlist* + *switchinnerflag* value = *0* or *1* + *0* = do not use inner switching function + *1* = use inner switching function + *sinner* values = *sinnerlist* *sinnerlist* = *ntypes* values of *Sinner* (distance units) + *dinner* values = *dinnerlist* *dinnerlist* = *ntypes* values of *Dinner* (distance units) Examples @@ -73,7 +77,7 @@ Examples compute vb all sna/atom 1.4 0.95 6 2.0 1.0 compute snap all snap 1.4 0.95 6 2.0 1.0 compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 chem 2 0 1 - compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 switchinnerflag 1.35 1.6 0.25 0.3 + compute snap all snap 1.0 0.99363 6 3.81 3.83 1.0 0.93 switchinnerflag 1 sinner 1.35 1.6 dinner 0.25 0.3 Description """"""""""" @@ -312,7 +316,8 @@ the resulting bispectrum rows are :math:`B_{i,k}` instead of just :math:`B_k`. In this case, the entries in the final column for these rows are set to zero. -The keyword *switchinnerflag* activates an additional radial switching +The keyword *switchinnerflag* with value 1 +activates an additional radial switching function similar to :math:`f_c(r)` above, but acting to switch off smoothly contributions from neighbor atoms at short separation distances. This is useful when SNAP is used in combination with a simple @@ -328,9 +333,9 @@ factor :math:`f_{inner}(r)` defined as follows: where the switching region is centered at :math:`S_{inner}` and it extends a distance :math:`D_{inner}` to the left and to the right of this. -The keyword is followed by the *ntypes* -values for :math:`S_{inner}` and the *ntypes* -values for :math:`D_{inner}`. +With this option, additional keywords *sinner* and *dinner* must be used, +each followed by *ntypes* +values for :math:`S_{inner}` and :math:`D_{inner}`, respectively. When the central atom and the neighbor atom have different types, the values of :math:`S_{inner}` and :math:`D_{inner}` are the arithmetic means of the values for both types. @@ -453,7 +458,7 @@ Default The optional keyword defaults are *rmin0* = 0, *switchflag* = 1, *bzeroflag* = 1, *quadraticflag* = 0, -*bnormflag* = 0, *wselfallflag* = 0, *switchinnerflag* = off, +*bnormflag* = 0, *wselfallflag* = 0, *switchinnerflag* = 0, ---------- diff --git a/doc/src/pair_snap.rst b/doc/src/pair_snap.rst index f7cc7ff731..ebedb288c1 100644 --- a/doc/src/pair_snap.rst +++ b/doc/src/pair_snap.rst @@ -136,7 +136,7 @@ keyword/value pair. The required keywords are *rcutfac* and *twojmax*\ . Optional keywords are *rfac0*, *rmin0*, *switchflag*, *bzeroflag*, *quadraticflag*, *chemflag*, *bnormflag*, *wselfallflag*, *switchinnerflag*, -*Sinner*, *Dinner*, *chunksize*, and *parallelthresh*\ . +*sinner*, *dinner*, *chunksize*, and *parallelthresh*\ . The default values for these keywords are @@ -197,7 +197,7 @@ pair_coeff command, to avoid ambiguity in the number of coefficients. The keyword *switchinnerflag* activates an additional switching function that smoothly turns off contributions to the SNAP potential from neighbor atoms at short separations. If *switchinnerflag* is set to 1 then -the additional keywords *Sinner* and *Dinner* must also be provided. +the additional keywords *sinner* and *dinner* must also be provided. Each of these is followed by *nelements* values, where *nelements* is the number of unique elements appearing in appearing in the LAMMPS pair_coeff command. The element order should correspond to the order diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h index 6cf40b31a5..6bc034a163 100644 --- a/src/KOKKOS/pair_snap_kokkos_impl.h +++ b/src/KOKKOS/pair_snap_kokkos_impl.h @@ -602,7 +602,7 @@ void PairSNAPKokkos::coeff(int narg, char Kokkos::deep_copy(d_map,h_map); snaKK = SNAKokkos(rfac0,twojmax, - rmin0,switchflag,bzeroflag,chemflag,bnormflag,wselfallflag,nelements); + rmin0,switchflag,bzeroflag,chemflag,bnormflag,wselfallflag,nelements,switchinnerflag); snaKK.grow_rij(0,0); snaKK.init(); } diff --git a/src/KOKKOS/sna_kokkos.h b/src/KOKKOS/sna_kokkos.h index c33f5e486e..706f2d2f4c 100644 --- a/src/KOKKOS/sna_kokkos.h +++ b/src/KOKKOS/sna_kokkos.h @@ -103,7 +103,7 @@ inline SNAKokkos(const SNAKokkos& sna, const typename Kokkos::TeamPolicy::member_type& team); inline - SNAKokkos(real_type, int, real_type, int, int, int, int, int, int); +SNAKokkos(real_type, int, real_type, int, int, int, int, int, int, int); KOKKOS_INLINE_FUNCTION ~SNAKokkos(); @@ -192,13 +192,13 @@ inline void compute_deidrj_cpu(const typename Kokkos::TeamPolicy::member_type& team, int, int); // ForceSNAP KOKKOS_INLINE_FUNCTION - real_type compute_sfac(real_type, real_type); // add_uarraytot, compute_duarray + real_type compute_sfac(real_type, real_type, real_type, real_type); // add_uarraytot, compute_duarray KOKKOS_INLINE_FUNCTION - real_type compute_dsfac(real_type, real_type); // compute_duarray + real_type compute_dsfac(real_type, real_type, real_type, real_type); // compute_duarray KOKKOS_INLINE_FUNCTION - void compute_s_dsfac(const real_type, const real_type, real_type&, real_type&); // compute_cayley_klein + void compute_s_dsfac(const real_type, const real_type, const real_type, const real_type, real_type&, real_type&); // compute_cayley_klein #ifdef TIMING_INFO double* timers; @@ -214,6 +214,8 @@ inline t_sna_2i inside; t_sna_2d wj; t_sna_2d rcutij; + t_sna_2d sinnerij; + t_sna_2d dinnerij; t_sna_2i element; t_sna_3d dedr; int natom, nmax; @@ -297,7 +299,7 @@ inline void init_rootpqarray(); // init() KOKKOS_INLINE_FUNCTION - void add_uarraytot(const typename Kokkos::TeamPolicy::member_type& team, int, int, const real_type&, const real_type&, const real_type&, int); // compute_ui + void add_uarraytot(const typename Kokkos::TeamPolicy::member_type& team, int, int, const real_type&, const real_type&, const real_type&, const real_type&, const real_type&, int); // compute_ui KOKKOS_INLINE_FUNCTION void compute_uarray_cpu(const typename Kokkos::TeamPolicy::member_type& team, int, int, @@ -313,13 +315,19 @@ inline KOKKOS_INLINE_FUNCTION void compute_duarray_cpu(const typename Kokkos::TeamPolicy::member_type& team, int, int, const real_type&, const real_type&, const real_type&, // compute_duidrj_cpu - const real_type&, const real_type&, const real_type&, const real_type&, const real_type&); + const real_type&, const real_type&, const real_type&, const real_type&, const real_type&, + const real_type&, const real_type&); // Sets the style for the switching function // 0 = none // 1 = cosine int switch_flag; + // Sets the style for the inner switching function + // 0 = none + // 1 = cosine + int switch_inner_flag; + // Chem snap flags int chem_flag; int bnorm_flag; diff --git a/src/KOKKOS/sna_kokkos_impl.h b/src/KOKKOS/sna_kokkos_impl.h index 550548df55..546c07ee11 100644 --- a/src/KOKKOS/sna_kokkos_impl.h +++ b/src/KOKKOS/sna_kokkos_impl.h @@ -25,12 +25,13 @@ namespace LAMMPS_NS { static const double MY_PI = 3.14159265358979323846; // pi +static const double MY_PI2 = 1.57079632679489661923; // pi/2 template inline SNAKokkos::SNAKokkos(real_type rfac0_in, int twojmax_in, real_type rmin0_in, int switch_flag_in, int bzero_flag_in, - int chem_flag_in, int bnorm_flag_in, int wselfall_flag_in, int nelements_in) + int chem_flag_in, int bnorm_flag_in, int wselfall_flag_in, int nelements_in, int switch_inner_flag_in) { LAMMPS_NS::ExecutionSpace execution_space = ExecutionSpaceFromDevice::space; host_flag = (execution_space == LAMMPS_NS::Host); @@ -40,6 +41,7 @@ SNAKokkos::SNAKokkos(real_type rfac0_in, rfac0 = rfac0_in; rmin0 = rmin0_in; switch_flag = switch_flag_in; + switch_inner_flag = switch_inner_flag_in; bzero_flag = bzero_flag_in; chem_flag = chem_flag_in; @@ -295,6 +297,8 @@ void SNAKokkos::grow_rij(int newnatom, int rij = t_sna_3d(Kokkos::NoInit("sna:rij"),natom,nmax,3); wj = t_sna_2d(Kokkos::NoInit("sna:wj"),natom,nmax); rcutij = t_sna_2d(Kokkos::NoInit("sna:rcutij"),natom,nmax); + sinnerij = t_sna_2d(Kokkos::NoInit("sna:sinnerij"),natom,nmax); + dinnerij = t_sna_2d(Kokkos::NoInit("sna:dinnerij"),natom,nmax); inside = t_sna_2i(Kokkos::NoInit("sna:inside"),natom,nmax); element = t_sna_2i(Kokkos::NoInit("sna:element"),natom,nmax); dedr = t_sna_3d(Kokkos::NoInit("sna:dedr"),natom,nmax,3); @@ -371,6 +375,8 @@ void SNAKokkos::compute_cayley_klein(const const real_type rsq = x * x + y * y + z * z; const real_type r = sqrt(rsq); const real_type rcut = rcutij(iatom, jnbor); + const real_type sinner = sinnerij(iatom, jnbor); + const real_type dinner = dinnerij(iatom, jnbor); const real_type rscale0 = rfac0 * static_cast(MY_PI) / (rcut - rmin0); const real_type theta0 = (r - rmin0) * rscale0; const real_type sn = sin(theta0); @@ -380,7 +386,7 @@ void SNAKokkos::compute_cayley_klein(const const real_type wj_local = wj(iatom, jnbor); real_type sfac, dsfac; - compute_s_dsfac(r, rcut, sfac, dsfac); + compute_s_dsfac(r, rcut, sinner, dinner, sfac, dsfac); sfac *= wj_local; dsfac *= wj_local; @@ -1238,7 +1244,7 @@ void SNAKokkos::compute_ui_cpu(const typen z0 = r / tan(theta0); compute_uarray_cpu(team, iatom, jnbor, x, y, z, z0, r); - add_uarraytot(team, iatom, jnbor, r, wj(iatom,jnbor), rcutij(iatom,jnbor), element(iatom, jnbor)); + add_uarraytot(team, iatom, jnbor, r, wj(iatom,jnbor), rcutij(iatom,jnbor), sinnerij(iatom,jnbor), dinnerij(iatom,jnbor), element(iatom, jnbor)); } /* ---------------------------------------------------------------------- @@ -1541,7 +1547,7 @@ void SNAKokkos::compute_duidrj_cpu(const t z0 = r * cs / sn; dz0dr = z0 / r - (r*rscale0) * (rsq + z0 * z0) / rsq; - compute_duarray_cpu(team, iatom, jnbor, x, y, z, z0, r, dz0dr, wj(iatom,jnbor), rcutij(iatom,jnbor)); + compute_duarray_cpu(team, iatom, jnbor, x, y, z, z0, r, dz0dr, wj(iatom,jnbor), rcutij(iatom,jnbor), sinnerij(iatom,jnbor), dinnerij(iatom,jnbor)); } @@ -1623,9 +1629,10 @@ void SNAKokkos::compute_deidrj_cpu(const t template KOKKOS_INLINE_FUNCTION void SNAKokkos::add_uarraytot(const typename Kokkos::TeamPolicy::member_type& team, int iatom, int jnbor, - const real_type& r, const real_type& wj, const real_type& rcut, int jelem) + const real_type& r, const real_type& wj, const real_type& rcut, + const real_type& sinner, const real_type& dinner, int jelem) { - const real_type sfac = compute_sfac(r, rcut) * wj; + const real_type sfac = compute_sfac(r, rcut, sinner, dinner) * wj; Kokkos::parallel_for(Kokkos::ThreadVectorRange(team,twojmax+1), [&] (const int& j) { @@ -1746,7 +1753,8 @@ KOKKOS_INLINE_FUNCTION void SNAKokkos::compute_duarray_cpu(const typename Kokkos::TeamPolicy::member_type& team, int iatom, int jnbor, const real_type& x, const real_type& y, const real_type& z, const real_type& z0, const real_type& r, const real_type& dz0dr, - const real_type& wj, const real_type& rcut) + const real_type& wj, const real_type& rcut, + const real_type& sinner, const real_type& dinner) { real_type r0inv; real_type a_r, a_i, b_r, b_i; @@ -1872,8 +1880,8 @@ void SNAKokkos::compute_duarray_cpu(const }); } - real_type sfac = compute_sfac(r, rcut); - real_type dsfac = compute_dsfac(r, rcut); + real_type sfac = compute_sfac(r, rcut, sinner, dinner); + real_type dsfac = compute_dsfac(r, rcut, sinner, dinner); sfac *= wj; dsfac *= wj; @@ -2224,7 +2232,7 @@ int SNAKokkos::compute_ncoeff() template KOKKOS_INLINE_FUNCTION -real_type SNAKokkos::compute_sfac(real_type r, real_type rcut) +real_type SNAKokkos::compute_sfac(real_type r, real_type rcut, real_type sinner, real_type dinner) { constexpr real_type one = static_cast(1.0); constexpr real_type zero = static_cast(0.0); @@ -2235,18 +2243,30 @@ real_type SNAKokkos::compute_sfac(real_typ else if (r > rcut) return zero; else { real_type rcutfac = static_cast(MY_PI) / (rcut - rmin0); - return onehalf * (cos((r - rmin0) * rcutfac) + one); + if (switch_inner_flag == 0) + return onehalf * (cos((r - rmin0) * rcutfac) + one); + if (switch_inner_flag == 1) { + if (r >= sinner + dinner) + return onehalf * (cos((r - rmin0) * rcutfac) + one); + else if (r > sinner - dinner) { + real_type rcutfacinner = static_cast(MY_PI2) / dinner; + return onehalf * (cos((r - rmin0) * rcutfac) + one) * + onehalf * (one - cos(static_cast(MY_PI2) + (r - sinner) * rcutfacinner)); + } else return zero; + } + return zero; // dummy return } } - return zero; + return zero; // dummy return } /* ---------------------------------------------------------------------- */ template KOKKOS_INLINE_FUNCTION -real_type SNAKokkos::compute_dsfac(real_type r, real_type rcut) +real_type SNAKokkos::compute_dsfac(real_type r, real_type rcut, real_type sinner, real_type dinner) { + constexpr real_type one = static_cast(1.0); constexpr real_type zero = static_cast(0.0); constexpr real_type onehalf = static_cast(0.5); if (switch_flag == 0) return zero; @@ -2258,12 +2278,12 @@ real_type SNAKokkos::compute_dsfac(real_ty return -onehalf * sin((r - rmin0) * rcutfac) * rcutfac; } } - return zero; + return zero; // dummy return } template KOKKOS_INLINE_FUNCTION -void SNAKokkos::compute_s_dsfac(const real_type r, const real_type rcut, real_type& sfac, real_type& dsfac) { +void SNAKokkos::compute_s_dsfac(const real_type r, const real_type rcut, const real_type sinner, const real_type dinner, real_type& sfac, real_type& dsfac) { constexpr real_type one = static_cast(1.0); constexpr real_type zero = static_cast(0.0); constexpr real_type onehalf = static_cast(0.5); @@ -2356,6 +2376,8 @@ double SNAKokkos::memory_usage() bytes += natom * nmax * sizeof(real_type); // inside bytes += natom * nmax * sizeof(real_type); // wj bytes += natom * nmax * sizeof(real_type); // rcutij + bytes += natom * nmax * sizeof(real_type); // sinnerij + bytes += natom * nmax * sizeof(real_type); // dinnerij return bytes; } diff --git a/src/ML-SNAP/compute_sna_atom.cpp b/src/ML-SNAP/compute_sna_atom.cpp index cec705770d..5f24ebeaa6 100644 --- a/src/ML-SNAP/compute_sna_atom.cpp +++ b/src/ML-SNAP/compute_sna_atom.cpp @@ -32,7 +32,7 @@ using namespace LAMMPS_NS; ComputeSNAAtom::ComputeSNAAtom(LAMMPS *lmp, int narg, char **arg) : Compute(lmp, narg, arg), cutsq(nullptr), list(nullptr), sna(nullptr), - radelem(nullptr), wjelem(nullptr), rinnerelem(nullptr), drinnerelem(nullptr) + radelem(nullptr), wjelem(nullptr), sinnerelem(nullptr), dinnerelem(nullptr) { double rmin0, rfac0; @@ -85,6 +85,11 @@ ComputeSNAAtom::ComputeSNAAtom(LAMMPS *lmp, int narg, char **arg) : } } + // set local input checks + + int sinnerflag = 0; + int dinnerflag = 0; + // process optional args int iarg = nargmin; @@ -134,21 +139,37 @@ ComputeSNAAtom::ComputeSNAAtom(LAMMPS *lmp, int narg, char **arg) : wselfallflag = atoi(arg[iarg+1]); iarg += 2; } else if (strcmp(arg[iarg],"switchinnerflag") == 0) { - if (iarg+1+2*ntypes > narg) + if (iarg+2 > narg) error->all(FLERR,"Illegal compute sna/atom command"); - switchinnerflag = 1; + switchinnerflag = atoi(arg[iarg+1]); + iarg += 2; + } else if (strcmp(arg[iarg],"sinner") == 0) { iarg++; - memory->create(rinnerelem,ntypes+1,"sna/atom:rinnerelem"); - memory->create(drinnerelem,ntypes+1,"sna/atom:drinnerelem"); + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute sna/atom command"); + memory->create(sinnerelem,ntypes+1,"sna/atom:sinnerelem"); for (int i = 0; i < ntypes; i++) - rinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerflag = 1; iarg += ntypes; + } else if (strcmp(arg[iarg],"dinner") == 0) { + iarg++; + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute sna/atom command"); + memory->create(dinnerelem,ntypes+1,"sna/atom:dinnerelem"); for (int i = 0; i < ntypes; i++) - drinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerflag = 1; iarg += ntypes; } else error->all(FLERR,"Illegal compute sna/atom command"); } + if (switchinnerflag && !(sinnerflag && dinnerflag)) + error->all(FLERR,"Illegal compute sna/atom command: switchinnerflag = 1, missing sinner/dinner keyword"); + + if (!switchinnerflag && (sinnerflag || dinnerflag)) + error->all(FLERR,"Illegal compute sna/atom command: switchinnerflag = 0, unexpected sinner/dinner keyword"); + snaptr = new SNA(lmp, rfac0, twojmax, rmin0, switchflag, bzeroflag, chemflag, bnormflag, wselfallflag, @@ -176,8 +197,8 @@ ComputeSNAAtom::~ComputeSNAAtom() if (chemflag) memory->destroy(map); if (switchinnerflag) { - memory->destroy(rinnerelem); - memory->destroy(drinnerelem); + memory->destroy(sinnerelem); + memory->destroy(dinnerelem); } } @@ -282,8 +303,8 @@ void ComputeSNAAtom::compute_peratom() snaptr->wj[ninside] = wjelem[jtype]; snaptr->rcutij[ninside] = (radi+radelem[jtype])*rcutfac; if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5*(rinnerelem[itype]+rinnerelem[jtype]); - snaptr->drinnerij[ninside] = 0.5*(drinnerelem[itype]+drinnerelem[jtype]); + snaptr->sinnerij[ninside] = 0.5*(sinnerelem[itype]+sinnerelem[jtype]); + snaptr->dinnerij[ninside] = 0.5*(dinnerelem[itype]+dinnerelem[jtype]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; diff --git a/src/ML-SNAP/compute_sna_atom.h b/src/ML-SNAP/compute_sna_atom.h index 5e62fe47d9..a7140c3949 100644 --- a/src/ML-SNAP/compute_sna_atom.h +++ b/src/ML-SNAP/compute_sna_atom.h @@ -45,8 +45,8 @@ class ComputeSNAAtom : public Compute { int *map; // map types to [0,nelements) int nelements, chemflag; int switchinnerflag; - double *rinnerelem; - double *drinnerelem; + double *sinnerelem; + double *dinnerelem; class SNA *snaptr; double cutmax; int quadraticflag; diff --git a/src/ML-SNAP/compute_snad_atom.cpp b/src/ML-SNAP/compute_snad_atom.cpp index 4874a53c9c..838d8a85c9 100644 --- a/src/ML-SNAP/compute_snad_atom.cpp +++ b/src/ML-SNAP/compute_snad_atom.cpp @@ -32,7 +32,7 @@ using namespace LAMMPS_NS; ComputeSNADAtom::ComputeSNADAtom(LAMMPS *lmp, int narg, char **arg) : Compute(lmp, narg, arg), cutsq(nullptr), list(nullptr), snad(nullptr), - radelem(nullptr), wjelem(nullptr), rinnerelem(nullptr), drinnerelem(nullptr) + radelem(nullptr), wjelem(nullptr), sinnerelem(nullptr), dinnerelem(nullptr) { double rfac0, rmin0; int twojmax, switchflag, bzeroflag, bnormflag, wselfallflag; @@ -57,8 +57,8 @@ ComputeSNADAtom::ComputeSNADAtom(LAMMPS *lmp, int narg, char **arg) : // process required arguments - memory->create(radelem,ntypes+1,"sna/atom:radelem"); // offset by 1 to match up with types - memory->create(wjelem,ntypes+1,"sna/atom:wjelem"); + memory->create(radelem,ntypes+1,"snad/atom:radelem"); // offset by 1 to match up with types + memory->create(wjelem,ntypes+1,"snad/atom:wjelem"); rcutfac = atof(arg[3]); rfac0 = atof(arg[4]); twojmax = atoi(arg[5]); @@ -71,7 +71,7 @@ ComputeSNADAtom::ComputeSNADAtom(LAMMPS *lmp, int narg, char **arg) : double cut; cutmax = 0.0; - memory->create(cutsq,ntypes+1,ntypes+1,"sna/atom:cutsq"); + memory->create(cutsq,ntypes+1,ntypes+1,"snad/atom:cutsq"); for (int i = 1; i <= ntypes; i++) { cut = 2.0*radelem[i]*rcutfac; if (cut > cutmax) cutmax = cut; @@ -82,6 +82,11 @@ ComputeSNADAtom::ComputeSNADAtom(LAMMPS *lmp, int narg, char **arg) : } } + // set local input checks + + int sinnerflag = 0; + int dinnerflag = 0; + // process optional args int iarg = nargmin; @@ -131,21 +136,38 @@ ComputeSNADAtom::ComputeSNADAtom(LAMMPS *lmp, int narg, char **arg) : wselfallflag = atoi(arg[iarg+1]); iarg += 2; } else if (strcmp(arg[iarg],"switchinnerflag") == 0) { - if (iarg+1+2*ntypes > narg) + if (iarg+2 > narg) error->all(FLERR,"Illegal compute snad/atom command"); - switchinnerflag = 1; + switchinnerflag = atoi(arg[iarg+1]); + iarg += 2; + } else if (strcmp(arg[iarg],"sinner") == 0) { iarg++; - memory->create(rinnerelem,ntypes+1,"snad/atom:rinnerelem"); - memory->create(drinnerelem,ntypes+1,"snad/atom:drinnerelem"); + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute snad/atom command"); + memory->create(sinnerelem,ntypes+1,"snad/atom:sinnerelem"); for (int i = 0; i < ntypes; i++) - rinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerflag = 1; iarg += ntypes; + } else if (strcmp(arg[iarg],"dinner") == 0) { + iarg++; + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute snad/atom command"); + memory->create(dinnerelem,ntypes+1,"snad/atom:dinnerelem"); for (int i = 0; i < ntypes; i++) - drinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerflag = 1; iarg += ntypes; } else error->all(FLERR,"Illegal compute snad/atom command"); } + if (switchinnerflag && !(sinnerflag && dinnerflag)) + error->all(FLERR,"Illegal compute snad/atom command: switchinnerflag = 1, missing sinner/dinner keyword"); + + if (!switchinnerflag && (sinnerflag || dinnerflag)) + error->all(FLERR,"Illegal compute snad/atom command: switchinnerflag = 0, unexpected sinner/dinner keyword"); + + snaptr = new SNA(lmp, rfac0, twojmax, rmin0, switchflag, bzeroflag, chemflag, bnormflag, wselfallflag, @@ -177,8 +199,8 @@ ComputeSNADAtom::~ComputeSNADAtom() if (chemflag) memory->destroy(map); if (switchinnerflag) { - memory->destroy(rinnerelem); - memory->destroy(drinnerelem); + memory->destroy(sinnerelem); + memory->destroy(dinnerelem); } } @@ -190,7 +212,7 @@ void ComputeSNADAtom::init() error->all(FLERR,"Compute snad/atom requires a pair style be defined"); if (cutmax > force->pair->cutforce) - error->all(FLERR,"Compute sna/atom cutoff is longer than pairwise cutoff"); + error->all(FLERR,"Compute snad/atom cutoff is longer than pairwise cutoff"); // need an occasional full neighbor list @@ -299,8 +321,8 @@ void ComputeSNADAtom::compute_peratom() snaptr->wj[ninside] = wjelem[jtype]; snaptr->rcutij[ninside] = (radi+radelem[jtype])*rcutfac; if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5*(rinnerelem[itype]+rinnerelem[jtype]); - snaptr->drinnerij[ninside] = 0.5*(drinnerelem[itype]+drinnerelem[jtype]); + snaptr->sinnerij[ninside] = 0.5*(sinnerelem[itype]+sinnerelem[jtype]); + snaptr->dinnerij[ninside] = 0.5*(dinnerelem[itype]+dinnerelem[jtype]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; diff --git a/src/ML-SNAP/compute_snad_atom.h b/src/ML-SNAP/compute_snad_atom.h index c03458446f..2e2842a0c1 100644 --- a/src/ML-SNAP/compute_snad_atom.h +++ b/src/ML-SNAP/compute_snad_atom.h @@ -47,8 +47,8 @@ class ComputeSNADAtom : public Compute { int *map; // map types to [0,nelements) int nelements, chemflag; int switchinnerflag; - double *rinnerelem; - double *drinnerelem; + double *sinnerelem; + double *dinnerelem; class SNA *snaptr; double cutmax; int quadraticflag; diff --git a/src/ML-SNAP/compute_snap.cpp b/src/ML-SNAP/compute_snap.cpp index 04dac2acba..83f56e338c 100644 --- a/src/ML-SNAP/compute_snap.cpp +++ b/src/ML-SNAP/compute_snap.cpp @@ -35,7 +35,7 @@ enum{SCALAR,VECTOR,ARRAY}; ComputeSnap::ComputeSnap(LAMMPS *lmp, int narg, char **arg) : Compute(lmp, narg, arg), cutsq(nullptr), list(nullptr), snap(nullptr), snapall(nullptr), snap_peratom(nullptr), radelem(nullptr), wjelem(nullptr), - rinnerelem(nullptr), drinnerelem(nullptr), snaptr(nullptr) + sinnerelem(nullptr), dinnerelem(nullptr), snaptr(nullptr) { array_flag = 1; @@ -89,6 +89,11 @@ ComputeSnap::ComputeSnap(LAMMPS *lmp, int narg, char **arg) : } } + // set local input checks + + int sinnerflag = 0; + int dinnerflag = 0; + // process optional args int iarg = nargmin; @@ -115,7 +120,7 @@ ComputeSnap::ComputeSnap(LAMMPS *lmp, int narg, char **arg) : quadraticflag = atoi(arg[iarg+1]); iarg += 2; } else if (strcmp(arg[iarg],"chem") == 0) { - if (iarg+2 > narg) + if (iarg+2+ntypes > narg) error->all(FLERR,"Illegal compute snap command"); chemflag = 1; memory->create(map,ntypes+1,"compute_snap:map"); @@ -143,21 +148,37 @@ ComputeSnap::ComputeSnap(LAMMPS *lmp, int narg, char **arg) : bikflag = atoi(arg[iarg+1]); iarg += 2; } else if (strcmp(arg[iarg],"switchinnerflag") == 0) { - if (iarg+1+2*ntypes > narg) + if (iarg+2 > narg) error->all(FLERR,"Illegal compute snap command"); - switchinnerflag = 1; + switchinnerflag = atoi(arg[iarg+1]); + iarg += 2; + } else if (strcmp(arg[iarg],"sinner") == 0) { iarg++; - memory->create(rinnerelem,ntypes+1,"snap:rinnerelem"); - memory->create(drinnerelem,ntypes+1,"snap:drinnerelem"); + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute snap command"); + memory->create(sinnerelem,ntypes+1,"snap:sinnerelem"); for (int i = 0; i < ntypes; i++) - rinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerflag = 1; iarg += ntypes; + } else if (strcmp(arg[iarg],"dinner") == 0) { + iarg++; + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute snap command"); + memory->create(dinnerelem,ntypes+1,"snap:dinnerelem"); for (int i = 0; i < ntypes; i++) - drinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerflag = 1; iarg += ntypes; } else error->all(FLERR,"Illegal compute snap command"); } + if (switchinnerflag && !(sinnerflag && dinnerflag)) + error->all(FLERR,"Illegal compute snap command: switchinnerflag = 1, missing sinner/dinner keyword"); + + if (!switchinnerflag && (sinnerflag || dinnerflag)) + error->all(FLERR,"Illegal compute snap command: switchinnerflag = 0, unexpected sinner/dinner keyword"); + snaptr = new SNA(lmp, rfac0, twojmax, rmin0, switchflag, bzeroflag, chemflag, bnormflag, wselfallflag, @@ -198,8 +219,8 @@ ComputeSnap::~ComputeSnap() if (chemflag) memory->destroy(map); if (switchinnerflag) { - memory->destroy(rinnerelem); - memory->destroy(drinnerelem); + memory->destroy(sinnerelem); + memory->destroy(dinnerelem); } } @@ -353,8 +374,8 @@ void ComputeSnap::compute_array() snaptr->wj[ninside] = wjelem[jtype]; snaptr->rcutij[ninside] = (radi+radelem[jtype])*rcutfac; if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5*(rinnerelem[itype]+rinnerelem[jtype]); - snaptr->drinnerij[ninside] = 0.5*(drinnerelem[itype]+drinnerelem[jtype]); + snaptr->sinnerij[ninside] = 0.5*(sinnerelem[itype]+sinnerelem[jtype]); + snaptr->dinnerij[ninside] = 0.5*(dinnerelem[itype]+dinnerelem[jtype]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; diff --git a/src/ML-SNAP/compute_snap.h b/src/ML-SNAP/compute_snap.h index fa31de747d..bf631539ba 100644 --- a/src/ML-SNAP/compute_snap.h +++ b/src/ML-SNAP/compute_snap.h @@ -47,8 +47,8 @@ class ComputeSnap : public Compute { int *map; // map types to [0,nelements) int nelements, chemflag; int switchinnerflag; - double *rinnerelem; - double *drinnerelem; + double *sinnerelem; + double *dinnerelem; class SNA *snaptr; double cutmax; int quadraticflag; diff --git a/src/ML-SNAP/compute_snav_atom.cpp b/src/ML-SNAP/compute_snav_atom.cpp index d3999381f2..0baa4bc110 100644 --- a/src/ML-SNAP/compute_snav_atom.cpp +++ b/src/ML-SNAP/compute_snav_atom.cpp @@ -31,7 +31,7 @@ using namespace LAMMPS_NS; ComputeSNAVAtom::ComputeSNAVAtom(LAMMPS *lmp, int narg, char **arg) : Compute(lmp, narg, arg), cutsq(nullptr), list(nullptr), snav(nullptr), - radelem(nullptr), wjelem(nullptr), rinnerelem(nullptr), drinnerelem(nullptr) + radelem(nullptr), wjelem(nullptr), sinnerelem(nullptr), dinnerelem(nullptr) { double rfac0, rmin0; int twojmax, switchflag, bzeroflag, bnormflag, wselfallflag; @@ -56,8 +56,8 @@ ComputeSNAVAtom::ComputeSNAVAtom(LAMMPS *lmp, int narg, char **arg) : // process required arguments - memory->create(radelem,ntypes+1,"sna/atom:radelem"); // offset by 1 to match up with types - memory->create(wjelem,ntypes+1,"sna/atom:wjelem"); + memory->create(radelem,ntypes+1,"snav/atom:radelem"); // offset by 1 to match up with types + memory->create(wjelem,ntypes+1,"snav/atom:wjelem"); rcutfac = atof(arg[3]); rfac0 = atof(arg[4]); twojmax = atoi(arg[5]); @@ -67,7 +67,7 @@ ComputeSNAVAtom::ComputeSNAVAtom(LAMMPS *lmp, int narg, char **arg) : wjelem[i+1] = atof(arg[6+ntypes+i]); // construct cutsq double cut; - memory->create(cutsq,ntypes+1,ntypes+1,"sna/atom:cutsq"); + memory->create(cutsq,ntypes+1,ntypes+1,"snav/atom:cutsq"); for (int i = 1; i <= ntypes; i++) { cut = 2.0*radelem[i]*rcutfac; cutsq[i][i] = cut*cut; @@ -77,6 +77,11 @@ ComputeSNAVAtom::ComputeSNAVAtom(LAMMPS *lmp, int narg, char **arg) : } } + // set local input checks + + int sinnerflag = 0; + int dinnerflag = 0; + // process optional args int iarg = nargmin; @@ -104,7 +109,7 @@ ComputeSNAVAtom::ComputeSNAVAtom(LAMMPS *lmp, int narg, char **arg) : iarg += 2; } else if (strcmp(arg[iarg],"chem") == 0) { if (iarg+2+ntypes > narg) - error->all(FLERR,"Illegal compute sna/atom command"); + error->all(FLERR,"Illegal compute snav/atom command"); chemflag = 1; memory->create(map,ntypes+1,"compute_sna_atom:map"); nelements = utils::inumeric(FLERR,arg[iarg+1],false,lmp); @@ -126,21 +131,37 @@ ComputeSNAVAtom::ComputeSNAVAtom(LAMMPS *lmp, int narg, char **arg) : wselfallflag = atoi(arg[iarg+1]); iarg += 2; } else if (strcmp(arg[iarg],"switchinnerflag") == 0) { - if (iarg+1+2*ntypes > narg) + if (iarg+2 > narg) error->all(FLERR,"Illegal compute snav/atom command"); - switchinnerflag = 1; + switchinnerflag = atoi(arg[iarg+1]); + iarg += 2; + } else if (strcmp(arg[iarg],"sinner") == 0) { iarg++; - memory->create(rinnerelem,ntypes+1,"snav/atom:rinnerelem"); - memory->create(drinnerelem,ntypes+1,"snav/atom:drinnerelem"); + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute snav/atom command"); + memory->create(sinnerelem,ntypes+1,"snav/atom:sinnerelem"); for (int i = 0; i < ntypes; i++) - rinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + sinnerflag = 1; iarg += ntypes; + } else if (strcmp(arg[iarg],"dinner") == 0) { + iarg++; + if (iarg+ntypes > narg) + error->all(FLERR,"Illegal compute snav/atom command"); + memory->create(dinnerelem,ntypes+1,"snav/atom:dinnerelem"); for (int i = 0; i < ntypes; i++) - drinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerelem[i+1] = utils::numeric(FLERR,arg[iarg+i],false,lmp); + dinnerflag = 1; iarg += ntypes; } else error->all(FLERR,"Illegal compute snav/atom command"); } + if (switchinnerflag && !(sinnerflag && dinnerflag)) + error->all(FLERR,"Illegal compute snav/atom command: switchinnerflag = 1, missing sinner/dinner keyword"); + + if (!switchinnerflag && (sinnerflag || dinnerflag)) + error->all(FLERR,"Illegal compute snav/atom command: switchinnerflag = 0, unexpected sinner/dinner keyword"); + snaptr = new SNA(lmp, rfac0, twojmax, rmin0, switchflag, bzeroflag, chemflag, bnormflag, wselfallflag, @@ -171,8 +192,8 @@ ComputeSNAVAtom::~ComputeSNAVAtom() if (chemflag) memory->destroy(map); if (switchinnerflag) { - memory->destroy(rinnerelem); - memory->destroy(drinnerelem); + memory->destroy(sinnerelem); + memory->destroy(dinnerelem); } } @@ -291,8 +312,8 @@ void ComputeSNAVAtom::compute_peratom() snaptr->wj[ninside] = wjelem[jtype]; snaptr->rcutij[ninside] = (radi+radelem[jtype])*rcutfac; if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5*(rinnerelem[itype]+rinnerelem[jtype]); - snaptr->drinnerij[ninside] = 0.5*(drinnerelem[itype]+drinnerelem[jtype]); + snaptr->sinnerij[ninside] = 0.5*(sinnerelem[itype]+sinnerelem[jtype]); + snaptr->dinnerij[ninside] = 0.5*(dinnerelem[itype]+dinnerelem[jtype]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; diff --git a/src/ML-SNAP/compute_snav_atom.h b/src/ML-SNAP/compute_snav_atom.h index be442f3a49..335f994323 100644 --- a/src/ML-SNAP/compute_snav_atom.h +++ b/src/ML-SNAP/compute_snav_atom.h @@ -47,8 +47,8 @@ class ComputeSNAVAtom : public Compute { int *map; // map types to [0,nelements) int nelements, chemflag; int switchinnerflag; - double *rinnerelem; - double *drinnerelem; + double *sinnerelem; + double *dinnerelem; class SNA *snaptr; int quadraticflag; }; diff --git a/src/ML-SNAP/pair_snap.cpp b/src/ML-SNAP/pair_snap.cpp index 3ba14c5872..9c076710da 100644 --- a/src/ML-SNAP/pair_snap.cpp +++ b/src/ML-SNAP/pair_snap.cpp @@ -45,8 +45,8 @@ PairSNAP::PairSNAP(LAMMPS *lmp) : Pair(lmp) radelem = nullptr; wjelem = nullptr; coeffelem = nullptr; - rinnerelem = nullptr; - drinnerelem = nullptr; + sinnerelem = nullptr; + dinnerelem = nullptr; beta_max = 0; beta = nullptr; @@ -65,8 +65,8 @@ PairSNAP::~PairSNAP() memory->destroy(coeffelem); if (switchinnerflag) { - memory->destroy(rinnerelem); - memory->destroy(drinnerelem); + memory->destroy(sinnerelem); + memory->destroy(dinnerelem); } memory->destroy(beta); @@ -158,8 +158,8 @@ void PairSNAP::compute(int eflag, int vflag) snaptr->wj[ninside] = wjelem[jelem]; snaptr->rcutij[ninside] = (radi + radelem[jelem])*rcutfac; if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5*(rinnerelem[ielem]+rinnerelem[jelem]); - snaptr->drinnerij[ninside] = 0.5*(drinnerelem[ielem]+drinnerelem[jelem]); + snaptr->sinnerij[ninside] = 0.5*(sinnerelem[ielem]+sinnerelem[jelem]); + snaptr->dinnerij[ninside] = 0.5*(dinnerelem[ielem]+dinnerelem[jelem]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; @@ -331,8 +331,8 @@ void PairSNAP::compute_bispectrum() snaptr->wj[ninside] = wjelem[jelem]; snaptr->rcutij[ninside] = (radi + radelem[jelem])*rcutfac; if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5*(rinnerelem[ielem]+rinnerelem[jelem]); - snaptr->drinnerij[ninside] = 0.5*(drinnerelem[ielem]+drinnerelem[jelem]); + snaptr->sinnerij[ninside] = 0.5*(sinnerelem[ielem]+sinnerelem[jelem]); + snaptr->dinnerij[ninside] = 0.5*(dinnerelem[ielem]+dinnerelem[jelem]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; @@ -516,13 +516,13 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) memory->destroy(radelem); memory->destroy(wjelem); memory->destroy(coeffelem); - memory->destroy(rinnerelem); - memory->destroy(drinnerelem); + memory->destroy(sinnerelem); + memory->destroy(dinnerelem); memory->create(radelem,nelements,"pair:radelem"); memory->create(wjelem,nelements,"pair:wjelem"); memory->create(coeffelem,nelements,ncoeffall,"pair:coeffelem"); - memory->create(rinnerelem,nelements,"pair:rinnerelem"); - memory->create(drinnerelem,nelements,"pair:drinnerelem"); + memory->create(sinnerelem,nelements,"pair:sinnerelem"); + memory->create(dinnerelem,nelements,"pair:dinnerelem"); // initialize checklist for all required nelements @@ -644,8 +644,8 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) // set local input checks - int rinnerflag = 0; - int drinnerflag = 0; + int sinnerflag = 0; + int dinnerflag = 0; // open SNAP parameter file on proc 0 @@ -690,7 +690,7 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) // check for keywords with one value per element - if (keywd == "rinner" || keywd == "drinner") { + if (keywd == "sinner" || keywd == "dinner") { if (words.size() != nelements+1) error->all(FLERR,"Incorrect SNAP parameter file"); @@ -700,20 +700,20 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) int iword = 1; - if (keywd == "rinner") { + if (keywd == "sinner") { keyval = words[iword]; for (int ielem = 0; ielem < nelements; ielem++) { - rinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); + sinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); iword++; } - rinnerflag = 1; - } else if (keywd == "drinner") { + sinnerflag = 1; + } else if (keywd == "dinner") { keyval = words[iword]; for (int ielem = 0; ielem < nelements; ielem++) { - drinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); + dinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); iword++; } - drinnerflag = 1; + dinnerflag = 1; } } else { @@ -765,10 +765,10 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) if (chemflag && nelemtmp != nelements) error->all(FLERR,"Incorrect SNAP parameter file"); - if (switchinnerflag && !(rinnerflag && drinnerflag)) + if (switchinnerflag && !(sinnerflag && dinnerflag)) error->all(FLERR,"Incorrect SNAP parameter file"); - if (!switchinnerflag && (rinnerflag || drinnerflag)) + if (!switchinnerflag && (sinnerflag || dinnerflag)) error->all(FLERR,"Incorrect SNAP parameter file"); } diff --git a/src/ML-SNAP/pair_snap.h b/src/ML-SNAP/pair_snap.h index 77ceeabd4f..f260704801 100644 --- a/src/ML-SNAP/pair_snap.h +++ b/src/ML-SNAP/pair_snap.h @@ -60,8 +60,8 @@ class PairSNAP : public Pair { int twojmax, switchflag, bzeroflag, bnormflag; int chemflag, wselfallflag; int switchinnerflag; // inner cutoff switch - double *rinnerelem; // element inner cutoff - double *drinnerelem; // element inner cutoff range + double *sinnerelem; // element inner cutoff midpoint + double *dinnerelem; // element inner cutoff half-width int chunksize,parallel_thresh; double rfac0, rmin0, wj1, wj2; int rcutfacflag, twojmaxflag; // flags for required parameters diff --git a/src/ML-SNAP/sna.cpp b/src/ML-SNAP/sna.cpp index 717d0014f5..0bd66d0bce 100644 --- a/src/ML-SNAP/sna.cpp +++ b/src/ML-SNAP/sna.cpp @@ -145,8 +145,8 @@ SNA::SNA(LAMMPS* lmp, double rfac0_in, int twojmax_in, inside = nullptr; wj = nullptr; rcutij = nullptr; - rinnerij = nullptr; - drinnerij = nullptr; + sinnerij = nullptr; + dinnerij = nullptr; element = nullptr; nmax = 0; idxz = nullptr; @@ -176,10 +176,8 @@ SNA::~SNA() memory->destroy(inside); memory->destroy(wj); memory->destroy(rcutij); - if (switch_inner_flag) { - memory->destroy(rinnerij); - memory->destroy(drinnerij); - } + memory->destroy(sinnerij); + memory->destroy(dinnerij); if (chem_flag) memory->destroy(element); memory->destroy(ulist_r_ij); memory->destroy(ulist_i_ij); @@ -327,10 +325,8 @@ void SNA::grow_rij(int newnmax) memory->destroy(inside); memory->destroy(wj); memory->destroy(rcutij); - if (switch_inner_flag) { - memory->destroy(rinnerij); - memory->destroy(drinnerij); - } + memory->destroy(sinnerij); + memory->destroy(dinnerij); if (chem_flag) memory->destroy(element); memory->destroy(ulist_r_ij); memory->destroy(ulist_i_ij); @@ -338,10 +334,8 @@ void SNA::grow_rij(int newnmax) memory->create(inside, nmax, "pair:inside"); memory->create(wj, nmax, "pair:wj"); memory->create(rcutij, nmax, "pair:rcutij"); - if (switch_inner_flag) { - memory->create(rinnerij, nmax, "pair:rinnerij"); - memory->create(drinnerij, nmax, "pair:drinnerij"); - } + memory->create(sinnerij, nmax, "pair:sinnerij"); + memory->create(dinnerij, nmax, "pair:dinnerij"); if (chem_flag) memory->create(element, nmax, "sna:element"); memory->create(ulist_r_ij, nmax, idxu_max, "sna:ulist_ij"); memory->create(ulist_i_ij, nmax, idxu_max, "sna:ulist_ij"); @@ -1020,12 +1014,9 @@ void SNA::add_uarraytot(double r, int jj) double sfac; int jelem; - sfac = compute_sfac(r, rcutij[jj]); + sfac = compute_sfac(r, rcutij[jj], sinnerij[jj], dinnerij[jj]); sfac *= wj[jj]; - if (switch_inner_flag) - sfac *= compute_sfac_inner(r, rinnerij[jj], drinnerij[jj]); - if (chem_flag) jelem = element[jj]; else jelem = 0; @@ -1268,14 +1259,8 @@ void SNA::compute_duarray(double x, double y, double z, } } - double sfac = compute_sfac(r, rcut); - double dsfac = compute_dsfac(r, rcut); - - if (switch_inner_flag) { - double sfac_inner = compute_sfac_inner(r, rinnerij[jj], drinnerij[jj]); - dsfac = dsfac*sfac_inner + sfac*compute_dsfac_inner(r, rinnerij[jj], drinnerij[jj]); - sfac *= sfac_inner; - } + double sfac = compute_sfac(r, rcut, sinnerij[jj], dinnerij[jj]); + double dsfac = compute_dsfac(r, rcut, sinnerij[jj], dinnerij[jj]); sfac *= wj; dsfac *= wj; @@ -1339,10 +1324,8 @@ double SNA::memory_usage() bytes += (double)nmax * sizeof(int); // inside bytes += (double)nmax * sizeof(double); // wj bytes += (double)nmax * sizeof(double); // rcutij - if (switch_inner_flag) { - bytes += (double)nmax * sizeof(double); // rinnerij - bytes += (double)nmax * sizeof(double); // drinnerij - } + bytes += (double)nmax * sizeof(double); // sinnerij + bytes += (double)nmax * sizeof(double); // dinnerij if (chem_flag) bytes += (double)nmax * sizeof(int); // element return bytes; @@ -1547,7 +1530,7 @@ void SNA::compute_ncoeff() /* ---------------------------------------------------------------------- */ -double SNA::compute_sfac(double r, double rcut) +double SNA::compute_sfac(double r, double rcut, double sinner, double dinner) { double sfac; if (switch_flag == 0) sfac = 1.0; @@ -1557,14 +1540,22 @@ double SNA::compute_sfac(double r, double rcut) double rcutfac = MY_PI / (rcut - rmin0); sfac = 0.5 * (cos((r - rmin0) * rcutfac) + 1.0); } + + if (switch_inner_flag == 1 && r < sinner + dinner) { + if (r > sinner - dinner) { + double rcutfac = MY_PI2 / dinner; + sfac *= 0.5 * (1.0 - cos(MY_PI2 + (r - sinner) * rcutfac)); + } else sfac = 0.0; + } + return sfac; } /* ---------------------------------------------------------------------- */ -double SNA::compute_dsfac(double r, double rcut) +double SNA::compute_dsfac(double r, double rcut, double sinner, double dinner) { - double dsfac; + double sfac, dsfac, sfac_inner, dsfac_inner; if (switch_flag == 0) dsfac = 0.0; else if (r <= rmin0) dsfac = 0.0; else if (r > rcut) dsfac = 0.0; @@ -1572,36 +1563,26 @@ double SNA::compute_dsfac(double r, double rcut) double rcutfac = MY_PI / (rcut - rmin0); dsfac = -0.5 * sin((r - rmin0) * rcutfac) * rcutfac; } - return dsfac; -} - -/* ---------------------------------------------------------------------- */ - -double SNA::compute_sfac_inner(double r, double rinner, double drinner) -{ - double sfac; - if (switch_inner_flag == 0) sfac = 1.0; - else if (r >= rinner + drinner) sfac = 1.0; - else if (r <= rinner - drinner) sfac = 0.0; - else { - double rcutfac = MY_PI2 / drinner; - sfac = 0.5 * (1.0 - cos(MY_PI2 + (r - rinner) * rcutfac)); - } - return sfac; -} - -/* ---------------------------------------------------------------------- */ - -double SNA::compute_dsfac_inner(double r, double rinner, double drinner) -{ - double dsfac; - if (switch_inner_flag == 0) dsfac = 0.0; - else if (r >= rinner + drinner) dsfac = 0.0; - else if (r <= rinner - drinner) dsfac = 0.0; - else { - double rcutfac = MY_PI2 / drinner; - dsfac = 0.5 * rcutfac * sin(MY_PI2 + (r - rinner) * rcutfac); - } + + // duplicated computation, but rarely visited + + if (switch_inner_flag == 1 && r < sinner + dinner) { + if (r > sinner - dinner) { + if (switch_flag == 0) sfac = 1.0; + else if (r <= rmin0) sfac = 1.0; + else if (r > rcut) sfac = 0.0; + else { + double rcutfac = MY_PI / (rcut - rmin0); + sfac = 0.5 * (cos((r - rmin0) * rcutfac) + 1.0); + } + + double rcutfac = MY_PI2 / dinner; + sfac_inner = 0.5 * (1.0 - cos(MY_PI2 + (r - sinner) * rcutfac)); + dsfac_inner = 0.5 * rcutfac * sin(MY_PI2 + (r - sinner) * rcutfac); + dsfac = dsfac*sfac_inner + sfac*dsfac_inner; + } else dsfac = 0.0; + } + return dsfac; } diff --git a/src/ML-SNAP/sna.h b/src/ML-SNAP/sna.h index 03d77d52b4..6eaf774941 100644 --- a/src/ML-SNAP/sna.h +++ b/src/ML-SNAP/sna.h @@ -56,11 +56,8 @@ class SNA : protected Pointers { void compute_duidrj(int); void compute_dbidrj(); void compute_deidrj(double *); - double compute_sfac(double, double); - double compute_dsfac(double, double); - - double compute_sfac_inner(double, double, double); - double compute_dsfac_inner(double, double, double); + double compute_sfac(double, double, double, double); + double compute_dsfac(double, double, double, double); // public bispectrum data @@ -80,8 +77,8 @@ class SNA : protected Pointers { // only allocated for switch_inner_flag=1 - double *rinnerij; // short inner cutoff list - double *drinnerij;// short inner range list + double *sinnerij; // short inner cutoff midpoint list + double *dinnerij; // short inner half-width list // only allocated for chem_flag=1 From 2c71d0eea25c3d4cf0bec407596e57e81d65c3bc Mon Sep 17 00:00:00 2001 From: Aidan Thompson Date: Fri, 22 Apr 2022 16:43:19 -0600 Subject: [PATCH 04/16] Ported to KOKKOS, untested --- src/KOKKOS/sna_kokkos_impl.h | 94 +++++++++++++++++++++++++++--------- src/ML-SNAP/sna.cpp | 34 ++++++++----- 2 files changed, 92 insertions(+), 36 deletions(-) diff --git a/src/KOKKOS/sna_kokkos_impl.h b/src/KOKKOS/sna_kokkos_impl.h index 546c07ee11..c9c1bfbb99 100644 --- a/src/KOKKOS/sna_kokkos_impl.h +++ b/src/KOKKOS/sna_kokkos_impl.h @@ -2234,29 +2234,30 @@ template KOKKOS_INLINE_FUNCTION real_type SNAKokkos::compute_sfac(real_type r, real_type rcut, real_type sinner, real_type dinner) { + real_type sfac_outer; constexpr real_type one = static_cast(1.0); constexpr real_type zero = static_cast(0.0); constexpr real_type onehalf = static_cast(0.5); - if (switch_flag == 0) return one; + if (switch_flag == 0) sfac_outer = one; if (switch_flag == 1) { - if (r <= rmin0) return one; + if (r <= rmin0) sfac_outer = one; else if (r > rcut) return zero; else { real_type rcutfac = static_cast(MY_PI) / (rcut - rmin0); - if (switch_inner_flag == 0) - return onehalf * (cos((r - rmin0) * rcutfac) + one); - if (switch_inner_flag == 1) { - if (r >= sinner + dinner) - return onehalf * (cos((r - rmin0) * rcutfac) + one); - else if (r > sinner - dinner) { - real_type rcutfacinner = static_cast(MY_PI2) / dinner; - return onehalf * (cos((r - rmin0) * rcutfac) + one) * - onehalf * (one - cos(static_cast(MY_PI2) + (r - sinner) * rcutfacinner)); - } else return zero; - } - return zero; // dummy return + sfac_outer = onehalf * (cos((r - rmin0) * rcutfac) + one); } } + + if (switch_inner_flag == 0) return sfac_outer; + if (switch_inner_flag == 1) { + if (r >= sinner + dinner) + return sfac_outer; + else if (r > sinner - dinner) { + real_type rcutfac = static_cast(MY_PI2) / dinner; + return sfac_outer * + onehalf * (one - cos(static_cast(MY_PI2) + (r - sinner) * rcutfac)); + } else return zero; + } return zero; // dummy return } @@ -2266,41 +2267,86 @@ template KOKKOS_INLINE_FUNCTION real_type SNAKokkos::compute_dsfac(real_type r, real_type rcut, real_type sinner, real_type dinner) { + real_type sfac_outer, dsfac_outer, sfac_inner, dsfac_inner; constexpr real_type one = static_cast(1.0); constexpr real_type zero = static_cast(0.0); constexpr real_type onehalf = static_cast(0.5); - if (switch_flag == 0) return zero; + if (switch_flag == 0) dsfac_outer = zero; if (switch_flag == 1) { - if (r <= rmin0) return zero; + if (r <= rmin0) dsfac_outer = zero; else if (r > rcut) return zero; else { real_type rcutfac = static_cast(MY_PI) / (rcut - rmin0); - return -onehalf * sin((r - rmin0) * rcutfac) * rcutfac; + dsfac_outer = -onehalf * sin((r - rmin0) * rcutfac) * rcutfac; } } + + if (switch_inner_flag == 0) return dsfac_outer; + if (switch_inner_flag == 1) { + if (r >= sinner + dinner) + return dsfac_outer; + else if (r > sinner - dinner) { + + // calculate sfac_outer + + if (switch_flag == 0) sfac_outer = one; + if (switch_flag == 1) { + if (r <= rmin0) sfac_outer = one; + else if (r > rcut) sfac_outer = zero; + else { + real_type rcutfac = static_cast(MY_PI) / (rcut - rmin0); + sfac_outer = onehalf * (cos((r - rmin0) * rcutfac) + one); + } + } + + // calculate sfac_inner + + real_type rcutfac = static_cast(MY_PI2) / dinner; + sfac_inner = onehalf * (one - cos(static_cast(MY_PI2) + (r - sinner) * rcutfac)); + dsfac_inner = onehalf * rcutfac * sin(static_cast(MY_PI2) + (r - sinner) * rcutfac); + return dsfac_outer * sfac_inner + sfac_outer * dsfac_inner; + + } else return zero; + } return zero; // dummy return } template KOKKOS_INLINE_FUNCTION void SNAKokkos::compute_s_dsfac(const real_type r, const real_type rcut, const real_type sinner, const real_type dinner, real_type& sfac, real_type& dsfac) { + + real_type sfac_outer, dsfac_outer, sfac_inner, dsfac_inner; constexpr real_type one = static_cast(1.0); constexpr real_type zero = static_cast(0.0); constexpr real_type onehalf = static_cast(0.5); - if (switch_flag == 0) { sfac = zero; dsfac = zero; } + + if (switch_flag == 0) { sfac_outer = zero; dsfac_outer = zero; } else if (switch_flag == 1) { - if (r <= rmin0) { sfac = one; dsfac = zero; } - else if (r > rcut) { sfac = zero; dsfac = zero; } + if (r <= rmin0) { sfac_outer = one; dsfac_outer = zero; } + else if (r > rcut) { sfac = zero; dsfac = zero; return; } else { const real_type rcutfac = static_cast(MY_PI) / (rcut - rmin0); const real_type theta0 = (r - rmin0) * rcutfac; const real_type sn = sin(theta0); const real_type cs = cos(theta0); - sfac = onehalf * (cs + one); - dsfac = -onehalf * sn * rcutfac; - + sfac_outer = onehalf * (cs + one); + dsfac_outer = -onehalf * sn * rcutfac; } - } else { sfac = zero; dsfac = zero; } + } else { sfac = zero; dsfac = zero; return; } // dummy return + + if (switch_inner_flag == 0) { sfac = sfac_outer; dsfac = dsfac_outer; return; } + else if (switch_inner_flag == 1) { + if (r >= sinner + dinner) { sfac = sfac_outer; dsfac = dsfac_outer; return; } + else if (r > sinner - dinner) { + real_type rcutfac = static_cast(MY_PI2) / dinner; + sfac_inner = onehalf * (one - cos(static_cast(MY_PI2) + (r - sinner) * rcutfac)); + dsfac_inner = onehalf * rcutfac * sin(static_cast(MY_PI2) + (r - sinner) * rcutfac); + sfac = sfac_outer * sfac_inner; + dsfac = dsfac_outer * sfac_inner + sfac_outer * dsfac_inner; + return; + } else { sfac = zero; dsfac = zero; return; } + } else { sfac = zero; dsfac = zero; return; } // dummy return + } /* ---------------------------------------------------------------------- diff --git a/src/ML-SNAP/sna.cpp b/src/ML-SNAP/sna.cpp index 0bd66d0bce..2b86774d5b 100644 --- a/src/ML-SNAP/sna.cpp +++ b/src/ML-SNAP/sna.cpp @@ -1533,6 +1533,9 @@ void SNA::compute_ncoeff() double SNA::compute_sfac(double r, double rcut, double sinner, double dinner) { double sfac; + + // calculate sfac = sfac_outer + if (switch_flag == 0) sfac = 1.0; else if (r <= rmin0) sfac = 1.0; else if (r > rcut) sfac = 0.0; @@ -1541,6 +1544,8 @@ double SNA::compute_sfac(double r, double rcut, double sinner, double dinner) sfac = 0.5 * (cos((r - rmin0) * rcutfac) + 1.0); } + // calculate sfac *= sfac_inner, rarely visited + if (switch_inner_flag == 1 && r < sinner + dinner) { if (r > sinner - dinner) { double rcutfac = MY_PI2 / dinner; @@ -1555,33 +1560,38 @@ double SNA::compute_sfac(double r, double rcut, double sinner, double dinner) double SNA::compute_dsfac(double r, double rcut, double sinner, double dinner) { - double sfac, dsfac, sfac_inner, dsfac_inner; - if (switch_flag == 0) dsfac = 0.0; - else if (r <= rmin0) dsfac = 0.0; - else if (r > rcut) dsfac = 0.0; + double dsfac, sfac_outer, dsfac_outer, sfac_inner, dsfac_inner; + if (switch_flag == 0) dsfac_outer = 0.0; + else if (r <= rmin0) dsfac_outer = 0.0; + else if (r > rcut) dsfac_outer = 0.0; else { double rcutfac = MY_PI / (rcut - rmin0); - dsfac = -0.5 * sin((r - rmin0) * rcutfac) * rcutfac; + dsfac_outer = -0.5 * sin((r - rmin0) * rcutfac) * rcutfac; } - // duplicated computation, but rarely visited + // some duplicated computation, but rarely visited if (switch_inner_flag == 1 && r < sinner + dinner) { if (r > sinner - dinner) { - if (switch_flag == 0) sfac = 1.0; - else if (r <= rmin0) sfac = 1.0; - else if (r > rcut) sfac = 0.0; + + // calculate sfac_outer + + if (switch_flag == 0) sfac_outer = 1.0; + else if (r <= rmin0) sfac_outer = 1.0; + else if (r > rcut) sfac_outer = 0.0; else { double rcutfac = MY_PI / (rcut - rmin0); - sfac = 0.5 * (cos((r - rmin0) * rcutfac) + 1.0); + sfac_outer = 0.5 * (cos((r - rmin0) * rcutfac) + 1.0); } + // calculate sfac_inner + double rcutfac = MY_PI2 / dinner; sfac_inner = 0.5 * (1.0 - cos(MY_PI2 + (r - sinner) * rcutfac)); dsfac_inner = 0.5 * rcutfac * sin(MY_PI2 + (r - sinner) * rcutfac); - dsfac = dsfac*sfac_inner + sfac*dsfac_inner; + dsfac = dsfac_outer*sfac_inner + sfac_outer*dsfac_inner; } else dsfac = 0.0; - } + } else dsfac = dsfac_outer; return dsfac; } From dc605d35eed20c2e2554680681da1bbe57cbd60d Mon Sep 17 00:00:00 2001 From: Aidan Thompson Date: Mon, 25 Apr 2022 13:58:41 -0600 Subject: [PATCH 05/16] Updated ML-IAP package --- src/ML-IAP/mliap_descriptor_snap.cpp | 53 +++++++++++++--------------- src/ML-IAP/mliap_descriptor_snap.h | 4 +-- 2 files changed, 27 insertions(+), 30 deletions(-) diff --git a/src/ML-IAP/mliap_descriptor_snap.cpp b/src/ML-IAP/mliap_descriptor_snap.cpp index 7badf7b641..195eaa07fa 100644 --- a/src/ML-IAP/mliap_descriptor_snap.cpp +++ b/src/ML-IAP/mliap_descriptor_snap.cpp @@ -41,8 +41,8 @@ MLIAPDescriptorSNAP::MLIAPDescriptorSNAP(LAMMPS *_lmp, char *paramfilename) : ML radelem = nullptr; wjelem = nullptr; snaptr = nullptr; - rinnerelem = nullptr; - drinnerelem = nullptr; + sinnerelem = nullptr; + dinnerelem = nullptr; read_paramfile(paramfilename); @@ -59,11 +59,8 @@ MLIAPDescriptorSNAP::~MLIAPDescriptorSNAP() memory->destroy(radelem); memory->destroy(wjelem); delete snaptr; - - if (switchinnerflag) { - memory->destroy(rinnerelem); - memory->destroy(drinnerelem); - } + memory->destroy(sinnerelem); + memory->destroy(dinnerelem); } /* ---------------------------------------------------------------------- @@ -94,8 +91,8 @@ void MLIAPDescriptorSNAP::compute_descriptors(class MLIAPData *data) snaptr->wj[ninside] = wjelem[jelem]; snaptr->rcutij[ninside] = sqrt(cutsq[ielem][jelem]); if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5 * (rinnerelem[ielem] + rinnerelem[jelem]); - snaptr->drinnerij[ninside] = 0.5 * (drinnerelem[ielem] + drinnerelem[jelem]); + snaptr->sinnerij[ninside] = 0.5 * (sinnerelem[ielem] + sinnerelem[jelem]); + snaptr->dinnerij[ninside] = 0.5 * (dinnerelem[ielem] + dinnerelem[jelem]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; @@ -150,8 +147,8 @@ void MLIAPDescriptorSNAP::compute_forces(class MLIAPData *data) snaptr->wj[ninside] = wjelem[jelem]; snaptr->rcutij[ninside] = sqrt(cutsq[ielem][jelem]); if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5 * (rinnerelem[ielem] + rinnerelem[jelem]); - snaptr->drinnerij[ninside] = 0.5 * (drinnerelem[ielem] + drinnerelem[jelem]); + snaptr->sinnerij[ninside] = 0.5 * (sinnerelem[ielem] + sinnerelem[jelem]); + snaptr->dinnerij[ninside] = 0.5 * (dinnerelem[ielem] + dinnerelem[jelem]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; @@ -222,8 +219,8 @@ void MLIAPDescriptorSNAP::compute_force_gradients(class MLIAPData *data) snaptr->wj[ninside] = wjelem[jelem]; snaptr->rcutij[ninside] = sqrt(cutsq[ielem][jelem]); if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5 * (rinnerelem[ielem] + rinnerelem[jelem]); - snaptr->drinnerij[ninside] = 0.5 * (drinnerelem[ielem] + drinnerelem[jelem]); + snaptr->sinnerij[ninside] = 0.5 * (sinnerelem[ielem] + sinnerelem[jelem]); + snaptr->dinnerij[ninside] = 0.5 * (dinnerelem[ielem] + dinnerelem[jelem]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; @@ -292,8 +289,8 @@ void MLIAPDescriptorSNAP::compute_descriptor_gradients(class MLIAPData *data) snaptr->wj[ninside] = wjelem[jelem]; snaptr->rcutij[ninside] = sqrt(cutsq[ielem][jelem]); if (switchinnerflag) { - snaptr->rinnerij[ninside] = 0.5 * (rinnerelem[ielem] + rinnerelem[jelem]); - snaptr->drinnerij[ninside] = 0.5 * (drinnerelem[ielem] + drinnerelem[jelem]); + snaptr->sinnerij[ninside] = 0.5 * (sinnerelem[ielem] + sinnerelem[jelem]); + snaptr->dinnerij[ninside] = 0.5 * (dinnerelem[ielem] + dinnerelem[jelem]); } if (chemflag) snaptr->element[ninside] = jelem; ninside++; @@ -364,8 +361,8 @@ void MLIAPDescriptorSNAP::read_paramfile(char *paramfilename) // set local input checks - int rinnerflag = 0; - int drinnerflag = 0; + int sinnerflag = 0; + int dinnerflag = 0; for (int i = 0; i < nelements; i++) delete[] elements[i]; delete[] elements; @@ -414,7 +411,7 @@ void MLIAPDescriptorSNAP::read_paramfile(char *paramfilename) // check for keywords with one value per element if ((keywd == "elems") || (keywd == "radelems") || (keywd == "welems") || - (keywd == "rinnerelems") || (keywd == "drinnerelems")) { + (keywd == "sinnerelems") || (keywd == "dinnerelems")) { if ((nelementsflag == 0) || ((int)words.count() != nelements + 1)) error->all(FLERR, "Incorrect SNAP parameter file"); @@ -433,14 +430,14 @@ void MLIAPDescriptorSNAP::read_paramfile(char *paramfilename) for (int ielem = 0; ielem < nelements; ielem++) wjelem[ielem] = utils::numeric(FLERR, words.next(), false, lmp); wjelemflag = 1; - } else if (keywd == "rinnerelems") { + } else if (keywd == "sinnerelems") { for (int ielem = 0; ielem < nelements; ielem++) - rinnerelem[ielem] = utils::numeric(FLERR, words.next(), false, lmp); - rinnerflag = 1; - } else if (keywd == "drinnerelems") { + sinnerelem[ielem] = utils::numeric(FLERR, words.next(), false, lmp); + sinnerflag = 1; + } else if (keywd == "dinnerelems") { for (int ielem = 0; ielem < nelements; ielem++) - drinnerelem[ielem] = utils::numeric(FLERR, words.next(), false, lmp); - rinnerflag = 1; + dinnerelem[ielem] = utils::numeric(FLERR, words.next(), false, lmp); + dinnerflag = 1; } } else { @@ -457,8 +454,8 @@ void MLIAPDescriptorSNAP::read_paramfile(char *paramfilename) elements = new char *[nelements]; memory->create(radelem, nelements, "mliap_snap_descriptor:radelem"); memory->create(wjelem, nelements, "mliap_snap_descriptor:wjelem"); - memory->create(rinnerelem, nelements, "mliap_snap_descriptor:rinner"); - memory->create(drinnerelem, nelements, "mliap_snap_descriptor:drinner"); + memory->create(sinnerelem, nelements, "mliap_snap_descriptor:sinner"); + memory->create(dinnerelem, nelements, "mliap_snap_descriptor:dinner"); nelementsflag = 1; } else if (keywd == "rcutfac") { rcutfac = utils::numeric(FLERR, keyval, false, lmp); @@ -491,10 +488,10 @@ void MLIAPDescriptorSNAP::read_paramfile(char *paramfilename) !wjelemflag) error->all(FLERR, "Incorrect SNAP parameter file"); - if (switchinnerflag && !(rinnerflag && drinnerflag)) + if (switchinnerflag && !(sinnerflag && dinnerflag)) error->all(FLERR, "Incorrect SNAP parameter file"); - if (!switchinnerflag && (rinnerflag || drinnerflag)) + if (!switchinnerflag && (sinnerflag || dinnerflag)) error->all(FLERR, "Incorrect SNAP parameter file"); // construct cutsq diff --git a/src/ML-IAP/mliap_descriptor_snap.h b/src/ML-IAP/mliap_descriptor_snap.h index 2136eefc5e..b49ff6d2cc 100644 --- a/src/ML-IAP/mliap_descriptor_snap.h +++ b/src/ML-IAP/mliap_descriptor_snap.h @@ -42,8 +42,8 @@ class MLIAPDescriptorSNAP : public MLIAPDescriptor { int switchinnerflag; double rfac0, rmin0; - double* rinnerelem; - double* drinnerelem; + double* sinnerelem; + double* dinnerelem; }; } // namespace LAMMPS_NS From 26278643d4660fa35f712705378be384914913fb Mon Sep 17 00:00:00 2001 From: Stan Moore Date: Tue, 26 Apr 2022 08:22:00 -0600 Subject: [PATCH 06/16] Port changes to Kokkos --- src/KOKKOS/pair_snap_kokkos.h | 2 ++ src/KOKKOS/pair_snap_kokkos_impl.h | 32 ++++++++++++++++++++++-------- 2 files changed, 26 insertions(+), 8 deletions(-) diff --git a/src/KOKKOS/pair_snap_kokkos.h b/src/KOKKOS/pair_snap_kokkos.h index bdad113c18..1c0bf38699 100644 --- a/src/KOKKOS/pair_snap_kokkos.h +++ b/src/KOKKOS/pair_snap_kokkos.h @@ -280,6 +280,8 @@ inline double dist2(double* x,double* y); Kokkos::View d_radelem; // element radii Kokkos::View d_wjelem; // elements weights Kokkos::View d_coeffelem; // element bispectrum coefficients + Kokkos::View d_sinnerelem; // element inner cutoff midpoint + Kokkos::View d_dinnerelem; // element inner cutoff half-width Kokkos::View d_map; // mapping from atom types to elements Kokkos::View d_ninside; // ninside for all atoms in list Kokkos::View d_beta; // betas for all atoms in list diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h index 6bc034a163..7c621176cd 100644 --- a/src/KOKKOS/pair_snap_kokkos_impl.h +++ b/src/KOKKOS/pair_snap_kokkos_impl.h @@ -578,15 +578,21 @@ void PairSNAPKokkos::coeff(int narg, char d_radelem = Kokkos::View("pair:radelem",nelements); d_wjelem = Kokkos::View("pair:wjelem",nelements); d_coeffelem = Kokkos::View("pair:coeffelem",nelements,ncoeffall); + d_sinnerelem = Kokkos::View("pair:sinnerelem",nelements); + d_dinnerelem = Kokkos::View("pair:dinnerelem",nelements); auto h_radelem = Kokkos::create_mirror_view(d_radelem); auto h_wjelem = Kokkos::create_mirror_view(d_wjelem); auto h_coeffelem = Kokkos::create_mirror_view(d_coeffelem); + auto h_sinnerelem = Kokkos::create_mirror_view(d_sinnerelem); + auto h_dinnerelem = Kokkos::create_mirror_view(d_dinnerelem); auto h_map = Kokkos::create_mirror_view(d_map); for (int ielem = 0; ielem < nelements; ielem++) { h_radelem(ielem) = radelem[ielem]; h_wjelem(ielem) = wjelem[ielem]; + h_sinnerelem(ielem) = sinnerelem[ielem]; + h_dinnerelem(ielem) = dinnerelem[ielem]; for (int jcoeff = 0; jcoeff < ncoeffall; jcoeff++) { h_coeffelem(ielem,jcoeff) = coeffelem[ielem][jcoeff]; } @@ -599,6 +605,8 @@ void PairSNAPKokkos::coeff(int narg, char Kokkos::deep_copy(d_radelem,h_radelem); Kokkos::deep_copy(d_wjelem,h_wjelem); Kokkos::deep_copy(d_coeffelem,h_coeffelem); + Kokkos::deep_copy(d_sinnerelem,h_sinnerelem); + Kokkos::deep_copy(d_dinnerelem,h_dinnerelem); Kokkos::deep_copy(d_map,h_map); snaKK = SNAKokkos(rfac0,twojmax, @@ -724,15 +732,19 @@ void PairSNAPKokkos::operator() (TagPairSN const F_FLOAT dx = x(j,0) - xtmp; const F_FLOAT dy = x(j,1) - ytmp; const F_FLOAT dz = x(j,2) - ztmp; - const int elem_j = d_map[jtype]; + const int jelem = d_map[jtype]; my_sna.rij(ii,offset,0) = static_cast(dx); my_sna.rij(ii,offset,1) = static_cast(dy); my_sna.rij(ii,offset,2) = static_cast(dz); - my_sna.wj(ii,offset) = static_cast(d_wjelem[elem_j]); - my_sna.rcutij(ii,offset) = static_cast((radi + d_radelem[elem_j])*rcutfac); + my_sna.wj(ii,offset) = static_cast(d_wjelem[jelem]); + my_sna.rcutij(ii,offset) = static_cast((radi + d_radelem[jelem])*rcutfac); my_sna.inside(ii,offset) = j; + if (switchinnerflag) { + my_sna.sinnerij(ii,ninside) = 0.5*(d_sinnerelem[ielem] + d_sinnerelem[jelem]); + my_sna.dinnerij(ii,ninside) = 0.5*(d_dinnerelem[ielem] + d_dinnerelem[jelem]); + } if (chemflag) - my_sna.element(ii,offset) = elem_j; + my_sna.element(ii,offset) = jelem; else my_sna.element(ii,offset) = 0; } @@ -1080,18 +1092,22 @@ void PairSNAPKokkos::operator() (TagPairSN const int jtype = type(j); const F_FLOAT rsq = dx*dx + dy*dy + dz*dz; - const int elem_j = d_map[jtype]; + const int jelem = d_map[jtype]; if (rsq < rnd_cutsq(itype,jtype)) { if (final) { my_sna.rij(ii,offset,0) = static_cast(dx); my_sna.rij(ii,offset,1) = static_cast(dy); my_sna.rij(ii,offset,2) = static_cast(dz); - my_sna.wj(ii,offset) = static_cast(d_wjelem[elem_j]); - my_sna.rcutij(ii,offset) = static_cast((radi + d_radelem[elem_j])*rcutfac); + my_sna.wj(ii,offset) = static_cast(d_wjelem[jelem]); + my_sna.rcutij(ii,offset) = static_cast((radi + d_radelem[jelem])*rcutfac); my_sna.inside(ii,offset) = j; + if (switchinnerflag) { + my_sna.sinnerij(ii,ninside) = 0.5*(d_sinnerelem[ielem] + d_sinnerelem[jelem]); + my_sna.dinnerij(ii,ninside) = 0.5*(d_dinnerelem[ielem] + d_dinnerelem[jelem]); + } if (chemflag) - my_sna.element(ii,offset) = elem_j; + my_sna.element(ii,offset) = jelem; else my_sna.element(ii,offset) = 0; } From a3c6baad4ca4ec173dedafbaf87f15d0abcf6b61 Mon Sep 17 00:00:00 2001 From: Aidan Thompson Date: Wed, 27 Apr 2022 10:58:22 -0600 Subject: [PATCH 07/16] Removed CPU specific options --- src/KOKKOS/pair_snap_kokkos_impl.h | 31 ++++++++++++++++-------------- 1 file changed, 17 insertions(+), 14 deletions(-) diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h index 7c621176cd..4d2edf037c 100644 --- a/src/KOKKOS/pair_snap_kokkos_impl.h +++ b/src/KOKKOS/pair_snap_kokkos_impl.h @@ -89,15 +89,17 @@ PairSNAPKokkos::~PairSNAPKokkos() template void PairSNAPKokkos::init_style() { - if (host_flag) { - if (lmp->kokkos->nthreads > 1) - if (comm->me == 0) - utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " - "CPU thread, even if more threads are requested\n"); - PairSNAP::init_style(); - return; - } + // Disabled for DEBUGGING on CPU + // if (host_flag) { + // if (lmp->kokkos->nthreads > 1) + // if (comm->me == 0) + // utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " + // "CPU thread, even if more threads are requested\n"); + + // PairSNAP::init_style(); + // return; + // } if (force->newton_pair == 0) error->all(FLERR,"Pair style SNAP requires newton pair on"); @@ -139,12 +141,13 @@ struct FindMaxNumNeighs { template void PairSNAPKokkos::compute(int eflag_in, int vflag_in) { - if (host_flag) { - atomKK->sync(Host,X_MASK|TYPE_MASK); - PairSNAP::compute(eflag_in,vflag_in); - atomKK->modified(Host,F_MASK); - return; - } + // Disabled for DEBUGGING on CPU + // if (host_flag) { + // atomKK->sync(Host,X_MASK|TYPE_MASK); + // PairSNAP::compute(eflag_in,vflag_in); + // atomKK->modified(Host,F_MASK); + // return; + // } eflag = eflag_in; vflag = vflag_in; From a3e0e1a6ebc3956364c72c23e793e58e9f843ab3 Mon Sep 17 00:00:00 2001 From: megmcca <90424440+megmcca@users.noreply.github.com> Date: Thu, 28 Apr 2022 15:57:53 -0600 Subject: [PATCH 08/16] fixed multielement assignment error --- src/ML-SNAP/pair_snap.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/ML-SNAP/pair_snap.cpp b/src/ML-SNAP/pair_snap.cpp index 6dcdaa67e9..089f594fff 100644 --- a/src/ML-SNAP/pair_snap.cpp +++ b/src/ML-SNAP/pair_snap.cpp @@ -701,15 +701,15 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) int iword = 1; if (keywd == "sinner") { - keyval = words[iword]; for (int ielem = 0; ielem < nelements; ielem++) { + keyval = words[iword]; sinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); iword++; } sinnerflag = 1; } else if (keywd == "dinner") { - keyval = words[iword]; for (int ielem = 0; ielem < nelements; ielem++) { + keyval = words[iword]; dinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); iword++; } From 2292caef3898ee5cfc0cb8a774cd8ec4b23fac5b Mon Sep 17 00:00:00 2001 From: megmcca <90424440+megmcca@users.noreply.github.com> Date: Fri, 29 Apr 2022 07:48:02 -0600 Subject: [PATCH 09/16] updated sinner/dinner log output now logmsg displays all element params for sinner, dinner on same line (replaces only first param + ...) --- src/ML-SNAP/pair_snap.cpp | 14 ++++++++++---- 1 file changed, 10 insertions(+), 4 deletions(-) diff --git a/src/ML-SNAP/pair_snap.cpp b/src/ML-SNAP/pair_snap.cpp index 089f594fff..56c3250c3f 100644 --- a/src/ML-SNAP/pair_snap.cpp +++ b/src/ML-SNAP/pair_snap.cpp @@ -688,16 +688,17 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) auto keywd = words[0]; auto keyval = words[1]; - // check for keywords with one value per element + // check for keywords with more than one value per element if (keywd == "sinner" || keywd == "dinner") { if ((int)words.size() != nelements+1) error->all(FLERR,"Incorrect SNAP parameter file"); + + // innerlogstr collects all values of sinner or dinner for log output below - if (comm->me == 0) - utils::logmesg(lmp,"SNAP keyword {} {} ... \n", keywd, keyval); - + std::string innerlogstr; + int iword = 1; if (keywd == "sinner") { @@ -705,6 +706,7 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) keyval = words[iword]; sinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); iword++; + innerlogstr += keyval + " "; } sinnerflag = 1; } else if (keywd == "dinner") { @@ -712,9 +714,13 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) keyval = words[iword]; dinnerelem[ielem] = utils::numeric(FLERR,keyval,false,lmp); iword++; + innerlogstr += keyval + " "; } dinnerflag = 1; } + + if (comm->me == 0) + utils::logmesg(lmp,"SNAP keyword {} {} ... \n", keywd, innerlogstr); } else { From c34efcc718111a91af12e35af72b21635d81830f Mon Sep 17 00:00:00 2001 From: megmcca <90424440+megmcca@users.noreply.github.com> Date: Fri, 29 Apr 2022 07:53:39 -0600 Subject: [PATCH 10/16] fixed multielement assignment error tested in CPU and GPU mode, gets same output for eneriges and forces as build without KOKKOS for various settings of sinner, dinner --- src/KOKKOS/pair_snap_kokkos_impl.h | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h index 4d2edf037c..65684d58ba 100644 --- a/src/KOKKOS/pair_snap_kokkos_impl.h +++ b/src/KOKKOS/pair_snap_kokkos_impl.h @@ -743,8 +743,8 @@ void PairSNAPKokkos::operator() (TagPairSN my_sna.rcutij(ii,offset) = static_cast((radi + d_radelem[jelem])*rcutfac); my_sna.inside(ii,offset) = j; if (switchinnerflag) { - my_sna.sinnerij(ii,ninside) = 0.5*(d_sinnerelem[ielem] + d_sinnerelem[jelem]); - my_sna.dinnerij(ii,ninside) = 0.5*(d_dinnerelem[ielem] + d_dinnerelem[jelem]); + my_sna.sinnerij(ii,offset) = 0.5*(d_sinnerelem[ielem] + d_sinnerelem[jelem]); + my_sna.dinnerij(ii,offset) = 0.5*(d_dinnerelem[ielem] + d_dinnerelem[jelem]); } if (chemflag) my_sna.element(ii,offset) = jelem; @@ -1106,8 +1106,8 @@ void PairSNAPKokkos::operator() (TagPairSN my_sna.rcutij(ii,offset) = static_cast((radi + d_radelem[jelem])*rcutfac); my_sna.inside(ii,offset) = j; if (switchinnerflag) { - my_sna.sinnerij(ii,ninside) = 0.5*(d_sinnerelem[ielem] + d_sinnerelem[jelem]); - my_sna.dinnerij(ii,ninside) = 0.5*(d_dinnerelem[ielem] + d_dinnerelem[jelem]); + my_sna.sinnerij(ii,offset) = 0.5*(d_sinnerelem[ielem] + d_sinnerelem[jelem]); + my_sna.dinnerij(ii,offset) = 0.5*(d_dinnerelem[ielem] + d_dinnerelem[jelem]); } if (chemflag) my_sna.element(ii,offset) = jelem; From 641e769bdea8a5494320462b6d55ca3d088e9fd4 Mon Sep 17 00:00:00 2001 From: Stan Moore Date: Fri, 29 Apr 2022 09:03:44 -0600 Subject: [PATCH 11/16] Revert a3c6baad4ca4ec173dedafbaf87f15d0abcf6b61 --- src/KOKKOS/pair_snap_kokkos_impl.h | 30 ++++++++++++++---------------- 1 file changed, 14 insertions(+), 16 deletions(-) diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h index 65684d58ba..2b49345371 100644 --- a/src/KOKKOS/pair_snap_kokkos_impl.h +++ b/src/KOKKOS/pair_snap_kokkos_impl.h @@ -90,16 +90,15 @@ template void PairSNAPKokkos::init_style() { - // Disabled for DEBUGGING on CPU - // if (host_flag) { - // if (lmp->kokkos->nthreads > 1) - // if (comm->me == 0) - // utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " - // "CPU thread, even if more threads are requested\n"); + if (host_flag) { + if (lmp->kokkos->nthreads > 1) + if (comm->me == 0) + utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " + "CPU thread, even if more threads are requested\n"); - // PairSNAP::init_style(); - // return; - // } + PairSNAP::init_style(); + return; + } if (force->newton_pair == 0) error->all(FLERR,"Pair style SNAP requires newton pair on"); @@ -141,13 +140,12 @@ struct FindMaxNumNeighs { template void PairSNAPKokkos::compute(int eflag_in, int vflag_in) { - // Disabled for DEBUGGING on CPU - // if (host_flag) { - // atomKK->sync(Host,X_MASK|TYPE_MASK); - // PairSNAP::compute(eflag_in,vflag_in); - // atomKK->modified(Host,F_MASK); - // return; - // } + if (host_flag) { + atomKK->sync(Host,X_MASK|TYPE_MASK); + PairSNAP::compute(eflag_in,vflag_in); + atomKK->modified(Host,F_MASK); + return; + } eflag = eflag_in; vflag = vflag_in; From 69a79ba082350f3efc7f440c051ee329bf1ef0c7 Mon Sep 17 00:00:00 2001 From: Stan Moore Date: Fri, 29 Apr 2022 09:17:24 -0600 Subject: [PATCH 12/16] whitespace --- src/KOKKOS/pair_snap_kokkos_impl.h | 29 ++++++++++++++--------------- 1 file changed, 14 insertions(+), 15 deletions(-) diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h index 2b49345371..e39f42058f 100644 --- a/src/KOKKOS/pair_snap_kokkos_impl.h +++ b/src/KOKKOS/pair_snap_kokkos_impl.h @@ -89,16 +89,15 @@ PairSNAPKokkos::~PairSNAPKokkos() template void PairSNAPKokkos::init_style() { + if (host_flag) { + if (lmp->kokkos->nthreads > 1) + if (comm->me == 0) + utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " + "CPU thread, even if more threads are requested\n"); - if (host_flag) { - if (lmp->kokkos->nthreads > 1) - if (comm->me == 0) - utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " - "CPU thread, even if more threads are requested\n"); - - PairSNAP::init_style(); - return; - } + PairSNAP::init_style(); + return; + } if (force->newton_pair == 0) error->all(FLERR,"Pair style SNAP requires newton pair on"); @@ -140,12 +139,12 @@ struct FindMaxNumNeighs { template void PairSNAPKokkos::compute(int eflag_in, int vflag_in) { - if (host_flag) { - atomKK->sync(Host,X_MASK|TYPE_MASK); - PairSNAP::compute(eflag_in,vflag_in); - atomKK->modified(Host,F_MASK); - return; - } + if (host_flag) { + atomKK->sync(Host,X_MASK|TYPE_MASK); + PairSNAP::compute(eflag_in,vflag_in); + atomKK->modified(Host,F_MASK); + return; + } eflag = eflag_in; vflag = vflag_in; From 57098f3df7e94ddeaaed1b4ab5628fcb3df66b6d Mon Sep 17 00:00:00 2001 From: Stan Moore Date: Fri, 29 Apr 2022 09:22:14 -0600 Subject: [PATCH 13/16] Remove extra files --- src/MAKE/OPTIONS/Makefile.kokkos_omp_db_asan | 126 ------------------- src/Si.tersoff | 18 --- src/Si.tersoff.mod | 27 ---- src/Si.tersoff.modc | 10 -- 4 files changed, 181 deletions(-) delete mode 100644 src/MAKE/OPTIONS/Makefile.kokkos_omp_db_asan delete mode 100644 src/Si.tersoff delete mode 100644 src/Si.tersoff.mod delete mode 100644 src/Si.tersoff.modc diff --git a/src/MAKE/OPTIONS/Makefile.kokkos_omp_db_asan b/src/MAKE/OPTIONS/Makefile.kokkos_omp_db_asan deleted file mode 100644 index 72427b9e5f..0000000000 --- a/src/MAKE/OPTIONS/Makefile.kokkos_omp_db_asan +++ /dev/null @@ -1,126 +0,0 @@ -# kokkos_omp = KOKKOS/OMP package, MPI with its default compiler - -SHELL = /bin/sh - -# --------------------------------------------------------------------- -# compiler/linker settings -# specify flags and libraries needed for your compiler - -CC = mpicxx -CCFLAGS = -g -O0 -fsanitize=address -SHFLAGS = -fPIC -DEPFLAGS = -M - -LINK = mpicxx -LINKFLAGS = -g -O0 -fsanitize=address -LIB = -SIZE = size - -ARCHIVE = ar -ARFLAGS = -rc -SHLIBFLAGS = -shared -KOKKOS_DEVICES = OpenMP -KOKKOS_DEBUG = yes - -# --------------------------------------------------------------------- -# LAMMPS-specific settings, all OPTIONAL -# specify settings for LAMMPS features you will use -# if you change any -D setting, do full re-compile after "make clean" - -# LAMMPS ifdef settings -# see possible settings in Section 3.5 of the manual - -LMP_INC = -DLAMMPS_GZIP - -# MPI library -# see discussion in Section 3.4 of the manual -# MPI wrapper compiler/linker can provide this info -# can point to dummy MPI library in src/STUBS as in Makefile.serial -# use -D MPICH and OMPI settings in INC to avoid C++ lib conflicts -# INC = path for mpi.h, MPI compiler settings -# PATH = path for MPI library -# LIB = name of MPI library - -MPI_INC = -DMPICH_SKIP_MPICXX -DOMPI_SKIP_MPICXX=1 -MPI_PATH = -MPI_LIB = - -# FFT library -# see discussion in Section 3.5.2 of manual -# can be left blank to use provided KISS FFT library -# INC = -DFFT setting, e.g. -DFFT_FFTW, FFT compiler settings -# PATH = path for FFT library -# LIB = name of FFT library - -FFT_INC = -FFT_PATH = -FFT_LIB = - -# JPEG and/or PNG library -# see discussion in Section 3.5.4 of manual -# only needed if -DLAMMPS_JPEG or -DLAMMPS_PNG listed with LMP_INC -# INC = path(s) for jpeglib.h and/or png.h -# PATH = path(s) for JPEG library and/or PNG library -# LIB = name(s) of JPEG library and/or PNG library - -JPG_INC = -JPG_PATH = -JPG_LIB = - -# library for loading shared objects (defaults to -ldl, should be empty on Windows) -# uncomment to change the default - -# override DYN_LIB = - -# --------------------------------------------------------------------- -# build rules and dependencies -# do not edit this section - -include Makefile.package.settings -include Makefile.package - -EXTRA_INC = $(LMP_INC) $(PKG_INC) $(MPI_INC) $(FFT_INC) $(JPG_INC) $(PKG_SYSINC) -EXTRA_PATH = $(PKG_PATH) $(MPI_PATH) $(FFT_PATH) $(JPG_PATH) $(PKG_SYSPATH) -EXTRA_LIB = $(PKG_LIB) $(MPI_LIB) $(FFT_LIB) $(JPG_LIB) $(PKG_SYSLIB) $(DYN_LIB) -EXTRA_CPP_DEPENDS = $(PKG_CPP_DEPENDS) -EXTRA_LINK_DEPENDS = $(PKG_LINK_DEPENDS) - -# Path to src files - -vpath %.cpp .. -vpath %.h .. - -# Link target - -$(EXE): main.o $(LMPLIB) $(EXTRA_LINK_DEPENDS) - $(LINK) $(LINKFLAGS) main.o $(EXTRA_PATH) $(LMPLINK) $(EXTRA_LIB) $(LIB) -o $@ - $(SIZE) $@ - -# Library targets - -$(ARLIB): $(OBJ) $(EXTRA_LINK_DEPENDS) - @rm -f ../$(ARLIB) - $(ARCHIVE) $(ARFLAGS) ../$(ARLIB) $(OBJ) - @rm -f $(ARLIB) - @ln -s ../$(ARLIB) $(ARLIB) - -$(SHLIB): $(OBJ) $(EXTRA_LINK_DEPENDS) - $(CC) $(CCFLAGS) $(SHFLAGS) $(SHLIBFLAGS) $(EXTRA_PATH) -o ../$(SHLIB) \ - $(OBJ) $(EXTRA_LIB) $(LIB) - @rm -f $(SHLIB) - @ln -s ../$(SHLIB) $(SHLIB) - -# Compilation rules - -%.o:%.cpp - $(CC) $(CCFLAGS) $(SHFLAGS) $(EXTRA_INC) -c $< - -# Individual dependencies - -depend : fastdep.exe $(SRC) - @./fastdep.exe $(EXTRA_INC) -- $^ > .depend || exit 1 - -fastdep.exe: ../DEPEND/fastdep.c - cc -O -o $@ $< - -sinclude .depend diff --git a/src/Si.tersoff b/src/Si.tersoff deleted file mode 100644 index fa565df42b..0000000000 --- a/src/Si.tersoff +++ /dev/null @@ -1,18 +0,0 @@ -# DATE: 2007-10-25 UNITS: metal CONTRIBUTOR: Aidan Thompson, athomps@sandia.gov CITATION: Tersoff, Phys Rev B, 37, 6991 (1988) - -# Tersoff parameters for various elements and mixtures -# multiple entries can be added to this file, LAMMPS reads the ones it needs -# these entries are in LAMMPS "metal" units: -# A,B = eV; lambda1,lambda2,lambda3 = 1/Angstroms; R,D = Angstroms -# other quantities are unitless - -# This is the Si parameterization from a particular Tersoff paper: -# J. Tersoff, PRB, 37, 6991 (1988) -# See the SiCGe.tersoff file for different Si variants. - -# format of a single entry (one or more lines): -# element 1, element 2, element 3, -# m, gamma, lambda3, c, d, costheta0, n, beta, lambda2, B, R, D, lambda1, A - -Si Si Si 3.0 1.0 1.3258 4.8381 2.0417 0.0000 22.956 - 0.33675 1.3258 95.373 3.0 0.2 3.2394 3264.7 diff --git a/src/Si.tersoff.mod b/src/Si.tersoff.mod deleted file mode 100644 index 016c247638..0000000000 --- a/src/Si.tersoff.mod +++ /dev/null @@ -1,27 +0,0 @@ -# DATE: 2013-07-26 UNITS: metal CONTRIBUTOR: Unknown CITATION: Kumagai, Izumi, Hara and Sakai, Comp Mat Sci, 39, 457 (2007) - -# Modified Tersoff potential (named MOD potential) parameters for various elements and mixtures -# multiple entries can be added to this file, LAMMPS reads the ones it needs -# these entries are in LAMMPS "metal" units: -# A,B = eV; lambda1,lambda2 = 1/Angstroms; R,D = Angstroms -# other quantities are unitless - -# MOD. This is the Si parameterization from the paper: -# [1] T. Kumagai, S. Izumi, S. Hara, S. Sakai, Comp. Mat. Sci., 39, 457 (2007) - -# format of a single entry (one or more lines): -# element 1, element 2, element 3, -# beta, alpha, h, eta, -# beta_ters, lambda2, B, R, D, lambda1, A, n -# c1, c2, c3, c4, c5 - -# Notes -# 1) beta_ters must be equal to 1.0 -# 2) R = (R1+R2)/2, where R1 and R2 determined at [1] (R1 = 2.7A, R2 = 3.3A) -# 3) D = (R2-R1)/2, where R1 and R2 determined at [1] (R1 = 2.7A, R2 = 3.3A) -# 4) n = 1.0/(2*delta), where delta determined at [1] (eta*delta = 0.53298909) - -#mod -Si Si Si 1.0 2.3890327 -.365 1.0 - 1.0 1.345797 121.00047 3.0 0.3 3.2300135 3281.5905 .93810551 - 0.20173476 730418.72 1000000.0 1.0 26.0 diff --git a/src/Si.tersoff.modc b/src/Si.tersoff.modc deleted file mode 100644 index ce184960bf..0000000000 --- a/src/Si.tersoff.modc +++ /dev/null @@ -1,10 +0,0 @@ -# DATE: 2016-11-09 UNITS: metal CONTRIBUTOR: Ganga P Purja Pun (George Mason University, Fairfax) CITATION: Unknown -# -# Format: -# element1 element2 element3 -# beta alpha h eta -# beta_ters lam2 B R D lam1 A -# n c1 c2 c3 c4 c5 C -Si Si Si 3.00000000 1.80536502 -0.38136087 2.16152496 -1 1.39343356 117.78072440 2.87478837 0.33090566 3.18011795 3198.51383127 -1.98633876 0.20123243 614230.04310619 996439.09714140 3.33560562 25.20963770 -0.00592042 From 894b7810b244bbba793db75e77255a907182b109 Mon Sep 17 00:00:00 2001 From: Stan Moore Date: Fri, 29 Apr 2022 09:32:24 -0600 Subject: [PATCH 14/16] whitespace --- src/KOKKOS/pair_snap_kokkos_impl.h | 4 ++-- src/KOKKOS/sna_kokkos.h | 4 ++-- src/KOKKOS/sna_kokkos_impl.h | 30 +++++++++++++++--------------- 3 files changed, 19 insertions(+), 19 deletions(-) diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h index e39f42058f..74264ca47e 100644 --- a/src/KOKKOS/pair_snap_kokkos_impl.h +++ b/src/KOKKOS/pair_snap_kokkos_impl.h @@ -92,8 +92,8 @@ void PairSNAPKokkos::init_style() if (host_flag) { if (lmp->kokkos->nthreads > 1) if (comm->me == 0) - utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " - "CPU thread, even if more threads are requested\n"); + utils::logmesg(lmp,"Pair style snap/kk currently only runs on a single " + "CPU thread, even if more threads are requested\n"); PairSNAP::init_style(); return; diff --git a/src/KOKKOS/sna_kokkos.h b/src/KOKKOS/sna_kokkos.h index 4f10a9ac0b..7ed2736947 100644 --- a/src/KOKKOS/sna_kokkos.h +++ b/src/KOKKOS/sna_kokkos.h @@ -315,8 +315,8 @@ inline KOKKOS_INLINE_FUNCTION void compute_duarray_cpu(const typename Kokkos::TeamPolicy::member_type& team, int, int, const real_type&, const real_type&, const real_type&, // compute_duidrj_cpu - const real_type&, const real_type&, const real_type&, const real_type&, const real_type&, - const real_type&, const real_type&); + const real_type&, const real_type&, const real_type&, const real_type&, const real_type&, + const real_type&, const real_type&); // Sets the style for the switching function // 0 = none diff --git a/src/KOKKOS/sna_kokkos_impl.h b/src/KOKKOS/sna_kokkos_impl.h index c9c1bfbb99..77130b9781 100644 --- a/src/KOKKOS/sna_kokkos_impl.h +++ b/src/KOKKOS/sna_kokkos_impl.h @@ -31,7 +31,7 @@ template inline SNAKokkos::SNAKokkos(real_type rfac0_in, int twojmax_in, real_type rmin0_in, int switch_flag_in, int bzero_flag_in, - int chem_flag_in, int bnorm_flag_in, int wselfall_flag_in, int nelements_in, int switch_inner_flag_in) + int chem_flag_in, int bnorm_flag_in, int wselfall_flag_in, int nelements_in, int switch_inner_flag_in) { LAMMPS_NS::ExecutionSpace execution_space = ExecutionSpaceFromDevice::space; host_flag = (execution_space == LAMMPS_NS::Host); @@ -1753,8 +1753,8 @@ KOKKOS_INLINE_FUNCTION void SNAKokkos::compute_duarray_cpu(const typename Kokkos::TeamPolicy::member_type& team, int iatom, int jnbor, const real_type& x, const real_type& y, const real_type& z, const real_type& z0, const real_type& r, const real_type& dz0dr, - const real_type& wj, const real_type& rcut, - const real_type& sinner, const real_type& dinner) + const real_type& wj, const real_type& rcut, + const real_type& sinner, const real_type& dinner) { real_type r0inv; real_type a_r, a_i, b_r, b_i; @@ -2247,15 +2247,15 @@ real_type SNAKokkos::compute_sfac(real_typ sfac_outer = onehalf * (cos((r - rmin0) * rcutfac) + one); } } - + if (switch_inner_flag == 0) return sfac_outer; if (switch_inner_flag == 1) { - if (r >= sinner + dinner) - return sfac_outer; + if (r >= sinner + dinner) + return sfac_outer; else if (r > sinner - dinner) { real_type rcutfac = static_cast(MY_PI2) / dinner; - return sfac_outer * - onehalf * (one - cos(static_cast(MY_PI2) + (r - sinner) * rcutfac)); + return sfac_outer * + onehalf * (one - cos(static_cast(MY_PI2) + (r - sinner) * rcutfac)); } else return zero; } return zero; // dummy return @@ -2283,7 +2283,7 @@ real_type SNAKokkos::compute_dsfac(real_ty if (switch_inner_flag == 0) return dsfac_outer; if (switch_inner_flag == 1) { - if (r >= sinner + dinner) + if (r >= sinner + dinner) return dsfac_outer; else if (r > sinner - dinner) { @@ -2291,12 +2291,12 @@ real_type SNAKokkos::compute_dsfac(real_ty if (switch_flag == 0) sfac_outer = one; if (switch_flag == 1) { - if (r <= rmin0) sfac_outer = one; - else if (r > rcut) sfac_outer = zero; - else { - real_type rcutfac = static_cast(MY_PI) / (rcut - rmin0); - sfac_outer = onehalf * (cos((r - rmin0) * rcutfac) + one); - } + if (r <= rmin0) sfac_outer = one; + else if (r > rcut) sfac_outer = zero; + else { + real_type rcutfac = static_cast(MY_PI) / (rcut - rmin0); + sfac_outer = onehalf * (cos((r - rmin0) * rcutfac) + one); + } } // calculate sfac_inner From 2f8e70818405ddca27ad478456fd495eeb336973 Mon Sep 17 00:00:00 2001 From: Stan Moore Date: Fri, 29 Apr 2022 09:37:27 -0600 Subject: [PATCH 15/16] more whitespace --- src/KOKKOS/pair_snap_kokkos.h | 16 ++++++++-------- src/KOKKOS/sna_kokkos.h | 25 +++++++++++++------------ 2 files changed, 21 insertions(+), 20 deletions(-) diff --git a/src/KOKKOS/pair_snap_kokkos.h b/src/KOKKOS/pair_snap_kokkos.h index 1c0bf38699..2cd0003478 100644 --- a/src/KOKKOS/pair_snap_kokkos.h +++ b/src/KOKKOS/pair_snap_kokkos.h @@ -73,7 +73,7 @@ struct TagPairSNAPComputeDeidrjCPU{}; template class PairSNAPKokkos : public PairSNAP { -public: + public: enum {EnabledNeighFlags=FULL|HALF|HALFTHREAD}; enum {COUL_FLAG=0}; typedef DeviceType device_type; @@ -228,7 +228,7 @@ public: const F_FLOAT &fx, const F_FLOAT &fy, const F_FLOAT &fz, const F_FLOAT &delx, const F_FLOAT &dely, const F_FLOAT &delz) const; -protected: + protected: typename AT::t_neighbors_2d d_neighbors; typename AT::t_int_1d_randomread d_ilist; typename AT::t_int_1d_randomread d_numneigh; @@ -252,9 +252,9 @@ protected: void allocate() override; //void read_files(char *, char *); /*template -inline int equal(double* x,double* y); + inline int equal(double* x,double* y); template -inline double dist2(double* x,double* y); + inline double dist2(double* x,double* y); double extra_cutoff(); void load_balance(); void set_sna_to_shared(int snaid,int i); @@ -331,10 +331,10 @@ inline double dist2(double* x,double* y); template class PairSNAPKokkosDevice : public PairSNAPKokkos { -private: + private: using Base = PairSNAPKokkos; -public: + public: PairSNAPKokkosDevice(class LAMMPS *); @@ -350,10 +350,10 @@ public: template class PairSNAPKokkosHost : public PairSNAPKokkos { -private: + private: using Base = PairSNAPKokkos; -public: + public: PairSNAPKokkosHost(class LAMMPS *); diff --git a/src/KOKKOS/sna_kokkos.h b/src/KOKKOS/sna_kokkos.h index 7ed2736947..fe70129660 100644 --- a/src/KOKKOS/sna_kokkos.h +++ b/src/KOKKOS/sna_kokkos.h @@ -64,7 +64,7 @@ struct alignas(8) FullHalfMapper { template class SNAKokkos { -public: + public: using real_type = real_type_; using complex = SNAComplex; static constexpr int vector_length = vector_length_; @@ -97,21 +97,22 @@ public: typedef Kokkos::View t_sna_3c3; typedef Kokkos::View t_sna_5c; -inline + inline SNAKokkos() {}; + KOKKOS_INLINE_FUNCTION SNAKokkos(const SNAKokkos& sna, const typename Kokkos::TeamPolicy::member_type& team); -inline -SNAKokkos(real_type, int, real_type, int, int, int, int, int, int, int); + inline + SNAKokkos(real_type, int, real_type, int, int, int, int, int, int, int); KOKKOS_INLINE_FUNCTION ~SNAKokkos(); -inline + inline void build_indexlist(); // SNAKokkos() -inline + inline void init(); // double memory_usage(); @@ -257,7 +258,7 @@ inline int ndoubles; int ntriples; -private: + private: real_type rmin0, rfac0; //use indexlist instead of loops, constructor generates these @@ -266,13 +267,13 @@ private: Kokkos::View idxb; Kokkos::View idxcg_block; -public: + public: Kokkos::View idxu_block; Kokkos::View idxu_half_block; Kokkos::View idxu_cache_block; Kokkos::View idxu_full_half; -private: + private: Kokkos::View idxz_block; Kokkos::View idxb_block; @@ -292,10 +293,10 @@ private: KOKKOS_INLINE_FUNCTION void create_thread_scratch_arrays(const typename Kokkos::TeamPolicy::member_type& team); // SNAKokkos() -inline + inline void init_clebsch_gordan(); // init() -inline + inline void init_rootpqarray(); // init() KOKKOS_INLINE_FUNCTION @@ -310,7 +311,7 @@ inline inline double deltacg(int, int, int); // init_clebsch_gordan -inline + inline int compute_ncoeff(); // SNAKokkos() KOKKOS_INLINE_FUNCTION void compute_duarray_cpu(const typename Kokkos::TeamPolicy::member_type& team, int, int, From 37d2d9e013f480afca3cf50ef9f7ee4a32453a31 Mon Sep 17 00:00:00 2001 From: Stan Moore Date: Fri, 29 Apr 2022 10:56:34 -0600 Subject: [PATCH 16/16] Remove trailing spaces --- doc/src/compute_sna_atom.rst | 4 ++-- src/ML-SNAP/pair_snap.cpp | 6 +++--- src/ML-SNAP/sna.cpp | 12 ++++++------ 3 files changed, 11 insertions(+), 11 deletions(-) diff --git a/doc/src/compute_sna_atom.rst b/doc/src/compute_sna_atom.rst index e34c9a8889..54a6df02a2 100644 --- a/doc/src/compute_sna_atom.rst +++ b/doc/src/compute_sna_atom.rst @@ -335,7 +335,7 @@ where the switching region is centered at :math:`S_{inner}` and it extends a dis to the left and to the right of this. With this option, additional keywords *sinner* and *dinner* must be used, each followed by *ntypes* -values for :math:`S_{inner}` and :math:`D_{inner}`, respectively. +values for :math:`S_{inner}` and :math:`D_{inner}`, respectively. When the central atom and the neighbor atom have different types, the values of :math:`S_{inner}` and :math:`D_{inner}` are the arithmetic means of the values for both types. @@ -458,7 +458,7 @@ Default The optional keyword defaults are *rmin0* = 0, *switchflag* = 1, *bzeroflag* = 1, *quadraticflag* = 0, -*bnormflag* = 0, *wselfallflag* = 0, *switchinnerflag* = 0, +*bnormflag* = 0, *wselfallflag* = 0, *switchinnerflag* = 0, ---------- diff --git a/src/ML-SNAP/pair_snap.cpp b/src/ML-SNAP/pair_snap.cpp index 56c3250c3f..07a8237ab5 100644 --- a/src/ML-SNAP/pair_snap.cpp +++ b/src/ML-SNAP/pair_snap.cpp @@ -694,11 +694,11 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) if ((int)words.size() != nelements+1) error->all(FLERR,"Incorrect SNAP parameter file"); - + // innerlogstr collects all values of sinner or dinner for log output below std::string innerlogstr; - + int iword = 1; if (keywd == "sinner") { @@ -718,7 +718,7 @@ void PairSNAP::read_files(char *coefffilename, char *paramfilename) } dinnerflag = 1; } - + if (comm->me == 0) utils::logmesg(lmp,"SNAP keyword {} {} ... \n", keywd, innerlogstr); diff --git a/src/ML-SNAP/sna.cpp b/src/ML-SNAP/sna.cpp index 2b86774d5b..33937b9c45 100644 --- a/src/ML-SNAP/sna.cpp +++ b/src/ML-SNAP/sna.cpp @@ -1552,7 +1552,7 @@ double SNA::compute_sfac(double r, double rcut, double sinner, double dinner) sfac *= 0.5 * (1.0 - cos(MY_PI2 + (r - sinner) * rcutfac)); } else sfac = 0.0; } - + return sfac; } @@ -1568,14 +1568,14 @@ double SNA::compute_dsfac(double r, double rcut, double sinner, double dinner) double rcutfac = MY_PI / (rcut - rmin0); dsfac_outer = -0.5 * sin((r - rmin0) * rcutfac) * rcutfac; } - - // some duplicated computation, but rarely visited + + // some duplicated computation, but rarely visited if (switch_inner_flag == 1 && r < sinner + dinner) { if (r > sinner - dinner) { // calculate sfac_outer - + if (switch_flag == 0) sfac_outer = 1.0; else if (r <= rmin0) sfac_outer = 1.0; else if (r > rcut) sfac_outer = 0.0; @@ -1585,14 +1585,14 @@ double SNA::compute_dsfac(double r, double rcut, double sinner, double dinner) } // calculate sfac_inner - + double rcutfac = MY_PI2 / dinner; sfac_inner = 0.5 * (1.0 - cos(MY_PI2 + (r - sinner) * rcutfac)); dsfac_inner = 0.5 * rcutfac * sin(MY_PI2 + (r - sinner) * rcutfac); dsfac = dsfac_outer*sfac_inner + sfac_outer*dsfac_inner; } else dsfac = 0.0; } else dsfac = dsfac_outer; - + return dsfac; }