Completed inner cutoff, KOKKOS still in progress
This commit is contained in:
Aidan Thompson
@ -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,
|
||||
|
||||
----------
|
||||
|
||||
|
||||
@ -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
|
||||
|
||||
@ -602,7 +602,7 @@ void PairSNAPKokkos<DeviceType, real_type, vector_length>::coeff(int narg, char
|
||||
Kokkos::deep_copy(d_map,h_map);
|
||||
|
||||
snaKK = SNAKokkos<DeviceType, real_type, vector_length>(rfac0,twojmax,
|
||||
rmin0,switchflag,bzeroflag,chemflag,bnormflag,wselfallflag,nelements);
|
||||
rmin0,switchflag,bzeroflag,chemflag,bnormflag,wselfallflag,nelements,switchinnerflag);
|
||||
snaKK.grow_rij(0,0);
|
||||
snaKK.init();
|
||||
}
|
||||
|
||||
@ -103,7 +103,7 @@ inline
|
||||
SNAKokkos(const SNAKokkos<DeviceType,real_type,vector_length>& sna, const typename Kokkos::TeamPolicy<DeviceType>::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<DeviceType>::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<DeviceType>::member_type& team, int, int, const real_type&, const real_type&, const real_type&, int); // compute_ui
|
||||
void add_uarraytot(const typename Kokkos::TeamPolicy<DeviceType>::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<DeviceType>::member_type& team, int, int,
|
||||
@ -313,13 +315,19 @@ inline
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void compute_duarray_cpu(const typename Kokkos::TeamPolicy<DeviceType>::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;
|
||||
|
||||
@ -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<class DeviceType, typename real_type, int vector_length>
|
||||
inline
|
||||
SNAKokkos<DeviceType, real_type, vector_length>::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<DeviceType>::space;
|
||||
host_flag = (execution_space == LAMMPS_NS::Host);
|
||||
@ -40,6 +41,7 @@ SNAKokkos<DeviceType, real_type, vector_length>::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<DeviceType, real_type, vector_length>::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<DeviceType, real_type, vector_length>::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<real_type>(MY_PI) / (rcut - rmin0);
|
||||
const real_type theta0 = (r - rmin0) * rscale0;
|
||||
const real_type sn = sin(theta0);
|
||||
@ -380,7 +386,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::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<DeviceType, real_type, vector_length>::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<DeviceType, real_type, vector_length>::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<DeviceType, real_type, vector_length>::compute_deidrj_cpu(const t
|
||||
template<class DeviceType, typename real_type, int vector_length>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void SNAKokkos<DeviceType, real_type, vector_length>::add_uarraytot(const typename Kokkos::TeamPolicy<DeviceType>::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<DeviceType, real_type, vector_length>::compute_duarray_cpu(const typename Kokkos::TeamPolicy<DeviceType>::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<DeviceType, real_type, vector_length>::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<DeviceType, real_type, vector_length>::compute_ncoeff()
|
||||
|
||||
template<class DeviceType, typename real_type, int vector_length>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_sfac(real_type r, real_type rcut)
|
||||
real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_sfac(real_type r, real_type rcut, real_type sinner, real_type dinner)
|
||||
{
|
||||
constexpr real_type one = static_cast<real_type>(1.0);
|
||||
constexpr real_type zero = static_cast<real_type>(0.0);
|
||||
@ -2235,18 +2243,30 @@ real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_sfac(real_typ
|
||||
else if (r > rcut) return zero;
|
||||
else {
|
||||
real_type rcutfac = static_cast<real_type>(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<real_type>(MY_PI2) / dinner;
|
||||
return onehalf * (cos((r - rmin0) * rcutfac) + one) *
|
||||
onehalf * (one - cos(static_cast<real_type>(MY_PI2) + (r - sinner) * rcutfacinner));
|
||||
} else return zero;
|
||||
}
|
||||
return zero; // dummy return
|
||||
}
|
||||
}
|
||||
return zero;
|
||||
return zero; // dummy return
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType, typename real_type, int vector_length>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_dsfac(real_type r, real_type rcut)
|
||||
real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_dsfac(real_type r, real_type rcut, real_type sinner, real_type dinner)
|
||||
{
|
||||
constexpr real_type one = static_cast<real_type>(1.0);
|
||||
constexpr real_type zero = static_cast<real_type>(0.0);
|
||||
constexpr real_type onehalf = static_cast<real_type>(0.5);
|
||||
if (switch_flag == 0) return zero;
|
||||
@ -2258,12 +2278,12 @@ real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_dsfac(real_ty
|
||||
return -onehalf * sin((r - rmin0) * rcutfac) * rcutfac;
|
||||
}
|
||||
}
|
||||
return zero;
|
||||
return zero; // dummy return
|
||||
}
|
||||
|
||||
template<class DeviceType, typename real_type, int vector_length>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void SNAKokkos<DeviceType, real_type, vector_length>::compute_s_dsfac(const real_type r, const real_type rcut, real_type& sfac, real_type& dsfac) {
|
||||
void SNAKokkos<DeviceType, real_type, vector_length>::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<real_type>(1.0);
|
||||
constexpr real_type zero = static_cast<real_type>(0.0);
|
||||
constexpr real_type onehalf = static_cast<real_type>(0.5);
|
||||
@ -2356,6 +2376,8 @@ double SNAKokkos<DeviceType, real_type, vector_length>::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;
|
||||
}
|
||||
|
||||
@ -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++;
|
||||
|
||||
@ -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;
|
||||
|
||||
@ -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++;
|
||||
|
||||
@ -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;
|
||||
|
||||
@ -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++;
|
||||
|
||||
@ -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;
|
||||
|
||||
@ -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++;
|
||||
|
||||
@ -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;
|
||||
};
|
||||
|
||||
@ -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");
|
||||
}
|
||||
|
||||
|
||||
@ -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
|
||||
|
||||
@ -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;
|
||||
}
|
||||
|
||||
|
||||
@ -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
|
||||
|
||||
|
||||
Reference in New Issue
Block a user