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; }