Adding Kokkos files
This commit is contained in:
Stan Moore
@ -113,6 +113,10 @@ if (test $1 = "USER-CG-CMM") then
|
||||
depend USER-OMP
|
||||
fi
|
||||
|
||||
if (test $1 = "USER-DPD") then
|
||||
depend KOKKOS
|
||||
fi
|
||||
|
||||
if (test $1 = "USER-FEP") then
|
||||
depend USER-OMP
|
||||
fi
|
||||
|
||||
@ -73,6 +73,8 @@ action domain_kokkos.cpp
|
||||
action domain_kokkos.h
|
||||
action fix_deform_kokkos.cpp
|
||||
action fix_deform_kokkos.h
|
||||
action fix_eos_table_rx_kokkos.cpp fix_eos_table_rx.cpp
|
||||
action fix_eos_table_rx_kokkos.h fix_eos_table_rx.h
|
||||
action fix_langevin_kokkos.cpp
|
||||
action fix_langevin_kokkos.h
|
||||
action fix_nh_kokkos.cpp
|
||||
@ -171,6 +173,8 @@ action pair_lj_gromacs_kokkos.cpp
|
||||
action pair_lj_gromacs_kokkos.h
|
||||
action pair_lj_sdk_kokkos.cpp pair_lj_sdk.cpp
|
||||
action pair_lj_sdk_kokkos.h pair_lj_sdk.h
|
||||
action pair_multi_lucy_rx_kokkos.cpp pair_multi_lucy_rx.cpp
|
||||
action pair_multi_lucy_rx_kokkos.h pair_multi_lucy_rx.h
|
||||
action pair_reax_c_kokkos.cpp pair_reax_c.cpp
|
||||
action pair_reax_c_kokkos.h pair_reax_c.h
|
||||
action pair_sw_kokkos.cpp pair_sw.cpp
|
||||
@ -179,6 +183,8 @@ action pair_vashishta_kokkos.cpp pair_vashishta.cpp
|
||||
action pair_vashishta_kokkos.h pair_vashishta.h
|
||||
action pair_table_kokkos.cpp
|
||||
action pair_table_kokkos.h
|
||||
action pair_table_rx_kokkos.cpp pair_table_rx.cpp
|
||||
action pair_table_rx_kokkos.h pair_table_rx.h
|
||||
action pair_tersoff_kokkos.cpp pair_tersoff.cpp
|
||||
action pair_tersoff_kokkos.h pair_tersoff.h
|
||||
action pair_tersoff_mod_kokkos.cpp pair_tersoff_mod.cpp
|
||||
|
||||
354
src/KOKKOS/fix_eos_table_rx_kokkos.cpp
Normal file
354
src/KOKKOS/fix_eos_table_rx_kokkos.cpp
Normal file
@ -0,0 +1,354 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Stan Moore (Sandia)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <stdlib.h>
|
||||
#include <string.h>
|
||||
#include "fix_eos_table_rx_kokkos.h"
|
||||
#include "atom_kokkos.h"
|
||||
#include "error.h"
|
||||
#include "force.h"
|
||||
#include "memory.h"
|
||||
#include "comm.h"
|
||||
#include <math.h>
|
||||
#include "modify.h"
|
||||
#include "atom_masks.h"
|
||||
|
||||
#define MAXLINE 1024
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
using namespace FixConst;
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
FixEOStableRXKokkos<DeviceType>::FixEOStableRXKokkos(LAMMPS *lmp, int narg, char **arg) :
|
||||
FixEOStableRX(lmp, narg, arg)
|
||||
{
|
||||
atomKK = (AtomKokkos *) atom;
|
||||
execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
|
||||
datamask_read = X_MASK | F_MASK | TYPE_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
FixEOStableRXKokkos<DeviceType>::~FixEOStableRXKokkos()
|
||||
{
|
||||
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void FixEOStableRXKokkos<DeviceType>::setup(int vflag)
|
||||
{
|
||||
int nlocal = atom->nlocal;
|
||||
mask = atomKK->k_mask.view<DeviceType>();
|
||||
uCond = atomKK->k_uCond.view<DeviceType>();
|
||||
uMech = atomKK->k_uMech.view<DeviceType>();
|
||||
uChem = atomKK->k_uChem.view<DeviceType>();
|
||||
dpdTheta= atomKK->k_dpdTheta.view<DeviceType>();
|
||||
uCG = atomKK->k_uCG.view<DeviceType>();
|
||||
uCGnew = atomKK->k_uCGnew.view<DeviceType>();
|
||||
double duChem;
|
||||
|
||||
for (int i = 0; i < nlocal; i++) // parallel_for
|
||||
if (mask[i] & groupbit){
|
||||
duChem = uCG[i] - uCGnew[i];
|
||||
uChem[i] += duChem;
|
||||
uCG[i] = 0.0;
|
||||
uCGnew[i] = 0.0;
|
||||
}
|
||||
|
||||
// Communicate the updated momenta and velocities to all nodes
|
||||
comm->forward_comm_fix(this);
|
||||
|
||||
for (int i = 0; i < nlocal; i++) // parallel_for
|
||||
if (mask[i] & groupbit)
|
||||
temperature_lookup(i,uCond[i]+uMech[i]+uChem[i],dpdTheta[i]);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void FixEOStableRXKokkos<DeviceType>::init()
|
||||
{
|
||||
int nlocal = atom->nlocal;
|
||||
mask = atomKK->k_mask.view<DeviceType>();
|
||||
uCond = atomKK->k_uCond.view<DeviceType>();
|
||||
uMech = atomKK->k_uMech.view<DeviceType>();
|
||||
uChem = atomKK->k_uChem.view<DeviceType>();
|
||||
dpdTheta= atomKK->k_dpdTheta.view<DeviceType>();
|
||||
double tmp;
|
||||
|
||||
if(this->restart_reset){
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
if (mask[i] & groupbit)
|
||||
temperature_lookup(i,uCond[i]+uMech[i]+uChem[i],dpdTheta[i]);
|
||||
} else {
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
if (mask[i] & groupbit) {
|
||||
if(dpdTheta[i] <= 0.0)
|
||||
error->one(FLERR,"Internal temperature <= zero");
|
||||
energy_lookup(i,dpdTheta[i],tmp);
|
||||
uCond[i] = tmp / 2.0;
|
||||
uMech[i] = tmp / 2.0;
|
||||
uChem[i] = 0.0;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void FixEOStableRXKokkos<DeviceType>::post_integrate()
|
||||
{
|
||||
int nlocal = atom->nlocal;
|
||||
mask = atomKK->k_mask.view<DeviceType>();
|
||||
uCond = atomKK->k_uCond.view<DeviceType>();
|
||||
uMech = atomKK->k_uMech.view<DeviceType>();
|
||||
uChem = atomKK->k_uChem.view<DeviceType>();
|
||||
dpdTheta= atomKK->k_dpdTheta.view<DeviceType>();
|
||||
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
if (mask[i] & groupbit){
|
||||
temperature_lookup(i,uCond[i]+uMech[i]+uChem[i],dpdTheta[i]);
|
||||
if(dpdTheta[i] <= 0.0)
|
||||
error->one(FLERR,"Internal temperature <= zero");
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void FixEOStableRXKokkos<DeviceType>::end_of_step()
|
||||
{
|
||||
int nlocal = atom->nlocal;
|
||||
mask = atomKK->k_mask.view<DeviceType>();
|
||||
uCond = atomKK->k_uCond.view<DeviceType>();
|
||||
uMech = atomKK->k_uMech.view<DeviceType>();
|
||||
uChem = atomKK->k_uChem.view<DeviceType>();
|
||||
dpdTheta= atomKK->k_dpdTheta.view<DeviceType>();
|
||||
uCG = atomKK->k_uCG.view<DeviceType>();
|
||||
uCGnew = atomKK->k_uCGnew.view<DeviceType>();
|
||||
double duChem;
|
||||
|
||||
// Communicate the ghost uCGnew
|
||||
comm->reverse_comm_fix(this);
|
||||
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
if (mask[i] & groupbit){
|
||||
duChem = uCG[i] - uCGnew[i];
|
||||
uChem[i] += duChem;
|
||||
uCG[i] = 0.0;
|
||||
uCGnew[i] = 0.0;
|
||||
}
|
||||
|
||||
// Communicate the updated momenta and velocities to all nodes
|
||||
comm->forward_comm_fix(this);
|
||||
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
if (mask[i] & groupbit){
|
||||
temperature_lookup(i,uCond[i]+uMech[i]+uChem[i],dpdTheta[i]);
|
||||
if(dpdTheta[i] <= 0.0)
|
||||
error->one(FLERR,"Internal temperature <= zero");
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
calculate potential ui at temperature thetai
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void FixEOStableRXKokkos<DeviceType>::energy_lookup(int id, double thetai, double &ui) const
|
||||
{
|
||||
int itable;
|
||||
double fraction, uTmp, nTotal;
|
||||
|
||||
ui = 0.0;
|
||||
nTotal = 0.0;
|
||||
for(int ispecies=0;ispecies<nspecies;ispecies++){
|
||||
Table *tb = &tables[ispecies];
|
||||
thetai = MAX(thetai,tb->lo);
|
||||
thetai = MIN(thetai,tb->hi);
|
||||
|
||||
if (tabstyle == LINEAR) {
|
||||
itable = static_cast<int> ((thetai - tb->lo) * tb->invdelta);
|
||||
fraction = (thetai - tb->r[itable]) * tb->invdelta;
|
||||
uTmp = tb->e[itable] + fraction*tb->de[itable];
|
||||
|
||||
uTmp += dHf[ispecies];
|
||||
// mol fraction form:
|
||||
ui += atom->dvector[ispecies][id]*uTmp;
|
||||
nTotal += atom->dvector[ispecies][id];
|
||||
}
|
||||
}
|
||||
ui = ui - double(nTotal+1.5)*force->boltz*thetai;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
calculate temperature thetai at energy ui
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void FixEOStableRXKokkos<DeviceType>::temperature_lookup(int id, double ui, double &thetai) const
|
||||
{
|
||||
Table *tb = &tables[0];
|
||||
|
||||
int it;
|
||||
double t1,t2,u1,u2,f1,f2;
|
||||
double maxit = 100;
|
||||
double temp;
|
||||
double delta = 0.001;
|
||||
|
||||
// Store the current thetai in t1
|
||||
t1 = MAX(thetai,tb->lo);
|
||||
t1 = MIN(t1,tb->hi);
|
||||
if(t1==tb->hi) delta = -delta;
|
||||
|
||||
// Compute u1 at thetai
|
||||
energy_lookup(id,t1,u1);
|
||||
|
||||
// Compute f1
|
||||
f1 = u1 - ui;
|
||||
|
||||
// Compute guess of t2
|
||||
t2 = (1.0 + delta)*t1;
|
||||
|
||||
// Compute u2 at t2
|
||||
energy_lookup(id,t2,u2);
|
||||
|
||||
// Compute f1
|
||||
f2 = u2 - ui;
|
||||
|
||||
// Apply the Secant Method
|
||||
for(it=0; it<maxit; it++){
|
||||
if(fabs(f2-f1)<1e-15){
|
||||
if(isnan(f1) || isnan(f2)) error->one(FLERR,"NaN detected in secant solver.");
|
||||
temp = t1;
|
||||
temp = MAX(temp,tb->lo);
|
||||
temp = MIN(temp,tb->hi);
|
||||
char str[256];
|
||||
sprintf(str,"Secant solver did not converge because table bounds were exceeded: it=%d id=%d ui=%lf thetai=%lf t1=%lf t2=%lf f1=%lf f2=%lf dpdTheta=%lf\n",it,id,ui,thetai,t1,t2,f1,f2,temp);
|
||||
error->warning(FLERR,str);
|
||||
break;
|
||||
}
|
||||
temp = t2 - f2*(t2-t1)/(f2-f1);
|
||||
if(fabs(temp-t2) < 1e-6) break;
|
||||
f1 = f2;
|
||||
t1 = t2;
|
||||
t2 = temp;
|
||||
energy_lookup(id,t2,u2);
|
||||
f2 = u2 - ui;
|
||||
}
|
||||
if(it==maxit){
|
||||
char str[256];
|
||||
sprintf(str,"Maxit exceeded in secant solver: id=%d ui=%lf thetai=%lf t1=%lf t2=%lf f1=%lf f2=%lf\n",id,ui,thetai,t1,t2,f1,f2);
|
||||
if(isnan(f1) || isnan(f2) || isnan(ui) || isnan(thetai) || isnan(t1) || isnan(t2))
|
||||
error->one(FLERR,"NaN detected in secant solver.");
|
||||
error->one(FLERR,str);
|
||||
}
|
||||
thetai = temp;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
int FixEOStableRXKokkos<DeviceType>::pack_forward_comm(int n, int *list, double *buf, int pbc_flag, int *pbc)
|
||||
{
|
||||
int ii,jj,m;
|
||||
uChem = atomKK->k_uChem.view<DeviceType>();
|
||||
uCG = atomKK->k_uCG.view<DeviceType>();
|
||||
uCGnew = atomKK->k_uCGnew.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
for (ii = 0; ii < n; ii++) {
|
||||
jj = list[ii];
|
||||
buf[m++] = uChem[jj];
|
||||
buf[m++] = uCG[jj];
|
||||
buf[m++] = uCGnew[jj];
|
||||
}
|
||||
return m;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void FixEOStableRXKokkos<DeviceType>::unpack_forward_comm(int n, int first, double *buf)
|
||||
{
|
||||
int ii,m,last;
|
||||
uChem = atomKK->k_uChem.view<DeviceType>();
|
||||
uCG = atomKK->k_uCG.view<DeviceType>();
|
||||
uCGnew = atomKK->k_uCGnew.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
last = first + n ;
|
||||
for (ii = first; ii < last; ii++){
|
||||
uChem[ii] = buf[m++];
|
||||
uCG[ii] = buf[m++];
|
||||
uCGnew[ii] = buf[m++];
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
int FixEOStableRXKokkos<DeviceType>::pack_reverse_comm(int n, int first, double *buf)
|
||||
{
|
||||
int i,m,last;
|
||||
uCG = atomKK->k_uCG.view<DeviceType>();
|
||||
uCGnew = atomKK->k_uCGnew.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
last = first + n;
|
||||
for (i = first; i < last; i++) {
|
||||
buf[m++] = uCG[i];
|
||||
buf[m++] = uCGnew[i];
|
||||
}
|
||||
return m;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void FixEOStableRXKokkos<DeviceType>::unpack_reverse_comm(int n, int *list, double *buf)
|
||||
{
|
||||
int i,j,m;
|
||||
uCG = atomKK->k_uCG.view<DeviceType>();
|
||||
uCGnew = atomKK->k_uCGnew.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
for (i = 0; i < n; i++) {
|
||||
j = list[i];
|
||||
|
||||
uCG[j] += buf[m++];
|
||||
uCGnew[j] += buf[m++];
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
template class FixEOStableRXKokkos<LMPDeviceType>;
|
||||
#ifdef KOKKOS_HAVE_CUDA
|
||||
template class FixEOStableRXKokkos<LMPHostType>;
|
||||
#endif
|
||||
}
|
||||
152
src/KOKKOS/fix_eos_table_rx_kokkos.h
Normal file
152
src/KOKKOS/fix_eos_table_rx_kokkos.h
Normal file
@ -0,0 +1,152 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef FIX_CLASS
|
||||
|
||||
FixStyle(eos/table/rx/kk,FixEOStableRXKokkos<LMPDeviceType>)
|
||||
FixStyle(eos/table/rx/kk/device,FixEOStableRXKokkos<LMPDeviceType>)
|
||||
FixStyle(eos/table/rx/kk/host,FixEOStableRXKokkos<LMPHostType>)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_FIX_EOS_TABLE_RX_KOKKOS_H
|
||||
#define LMP_FIX_EOS_TABLE_RX_KOKKOS_H
|
||||
|
||||
#include "fix_eos_table_rx.h"
|
||||
#include "kokkos_type.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
template<class DeviceType>
|
||||
class FixEOStableRXKokkos : public FixEOStableRX {
|
||||
public:
|
||||
typedef DeviceType device_type;
|
||||
typedef ArrayTypes<DeviceType> AT;
|
||||
typedef EV_FLOAT value_type;
|
||||
|
||||
FixEOStableRXKokkos(class LAMMPS *, int, char **);
|
||||
virtual ~FixEOStableRXKokkos();
|
||||
void setup(int);
|
||||
void init();
|
||||
void post_integrate();
|
||||
void end_of_step();
|
||||
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void energy_lookup(int, double, double &) const;
|
||||
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void temperature_lookup(int, double, double &) const;
|
||||
|
||||
protected:
|
||||
//struct Table {
|
||||
// int ninput;
|
||||
// double lo,hi;
|
||||
// double *rfile,*efile;
|
||||
// double *e2file;
|
||||
// double delta,invdelta,deltasq6;
|
||||
// double *r,*e,*de,*e2;
|
||||
//};
|
||||
//Table *tables, *tables2;
|
||||
|
||||
void allocate();
|
||||
|
||||
//double *dHf;
|
||||
|
||||
typename AT::t_int_1d mask;
|
||||
typename AT::t_efloat_1d uCond,uMech,uChem,uCG,uCGnew,rho,dpdTheta,duChem;
|
||||
|
||||
int pack_reverse_comm(int, int, double *);
|
||||
void unpack_reverse_comm(int, int *, double *);
|
||||
int pack_forward_comm(int , int *, double *, int, int *);
|
||||
void unpack_forward_comm(int , int , double *);
|
||||
|
||||
//int *eosSpecies;
|
||||
};
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Illegal ... command
|
||||
|
||||
Self-explanatory. Check the input script syntax and compare to the
|
||||
documentation for the command. You can use -echo screen as a
|
||||
command-line option when running LAMMPS to see the offending line.
|
||||
|
||||
E: FixEOStableRXKokkos requires a fix rx command.
|
||||
|
||||
The fix rx command must come before the pair style command in the input file
|
||||
|
||||
E: There are no rx species specified
|
||||
|
||||
There must be at least one species specified through the fix rx command
|
||||
|
||||
E: Invalid eos/table/rx length
|
||||
|
||||
The eos/table/rx table must have more than one entry.
|
||||
|
||||
E: eos/table/rx values are not increasing
|
||||
|
||||
The equation-of-state must an increasing function
|
||||
|
||||
E: Internal temperature <= zero.
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Cannot open eos table/rx potential file %s
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Incorrect format in eos table/rx file
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Cannot open file %s
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Did not find keyword in table file
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Illegal fix eos/table/rx command
|
||||
|
||||
Incorrect number of arguments specified for the fix eos/table/rx command.
|
||||
|
||||
E: Invalid keyword in fix eos/table/rx parameters
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: The number of columns in fix eos/table/rx does not match the number of species.
|
||||
|
||||
Self-explanatory. Check format for fix eos/table/rx file.
|
||||
|
||||
E: fix eos/table/rx parameters did not set N
|
||||
|
||||
The number of table entries was not set in the eos/table/rx file
|
||||
|
||||
W: Secant solver did not converge because table bounds were exceeded
|
||||
|
||||
The secant solver failed to converge, resulting in the lower or upper table bound temperature to be returned
|
||||
|
||||
E: NaN detected in secant solver.
|
||||
|
||||
Self-explanatory.
|
||||
|
||||
E: Maxit exceeded in secant solver
|
||||
|
||||
The maximum number of interations was exceeded in the secant solver
|
||||
|
||||
*/
|
||||
@ -31,6 +31,8 @@
|
||||
#include "modify.h"
|
||||
#include "fix.h"
|
||||
#include <float.h>
|
||||
#include "atom_masks.h"
|
||||
#include "neigh_request.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
using namespace MathConst;
|
||||
@ -50,7 +52,10 @@ using namespace MathSpecial;
|
||||
template<class DeviceType>
|
||||
PairExp6rxKokkos<DeviceType>::PairExp6rxKokkos(LAMMPS *lmp) : PairExp6rx(lmp)
|
||||
{
|
||||
|
||||
atomKK = (AtomKokkos *) atom;
|
||||
execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
|
||||
datamask_read = X_MASK | F_MASK | TYPE_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
@ -63,6 +68,39 @@ PairExp6rxKokkos<DeviceType>::~PairExp6rxKokkos()
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairExp6rxKokkos<DeviceType>::init_style()
|
||||
{
|
||||
PairExp6rxKokkos::init_style();
|
||||
|
||||
// irequest = neigh request made by parent class
|
||||
|
||||
neighflag = lmp->kokkos->neighflag;
|
||||
int irequest = neighbor->nrequest - 1;
|
||||
|
||||
neighbor->requests[irequest]->
|
||||
kokkos_host = Kokkos::Impl::is_same<DeviceType,LMPHostType>::value &&
|
||||
!Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
|
||||
neighbor->requests[irequest]->
|
||||
kokkos_device = Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
|
||||
|
||||
if (neighflag == FULL) {
|
||||
neighbor->requests[irequest]->full = 1;
|
||||
neighbor->requests[irequest]->half = 0;
|
||||
neighbor->requests[irequest]->full_cluster = 0;
|
||||
neighbor->requests[irequest]->ghost = 1;
|
||||
} else if (neighflag == HALF || neighflag == HALFTHREAD) {
|
||||
neighbor->requests[irequest]->full = 0;
|
||||
neighbor->requests[irequest]->half = 1;
|
||||
neighbor->requests[irequest]->full_cluster = 0;
|
||||
neighbor->requests[irequest]->ghost = 1;
|
||||
} else {
|
||||
error->all(FLERR,"Cannot use chosen neighbor list style with reax/c/kk");
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairExp6rxKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
|
||||
{
|
||||
@ -270,14 +308,14 @@ void PairExp6rxKokkos<DeviceType>::operator()(TagPairExp6rxCompute<NEIGHFLAG,NEW
|
||||
rsq = delx*delx + dely*dely + delz*delz;
|
||||
jtype = type[j];
|
||||
|
||||
if (rsq < cutsq[itype][jtype]) {
|
||||
if (rsq < d_cutsq(itype,jtype)) { // optimize
|
||||
r2inv = 1.0/rsq;
|
||||
r6inv = r2inv*r2inv*r2inv;
|
||||
|
||||
r = sqrt(rsq);
|
||||
rCut2inv = 1.0/cutsq[itype][jtype];
|
||||
rCut2inv = 1.0/d_cutsq(itype,jtype);
|
||||
rCut6inv = rCut2inv*rCut2inv*rCut2inv;
|
||||
rCut = sqrt(cutsq[itype][jtype]);
|
||||
rCut = sqrt(d_cutsq(itype,jtype));
|
||||
rCutInv = 1.0/rCut;
|
||||
|
||||
//
|
||||
|
||||
@ -63,7 +63,11 @@ class PairExp6rxKokkos : public PairExp6rx {
|
||||
|
||||
PairExp6rxKokkos(class LAMMPS *);
|
||||
virtual ~PairExp6rxKokkos();
|
||||
virtual void compute(int, int);
|
||||
void compute(int, int);
|
||||
void init_style();
|
||||
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairExp6rxgetParamsEXP6, const int&) const;
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR, int EVFLAG>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
@ -73,9 +77,6 @@ class PairExp6rxKokkos : public PairExp6rx {
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairExp6rxCompute<NEIGHFLAG,NEWTON_PAIR,EVFLAG>, const int&) const;
|
||||
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairExp6rxgetParamsEXP6, const int&) const;
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void ev_tally(EV_FLOAT &ev, const int &i, const int &j,
|
||||
|
||||
791
src/KOKKOS/pair_multi_lucy_rx_kokkos.cpp
Normal file
791
src/KOKKOS/pair_multi_lucy_rx_kokkos.cpp
Normal file
@ -0,0 +1,791 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------------------------
|
||||
Contributing authors:
|
||||
Stan Moore (Sandia)
|
||||
|
||||
Please cite the related publications:
|
||||
J.D. Moore, B.C. Barnes, S. Izvekov, M. Lisal, M.S. Sellers, D.E. Taylor & J.K. Brennan
|
||||
"A coarse-grain force field for RDX: Density dependent and energy conserving"
|
||||
The Journal of Chemical Physics, 2016, 144, 104501.
|
||||
------------------------------------------------------------------------------------------- */
|
||||
|
||||
#include <mpi.h>
|
||||
#include <math.h>
|
||||
#include "math_const.h"
|
||||
#include <stdlib.h>
|
||||
#include <string.h>
|
||||
#include "pair_multi_lucy_rx_kokkos.h"
|
||||
#include "atom_kokkos.h"
|
||||
#include "force.h"
|
||||
#include "comm.h"
|
||||
#include "neigh_list.h"
|
||||
#include "memory.h"
|
||||
#include "error.h"
|
||||
#include "citeme.h"
|
||||
#include "modify.h"
|
||||
#include "fix.h"
|
||||
#include "atom_masks.h"
|
||||
#include "neigh_request.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
enum{NONE,RLINEAR,RSQ};
|
||||
|
||||
#define MAXLINE 1024
|
||||
|
||||
#define oneFluidParameter (-1)
|
||||
#define isOneFluid(_site) ( (_site) == oneFluidParameter )
|
||||
|
||||
static const char cite_pair_multi_lucy_rx[] =
|
||||
"pair_style multi/lucy/rx command:\n\n"
|
||||
"@Article{Moore16,\n"
|
||||
" author = {J.D. Moore, B.C. Barnes, S. Izvekov, M. Lisal, M.S. Sellers, D.E. Taylor and J. K. Brennan},\n"
|
||||
" title = {A coarse-grain force field for RDX: Density dependent and energy conserving},\n"
|
||||
" journal = {J. Chem. Phys.},\n"
|
||||
" year = 2016,\n"
|
||||
" volume = 144\n"
|
||||
" pages = {104501}\n"
|
||||
"}\n\n";
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
PairMultiLucyRXKokkos<DeviceType>::PairMultiLucyRXKokkos(LAMMPS *lmp) : PairMultiLucyRX(lmp)
|
||||
{
|
||||
respa_enable = 0;
|
||||
|
||||
atomKK = (AtomKokkos *) atom;
|
||||
execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
|
||||
datamask_read = X_MASK | F_MASK | TYPE_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
PairMultiLucyRXKokkos<DeviceType>::~PairMultiLucyRXKokkos()
|
||||
{
|
||||
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairMultiLucyRXKokkos<DeviceType>::init_style()
|
||||
{
|
||||
PairMultiLucyRX::init_style();
|
||||
|
||||
// irequest = neigh request made by parent class
|
||||
|
||||
neighflag = lmp->kokkos->neighflag;
|
||||
int irequest = neighbor->nrequest - 1;
|
||||
|
||||
neighbor->requests[irequest]->
|
||||
kokkos_host = Kokkos::Impl::is_same<DeviceType,LMPHostType>::value &&
|
||||
!Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
|
||||
neighbor->requests[irequest]->
|
||||
kokkos_device = Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
|
||||
|
||||
if (neighflag == FULL) {
|
||||
neighbor->requests[irequest]->full = 1;
|
||||
neighbor->requests[irequest]->half = 0;
|
||||
neighbor->requests[irequest]->full_cluster = 0;
|
||||
neighbor->requests[irequest]->ghost = 1;
|
||||
} else if (neighflag == HALF || neighflag == HALFTHREAD) {
|
||||
neighbor->requests[irequest]->full = 0;
|
||||
neighbor->requests[irequest]->half = 1;
|
||||
neighbor->requests[irequest]->full_cluster = 0;
|
||||
neighbor->requests[irequest]->ghost = 1;
|
||||
} else {
|
||||
error->all(FLERR,"Cannot use chosen neighbor list style with reax/c/kk");
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairMultiLucyRXKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
|
||||
{
|
||||
eflag = eflag_in;
|
||||
vflag = vflag_in;
|
||||
|
||||
double evdwl,evdwlOld;
|
||||
|
||||
evdwlOld = 0.0;
|
||||
evdwl = 0.0;
|
||||
if (neighflag == FULL) no_virial_fdotr_compute = 1;
|
||||
if (eflag || vflag) ev_setup(eflag,vflag);
|
||||
else evflag = vflag_fdotr = 0;
|
||||
|
||||
// reallocate per-atom arrays if necessary
|
||||
|
||||
if (eflag_atom) {
|
||||
memory->destroy_kokkos(k_eatom,eatom);
|
||||
memory->create_kokkos(k_eatom,eatom,maxeatom,"pair:eatom");
|
||||
d_eatom = k_eatom.d_view;
|
||||
}
|
||||
if (vflag_atom) {
|
||||
memory->destroy_kokkos(k_vatom,vatom);
|
||||
memory->create_kokkos(k_vatom,vatom,maxvatom,6,"pair:vatom");
|
||||
d_vatom = k_vatom.d_view;
|
||||
}
|
||||
|
||||
x = atomKK->k_x.view<DeviceType>();
|
||||
f = atomKK->k_f.view<DeviceType>();
|
||||
type = atomKK->k_type.view<DeviceType>();
|
||||
uCG = atomKK->k_uCG.view<DeviceType>();
|
||||
uCGnew = atomKK->k_uCGnew.view<DeviceType>();
|
||||
dvector = atomKK->k_dvector.view<DeviceType>();
|
||||
rho = atomKK->k_rho.view<DeviceType>();
|
||||
|
||||
nlocal = atom->nlocal;
|
||||
int nghost = atom->nghost;
|
||||
int newton_pair = force->newton_pair;
|
||||
|
||||
{
|
||||
const int ntotal = nlocal + nghost;
|
||||
d_fractionOld1 = typename AT::t_float_1d("PairMultiLucyRX::fractionOld1",ntotal);
|
||||
d_fractionOld2 = typename AT::t_float_1d("PairMultiLucyRX::fractionOld2",ntotal);
|
||||
d_fraction1 = typename AT::t_float_1d("PairMultiLucyRX::fraction1",ntotal);
|
||||
d_fraction2 = typename AT::t_float_1d("PairMultiLucyRX::fraction2",ntotal);
|
||||
|
||||
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairMultiLucyRXgetParams>(0,ntotal),*this);
|
||||
}
|
||||
|
||||
const int inum = list->inum;
|
||||
NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list);
|
||||
d_numneigh = k_list->d_numneigh;
|
||||
d_neighbors = k_list->d_neighbors;
|
||||
d_ilist = k_list->d_ilist;
|
||||
|
||||
computeLocalDensity();
|
||||
|
||||
// loop over neighbors of my atoms
|
||||
|
||||
EV_FLOAT ev;
|
||||
|
||||
if (evflag) {
|
||||
Kokkos::parallel_reduce(Kokkos::RangePolicy<DeviceType, TagPairMultiLucyRXCompute<HALF,1,1> >(0,inum),*this,ev);
|
||||
} else {
|
||||
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairMultiLucyRXCompute<HALF,1,0> >(0,inum),*this);
|
||||
}
|
||||
|
||||
if (eflag_global) eng_vdwl += ev.evdwl;
|
||||
if (vflag_global) {
|
||||
virial[0] += ev.v[0];
|
||||
virial[1] += ev.v[1];
|
||||
virial[2] += ev.v[2];
|
||||
virial[3] += ev.v[3];
|
||||
virial[4] += ev.v[4];
|
||||
virial[5] += ev.v[5];
|
||||
}
|
||||
|
||||
if (vflag_fdotr) pair_virial_fdotr_compute(this);
|
||||
|
||||
if (eflag_atom) {
|
||||
k_eatom.template modify<DeviceType>();
|
||||
k_eatom.template sync<LMPHostType>();
|
||||
}
|
||||
|
||||
if (vflag_atom) {
|
||||
k_vatom.template modify<DeviceType>();
|
||||
k_vatom.template sync<LMPHostType>();
|
||||
}
|
||||
|
||||
copymode = 0;
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairMultiLucyRXKokkos<DeviceType>::operator()(TagPairMultiLucyRXgetParams, const int &i) const {
|
||||
getParams(i, d_fractionOld1[i], d_fractionOld2[i], d_fraction1[i], d_fraction2[i]);
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR, int EVFLAG>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairMultiLucyRXKokkos<DeviceType>::operator()(TagPairMultiLucyRXCompute<NEIGHFLAG,NEWTON_PAIR,EVFLAG>, const int &ii, EV_FLOAT& ev) const {
|
||||
int i,j,jj,inum,jnum,itype,jtype,itable;
|
||||
double xtmp,ytmp,ztmp,delx,dely,delz,evdwl,evdwlOld,fpair;
|
||||
double rsq;
|
||||
|
||||
double fractionOld1_i,fractionOld1_j;
|
||||
double fractionOld2_i,fractionOld2_j;
|
||||
double fraction1_i;
|
||||
|
||||
double pi = MathConst::MY_PI;
|
||||
double A_i, A_j;
|
||||
double fraction_i,fraction_j;
|
||||
int jtable;
|
||||
|
||||
Table *tb;
|
||||
|
||||
int tlm1 = tablength - 1;
|
||||
|
||||
i = d_ilist[ii];
|
||||
xtmp = x(i,0);
|
||||
ytmp = x(i,1);
|
||||
ztmp = x(i,2);
|
||||
itype = type[i];
|
||||
jnum = d_numneigh[i];
|
||||
|
||||
double fx_i = 0.0;
|
||||
double fy_i = 0.0;
|
||||
double fz_i = 0.0;
|
||||
|
||||
fractionOld1_i = d_fractionOld1[i];
|
||||
fractionOld2_i = d_fractionOld2[i];
|
||||
fraction1_i = d_fraction1[i];
|
||||
|
||||
for (jj = 0; jj < jnum; jj++) {
|
||||
int j = d_neighbors(i,jj);
|
||||
j &= NEIGHMASK;
|
||||
|
||||
delx = xtmp - x(j,0);
|
||||
dely = ytmp - x(j,1);
|
||||
delz = ztmp - x(j,2);
|
||||
rsq = delx*delx + dely*dely + delz*delz;
|
||||
jtype = type[j];
|
||||
|
||||
if (rsq < d_cutsq(itype,jtype)) { // optimize
|
||||
fpair = 0.0;
|
||||
|
||||
fractionOld1_j = d_fractionOld1[j];
|
||||
fractionOld2_j = d_fractionOld2[j];
|
||||
|
||||
tb = &tables[tabindex[itype][jtype]];
|
||||
if (rho[i]*rho[i] < tb->innersq || rho[j]*rho[j] < tb->innersq){
|
||||
//printf("Table inner cutoff = %lf\n",sqrt(tb->innersq));
|
||||
//printf("rho[%d]=%lf\n",i,rho[i]);
|
||||
//printf("rho[%d]=%lf\n",j,rho[j]);
|
||||
error->one(FLERR,"Density < table inner cutoff");
|
||||
}
|
||||
if (tabstyle == LOOKUP) {
|
||||
itable = static_cast<int> (((rho[i]*rho[i]) - tb->innersq) * tb->invdelta);
|
||||
jtable = static_cast<int> (((rho[j]*rho[j]) - tb->innersq) * tb->invdelta);
|
||||
if (itable >= tlm1 || jtable >= tlm1){
|
||||
//printf("Table outer index = %d\n",tlm1);
|
||||
//printf("itableIndex=%d rho[%d]=%lf\n",itable,i,rho[i]);
|
||||
//printf("jtableIndex=%d rho[%d]=%lf\n",jtable,j,rho[j]);
|
||||
error->one(FLERR,"Density > table outer cutoff");
|
||||
}
|
||||
A_i = tb->f[itable];
|
||||
A_j = tb->f[jtable];
|
||||
|
||||
const double rfactor = 1.0-sqrt(rsq/d_cutsq(itype,jtype));
|
||||
fpair = 0.5*(A_i + A_j)*(4.0-3.0*rfactor)*rfactor*rfactor*rfactor;
|
||||
fpair /= sqrt(rsq);
|
||||
|
||||
} else if (tabstyle == LINEAR) {
|
||||
itable = static_cast<int> ((rho[i]*rho[i] - tb->innersq) * tb->invdelta);
|
||||
jtable = static_cast<int> (((rho[j]*rho[j]) - tb->innersq) * tb->invdelta);
|
||||
if (itable >= tlm1 || jtable >= tlm1){
|
||||
//printf("Table outer index = %d\n",tlm1);
|
||||
//printf("itableIndex=%d rho[%d]=%lf\n",itable,i,rho[i]);
|
||||
//printf("jtableIndex=%d rho[%d]=%lf\n",jtable,j,rho[j]);
|
||||
error->one(FLERR,"Density > table outer cutoff");
|
||||
}
|
||||
if(itable<0) itable=0;
|
||||
if(itable>=tlm1) itable=tlm1;
|
||||
if(jtable<0) jtable=0;
|
||||
if(jtable>=tlm1)jtable=tlm1;
|
||||
|
||||
fraction_i = (((rho[i]*rho[i]) - tb->rsq[itable]) * tb->invdelta);
|
||||
fraction_j = (((rho[j]*rho[j]) - tb->rsq[jtable]) * tb->invdelta);
|
||||
if(itable==0) fraction_i=0.0;
|
||||
if(itable==tlm1) fraction_i=0.0;
|
||||
if(jtable==0) fraction_j=0.0;
|
||||
if(jtable==tlm1) fraction_j=0.0;
|
||||
|
||||
A_i = tb->f[itable] + fraction_i*tb->df[itable];
|
||||
A_j = tb->f[jtable] + fraction_j*tb->df[jtable];
|
||||
|
||||
const double rfactor = 1.0-sqrt(rsq/d_cutsq(itype,jtype));
|
||||
fpair = 0.5*(A_i + A_j)*(4.0-3.0*rfactor)*rfactor*rfactor*rfactor;
|
||||
fpair /= sqrt(rsq);
|
||||
|
||||
} else error->one(FLERR,"Only LOOKUP and LINEAR table styles have been implemented for pair multi/lucy/rx");
|
||||
|
||||
if (isite1 == isite2) fpair = sqrt(fractionOld1_i*fractionOld2_j)*fpair;
|
||||
else fpair = (sqrt(fractionOld1_i*fractionOld2_j) + sqrt(fractionOld2_i*fractionOld1_j))*fpair;
|
||||
|
||||
fx_i += delx*fpair;
|
||||
fy_i += dely*fpair;
|
||||
fz_i += delz*fpair;
|
||||
if (NEWTON_PAIR || j < nlocal) {
|
||||
f(j,0) -= delx*fpair;
|
||||
f(j,1) -= dely*fpair;
|
||||
f(j,2) -= delz*fpair;
|
||||
}
|
||||
//if (evflag) ev_tally(i,j,nlocal,newton_pair,0.0,0.0,fpair,delx,dely,delz);
|
||||
if (EVFLAG) this->template ev_tally<NEIGHFLAG,NEWTON_PAIR>(ev,i,j,0.0,fpair,delx,dely,delz);
|
||||
}
|
||||
}
|
||||
|
||||
f(i,0) += fx_i;
|
||||
f(i,1) += fy_i;
|
||||
f(i,2) += fz_i;
|
||||
|
||||
tb = &tables[tabindex[itype][itype]];
|
||||
itable = static_cast<int> (((rho[i]*rho[i]) - tb->innersq) * tb->invdelta);
|
||||
if (tabstyle == LOOKUP) evdwl = tb->e[itable];
|
||||
else if (tabstyle == LINEAR){
|
||||
if (itable >= tlm1){
|
||||
//printf("itableIndex=%d rho[%d]=%lf\n",itable,i,rho[i]);
|
||||
error->one(FLERR,"Density > table outer cutoff");
|
||||
}
|
||||
if(itable==0) fraction_i=0.0;
|
||||
else fraction_i = (((rho[i]*rho[i]) - tb->rsq[itable]) * tb->invdelta);
|
||||
evdwl = tb->e[itable] + fraction_i*tb->de[itable];
|
||||
} else error->one(FLERR,"Only LOOKUP and LINEAR table styles have been implemented for pair multi/lucy/rx");
|
||||
|
||||
evdwl *=(pi*d_cutsq(itype,itype)*d_cutsq(itype,itype))/84.0;
|
||||
evdwlOld = fractionOld1_i*evdwl;
|
||||
evdwl = fraction1_i*evdwl;
|
||||
|
||||
uCG[i] += evdwlOld;
|
||||
uCGnew[i] += evdwl;
|
||||
|
||||
evdwl = evdwlOld;
|
||||
|
||||
//if (evflag) ev_tally(0,0,nlocal,newton_pair,evdwl,0.0,0.0,0.0,0.0,0.0);
|
||||
if (EVFLAG) ev.evdwl += evdwl;
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR, int EVFLAG>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairMultiLucyRXKokkos<DeviceType>::operator()(TagPairMultiLucyRXCompute<NEIGHFLAG,NEWTON_PAIR,EVFLAG>, const int &ii) const {
|
||||
EV_FLOAT ev;
|
||||
this->template operator()<NEIGHFLAG,NEWTON_PAIR,EVFLAG>(TagPairMultiLucyRXCompute<NEIGHFLAG,NEWTON_PAIR,EVFLAG>(), ii, ev);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
set coeffs for one or more type pairs
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairMultiLucyRXKokkos<DeviceType>::coeff(int narg, char **arg)
|
||||
{
|
||||
if (narg != 6 && narg != 7) error->all(FLERR,"Illegal pair_coeff command");
|
||||
|
||||
bool rx_flag = false;
|
||||
for (int i = 0; i < modify->nfix; i++)
|
||||
if (strncmp(modify->fix[i]->style,"rx",2) == 0) rx_flag = true;
|
||||
if (!rx_flag) error->all(FLERR,"PairMultiLucyRXKokkos<DeviceType> requires a fix rx command.");
|
||||
|
||||
if (!allocated) allocate();
|
||||
|
||||
int ilo,ihi,jlo,jhi;
|
||||
force->bounds(FLERR,arg[0],atom->ntypes,ilo,ihi);
|
||||
force->bounds(FLERR,arg[1],atom->ntypes,jlo,jhi);
|
||||
|
||||
int me;
|
||||
MPI_Comm_rank(world,&me);
|
||||
tables = (Table *)
|
||||
memory->srealloc(tables,(ntables+1)*sizeof(Table),"pair:tables");
|
||||
Table *tb = &tables[ntables];
|
||||
null_table(tb);
|
||||
if (me == 0) read_table(tb,arg[2],arg[3]);
|
||||
bcast_table(tb);
|
||||
|
||||
nspecies = atom->nspecies_dpd;
|
||||
int n;
|
||||
n = strlen(arg[3]) + 1;
|
||||
site1 = new char[n];
|
||||
strcpy(site1,arg[4]);
|
||||
|
||||
n = strlen(arg[4]) + 1;
|
||||
site2 = new char[n];
|
||||
strcpy(site2,arg[5]);
|
||||
|
||||
// set table cutoff
|
||||
|
||||
if (narg == 7) tb->cut = force->numeric(FLERR,arg[6]);
|
||||
else if (tb->rflag) tb->cut = tb->rhi;
|
||||
else tb->cut = tb->rfile[tb->ninput-1];
|
||||
|
||||
// error check on table parameters
|
||||
// insure cutoff is within table
|
||||
|
||||
if (tb->ninput <= 1) error->one(FLERR,"Invalid pair table length");
|
||||
if (tb->rflag == 0) {
|
||||
rho_0 = tb->rfile[0];
|
||||
} else {
|
||||
rho_0 = tb->rlo;
|
||||
}
|
||||
|
||||
tb->match = 0;
|
||||
if (tabstyle == LINEAR && tb->ninput == tablength &&
|
||||
tb->rflag == RSQ) tb->match = 1;
|
||||
|
||||
// spline read-in values and compute r,e,f vectors within table
|
||||
|
||||
if (tb->match == 0) spline_table(tb);
|
||||
compute_table(tb);
|
||||
|
||||
// store ptr to table in tabindex
|
||||
|
||||
int count = 0;
|
||||
for (int i = ilo; i <= ihi; i++) {
|
||||
for (int j = MAX(jlo,i); j <= jhi; j++) {
|
||||
tabindex[i][j] = ntables;
|
||||
setflag[i][j] = 1;
|
||||
count++;
|
||||
}
|
||||
}
|
||||
|
||||
if (count == 0) error->all(FLERR,"Illegal pair_coeff command");
|
||||
ntables++;
|
||||
|
||||
// Match site* to isite values.
|
||||
|
||||
if (strcmp(site1, "1fluid") == 0)
|
||||
isite1 = oneFluidParameter;
|
||||
else {
|
||||
isite1 = nspecies;
|
||||
for (int ispecies = 0; ispecies < nspecies; ++ispecies)
|
||||
if (strcmp(site1, atom->dname[ispecies]) == 0){
|
||||
isite1 = ispecies;
|
||||
break;
|
||||
}
|
||||
|
||||
if (isite1 == nspecies)
|
||||
error->all(FLERR,"Pair_multi_lucy_rx site1 is invalid.");
|
||||
}
|
||||
|
||||
if (strcmp(site2, "1fluid") == 0)
|
||||
isite2 = oneFluidParameter;
|
||||
else {
|
||||
isite2 = nspecies;
|
||||
for (int ispecies = 0; ispecies < nspecies; ++ispecies)
|
||||
if (strcmp(site2, atom->dname[ispecies]) == 0){
|
||||
isite2 = ispecies;
|
||||
break;
|
||||
}
|
||||
|
||||
if (isite2 == nspecies)
|
||||
error->all(FLERR,"Pair_multi_lucy_rx site2 is invalid.");
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairMultiLucyRXKokkos<DeviceType>::computeLocalDensity()
|
||||
{
|
||||
x = atomKK->k_x.view<DeviceType>();
|
||||
type = atomKK->k_type.view<DeviceType>();
|
||||
rho = atomKK->k_rho.view<DeviceType>();
|
||||
nlocal = atom->nlocal;
|
||||
|
||||
//sync
|
||||
|
||||
const int inum = list->inum;
|
||||
NeighListKokkos<DeviceType>* k_list = static_cast<NeighListKokkos<DeviceType>*>(list);
|
||||
d_numneigh = k_list->d_numneigh;
|
||||
d_neighbors = k_list->d_neighbors;
|
||||
d_ilist = k_list->d_ilist;
|
||||
|
||||
const double pi = MathConst::MY_PI;
|
||||
|
||||
const bool newton_pair = force->newton_pair;
|
||||
one_type = (atom->ntypes == 1);
|
||||
|
||||
// Special cut-off values for when there's only one type.
|
||||
cutsq_type11 = cutsq[1][1];
|
||||
rcut_type11 = sqrt(cutsq_type11);
|
||||
factor_type11 = 84.0/(5.0*pi*rcut_type11*rcut_type11*rcut_type11);
|
||||
|
||||
// zero out density
|
||||
int m = nlocal;
|
||||
if (newton_pair) m += atom->nghost;
|
||||
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairMultiLucyRXZero>(0,m),*this);
|
||||
|
||||
// rho = density at each atom
|
||||
// loop over neighbors of my atoms
|
||||
if (newton_pair)
|
||||
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairMultiLucyRXComputeLocalDensity<HALF,1> >(0,inum),*this);
|
||||
else
|
||||
Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType, TagPairMultiLucyRXComputeLocalDensity<HALF,0> >(0,inum),*this);
|
||||
|
||||
if (newton_pair) comm->reverse_comm_pair(this);
|
||||
|
||||
comm->forward_comm_pair(this);
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairMultiLucyRXKokkos<DeviceType>::operator()(TagPairMultiLucyRXZero, const int &i) const {
|
||||
rho[i] = 0.0;
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairMultiLucyRXKokkos<DeviceType>::operator()(TagPairMultiLucyRXComputeLocalDensity<NEIGHFLAG,NEWTON_PAIR>, const int &ii) const {
|
||||
const int i = d_ilist[ii];
|
||||
|
||||
const double xtmp = x(i,0);
|
||||
const double ytmp = x(i,1);
|
||||
const double ztmp = x(i,2);
|
||||
|
||||
double rho_i = rho[i];
|
||||
|
||||
const int itype = type[i];
|
||||
const int jnum = d_numneigh[i];
|
||||
|
||||
const double pi = MathConst::MY_PI;
|
||||
|
||||
for (int jj = 0; jj < jnum; jj++){
|
||||
const int j = (d_neighbors(i,jj) & NEIGHMASK);
|
||||
const int jtype = type[j];
|
||||
|
||||
const double delx = xtmp - x(j,0);
|
||||
const double dely = ytmp - x(j,1);
|
||||
const double delz = ztmp - x(j,2);
|
||||
const double rsq = delx*delx + dely*dely + delz*delz;
|
||||
|
||||
if (one_type) {
|
||||
if (rsq < cutsq_type11) {
|
||||
const double rcut = rcut_type11;
|
||||
const double r_over_rcut = sqrt(rsq) / rcut;
|
||||
const double tmpFactor = 1.0 - r_over_rcut;
|
||||
const double tmpFactor4 = tmpFactor*tmpFactor*tmpFactor*tmpFactor;
|
||||
const double factor = factor_type11*(1.0 + 1.5*r_over_rcut)*tmpFactor4;
|
||||
rho_i += factor;
|
||||
if (NEWTON_PAIR || j < nlocal)
|
||||
rho[j] += factor;
|
||||
} else if (rsq < d_cutsq(itype,jtype)) {
|
||||
const double rcut = sqrt(d_cutsq(itype,jtype));
|
||||
const double tmpFactor = 1.0-sqrt(rsq)/rcut;
|
||||
const double tmpFactor4 = tmpFactor*tmpFactor*tmpFactor*tmpFactor;
|
||||
const double factor = (84.0/(5.0*pi*rcut*rcut*rcut))*(1.0+3.0*sqrt(rsq)/(2.0*rcut))*tmpFactor4;
|
||||
rho_i += factor;
|
||||
if (NEWTON_PAIR || j < nlocal)
|
||||
rho[j] += factor;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
rho[i] = rho_i;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairMultiLucyRXKokkos<DeviceType>::getParams(int id, double &fractionOld1, double &fractionOld2, double &fraction1, double &fraction2) const
|
||||
{
|
||||
double fractionOld, fraction;
|
||||
double nTotal, nTotalOld;
|
||||
|
||||
nTotal = 0.0;
|
||||
nTotalOld = 0.0;
|
||||
for (int ispecies = 0; ispecies < nspecies; ispecies++){
|
||||
nTotal += dvector(ispecies,id);
|
||||
nTotalOld += dvector(ispecies+nspecies,id);
|
||||
}
|
||||
|
||||
if (isOneFluid(isite1) == false){
|
||||
fractionOld1 = dvector(isite1+nspecies,id)/nTotalOld;
|
||||
fraction1 = dvector(isite1,id)/nTotal;
|
||||
}
|
||||
if (isOneFluid(isite2) == false){
|
||||
fractionOld2 = dvector(isite2+nspecies,id)/nTotalOld;
|
||||
fraction2 = dvector(isite2,id)/nTotal;
|
||||
}
|
||||
|
||||
if (isOneFluid(isite1) || isOneFluid(isite2)){
|
||||
fractionOld = 0.0;
|
||||
fraction = 0.0;
|
||||
|
||||
for (int ispecies = 0; ispecies < nspecies; ispecies++){
|
||||
if (isite1 == ispecies || isite2 == ispecies) continue;
|
||||
fractionOld += dvector(ispecies+nspecies,id) / nTotalOld;
|
||||
fraction += dvector(ispecies,id) / nTotal;
|
||||
}
|
||||
if (isOneFluid(isite1)){
|
||||
fractionOld1 = fractionOld;
|
||||
fraction1 = fraction;
|
||||
}
|
||||
if (isOneFluid(isite2)){
|
||||
fractionOld2 = fractionOld;
|
||||
fraction2 = fraction;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
int PairMultiLucyRXKokkos<DeviceType>::pack_forward_comm(int n, int *list, double *buf, int pbc_flag, int *pbc)
|
||||
{
|
||||
int i,j,m;
|
||||
rho = atomKK->k_rho.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
for (i = 0; i < n; i++) {
|
||||
j = list[i];
|
||||
buf[m++] = rho[j];
|
||||
}
|
||||
return m;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairMultiLucyRXKokkos<DeviceType>::unpack_forward_comm(int n, int first, double *buf)
|
||||
{
|
||||
int i,m,last;
|
||||
rho = atomKK->k_rho.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
last = first + n;
|
||||
for (i = first; i < last; i++) rho[i] = buf[m++];
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
int PairMultiLucyRXKokkos<DeviceType>::pack_reverse_comm(int n, int first, double *buf)
|
||||
{
|
||||
int i,m,last;
|
||||
rho = atomKK->k_rho.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
last = first + n;
|
||||
for (i = first; i < last; i++) buf[m++] = rho[i];
|
||||
return m;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairMultiLucyRXKokkos<DeviceType>::unpack_reverse_comm(int n, int *list, double *buf)
|
||||
{
|
||||
int i,j,m;
|
||||
rho = atomKK->k_rho.view<DeviceType>();
|
||||
|
||||
m = 0;
|
||||
for (i = 0; i < n; i++) {
|
||||
j = list[i];
|
||||
rho[j] += buf[m++];
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairMultiLucyRXKokkos<DeviceType>::ev_tally(EV_FLOAT &ev, const int &i, const int &j,
|
||||
const F_FLOAT &epair, const F_FLOAT &fpair, const F_FLOAT &delx,
|
||||
const F_FLOAT &dely, const F_FLOAT &delz) const
|
||||
{
|
||||
const int EFLAG = eflag;
|
||||
const int VFLAG = vflag_either;
|
||||
|
||||
// The eatom and vatom arrays are atomic for Half/Thread neighbor style
|
||||
Kokkos::View<E_FLOAT*, typename DAT::t_efloat_1d::array_layout,DeviceType,Kokkos::MemoryTraits<AtomicF<NEIGHFLAG>::value> > v_eatom = k_eatom.view<DeviceType>();
|
||||
Kokkos::View<F_FLOAT*[6], typename DAT::t_virial_array::array_layout,DeviceType,Kokkos::MemoryTraits<AtomicF<NEIGHFLAG>::value> > v_vatom = k_vatom.view<DeviceType>();
|
||||
|
||||
if (EFLAG) {
|
||||
if (eflag_atom) {
|
||||
const E_FLOAT epairhalf = 0.5 * epair;
|
||||
if (NEIGHFLAG!=FULL) {
|
||||
if (NEWTON_PAIR || i < nlocal) v_eatom[i] += epairhalf;
|
||||
if (NEWTON_PAIR || j < nlocal) v_eatom[j] += epairhalf;
|
||||
} else {
|
||||
v_eatom[i] += epairhalf;
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
if (VFLAG) {
|
||||
const E_FLOAT v0 = delx*delx*fpair;
|
||||
const E_FLOAT v1 = dely*dely*fpair;
|
||||
const E_FLOAT v2 = delz*delz*fpair;
|
||||
const E_FLOAT v3 = delx*dely*fpair;
|
||||
const E_FLOAT v4 = delx*delz*fpair;
|
||||
const E_FLOAT v5 = dely*delz*fpair;
|
||||
|
||||
if (vflag_global) {
|
||||
if (NEIGHFLAG!=FULL) {
|
||||
if (NEWTON_PAIR || i < nlocal) {
|
||||
ev.v[0] += 0.5*v0;
|
||||
ev.v[1] += 0.5*v1;
|
||||
ev.v[2] += 0.5*v2;
|
||||
ev.v[3] += 0.5*v3;
|
||||
ev.v[4] += 0.5*v4;
|
||||
ev.v[5] += 0.5*v5;
|
||||
}
|
||||
if (NEWTON_PAIR || j < nlocal) {
|
||||
ev.v[0] += 0.5*v0;
|
||||
ev.v[1] += 0.5*v1;
|
||||
ev.v[2] += 0.5*v2;
|
||||
ev.v[3] += 0.5*v3;
|
||||
ev.v[4] += 0.5*v4;
|
||||
ev.v[5] += 0.5*v5;
|
||||
}
|
||||
} else {
|
||||
ev.v[0] += 0.5*v0;
|
||||
ev.v[1] += 0.5*v1;
|
||||
ev.v[2] += 0.5*v2;
|
||||
ev.v[3] += 0.5*v3;
|
||||
ev.v[4] += 0.5*v4;
|
||||
ev.v[5] += 0.5*v5;
|
||||
}
|
||||
}
|
||||
|
||||
if (vflag_atom) {
|
||||
if (NEIGHFLAG!=FULL) {
|
||||
if (NEWTON_PAIR || i < nlocal) {
|
||||
v_vatom(i,0) += 0.5*v0;
|
||||
v_vatom(i,1) += 0.5*v1;
|
||||
v_vatom(i,2) += 0.5*v2;
|
||||
v_vatom(i,3) += 0.5*v3;
|
||||
v_vatom(i,4) += 0.5*v4;
|
||||
v_vatom(i,5) += 0.5*v5;
|
||||
}
|
||||
if (NEWTON_PAIR || j < nlocal) {
|
||||
v_vatom(j,0) += 0.5*v0;
|
||||
v_vatom(j,1) += 0.5*v1;
|
||||
v_vatom(j,2) += 0.5*v2;
|
||||
v_vatom(j,3) += 0.5*v3;
|
||||
v_vatom(j,4) += 0.5*v4;
|
||||
v_vatom(j,5) += 0.5*v5;
|
||||
}
|
||||
} else {
|
||||
v_vatom(i,0) += 0.5*v0;
|
||||
v_vatom(i,1) += 0.5*v1;
|
||||
v_vatom(i,2) += 0.5*v2;
|
||||
v_vatom(i,3) += 0.5*v3;
|
||||
v_vatom(i,4) += 0.5*v4;
|
||||
v_vatom(i,5) += 0.5*v5;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
template class PairMultiLucyRXKokkos<LMPDeviceType>;
|
||||
#ifdef KOKKOS_HAVE_CUDA
|
||||
template class PairMultiLucyRXKokkos<LMPHostType>;
|
||||
#endif
|
||||
}
|
||||
215
src/KOKKOS/pair_multi_lucy_rx_kokkos.h
Normal file
215
src/KOKKOS/pair_multi_lucy_rx_kokkos.h
Normal file
@ -0,0 +1,215 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef PAIR_CLASS
|
||||
|
||||
PairStyle(multi/lucy/rx/kk,PairMultiLucyRXKokkos<LMPDeviceType>)
|
||||
PairStyle(multi/lucy/rx/kk/device,PairMultiLucyRXKokkos<LMPDeviceType>)
|
||||
PairStyle(multi/lucy/rx/kk/host,PairMultiLucyRXKokkos<LMPHostType>)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_PAIR_MULTI_LUCY_RX_KOKKOS_H
|
||||
#define LMP_PAIR_MULTI_LUCY_RX_KOKKOS_H
|
||||
|
||||
|
||||
#include "pair_multi_lucy_rx.h"
|
||||
#include "pair_kokkos.h"
|
||||
#include "kokkos_type.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
struct TagPairMultiLucyRXgetParams{};
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR, int EVFLAG>
|
||||
struct TagPairMultiLucyRXCompute{};
|
||||
|
||||
struct TagPairMultiLucyRXZero{};
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR>
|
||||
struct TagPairMultiLucyRXComputeLocalDensity{};
|
||||
|
||||
template<class DeviceType>
|
||||
class PairMultiLucyRXKokkos : public PairMultiLucyRX {
|
||||
public:
|
||||
typedef DeviceType device_type;
|
||||
typedef ArrayTypes<DeviceType> AT;
|
||||
typedef EV_FLOAT value_type;
|
||||
|
||||
PairMultiLucyRXKokkos(class LAMMPS *);
|
||||
virtual ~PairMultiLucyRXKokkos();
|
||||
|
||||
void compute(int, int);
|
||||
void init_style();
|
||||
void coeff(int, char **);
|
||||
int pack_forward_comm(int, int *, double *, int, int *);
|
||||
void unpack_forward_comm(int, int, double *);
|
||||
int pack_reverse_comm(int, int, double *);
|
||||
void unpack_reverse_comm(int, int *, double *);
|
||||
void computeLocalDensity();
|
||||
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairMultiLucyRXgetParams, const int&) const;
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR, int EVFLAG>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairMultiLucyRXCompute<NEIGHFLAG,NEWTON_PAIR,EVFLAG>, const int&, EV_FLOAT&) const;
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR, int EVFLAG>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairMultiLucyRXCompute<NEIGHFLAG,NEWTON_PAIR,EVFLAG>, const int&) const;
|
||||
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairMultiLucyRXZero, const int&) const;
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void operator()(TagPairMultiLucyRXComputeLocalDensity<NEIGHFLAG,NEWTON_PAIR>, const int&) const;
|
||||
|
||||
template<int NEIGHFLAG, int NEWTON_PAIR>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void ev_tally(EV_FLOAT &ev, const int &i, const int &j,
|
||||
const F_FLOAT &epair, const F_FLOAT &fpair, const F_FLOAT &delx,
|
||||
const F_FLOAT &dely, const F_FLOAT &delz) const;
|
||||
|
||||
private:
|
||||
int nlocal;
|
||||
int neighflag;
|
||||
int eflag,vflag;
|
||||
|
||||
bool one_type;
|
||||
double cutsq_type11;
|
||||
double rcut_type11;
|
||||
double factor_type11;
|
||||
|
||||
//struct Table {
|
||||
// int ninput,rflag,fpflag,match;
|
||||
// double rlo,rhi,fplo,fphi,cut;
|
||||
// double *rfile,*efile,*ffile;
|
||||
// double *e2file,*f2file;
|
||||
// double innersq,delta,invdelta,deltasq6;
|
||||
// double *rsq,*drsq,*e,*de,*f,*df,*e2,*f2;
|
||||
//};
|
||||
//Table *tables;
|
||||
|
||||
int **tabindex;
|
||||
|
||||
//void read_table(Table *, char *, char *);
|
||||
//void param_extract(Table *, char *);
|
||||
|
||||
char *site1, *site2;
|
||||
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void getParams(int, double &, double &, double &, double &) const;
|
||||
|
||||
typename AT::t_float_1d d_fractionOld1,d_fractionOld2,d_fraction1,d_fraction2;
|
||||
|
||||
typename AT::t_x_array_randomread x;
|
||||
typename AT::t_f_array f;
|
||||
typename AT::t_int_1d_randomread type;
|
||||
typename AT::t_efloat_1d rho;
|
||||
typename AT::t_efloat_1d uCG, uCGnew;
|
||||
typename AT::t_float_2d dvector;
|
||||
|
||||
DAT::tdual_efloat_1d k_eatom;
|
||||
DAT::tdual_virial_array k_vatom;
|
||||
DAT::t_efloat_1d d_eatom;
|
||||
DAT::t_virial_array d_vatom;
|
||||
|
||||
typename AT::t_neighbors_2d d_neighbors;
|
||||
typename AT::t_int_1d_randomread d_ilist;
|
||||
typename AT::t_int_1d_randomread d_numneigh;
|
||||
|
||||
typename AT::tdual_ffloat_2d k_cutsq;
|
||||
typename AT::t_ffloat_2d d_cutsq;
|
||||
|
||||
friend void pair_virial_fdotr_compute<PairMultiLucyRXKokkos>(PairMultiLucyRXKokkos*);
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Pair multi/lucy/rx command requires atom_style with density (e.g. dpd, meso)
|
||||
|
||||
Self-explanatory
|
||||
|
||||
E: Density < table inner cutoff
|
||||
|
||||
The local density inner is smaller than the inner cutoff
|
||||
|
||||
E: Density > table inner cutoff
|
||||
|
||||
The local density inner is greater than the inner cutoff
|
||||
|
||||
E: Only LOOKUP and LINEAR table styles have been implemented for pair multi/lucy/rx
|
||||
|
||||
Self-explanatory
|
||||
|
||||
E: Illegal ... command
|
||||
|
||||
Self-explanatory. Check the input script syntax and compare to the
|
||||
documentation for the command. You can use -echo screen as a
|
||||
command-line option when running LAMMPS to see the offending line.
|
||||
|
||||
E: Unknown table style in pair_style command
|
||||
|
||||
Self-explanatory
|
||||
|
||||
E: Illegal number of pair table entries
|
||||
|
||||
There must be at least 2 table entries.
|
||||
|
||||
E: Illegal pair_coeff command
|
||||
|
||||
All pair coefficients must be set in the data file or by the
|
||||
pair_coeff command before running a simulation.
|
||||
|
||||
E: PairMultiLucyRXKokkos requires a fix rx command
|
||||
|
||||
The fix rx command must come before the pair style command in the input file
|
||||
|
||||
E: There are no rx species specified
|
||||
|
||||
There must be at least one species specified through the fix rx command
|
||||
|
||||
E: Invalid pair table length
|
||||
|
||||
Length of read-in pair table is invalid
|
||||
|
||||
E: All pair coeffs are not set
|
||||
|
||||
All pair coefficients must be set in the data file or by the
|
||||
pair_coeff command before running a simulation.
|
||||
|
||||
E: Cannot open file %s
|
||||
|
||||
The specified file cannot be opened. Check that the path and name are
|
||||
correct.
|
||||
|
||||
E: Did not find keyword in table file
|
||||
|
||||
Keyword used in pair_coeff command was not found in table file.
|
||||
|
||||
E: Invalid keyword in pair table parameters
|
||||
|
||||
Keyword used in list of table parameters is not recognized.
|
||||
|
||||
E: Pair table parameters did not set N
|
||||
|
||||
List of pair table parameters must include N setting.
|
||||
|
||||
*/
|
||||
@ -12,7 +12,7 @@
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Paul Crozier (SNL)
|
||||
Contributing author: Christian Trott (SNL)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <mpi.h>
|
||||
@ -41,7 +41,7 @@ enum{FULL,HALFTHREAD,HALF};
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
PairTableKokkos<DeviceType>::PairTableKokkos(LAMMPS *lmp) : Pair(lmp)
|
||||
PairTableKokkos<DeviceType>::PairTableKokkos(LAMMPS *lmp) : PairTable(lmp)
|
||||
{
|
||||
update_table = 0;
|
||||
atomKK = (AtomKokkos *) atom;
|
||||
@ -98,6 +98,7 @@ void PairTableKokkos<DeviceType>::compute_style(int eflag_in, int vflag_in)
|
||||
|
||||
if (neighflag == FULL || neighflag == FULLCLUSTER) no_virial_fdotr_compute = 1;
|
||||
|
||||
|
||||
if (eflag || vflag) ev_setup(eflag,vflag);
|
||||
else evflag = vflag_fdotr = 0;
|
||||
|
||||
@ -221,6 +222,7 @@ compute_fpair(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, c
|
||||
//if (rsq < d_table_const.innersq(tidx))
|
||||
// error->one(FLERR,"Pair distance < table inner cutoff");
|
||||
|
||||
|
||||
if (Specialisation::TabStyle == LOOKUP) {
|
||||
const int itable = static_cast<int> ((rsq - d_table_const.innersq(tidx)) * d_table_const.invdelta(tidx));
|
||||
//if (itable >= tlm1)
|
||||
@ -338,8 +340,6 @@ void PairTableKokkos<DeviceType>::create_kokkos_tables()
|
||||
memory->create_kokkos(d_table->drsq,h_table->drsq,ntables,ntable,"Table::drsq");
|
||||
}
|
||||
|
||||
|
||||
|
||||
for(int i=0; i < ntables; i++) {
|
||||
Table* tb = &tables[i];
|
||||
|
||||
@ -477,85 +477,6 @@ void PairTableKokkos<DeviceType>::settings(int narg, char **arg)
|
||||
tables = NULL;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
set coeffs for one or more type pairs
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::coeff(int narg, char **arg)
|
||||
{
|
||||
if (narg != 4 && narg != 5) error->all(FLERR,"Illegal pair_coeff command");
|
||||
if (!allocated) allocate();
|
||||
|
||||
int ilo,ihi,jlo,jhi;
|
||||
force->bounds(FLERR,arg[0],atom->ntypes,ilo,ihi);
|
||||
force->bounds(FLERR,arg[1],atom->ntypes,jlo,jhi);
|
||||
|
||||
int me;
|
||||
MPI_Comm_rank(world,&me);
|
||||
tables = (Table *)
|
||||
memory->srealloc(tables,(ntables+1)*sizeof(Table),"pair:tables");
|
||||
Table *tb = &tables[ntables];
|
||||
null_table(tb);
|
||||
if (me == 0) read_table(tb,arg[2],arg[3]);
|
||||
bcast_table(tb);
|
||||
|
||||
// set table cutoff
|
||||
|
||||
if (narg == 5) tb->cut = force->numeric(FLERR,arg[4]);
|
||||
else if (tb->rflag) tb->cut = tb->rhi;
|
||||
else tb->cut = tb->rfile[tb->ninput-1];
|
||||
|
||||
// error check on table parameters
|
||||
// insure cutoff is within table
|
||||
// for BITMAP tables, file values can be in non-ascending order
|
||||
|
||||
if (tb->ninput <= 1) error->one(FLERR,"Invalid pair table length");
|
||||
double rlo,rhi;
|
||||
if (tb->rflag == 0) {
|
||||
rlo = tb->rfile[0];
|
||||
rhi = tb->rfile[tb->ninput-1];
|
||||
} else {
|
||||
rlo = tb->rlo;
|
||||
rhi = tb->rhi;
|
||||
}
|
||||
if (tb->cut <= rlo || tb->cut > rhi)
|
||||
error->all(FLERR,"Invalid pair table cutoff");
|
||||
if (rlo <= 0.0) error->all(FLERR,"Invalid pair table cutoff");
|
||||
|
||||
// match = 1 if don't need to spline read-in tables
|
||||
// this is only the case if r values needed by final tables
|
||||
// exactly match r values read from file
|
||||
// for tabstyle SPLINE, always need to build spline tables
|
||||
|
||||
tb->match = 0;
|
||||
if (tabstyle == LINEAR && tb->ninput == tablength &&
|
||||
tb->rflag == RSQ && tb->rhi == tb->cut) tb->match = 1;
|
||||
if (tabstyle == BITMAP && tb->ninput == 1 << tablength &&
|
||||
tb->rflag == BMP && tb->rhi == tb->cut) tb->match = 1;
|
||||
if (tb->rflag == BMP && tb->match == 0)
|
||||
error->all(FLERR,"Bitmapped table in file does not match requested table");
|
||||
|
||||
// spline read-in values and compute r,e,f vectors within table
|
||||
|
||||
if (tb->match == 0) spline_table(tb);
|
||||
compute_table(tb);
|
||||
|
||||
// store ptr to table in tabindex
|
||||
|
||||
int count = 0;
|
||||
for (int i = ilo; i <= ihi; i++) {
|
||||
for (int j = MAX(jlo,i); j <= jhi; j++) {
|
||||
tabindex[i][j] = ntables;
|
||||
setflag[i][j] = 1;
|
||||
count++;
|
||||
}
|
||||
}
|
||||
|
||||
if (count == 0) error->all(FLERR,"Illegal pair_coeff command");
|
||||
ntables++;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
init for one type pair i,j and corresponding j,i
|
||||
------------------------------------------------------------------------- */
|
||||
@ -574,677 +495,6 @@ double PairTableKokkos<DeviceType>::init_one(int i, int j)
|
||||
return tables[tabindex[i][j]].cut;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
read a table section from a tabulated potential file
|
||||
only called by proc 0
|
||||
this function sets these values in Table:
|
||||
ninput,rfile,efile,ffile,rflag,rlo,rhi,fpflag,fplo,fphi,ntablebits
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::read_table(Table *tb, char *file, char *keyword)
|
||||
{
|
||||
char line[MAXLINE];
|
||||
|
||||
// open file
|
||||
|
||||
FILE *fp = force->open_potential(file);
|
||||
if (fp == NULL) {
|
||||
char str[128];
|
||||
sprintf(str,"Cannot open file %s",file);
|
||||
error->one(FLERR,str);
|
||||
}
|
||||
|
||||
// loop until section found with matching keyword
|
||||
|
||||
while (1) {
|
||||
if (fgets(line,MAXLINE,fp) == NULL)
|
||||
error->one(FLERR,"Did not find keyword in table file");
|
||||
if (strspn(line," \t\n\r") == strlen(line)) continue; // blank line
|
||||
if (line[0] == '#') continue; // comment
|
||||
char *word = strtok(line," \t\n\r");
|
||||
if (strcmp(word,keyword) == 0) break; // matching keyword
|
||||
fgets(line,MAXLINE,fp); // no match, skip section
|
||||
param_extract(tb,line);
|
||||
fgets(line,MAXLINE,fp);
|
||||
for (int i = 0; i < tb->ninput; i++) fgets(line,MAXLINE,fp);
|
||||
}
|
||||
|
||||
// read args on 2nd line of section
|
||||
// allocate table arrays for file values
|
||||
|
||||
fgets(line,MAXLINE,fp);
|
||||
param_extract(tb,line);
|
||||
memory->create(tb->rfile,tb->ninput,"pair:rfile");
|
||||
memory->create(tb->efile,tb->ninput,"pair:efile");
|
||||
memory->create(tb->ffile,tb->ninput,"pair:ffile");
|
||||
|
||||
// setup bitmap parameters for table to read in
|
||||
|
||||
tb->ntablebits = 0;
|
||||
int masklo,maskhi,nmask,nshiftbits;
|
||||
if (tb->rflag == BMP) {
|
||||
while (1 << tb->ntablebits < tb->ninput) tb->ntablebits++;
|
||||
if (1 << tb->ntablebits != tb->ninput)
|
||||
error->one(FLERR,"Bitmapped table is incorrect length in table file");
|
||||
init_bitmap(tb->rlo,tb->rhi,tb->ntablebits,masklo,maskhi,nmask,nshiftbits);
|
||||
}
|
||||
|
||||
// read r,e,f table values from file
|
||||
// if rflag set, compute r
|
||||
// if rflag not set, use r from file
|
||||
|
||||
int itmp;
|
||||
double rtmp;
|
||||
union_int_float_t rsq_lookup;
|
||||
|
||||
fgets(line,MAXLINE,fp);
|
||||
for (int i = 0; i < tb->ninput; i++) {
|
||||
fgets(line,MAXLINE,fp);
|
||||
sscanf(line,"%d %lg %lg %lg",&itmp,&rtmp,&tb->efile[i],&tb->ffile[i]);
|
||||
|
||||
if (tb->rflag == RLINEAR)
|
||||
rtmp = tb->rlo + (tb->rhi - tb->rlo)*i/(tb->ninput-1);
|
||||
else if (tb->rflag == RSQ) {
|
||||
rtmp = tb->rlo*tb->rlo +
|
||||
(tb->rhi*tb->rhi - tb->rlo*tb->rlo)*i/(tb->ninput-1);
|
||||
rtmp = sqrt(rtmp);
|
||||
} else if (tb->rflag == BMP) {
|
||||
rsq_lookup.i = i << nshiftbits;
|
||||
rsq_lookup.i |= masklo;
|
||||
if (rsq_lookup.f < tb->rlo*tb->rlo) {
|
||||
rsq_lookup.i = i << nshiftbits;
|
||||
rsq_lookup.i |= maskhi;
|
||||
}
|
||||
rtmp = sqrtf(rsq_lookup.f);
|
||||
}
|
||||
|
||||
tb->rfile[i] = rtmp;
|
||||
}
|
||||
|
||||
// close file
|
||||
|
||||
fclose(fp);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
broadcast read-in table info from proc 0 to other procs
|
||||
this function communicates these values in Table:
|
||||
ninput,rfile,efile,ffile,rflag,rlo,rhi,fpflag,fplo,fphi
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::bcast_table(Table *tb)
|
||||
{
|
||||
MPI_Bcast(&tb->ninput,1,MPI_INT,0,world);
|
||||
|
||||
int me;
|
||||
MPI_Comm_rank(world,&me);
|
||||
if (me > 0) {
|
||||
memory->create(tb->rfile,tb->ninput,"pair:rfile");
|
||||
memory->create(tb->efile,tb->ninput,"pair:efile");
|
||||
memory->create(tb->ffile,tb->ninput,"pair:ffile");
|
||||
}
|
||||
|
||||
MPI_Bcast(tb->rfile,tb->ninput,MPI_DOUBLE,0,world);
|
||||
MPI_Bcast(tb->efile,tb->ninput,MPI_DOUBLE,0,world);
|
||||
MPI_Bcast(tb->ffile,tb->ninput,MPI_DOUBLE,0,world);
|
||||
|
||||
MPI_Bcast(&tb->rflag,1,MPI_INT,0,world);
|
||||
if (tb->rflag) {
|
||||
MPI_Bcast(&tb->rlo,1,MPI_DOUBLE,0,world);
|
||||
MPI_Bcast(&tb->rhi,1,MPI_DOUBLE,0,world);
|
||||
}
|
||||
MPI_Bcast(&tb->fpflag,1,MPI_INT,0,world);
|
||||
if (tb->fpflag) {
|
||||
MPI_Bcast(&tb->fplo,1,MPI_DOUBLE,0,world);
|
||||
MPI_Bcast(&tb->fphi,1,MPI_DOUBLE,0,world);
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
build spline representation of e,f over entire range of read-in table
|
||||
this function sets these values in Table: e2file,f2file
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::spline_table(Table *tb)
|
||||
{
|
||||
memory->create(tb->e2file,tb->ninput,"pair:e2file");
|
||||
memory->create(tb->f2file,tb->ninput,"pair:f2file");
|
||||
|
||||
double ep0 = - tb->ffile[0];
|
||||
double epn = - tb->ffile[tb->ninput-1];
|
||||
spline(tb->rfile,tb->efile,tb->ninput,ep0,epn,tb->e2file);
|
||||
|
||||
if (tb->fpflag == 0) {
|
||||
tb->fplo = (tb->ffile[1] - tb->ffile[0]) / (tb->rfile[1] - tb->rfile[0]);
|
||||
tb->fphi = (tb->ffile[tb->ninput-1] - tb->ffile[tb->ninput-2]) /
|
||||
(tb->rfile[tb->ninput-1] - tb->rfile[tb->ninput-2]);
|
||||
}
|
||||
|
||||
double fp0 = tb->fplo;
|
||||
double fpn = tb->fphi;
|
||||
spline(tb->rfile,tb->ffile,tb->ninput,fp0,fpn,tb->f2file);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
extract attributes from parameter line in table section
|
||||
format of line: N value R/RSQ/BITMAP lo hi FP fplo fphi
|
||||
N is required, other params are optional
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::param_extract(Table *tb, char *line)
|
||||
{
|
||||
tb->ninput = 0;
|
||||
tb->rflag = NONE;
|
||||
tb->fpflag = 0;
|
||||
|
||||
char *word = strtok(line," \t\n\r\f");
|
||||
while (word) {
|
||||
if (strcmp(word,"N") == 0) {
|
||||
word = strtok(NULL," \t\n\r\f");
|
||||
tb->ninput = atoi(word);
|
||||
} else if (strcmp(word,"R") == 0 || strcmp(word,"RSQ") == 0 ||
|
||||
strcmp(word,"BITMAP") == 0) {
|
||||
if (strcmp(word,"R") == 0) tb->rflag = RLINEAR;
|
||||
else if (strcmp(word,"RSQ") == 0) tb->rflag = RSQ;
|
||||
else if (strcmp(word,"BITMAP") == 0) tb->rflag = BMP;
|
||||
word = strtok(NULL," \t\n\r\f");
|
||||
tb->rlo = atof(word);
|
||||
word = strtok(NULL," \t\n\r\f");
|
||||
tb->rhi = atof(word);
|
||||
} else if (strcmp(word,"FP") == 0) {
|
||||
tb->fpflag = 1;
|
||||
word = strtok(NULL," \t\n\r\f");
|
||||
tb->fplo = atof(word);
|
||||
word = strtok(NULL," \t\n\r\f");
|
||||
tb->fphi = atof(word);
|
||||
} else {
|
||||
error->one(FLERR,"Invalid keyword in pair table parameters");
|
||||
}
|
||||
word = strtok(NULL," \t\n\r\f");
|
||||
}
|
||||
|
||||
if (tb->ninput == 0) error->one(FLERR,"Pair table parameters did not set N");
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
compute r,e,f vectors from splined values
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::compute_table(Table *tb)
|
||||
{
|
||||
update_table = 1;
|
||||
int tlm1 = tablength-1;
|
||||
|
||||
// inner = inner table bound
|
||||
// cut = outer table bound
|
||||
// delta = table spacing in rsq for N-1 bins
|
||||
|
||||
double inner;
|
||||
if (tb->rflag) inner = tb->rlo;
|
||||
else inner = tb->rfile[0];
|
||||
tb->innersq = inner*inner;
|
||||
tb->delta = (tb->cut*tb->cut - tb->innersq) / tlm1;
|
||||
tb->invdelta = 1.0/tb->delta;
|
||||
|
||||
// direct lookup tables
|
||||
// N-1 evenly spaced bins in rsq from inner to cut
|
||||
// e,f = value at midpt of bin
|
||||
// e,f are N-1 in length since store 1 value at bin midpt
|
||||
// f is converted to f/r when stored in f[i]
|
||||
// e,f are never a match to read-in values, always computed via spline interp
|
||||
|
||||
if (tabstyle == LOOKUP) {
|
||||
memory->create(tb->e,tlm1,"pair:e");
|
||||
memory->create(tb->f,tlm1,"pair:f");
|
||||
|
||||
double r,rsq;
|
||||
for (int i = 0; i < tlm1; i++) {
|
||||
rsq = tb->innersq + (i+0.5)*tb->delta;
|
||||
r = sqrt(rsq);
|
||||
tb->e[i] = splint(tb->rfile,tb->efile,tb->e2file,tb->ninput,r);
|
||||
tb->f[i] = splint(tb->rfile,tb->ffile,tb->f2file,tb->ninput,r)/r;
|
||||
}
|
||||
}
|
||||
|
||||
// linear tables
|
||||
// N-1 evenly spaced bins in rsq from inner to cut
|
||||
// rsq,e,f = value at lower edge of bin
|
||||
// de,df values = delta from lower edge to upper edge of bin
|
||||
// rsq,e,f are N in length so de,df arrays can compute difference
|
||||
// f is converted to f/r when stored in f[i]
|
||||
// e,f can match read-in values, else compute via spline interp
|
||||
|
||||
if (tabstyle == LINEAR) {
|
||||
memory->create(tb->rsq,tablength,"pair:rsq");
|
||||
memory->create(tb->e,tablength,"pair:e");
|
||||
memory->create(tb->f,tablength,"pair:f");
|
||||
memory->create(tb->de,tlm1,"pair:de");
|
||||
memory->create(tb->df,tlm1,"pair:df");
|
||||
|
||||
double r,rsq;
|
||||
for (int i = 0; i < tablength; i++) {
|
||||
rsq = tb->innersq + i*tb->delta;
|
||||
r = sqrt(rsq);
|
||||
tb->rsq[i] = rsq;
|
||||
if (tb->match) {
|
||||
tb->e[i] = tb->efile[i];
|
||||
tb->f[i] = tb->ffile[i]/r;
|
||||
} else {
|
||||
tb->e[i] = splint(tb->rfile,tb->efile,tb->e2file,tb->ninput,r);
|
||||
tb->f[i] = splint(tb->rfile,tb->ffile,tb->f2file,tb->ninput,r)/r;
|
||||
}
|
||||
}
|
||||
|
||||
for (int i = 0; i < tlm1; i++) {
|
||||
tb->de[i] = tb->e[i+1] - tb->e[i];
|
||||
tb->df[i] = tb->f[i+1] - tb->f[i];
|
||||
}
|
||||
}
|
||||
|
||||
// cubic spline tables
|
||||
// N-1 evenly spaced bins in rsq from inner to cut
|
||||
// rsq,e,f = value at lower edge of bin
|
||||
// e2,f2 = spline coefficient for each bin
|
||||
// rsq,e,f,e2,f2 are N in length so have N-1 spline bins
|
||||
// f is converted to f/r after e is splined
|
||||
// e,f can match read-in values, else compute via spline interp
|
||||
|
||||
if (tabstyle == SPLINE) {
|
||||
memory->create(tb->rsq,tablength,"pair:rsq");
|
||||
memory->create(tb->e,tablength,"pair:e");
|
||||
memory->create(tb->f,tablength,"pair:f");
|
||||
memory->create(tb->e2,tablength,"pair:e2");
|
||||
memory->create(tb->f2,tablength,"pair:f2");
|
||||
|
||||
tb->deltasq6 = tb->delta*tb->delta / 6.0;
|
||||
|
||||
double r,rsq;
|
||||
for (int i = 0; i < tablength; i++) {
|
||||
rsq = tb->innersq + i*tb->delta;
|
||||
r = sqrt(rsq);
|
||||
tb->rsq[i] = rsq;
|
||||
if (tb->match) {
|
||||
tb->e[i] = tb->efile[i];
|
||||
tb->f[i] = tb->ffile[i]/r;
|
||||
} else {
|
||||
tb->e[i] = splint(tb->rfile,tb->efile,tb->e2file,tb->ninput,r);
|
||||
tb->f[i] = splint(tb->rfile,tb->ffile,tb->f2file,tb->ninput,r);
|
||||
}
|
||||
}
|
||||
|
||||
// ep0,epn = dh/dg at inner and at cut
|
||||
// h(r) = e(r) and g(r) = r^2
|
||||
// dh/dg = (de/dr) / 2r = -f/2r
|
||||
|
||||
double ep0 = - tb->f[0] / (2.0 * sqrt(tb->innersq));
|
||||
double epn = - tb->f[tlm1] / (2.0 * tb->cut);
|
||||
spline(tb->rsq,tb->e,tablength,ep0,epn,tb->e2);
|
||||
|
||||
// fp0,fpn = dh/dg at inner and at cut
|
||||
// h(r) = f(r)/r and g(r) = r^2
|
||||
// dh/dg = (1/r df/dr - f/r^2) / 2r
|
||||
// dh/dg in secant approx = (f(r2)/r2 - f(r1)/r1) / (g(r2) - g(r1))
|
||||
|
||||
double fp0,fpn;
|
||||
double secant_factor = 0.1;
|
||||
if (tb->fpflag) fp0 = (tb->fplo/sqrt(tb->innersq) - tb->f[0]/tb->innersq) /
|
||||
(2.0 * sqrt(tb->innersq));
|
||||
else {
|
||||
double rsq1 = tb->innersq;
|
||||
double rsq2 = rsq1 + secant_factor*tb->delta;
|
||||
fp0 = (splint(tb->rfile,tb->ffile,tb->f2file,tb->ninput,sqrt(rsq2)) /
|
||||
sqrt(rsq2) - tb->f[0] / sqrt(rsq1)) / (secant_factor*tb->delta);
|
||||
}
|
||||
|
||||
if (tb->fpflag && tb->cut == tb->rfile[tb->ninput-1]) fpn =
|
||||
(tb->fphi/tb->cut - tb->f[tlm1]/(tb->cut*tb->cut)) / (2.0 * tb->cut);
|
||||
else {
|
||||
double rsq2 = tb->cut * tb->cut;
|
||||
double rsq1 = rsq2 - secant_factor*tb->delta;
|
||||
fpn = (tb->f[tlm1] / sqrt(rsq2) -
|
||||
splint(tb->rfile,tb->ffile,tb->f2file,tb->ninput,sqrt(rsq1)) /
|
||||
sqrt(rsq1)) / (secant_factor*tb->delta);
|
||||
}
|
||||
|
||||
for (int i = 0; i < tablength; i++) tb->f[i] /= sqrt(tb->rsq[i]);
|
||||
spline(tb->rsq,tb->f,tablength,fp0,fpn,tb->f2);
|
||||
}
|
||||
|
||||
// bitmapped linear tables
|
||||
// 2^N bins from inner to cut, spaced in bitmapped manner
|
||||
// f is converted to f/r when stored in f[i]
|
||||
// e,f can match read-in values, else compute via spline interp
|
||||
|
||||
if (tabstyle == BITMAP) {
|
||||
double r;
|
||||
union_int_float_t rsq_lookup;
|
||||
int masklo,maskhi;
|
||||
|
||||
// linear lookup tables of length ntable = 2^n
|
||||
// stored value = value at lower edge of bin
|
||||
|
||||
init_bitmap(inner,tb->cut,tablength,masklo,maskhi,tb->nmask,tb->nshiftbits);
|
||||
int ntable = 1 << tablength;
|
||||
int ntablem1 = ntable - 1;
|
||||
|
||||
memory->create(tb->rsq,ntable,"pair:rsq");
|
||||
memory->create(tb->e,ntable,"pair:e");
|
||||
memory->create(tb->f,ntable,"pair:f");
|
||||
memory->create(tb->de,ntable,"pair:de");
|
||||
memory->create(tb->df,ntable,"pair:df");
|
||||
memory->create(tb->drsq,ntable,"pair:drsq");
|
||||
|
||||
union_int_float_t minrsq_lookup;
|
||||
minrsq_lookup.i = 0 << tb->nshiftbits;
|
||||
minrsq_lookup.i |= maskhi;
|
||||
|
||||
for (int i = 0; i < ntable; i++) {
|
||||
rsq_lookup.i = i << tb->nshiftbits;
|
||||
rsq_lookup.i |= masklo;
|
||||
if (rsq_lookup.f < tb->innersq) {
|
||||
rsq_lookup.i = i << tb->nshiftbits;
|
||||
rsq_lookup.i |= maskhi;
|
||||
}
|
||||
r = sqrtf(rsq_lookup.f);
|
||||
tb->rsq[i] = rsq_lookup.f;
|
||||
if (tb->match) {
|
||||
tb->e[i] = tb->efile[i];
|
||||
tb->f[i] = tb->ffile[i]/r;
|
||||
} else {
|
||||
tb->e[i] = splint(tb->rfile,tb->efile,tb->e2file,tb->ninput,r);
|
||||
tb->f[i] = splint(tb->rfile,tb->ffile,tb->f2file,tb->ninput,r)/r;
|
||||
}
|
||||
minrsq_lookup.f = MIN(minrsq_lookup.f,rsq_lookup.f);
|
||||
}
|
||||
|
||||
tb->innersq = minrsq_lookup.f;
|
||||
|
||||
for (int i = 0; i < ntablem1; i++) {
|
||||
tb->de[i] = tb->e[i+1] - tb->e[i];
|
||||
tb->df[i] = tb->f[i+1] - tb->f[i];
|
||||
tb->drsq[i] = 1.0/(tb->rsq[i+1] - tb->rsq[i]);
|
||||
}
|
||||
|
||||
// get the delta values for the last table entries
|
||||
// tables are connected periodically between 0 and ntablem1
|
||||
|
||||
tb->de[ntablem1] = tb->e[0] - tb->e[ntablem1];
|
||||
tb->df[ntablem1] = tb->f[0] - tb->f[ntablem1];
|
||||
tb->drsq[ntablem1] = 1.0/(tb->rsq[0] - tb->rsq[ntablem1]);
|
||||
|
||||
// get the correct delta values at itablemax
|
||||
// smallest r is in bin itablemin
|
||||
// largest r is in bin itablemax, which is itablemin-1,
|
||||
// or ntablem1 if itablemin=0
|
||||
|
||||
// deltas at itablemax only needed if corresponding rsq < cut*cut
|
||||
// if so, compute deltas between rsq and cut*cut
|
||||
// if tb->match, data at cut*cut is unavailable, so we'll take
|
||||
// deltas at itablemax-1 as a good approximation
|
||||
|
||||
double e_tmp,f_tmp;
|
||||
int itablemin = minrsq_lookup.i & tb->nmask;
|
||||
itablemin >>= tb->nshiftbits;
|
||||
int itablemax = itablemin - 1;
|
||||
if (itablemin == 0) itablemax = ntablem1;
|
||||
int itablemaxm1 = itablemax - 1;
|
||||
if (itablemax == 0) itablemaxm1 = ntablem1;
|
||||
rsq_lookup.i = itablemax << tb->nshiftbits;
|
||||
rsq_lookup.i |= maskhi;
|
||||
if (rsq_lookup.f < tb->cut*tb->cut) {
|
||||
if (tb->match) {
|
||||
tb->de[itablemax] = tb->de[itablemaxm1];
|
||||
tb->df[itablemax] = tb->df[itablemaxm1];
|
||||
tb->drsq[itablemax] = tb->drsq[itablemaxm1];
|
||||
} else {
|
||||
rsq_lookup.f = tb->cut*tb->cut;
|
||||
r = sqrtf(rsq_lookup.f);
|
||||
e_tmp = splint(tb->rfile,tb->efile,tb->e2file,tb->ninput,r);
|
||||
f_tmp = splint(tb->rfile,tb->ffile,tb->f2file,tb->ninput,r)/r;
|
||||
tb->de[itablemax] = e_tmp - tb->e[itablemax];
|
||||
tb->df[itablemax] = f_tmp - tb->f[itablemax];
|
||||
tb->drsq[itablemax] = 1.0/(rsq_lookup.f - tb->rsq[itablemax]);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
set all ptrs in a table to NULL, so can be freed safely
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::null_table(Table *tb)
|
||||
{
|
||||
tb->rfile = tb->efile = tb->ffile = NULL;
|
||||
tb->e2file = tb->f2file = NULL;
|
||||
tb->rsq = tb->drsq = tb->e = tb->de = NULL;
|
||||
tb->f = tb->df = tb->e2 = tb->f2 = NULL;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
free all arrays in a table
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::free_table(Table *tb)
|
||||
{
|
||||
memory->destroy(tb->rfile);
|
||||
memory->destroy(tb->efile);
|
||||
memory->destroy(tb->ffile);
|
||||
memory->destroy(tb->e2file);
|
||||
memory->destroy(tb->f2file);
|
||||
|
||||
memory->destroy(tb->rsq);
|
||||
memory->destroy(tb->drsq);
|
||||
memory->destroy(tb->e);
|
||||
memory->destroy(tb->de);
|
||||
memory->destroy(tb->f);
|
||||
memory->destroy(tb->df);
|
||||
memory->destroy(tb->e2);
|
||||
memory->destroy(tb->f2);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
spline and splint routines modified from Numerical Recipes
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::spline(double *x, double *y, int n,
|
||||
double yp1, double ypn, double *y2)
|
||||
{
|
||||
int i,k;
|
||||
double p,qn,sig,un;
|
||||
double *u = new double[n];
|
||||
|
||||
if (yp1 > 0.99e30) y2[0] = u[0] = 0.0;
|
||||
else {
|
||||
y2[0] = -0.5;
|
||||
u[0] = (3.0/(x[1]-x[0])) * ((y[1]-y[0]) / (x[1]-x[0]) - yp1);
|
||||
}
|
||||
for (i = 1; i < n-1; i++) {
|
||||
sig = (x[i]-x[i-1]) / (x[i+1]-x[i-1]);
|
||||
p = sig*y2[i-1] + 2.0;
|
||||
y2[i] = (sig-1.0) / p;
|
||||
u[i] = (y[i+1]-y[i]) / (x[i+1]-x[i]) - (y[i]-y[i-1]) / (x[i]-x[i-1]);
|
||||
u[i] = (6.0*u[i] / (x[i+1]-x[i-1]) - sig*u[i-1]) / p;
|
||||
}
|
||||
if (ypn > 0.99e30) qn = un = 0.0;
|
||||
else {
|
||||
qn = 0.5;
|
||||
un = (3.0/(x[n-1]-x[n-2])) * (ypn - (y[n-1]-y[n-2]) / (x[n-1]-x[n-2]));
|
||||
}
|
||||
y2[n-1] = (un-qn*u[n-2]) / (qn*y2[n-2] + 1.0);
|
||||
for (k = n-2; k >= 0; k--) y2[k] = y2[k]*y2[k+1] + u[k];
|
||||
|
||||
delete [] u;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
double PairTableKokkos<DeviceType>::splint(double *xa, double *ya, double *y2a, int n, double x)
|
||||
{
|
||||
int klo,khi,k;
|
||||
double h,b,a,y;
|
||||
|
||||
klo = 0;
|
||||
khi = n-1;
|
||||
while (khi-klo > 1) {
|
||||
k = (khi+klo) >> 1;
|
||||
if (xa[k] > x) khi = k;
|
||||
else klo = k;
|
||||
}
|
||||
h = xa[khi]-xa[klo];
|
||||
a = (xa[khi]-x) / h;
|
||||
b = (x-xa[klo]) / h;
|
||||
y = a*ya[klo] + b*ya[khi] +
|
||||
((a*a*a-a)*y2a[klo] + (b*b*b-b)*y2a[khi]) * (h*h)/6.0;
|
||||
return y;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 writes to restart file
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::write_restart(FILE *fp)
|
||||
{
|
||||
write_restart_settings(fp);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 reads from restart file, bcasts
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::read_restart(FILE *fp)
|
||||
{
|
||||
read_restart_settings(fp);
|
||||
allocate();
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 writes to restart file
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::write_restart_settings(FILE *fp)
|
||||
{
|
||||
fwrite(&tabstyle,sizeof(int),1,fp);
|
||||
fwrite(&tablength,sizeof(int),1,fp);
|
||||
fwrite(&ewaldflag,sizeof(int),1,fp);
|
||||
fwrite(&pppmflag,sizeof(int),1,fp);
|
||||
fwrite(&msmflag,sizeof(int),1,fp);
|
||||
fwrite(&dispersionflag,sizeof(int),1,fp);
|
||||
fwrite(&tip4pflag,sizeof(int),1,fp);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 reads from restart file, bcasts
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::read_restart_settings(FILE *fp)
|
||||
{
|
||||
if (comm->me == 0) {
|
||||
fread(&tabstyle,sizeof(int),1,fp);
|
||||
fread(&tablength,sizeof(int),1,fp);
|
||||
fread(&ewaldflag,sizeof(int),1,fp);
|
||||
fread(&pppmflag,sizeof(int),1,fp);
|
||||
fread(&msmflag,sizeof(int),1,fp);
|
||||
fread(&dispersionflag,sizeof(int),1,fp);
|
||||
fread(&tip4pflag,sizeof(int),1,fp);
|
||||
}
|
||||
MPI_Bcast(&tabstyle,1,MPI_INT,0,world);
|
||||
MPI_Bcast(&tablength,1,MPI_INT,0,world);
|
||||
MPI_Bcast(&ewaldflag,1,MPI_INT,0,world);
|
||||
MPI_Bcast(&pppmflag,1,MPI_INT,0,world);
|
||||
MPI_Bcast(&msmflag,1,MPI_INT,0,world);
|
||||
MPI_Bcast(&dispersionflag,1,MPI_INT,0,world);
|
||||
MPI_Bcast(&tip4pflag,1,MPI_INT,0,world);
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
double PairTableKokkos<DeviceType>::single(int i, int j, int itype, int jtype, double rsq,
|
||||
double factor_coul, double factor_lj,
|
||||
double &fforce)
|
||||
{
|
||||
int itable;
|
||||
double fraction,value,a,b,phi;
|
||||
int tlm1 = tablength - 1;
|
||||
|
||||
Table *tb = &tables[tabindex[itype][jtype]];
|
||||
if (rsq < tb->innersq) error->one(FLERR,"Pair distance < table inner cutoff");
|
||||
|
||||
if (tabstyle == LOOKUP) {
|
||||
itable = static_cast<int> ((rsq-tb->innersq) * tb->invdelta);
|
||||
if (itable >= tlm1) error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
fforce = factor_lj * tb->f[itable];
|
||||
} else if (tabstyle == LINEAR) {
|
||||
itable = static_cast<int> ((rsq-tb->innersq) * tb->invdelta);
|
||||
if (itable >= tlm1) error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
fraction = (rsq - tb->rsq[itable]) * tb->invdelta;
|
||||
value = tb->f[itable] + fraction*tb->df[itable];
|
||||
fforce = factor_lj * value;
|
||||
} else if (tabstyle == SPLINE) {
|
||||
itable = static_cast<int> ((rsq-tb->innersq) * tb->invdelta);
|
||||
if (itable >= tlm1) error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
b = (rsq - tb->rsq[itable]) * tb->invdelta;
|
||||
a = 1.0 - b;
|
||||
value = a * tb->f[itable] + b * tb->f[itable+1] +
|
||||
((a*a*a-a)*tb->f2[itable] + (b*b*b-b)*tb->f2[itable+1]) *
|
||||
tb->deltasq6;
|
||||
fforce = factor_lj * value;
|
||||
} else {
|
||||
union_int_float_t rsq_lookup;
|
||||
rsq_lookup.f = rsq;
|
||||
itable = rsq_lookup.i & tb->nmask;
|
||||
itable >>= tb->nshiftbits;
|
||||
fraction = (rsq_lookup.f - tb->rsq[itable]) * tb->drsq[itable];
|
||||
value = tb->f[itable] + fraction*tb->df[itable];
|
||||
fforce = factor_lj * value;
|
||||
}
|
||||
|
||||
if (tabstyle == LOOKUP)
|
||||
phi = tb->e[itable];
|
||||
else if (tabstyle == LINEAR || tabstyle == BITMAP)
|
||||
phi = tb->e[itable] + fraction*tb->de[itable];
|
||||
else
|
||||
phi = a * tb->e[itable] + b * tb->e[itable+1] +
|
||||
((a*a*a-a)*tb->e2[itable] + (b*b*b-b)*tb->e2[itable+1]) * tb->deltasq6;
|
||||
return factor_lj*phi;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
return the Coulomb cutoff for tabled potentials
|
||||
called by KSpace solvers which require that all pairwise cutoffs be the same
|
||||
loop over all tables not just those indexed by tabindex[i][j] since
|
||||
no way to know which tables are active since pair::init() not yet called
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void *PairTableKokkos<DeviceType>::extract(const char *str, int &dim)
|
||||
{
|
||||
if (strcmp(str,"cut_coul") != 0) return NULL;
|
||||
if (ntables == 0) error->all(FLERR,"All pair coeffs are not set");
|
||||
|
||||
double cut_coul = tables[0].cut;
|
||||
for (int m = 1; m < ntables; m++)
|
||||
if (tables[m].cut != cut_coul)
|
||||
error->all(FLERR,
|
||||
"Pair table cutoffs must all be equal to use with KSpace");
|
||||
dim = 0;
|
||||
return &tables[0].cut;
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableKokkos<DeviceType>::init_style()
|
||||
{
|
||||
|
||||
@ -22,7 +22,7 @@ PairStyle(table/kk/host,PairTableKokkos<LMPHostType>)
|
||||
#ifndef LMP_PAIR_TABLE_KOKKOS_H
|
||||
#define LMP_PAIR_TABLE_KOKKOS_H
|
||||
|
||||
#include "pair.h"
|
||||
#include "pair_table.h"
|
||||
#include "pair_kokkos.h"
|
||||
#include "neigh_list_kokkos.h"
|
||||
#include "atom_kokkos.h"
|
||||
@ -38,7 +38,7 @@ template <class DeviceType, int NEIGHFLAG, int TABSTYLE>
|
||||
class PairTableComputeFunctor;
|
||||
|
||||
template<class DeviceType>
|
||||
class PairTableKokkos : public Pair {
|
||||
class PairTableKokkos : public PairTable {
|
||||
public:
|
||||
|
||||
enum {EnabledNeighFlags=FULL|HALFTHREAD|HALF|N2|FULLCLUSTER};
|
||||
@ -59,18 +59,9 @@ class PairTableKokkos : public Pair {
|
||||
const NeighListKokkos<DeviceType> &list) const;
|
||||
*/
|
||||
void settings(int, char **);
|
||||
void coeff(int, char **);
|
||||
double init_one(int, int);
|
||||
void write_restart(FILE *);
|
||||
void read_restart(FILE *);
|
||||
void write_restart_settings(FILE *);
|
||||
void read_restart_settings(FILE *);
|
||||
double single(int, int, int, int, double, double, double, double &);
|
||||
void *extract(const char *, int &);
|
||||
|
||||
void init_style();
|
||||
|
||||
|
||||
protected:
|
||||
enum{LOOKUP,LINEAR,SPLINE,BITMAP};
|
||||
|
||||
@ -107,17 +98,6 @@ class PairTableKokkos : public Pair {
|
||||
typename ArrayTypes<LMPHostType>::t_ffloat_2d rsq,drsq,e,de,f,df,e2,f2;
|
||||
};
|
||||
|
||||
struct Table {
|
||||
int ninput,rflag,fpflag,match,ntablebits;
|
||||
int nshiftbits,nmask;
|
||||
double rlo,rhi,fplo,fphi,cut;
|
||||
double *rfile,*efile,*ffile;
|
||||
double *e2file,*f2file;
|
||||
double innersq,delta,invdelta,deltasq6;
|
||||
double *rsq,*drsq,*e,*de,*f,*df,*e2,*f2;
|
||||
};
|
||||
int ntables;
|
||||
Table *tables;
|
||||
TableDeviceConst d_table_const;
|
||||
TableDevice* d_table;
|
||||
TableHost* h_table;
|
||||
@ -128,15 +108,6 @@ class PairTableKokkos : public Pair {
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d d_cutsq;
|
||||
|
||||
void allocate();
|
||||
void read_table(Table *, char *, char *);
|
||||
void param_extract(Table *, char *);
|
||||
void bcast_table(Table *);
|
||||
void spline_table(Table *);
|
||||
void compute_table(Table *);
|
||||
void null_table(Table *);
|
||||
void free_table(Table *);
|
||||
void spline(double *, double *, int, double, double, double *);
|
||||
double splint(double *, double *, double *, int, double);
|
||||
|
||||
typename ArrayTypes<DeviceType>::t_x_array_randomread x;
|
||||
typename ArrayTypes<DeviceType>::t_x_array_const c_x;
|
||||
@ -213,11 +184,6 @@ class PairTableKokkos : public Pair {
|
||||
friend void pair_virial_fdotr_compute<PairTableKokkos>(PairTableKokkos*);
|
||||
};
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
@ -297,4 +263,10 @@ E: Cannot use chosen neighbor list style with lj/cut/kk
|
||||
|
||||
That style is not supported by Kokkos.
|
||||
|
||||
|
||||
|
||||
|
||||
*/
|
||||
|
||||
|
||||
|
||||
|
||||
634
src/KOKKOS/pair_table_rx_kokkos.cpp
Normal file
634
src/KOKKOS/pair_table_rx_kokkos.cpp
Normal file
@ -0,0 +1,634 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Stan Moore (SNL)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include <mpi.h>
|
||||
#include <math.h>
|
||||
#include <stdlib.h>
|
||||
#include <string.h>
|
||||
#include "pair_table_rx_kokkos.h"
|
||||
#include "kokkos.h"
|
||||
#include "atom.h"
|
||||
#include "force.h"
|
||||
#include "comm.h"
|
||||
#include "neighbor.h"
|
||||
#include "neigh_list.h"
|
||||
#include "neigh_request.h"
|
||||
#include "memory.h"
|
||||
#include "error.h"
|
||||
#include "atom_masks.h"
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
enum{NONE,RLINEAR,RSQ,BMP};
|
||||
enum{FULL,HALFTHREAD,HALF};
|
||||
|
||||
#define MAXLINE 1024
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
PairTableRXKokkos<DeviceType>::PairTableRXKokkos(LAMMPS *lmp) : PairTableRX(lmp)
|
||||
{
|
||||
update_table = 0;
|
||||
atomKK = (AtomKokkos *) atom;
|
||||
ntables = 0;
|
||||
tables = NULL;
|
||||
execution_space = ExecutionSpaceFromDevice<DeviceType>::space;
|
||||
datamask_read = X_MASK | F_MASK | TYPE_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
datamask_modify = F_MASK | ENERGY_MASK | VIRIAL_MASK;
|
||||
h_table = new TableHost();
|
||||
d_table = new TableDevice();
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
PairTableRXKokkos<DeviceType>::~PairTableRXKokkos()
|
||||
{
|
||||
/* for (int m = 0; m < ntables; m++) free_table(&tables[m]);
|
||||
memory->sfree(tables);
|
||||
|
||||
if (allocated) {
|
||||
memory->destroy(setflag);
|
||||
memory->destroy(cutsq);
|
||||
memory->destroy(tabindex);
|
||||
}*/
|
||||
delete h_table;
|
||||
delete d_table;
|
||||
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableRXKokkos<DeviceType>::compute(int eflag_in, int vflag_in)
|
||||
{
|
||||
if(update_table)
|
||||
create_kokkos_tables();
|
||||
if(tabstyle == LOOKUP)
|
||||
compute_style<LOOKUP>(eflag_in,vflag_in);
|
||||
if(tabstyle == LINEAR)
|
||||
compute_style<LINEAR>(eflag_in,vflag_in);
|
||||
if(tabstyle == SPLINE)
|
||||
compute_style<SPLINE>(eflag_in,vflag_in);
|
||||
if(tabstyle == BITMAP)
|
||||
compute_style<BITMAP>(eflag_in,vflag_in);
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
template<int TABSTYLE>
|
||||
void PairTableRXKokkos<DeviceType>::compute_style(int eflag_in, int vflag_in)
|
||||
{
|
||||
eflag = eflag_in;
|
||||
vflag = vflag_in;
|
||||
|
||||
if (neighflag == FULL || neighflag == FULLCLUSTER) no_virial_fdotr_compute = 1;
|
||||
|
||||
|
||||
if (eflag || vflag) ev_setup(eflag,vflag);
|
||||
else evflag = vflag_fdotr = 0;
|
||||
|
||||
atomKK->sync(execution_space,datamask_read);
|
||||
//k_cutsq.template sync<DeviceType>();
|
||||
//k_params.template sync<DeviceType>();
|
||||
if (eflag || vflag) atomKK->modified(execution_space,datamask_modify);
|
||||
else atomKK->modified(execution_space,F_MASK);
|
||||
|
||||
x = c_x = atomKK->k_x.view<DeviceType>();
|
||||
f = atomKK->k_f.view<DeviceType>();
|
||||
type = atomKK->k_type.view<DeviceType>();
|
||||
nlocal = atom->nlocal;
|
||||
nall = atom->nlocal + atom->nghost;
|
||||
special_lj[0] = force->special_lj[0];
|
||||
special_lj[1] = force->special_lj[1];
|
||||
special_lj[2] = force->special_lj[2];
|
||||
special_lj[3] = force->special_lj[3];
|
||||
newton_pair = force->newton_pair;
|
||||
d_cutsq = d_table->cutsq;
|
||||
// loop over neighbors of my atoms
|
||||
|
||||
EV_FLOAT ev;
|
||||
if(atom->ntypes > MAX_TYPES_STACKPARAMS) {
|
||||
if (neighflag == FULL) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,FULL,false,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
ff(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(list->inum,ff,ev);
|
||||
else Kokkos::parallel_for(list->inum,ff);
|
||||
} else if (neighflag == HALFTHREAD) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,HALFTHREAD,false,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
ff(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(list->inum,ff,ev);
|
||||
else Kokkos::parallel_for(list->inum,ff);
|
||||
} else if (neighflag == HALF) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,HALF,false,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(list->inum,f,ev);
|
||||
else Kokkos::parallel_for(list->inum,f);
|
||||
} else if (neighflag == N2) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,N2,false,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(nlocal,f,ev);
|
||||
else Kokkos::parallel_for(nlocal,f);
|
||||
} else if (neighflag == FULLCLUSTER) {
|
||||
typedef PairComputeFunctor<PairTableRXKokkos<DeviceType>,FULLCLUSTER,false,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f_type;
|
||||
f_type f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
#ifdef KOKKOS_HAVE_CUDA
|
||||
const int teamsize = Kokkos::Impl::is_same<DeviceType, Kokkos::Cuda>::value ? 32 : 1;
|
||||
#else
|
||||
const int teamsize = 1;
|
||||
#endif
|
||||
const int nteams = (list->inum*+teamsize-1)/teamsize;
|
||||
Kokkos::TeamPolicy<DeviceType> config(nteams,teamsize,NeighClusterSize);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(config,f,ev);
|
||||
else Kokkos::parallel_for(config,f);
|
||||
}
|
||||
} else {
|
||||
if (neighflag == FULL) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,FULL,true,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(list->inum,f,ev);
|
||||
else Kokkos::parallel_for(list->inum,f);
|
||||
} else if (neighflag == HALFTHREAD) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,HALFTHREAD,true,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(list->inum,f,ev);
|
||||
else Kokkos::parallel_for(list->inum,f);
|
||||
} else if (neighflag == HALF) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,HALF,true,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(list->inum,f,ev);
|
||||
else Kokkos::parallel_for(list->inum,f);
|
||||
} else if (neighflag == N2) {
|
||||
PairComputeFunctor<PairTableRXKokkos<DeviceType>,N2,true,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(nlocal,f,ev);
|
||||
else Kokkos::parallel_for(nlocal,f);
|
||||
} else if (neighflag == FULLCLUSTER) {
|
||||
typedef PairComputeFunctor<PairTableRXKokkos<DeviceType>,FULLCLUSTER,true,S_TableRXCompute<DeviceType,TABSTYLE> >
|
||||
f_type;
|
||||
f_type f(this,(NeighListKokkos<DeviceType>*) list);
|
||||
#ifdef KOKKOS_HAVE_CUDA
|
||||
const int teamsize = Kokkos::Impl::is_same<DeviceType, Kokkos::Cuda>::value ? 32 : 1;
|
||||
#else
|
||||
const int teamsize = 1;
|
||||
#endif
|
||||
const int nteams = (list->inum*+teamsize-1)/teamsize;
|
||||
Kokkos::TeamPolicy<DeviceType> config(nteams,teamsize,NeighClusterSize);
|
||||
if (eflag || vflag) Kokkos::parallel_reduce(config,f,ev);
|
||||
else Kokkos::parallel_for(config,f);
|
||||
}
|
||||
}
|
||||
|
||||
if (eflag) eng_vdwl += ev.evdwl;
|
||||
if (vflag_global) {
|
||||
virial[0] += ev.v[0];
|
||||
virial[1] += ev.v[1];
|
||||
virial[2] += ev.v[2];
|
||||
virial[3] += ev.v[3];
|
||||
virial[4] += ev.v[4];
|
||||
virial[5] += ev.v[5];
|
||||
}
|
||||
|
||||
if (vflag_fdotr) pair_virial_fdotr_compute(this);
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
template<bool STACKPARAMS, class Specialisation>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
F_FLOAT PairTableRXKokkos<DeviceType>::
|
||||
compute_fpair(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
|
||||
(void) i;
|
||||
(void) j;
|
||||
union_int_float_t rsq_lookup;
|
||||
double fpair;
|
||||
const int tidx = d_table_const.tabindex(itype,jtype);
|
||||
//const Table* const tb = &tables[tabindex[itype][jtype]];
|
||||
|
||||
//if (rsq < d_table_const.innersq(tidx))
|
||||
// error->one(FLERR,"Pair distance < table inner cutoff");
|
||||
|
||||
|
||||
if (Specialisation::TabStyle == LOOKUP) {
|
||||
const int itable = static_cast<int> ((rsq - d_table_const.innersq(tidx)) * d_table_const.invdelta(tidx));
|
||||
//if (itable >= tlm1)
|
||||
// error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
fpair = d_table_const.f(tidx,itable);
|
||||
} else if (Specialisation::TabStyle == LINEAR) {
|
||||
const int itable = static_cast<int> ((rsq - d_table_const.innersq(tidx)) * d_table_const.invdelta(tidx));
|
||||
//if (itable >= tlm1)
|
||||
// error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
const double fraction = (rsq - d_table_const.rsq(tidx,itable)) * d_table_const.invdelta(tidx);
|
||||
fpair = d_table_const.f(tidx,itable) + fraction*d_table_const.df(tidx,itable);
|
||||
} else if (Specialisation::TabStyle == SPLINE) {
|
||||
const int itable = static_cast<int> ((rsq - d_table_const.innersq(tidx)) * d_table_const.invdelta(tidx));
|
||||
//if (itable >= tlm1)
|
||||
// error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
const double b = (rsq - d_table_const.rsq(tidx,itable)) * d_table_const.invdelta(tidx);
|
||||
const double a = 1.0 - b;
|
||||
fpair = a * d_table_const.f(tidx,itable) + b * d_table_const.f(tidx,itable+1) +
|
||||
((a*a*a-a)*d_table_const.f2(tidx,itable) + (b*b*b-b)*d_table_const.f2(tidx,itable+1)) *
|
||||
d_table_const.deltasq6(tidx);
|
||||
} else {
|
||||
rsq_lookup.f = rsq;
|
||||
int itable = rsq_lookup.i & d_table_const.nmask(tidx);
|
||||
itable >>= d_table_const.nshiftbits(tidx);
|
||||
const double fraction = (rsq_lookup.f - d_table_const.rsq(tidx,itable)) * d_table_const.drsq(tidx,itable);
|
||||
fpair = d_table_const.f(tidx,itable) + fraction*d_table_const.df(tidx,itable);
|
||||
}
|
||||
return fpair;
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
template<bool STACKPARAMS, class Specialisation>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
F_FLOAT PairTableRXKokkos<DeviceType>::
|
||||
compute_evdwl(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
|
||||
(void) i;
|
||||
(void) j;
|
||||
double evdwl;
|
||||
union_int_float_t rsq_lookup;
|
||||
const int tidx = d_table_const.tabindex(itype,jtype);
|
||||
//const Table* const tb = &tables[tabindex[itype][jtype]];
|
||||
|
||||
//if (rsq < d_table_const.innersq(tidx))
|
||||
// error->one(FLERR,"Pair distance < table inner cutoff");
|
||||
|
||||
if (Specialisation::TabStyle == LOOKUP) {
|
||||
const int itable = static_cast<int> ((rsq - d_table_const.innersq(tidx)) * d_table_const.invdelta(tidx));
|
||||
//if (itable >= tlm1)
|
||||
// error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
evdwl = d_table_const.e(tidx,itable);
|
||||
} else if (Specialisation::TabStyle == LINEAR) {
|
||||
const int itable = static_cast<int> ((rsq - d_table_const.innersq(tidx)) * d_table_const.invdelta(tidx));
|
||||
//if (itable >= tlm1)
|
||||
// error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
const double fraction = (rsq - d_table_const.rsq(tidx,itable)) * d_table_const.invdelta(tidx);
|
||||
evdwl = d_table_const.e(tidx,itable) + fraction*d_table_const.de(tidx,itable);
|
||||
} else if (Specialisation::TabStyle == SPLINE) {
|
||||
const int itable = static_cast<int> ((rsq - d_table_const.innersq(tidx)) * d_table_const.invdelta(tidx));
|
||||
//if (itable >= tlm1)
|
||||
// error->one(FLERR,"Pair distance > table outer cutoff");
|
||||
const double b = (rsq - d_table_const.rsq(tidx,itable)) * d_table_const.invdelta(tidx);
|
||||
const double a = 1.0 - b;
|
||||
evdwl = a * d_table_const.e(tidx,itable) + b * d_table_const.e(tidx,itable+1) +
|
||||
((a*a*a-a)*d_table_const.e2(tidx,itable) + (b*b*b-b)*d_table_const.e2(tidx,itable+1)) *
|
||||
d_table_const.deltasq6(tidx);
|
||||
} else {
|
||||
rsq_lookup.f = rsq;
|
||||
int itable = rsq_lookup.i & d_table_const.nmask(tidx);
|
||||
itable >>= d_table_const.nshiftbits(tidx);
|
||||
const double fraction = (rsq_lookup.f - d_table_const.rsq(tidx,itable)) * d_table_const.drsq(tidx,itable);
|
||||
evdwl = d_table_const.e(tidx,itable) + fraction*d_table_const.de(tidx,itable);
|
||||
}
|
||||
return evdwl;
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableRXKokkos<DeviceType>::create_kokkos_tables()
|
||||
{
|
||||
const int tlm1 = tablength-1;
|
||||
|
||||
memory->create_kokkos(d_table->nshiftbits,h_table->nshiftbits,ntables,"Table::nshiftbits");
|
||||
memory->create_kokkos(d_table->nmask,h_table->nmask,ntables,"Table::nmask");
|
||||
memory->create_kokkos(d_table->innersq,h_table->innersq,ntables,"Table::innersq");
|
||||
memory->create_kokkos(d_table->invdelta,h_table->invdelta,ntables,"Table::invdelta");
|
||||
memory->create_kokkos(d_table->deltasq6,h_table->deltasq6,ntables,"Table::deltasq6");
|
||||
|
||||
if(tabstyle == LOOKUP) {
|
||||
memory->create_kokkos(d_table->e,h_table->e,ntables,tlm1,"Table::e");
|
||||
memory->create_kokkos(d_table->f,h_table->f,ntables,tlm1,"Table::f");
|
||||
}
|
||||
|
||||
if(tabstyle == LINEAR) {
|
||||
memory->create_kokkos(d_table->rsq,h_table->rsq,ntables,tablength,"Table::rsq");
|
||||
memory->create_kokkos(d_table->e,h_table->e,ntables,tablength,"Table::e");
|
||||
memory->create_kokkos(d_table->f,h_table->f,ntables,tablength,"Table::f");
|
||||
memory->create_kokkos(d_table->de,h_table->de,ntables,tlm1,"Table::de");
|
||||
memory->create_kokkos(d_table->df,h_table->df,ntables,tlm1,"Table::df");
|
||||
}
|
||||
|
||||
if(tabstyle == SPLINE) {
|
||||
memory->create_kokkos(d_table->rsq,h_table->rsq,ntables,tablength,"Table::rsq");
|
||||
memory->create_kokkos(d_table->e,h_table->e,ntables,tablength,"Table::e");
|
||||
memory->create_kokkos(d_table->f,h_table->f,ntables,tablength,"Table::f");
|
||||
memory->create_kokkos(d_table->e2,h_table->e2,ntables,tablength,"Table::e2");
|
||||
memory->create_kokkos(d_table->f2,h_table->f2,ntables,tablength,"Table::f2");
|
||||
}
|
||||
|
||||
if(tabstyle == BITMAP) {
|
||||
int ntable = 1 << tablength;
|
||||
memory->create_kokkos(d_table->rsq,h_table->rsq,ntables,ntable,"Table::rsq");
|
||||
memory->create_kokkos(d_table->e,h_table->e,ntables,ntable,"Table::e");
|
||||
memory->create_kokkos(d_table->f,h_table->f,ntables,ntable,"Table::f");
|
||||
memory->create_kokkos(d_table->de,h_table->de,ntables,ntable,"Table::de");
|
||||
memory->create_kokkos(d_table->df,h_table->df,ntables,ntable,"Table::df");
|
||||
memory->create_kokkos(d_table->drsq,h_table->drsq,ntables,ntable,"Table::drsq");
|
||||
}
|
||||
|
||||
for(int i=0; i < ntables; i++) {
|
||||
Table* tb = &tables[i];
|
||||
|
||||
h_table->nshiftbits[i] = tb->nshiftbits;
|
||||
h_table->nmask[i] = tb->nmask;
|
||||
h_table->innersq[i] = tb->innersq;
|
||||
h_table->invdelta[i] = tb->invdelta;
|
||||
h_table->deltasq6[i] = tb->deltasq6;
|
||||
|
||||
for(int j = 0; j<h_table->rsq.dimension_1(); j++)
|
||||
h_table->rsq(i,j) = tb->rsq[j];
|
||||
for(int j = 0; j<h_table->drsq.dimension_1(); j++)
|
||||
h_table->drsq(i,j) = tb->drsq[j];
|
||||
for(int j = 0; j<h_table->e.dimension_1(); j++)
|
||||
h_table->e(i,j) = tb->e[j];
|
||||
for(int j = 0; j<h_table->de.dimension_1(); j++)
|
||||
h_table->de(i,j) = tb->de[j];
|
||||
for(int j = 0; j<h_table->f.dimension_1(); j++)
|
||||
h_table->f(i,j) = tb->f[j];
|
||||
for(int j = 0; j<h_table->df.dimension_1(); j++)
|
||||
h_table->df(i,j) = tb->df[j];
|
||||
for(int j = 0; j<h_table->e2.dimension_1(); j++)
|
||||
h_table->e2(i,j) = tb->e2[j];
|
||||
for(int j = 0; j<h_table->f2.dimension_1(); j++)
|
||||
h_table->f2(i,j) = tb->f2[j];
|
||||
}
|
||||
|
||||
|
||||
Kokkos::deep_copy(d_table->nshiftbits,h_table->nshiftbits);
|
||||
Kokkos::deep_copy(d_table->nmask,h_table->nmask);
|
||||
Kokkos::deep_copy(d_table->innersq,h_table->innersq);
|
||||
Kokkos::deep_copy(d_table->invdelta,h_table->invdelta);
|
||||
Kokkos::deep_copy(d_table->deltasq6,h_table->deltasq6);
|
||||
Kokkos::deep_copy(d_table->rsq,h_table->rsq);
|
||||
Kokkos::deep_copy(d_table->drsq,h_table->drsq);
|
||||
Kokkos::deep_copy(d_table->e,h_table->e);
|
||||
Kokkos::deep_copy(d_table->de,h_table->de);
|
||||
Kokkos::deep_copy(d_table->f,h_table->f);
|
||||
Kokkos::deep_copy(d_table->df,h_table->df);
|
||||
Kokkos::deep_copy(d_table->e2,h_table->e2);
|
||||
Kokkos::deep_copy(d_table->f2,h_table->f2);
|
||||
Kokkos::deep_copy(d_table->tabindex,h_table->tabindex);
|
||||
|
||||
d_table_const.nshiftbits = d_table->nshiftbits;
|
||||
d_table_const.nmask = d_table->nmask;
|
||||
d_table_const.innersq = d_table->innersq;
|
||||
d_table_const.invdelta = d_table->invdelta;
|
||||
d_table_const.deltasq6 = d_table->deltasq6;
|
||||
d_table_const.rsq = d_table->rsq;
|
||||
d_table_const.drsq = d_table->drsq;
|
||||
d_table_const.e = d_table->e;
|
||||
d_table_const.de = d_table->de;
|
||||
d_table_const.f = d_table->f;
|
||||
d_table_const.df = d_table->df;
|
||||
d_table_const.e2 = d_table->e2;
|
||||
d_table_const.f2 = d_table->f2;
|
||||
|
||||
|
||||
Kokkos::deep_copy(d_table->cutsq,h_table->cutsq);
|
||||
update_table = 0;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
allocate all arrays
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableRXKokkos<DeviceType>::allocate()
|
||||
{
|
||||
allocated = 1;
|
||||
const int nt = atom->ntypes + 1;
|
||||
|
||||
memory->create(setflag,nt,nt,"pair:setflag");
|
||||
memory->create_kokkos(d_table->cutsq,h_table->cutsq,cutsq,nt,nt,"pair:cutsq");
|
||||
memory->create_kokkos(d_table->tabindex,h_table->tabindex,tabindex,nt,nt,"pair:tabindex");
|
||||
|
||||
d_table_const.cutsq = d_table->cutsq;
|
||||
d_table_const.tabindex = d_table->tabindex;
|
||||
memset(&setflag[0][0],0,nt*nt*sizeof(int));
|
||||
memset(&cutsq[0][0],0,nt*nt*sizeof(double));
|
||||
memset(&tabindex[0][0],0,nt*nt*sizeof(int));
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
global settings
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableRXKokkos<DeviceType>::settings(int narg, char **arg)
|
||||
{
|
||||
if (narg < 2) error->all(FLERR,"Illegal pair_style command");
|
||||
|
||||
// new settings
|
||||
|
||||
if (strcmp(arg[0],"lookup") == 0) tabstyle = LOOKUP;
|
||||
else if (strcmp(arg[0],"linear") == 0) tabstyle = LINEAR;
|
||||
else if (strcmp(arg[0],"spline") == 0) tabstyle = SPLINE;
|
||||
else if (strcmp(arg[0],"bitmap") == 0) tabstyle = BITMAP;
|
||||
else error->all(FLERR,"Unknown table style in pair_style command");
|
||||
|
||||
tablength = force->inumeric(FLERR,arg[1]);
|
||||
if (tablength < 2) error->all(FLERR,"Illegal number of pair table entries");
|
||||
|
||||
// optional keywords
|
||||
// assert the tabulation is compatible with a specific long-range solver
|
||||
|
||||
int iarg = 2;
|
||||
while (iarg < narg) {
|
||||
if (strcmp(arg[iarg],"ewald") == 0) ewaldflag = 1;
|
||||
else if (strcmp(arg[iarg],"pppm") == 0) pppmflag = 1;
|
||||
else if (strcmp(arg[iarg],"msm") == 0) msmflag = 1;
|
||||
else if (strcmp(arg[iarg],"dispersion") == 0) dispersionflag = 1;
|
||||
else if (strcmp(arg[iarg],"tip4p") == 0) tip4pflag = 1;
|
||||
else error->all(FLERR,"Illegal pair_style command");
|
||||
iarg++;
|
||||
}
|
||||
|
||||
// delete old tables, since cannot just change settings
|
||||
|
||||
for (int m = 0; m < ntables; m++) free_table(&tables[m]);
|
||||
memory->sfree(tables);
|
||||
|
||||
if (allocated) {
|
||||
memory->destroy(setflag);
|
||||
|
||||
d_table_const.tabindex = d_table->tabindex = typename ArrayTypes<DeviceType>::t_int_2d();
|
||||
h_table->tabindex = typename ArrayTypes<LMPHostType>::t_int_2d();
|
||||
|
||||
d_table_const.cutsq = d_table->cutsq = typename ArrayTypes<DeviceType>::t_ffloat_2d();
|
||||
h_table->cutsq = typename ArrayTypes<LMPHostType>::t_ffloat_2d();
|
||||
}
|
||||
allocated = 0;
|
||||
|
||||
ntables = 0;
|
||||
tables = NULL;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
init for one type pair i,j and corresponding j,i
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
template<class DeviceType>
|
||||
double PairTableRXKokkos<DeviceType>::init_one(int i, int j)
|
||||
{
|
||||
if (setflag[i][j] == 0) error->all(FLERR,"All pair coeffs are not set");
|
||||
|
||||
tabindex[j][i] = tabindex[i][j];
|
||||
|
||||
if(i<MAX_TYPES_STACKPARAMS+1 && j<MAX_TYPES_STACKPARAMS+1) {
|
||||
m_cutsq[j][i] = m_cutsq[i][j] = tables[tabindex[i][j]].cut*tables[tabindex[i][j]].cut;
|
||||
}
|
||||
|
||||
return tables[tabindex[i][j]].cut;
|
||||
}
|
||||
|
||||
template<class DeviceType>
|
||||
void PairTableRXKokkos<DeviceType>::init_style()
|
||||
{
|
||||
neighbor->request(this,instance_me);
|
||||
neighflag = lmp->kokkos->neighflag;
|
||||
int irequest = neighbor->nrequest - 1;
|
||||
|
||||
neighbor->requests[irequest]->
|
||||
kokkos_host = Kokkos::Impl::is_same<DeviceType,LMPHostType>::value &&
|
||||
!Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
|
||||
neighbor->requests[irequest]->
|
||||
kokkos_device = Kokkos::Impl::is_same<DeviceType,LMPDeviceType>::value;
|
||||
|
||||
if (neighflag == FULL) {
|
||||
neighbor->requests[irequest]->full = 1;
|
||||
neighbor->requests[irequest]->half = 0;
|
||||
neighbor->requests[irequest]->full_cluster = 0;
|
||||
} else if (neighflag == HALF || neighflag == HALFTHREAD) {
|
||||
neighbor->requests[irequest]->full = 0;
|
||||
neighbor->requests[irequest]->half = 1;
|
||||
neighbor->requests[irequest]->full_cluster = 0;
|
||||
} else if (neighflag == N2) {
|
||||
neighbor->requests[irequest]->full = 0;
|
||||
neighbor->requests[irequest]->half = 0;
|
||||
neighbor->requests[irequest]->full_cluster = 0;
|
||||
} else if (neighflag == FULLCLUSTER) {
|
||||
neighbor->requests[irequest]->full_cluster = 1;
|
||||
neighbor->requests[irequest]->full = 1;
|
||||
neighbor->requests[irequest]->half = 0;
|
||||
} else {
|
||||
error->all(FLERR,"Cannot use chosen neighbor list style with lj/cut/kk");
|
||||
}
|
||||
}
|
||||
|
||||
/*
|
||||
template <class DeviceType> template<int NEIGHFLAG>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
void PairTableRXKokkos<DeviceType>::
|
||||
ev_tally(EV_FLOAT &ev, const int &i, const int &j, const F_FLOAT &fpair,
|
||||
const F_FLOAT &delx, const F_FLOAT &dely, const F_FLOAT &delz) const
|
||||
{
|
||||
const int EFLAG = eflag;
|
||||
const int NEWTON_PAIR = newton_pair;
|
||||
const int VFLAG = vflag_either;
|
||||
|
||||
if (EFLAG) {
|
||||
if (eflag_atom) {
|
||||
E_FLOAT epairhalf = 0.5 * (ev.evdwl + ev.ecoul);
|
||||
if (NEWTON_PAIR || i < nlocal) eatom[i] += epairhalf;
|
||||
if (NEWTON_PAIR || j < nlocal) eatom[j] += epairhalf;
|
||||
}
|
||||
}
|
||||
|
||||
if (VFLAG) {
|
||||
const E_FLOAT v0 = delx*delx*fpair;
|
||||
const E_FLOAT v1 = dely*dely*fpair;
|
||||
const E_FLOAT v2 = delz*delz*fpair;
|
||||
const E_FLOAT v3 = delx*dely*fpair;
|
||||
const E_FLOAT v4 = delx*delz*fpair;
|
||||
const E_FLOAT v5 = dely*delz*fpair;
|
||||
|
||||
if (vflag_global) {
|
||||
if (NEIGHFLAG) {
|
||||
if (NEWTON_PAIR) {
|
||||
ev.v[0] += v0;
|
||||
ev.v[1] += v1;
|
||||
ev.v[2] += v2;
|
||||
ev.v[3] += v3;
|
||||
ev.v[4] += v4;
|
||||
ev.v[5] += v5;
|
||||
} else {
|
||||
if (i < nlocal) {
|
||||
ev.v[0] += 0.5*v0;
|
||||
ev.v[1] += 0.5*v1;
|
||||
ev.v[2] += 0.5*v2;
|
||||
ev.v[3] += 0.5*v3;
|
||||
ev.v[4] += 0.5*v4;
|
||||
ev.v[5] += 0.5*v5;
|
||||
}
|
||||
if (j < nlocal) {
|
||||
ev.v[0] += 0.5*v0;
|
||||
ev.v[1] += 0.5*v1;
|
||||
ev.v[2] += 0.5*v2;
|
||||
ev.v[3] += 0.5*v3;
|
||||
ev.v[4] += 0.5*v4;
|
||||
ev.v[5] += 0.5*v5;
|
||||
}
|
||||
}
|
||||
} else {
|
||||
ev.v[0] += 0.5*v0;
|
||||
ev.v[1] += 0.5*v1;
|
||||
ev.v[2] += 0.5*v2;
|
||||
ev.v[3] += 0.5*v3;
|
||||
ev.v[4] += 0.5*v4;
|
||||
ev.v[5] += 0.5*v5;
|
||||
}
|
||||
}
|
||||
|
||||
if (vflag_atom) {
|
||||
if (NEWTON_PAIR || i < nlocal) {
|
||||
d_vatom(i,0) += 0.5*v0;
|
||||
d_vatom(i,1) += 0.5*v1;
|
||||
d_vatom(i,2) += 0.5*v2;
|
||||
d_vatom(i,3) += 0.5*v3;
|
||||
d_vatom(i,4) += 0.5*v4;
|
||||
d_vatom(i,5) += 0.5*v5;
|
||||
}
|
||||
if (NEWTON_PAIR || (NEIGHFLAG && j < nlocal)) {
|
||||
d_vatom(j,0) += 0.5*v0;
|
||||
d_vatom(j,1) += 0.5*v1;
|
||||
d_vatom(j,2) += 0.5*v2;
|
||||
d_vatom(j,3) += 0.5*v3;
|
||||
d_vatom(j,4) += 0.5*v4;
|
||||
d_vatom(j,5) += 0.5*v5;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
*/
|
||||
template<class DeviceType>
|
||||
void PairTableRXKokkos<DeviceType>::cleanup_copy() {
|
||||
// WHY needed: this prevents parent copy from deallocating any arrays
|
||||
allocated = 0;
|
||||
cutsq = NULL;
|
||||
eatom = NULL;
|
||||
vatom = NULL;
|
||||
h_table=NULL; d_table=NULL;
|
||||
}
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
template class PairTableRXKokkos<LMPDeviceType>;
|
||||
#ifdef KOKKOS_HAVE_CUDA
|
||||
template class PairTableRXKokkos<LMPHostType>;
|
||||
#endif
|
||||
|
||||
}
|
||||
|
||||
269
src/KOKKOS/pair_table_rx_kokkos.h
Normal file
269
src/KOKKOS/pair_table_rx_kokkos.h
Normal file
@ -0,0 +1,269 @@
|
||||
/* -*- c++ -*- ----------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, Sandia National Laboratories
|
||||
Steve Plimpton, sjplimp@sandia.gov
|
||||
|
||||
Copyright (2003) Sandia Corporation. Under the terms of Contract
|
||||
DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
|
||||
certain rights in this software. This software is distributed under
|
||||
the GNU General Public License.
|
||||
|
||||
See the README file in the top-level LAMMPS directory.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef PAIR_CLASS
|
||||
|
||||
PairStyle(table/rx/kk,PairTableRXKokkos<LMPDeviceType>)
|
||||
PairStyle(table/rx/kk/device,PairTableRXKokkos<LMPDeviceType>)
|
||||
PairStyle(table/rx/kk/host,PairTableRXKokkos<LMPHostType>)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_PAIR_TABLE_RX_KOKKOS_H
|
||||
#define LMP_PAIR_TABLE_RX_KOKKOS_H
|
||||
|
||||
#include "pair_table_rx.h"
|
||||
#include "pair_kokkos.h"
|
||||
#include "neigh_list_kokkos.h"
|
||||
#include "atom_kokkos.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
template<class Device,int TABSTYLE>
|
||||
struct S_TableRXCompute {
|
||||
enum {TabStyle = TABSTYLE};
|
||||
};
|
||||
|
||||
template <class DeviceType, int NEIGHFLAG, int TABSTYLE>
|
||||
class PairTableRXComputeFunctor;
|
||||
|
||||
template<class DeviceType>
|
||||
class PairTableRXKokkos : public PairTableRX {
|
||||
public:
|
||||
|
||||
enum {EnabledNeighFlags=FULL|HALFTHREAD|HALF|N2|FULLCLUSTER};
|
||||
enum {COUL_FLAG=0};
|
||||
typedef DeviceType device_type;
|
||||
|
||||
PairTableRXKokkos(class LAMMPS *);
|
||||
virtual ~PairTableRXKokkos();
|
||||
|
||||
virtual void compute(int, int);
|
||||
|
||||
template<int TABSTYLE>
|
||||
void compute_style(int, int);
|
||||
|
||||
/*template<int EVFLAG, int NEIGHFLAG, int NEWTON_PAIR, int TABSTYLE>
|
||||
KOKKOS_FUNCTION
|
||||
EV_FLOAT compute_item(const int& i,
|
||||
const NeighListKokkos<DeviceType> &list) const;
|
||||
*/
|
||||
void settings(int, char **);
|
||||
double init_one(int, int);
|
||||
void init_style();
|
||||
|
||||
protected:
|
||||
enum{LOOKUP,LINEAR,SPLINE,BITMAP};
|
||||
|
||||
int tabstyle,tablength;
|
||||
/*struct TableDeviceConst {
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d_randomread cutsq;
|
||||
typename ArrayTypes<DeviceType>::t_int_2d_randomread tabindex;
|
||||
typename ArrayTypes<DeviceType>::t_int_1d_randomread nshiftbits,nmask;
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_1d_randomread innersq,invdelta,deltasq6;
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d_randomread rsq,drsq,e,de,f,df,e2,f2;
|
||||
};*/
|
||||
//Its faster not to use texture fetch if the number of tables is less than 32!
|
||||
struct TableDeviceConst {
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d cutsq;
|
||||
typename ArrayTypes<DeviceType>::t_int_2d tabindex;
|
||||
typename ArrayTypes<DeviceType>::t_int_1d nshiftbits,nmask;
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_1d innersq,invdelta,deltasq6;
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d_randomread rsq,drsq,e,de,f,df,e2,f2;
|
||||
};
|
||||
|
||||
struct TableDevice {
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d cutsq;
|
||||
typename ArrayTypes<DeviceType>::t_int_2d tabindex;
|
||||
typename ArrayTypes<DeviceType>::t_int_1d nshiftbits,nmask;
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_1d innersq,invdelta,deltasq6;
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d rsq,drsq,e,de,f,df,e2,f2;
|
||||
};
|
||||
|
||||
struct TableHost {
|
||||
typename ArrayTypes<LMPHostType>::t_ffloat_2d cutsq;
|
||||
typename ArrayTypes<LMPHostType>::t_int_2d tabindex;
|
||||
typename ArrayTypes<LMPHostType>::t_int_1d nshiftbits,nmask;
|
||||
typename ArrayTypes<LMPHostType>::t_ffloat_1d innersq,invdelta,deltasq6;
|
||||
typename ArrayTypes<LMPHostType>::t_ffloat_2d rsq,drsq,e,de,f,df,e2,f2;
|
||||
};
|
||||
|
||||
TableDeviceConst d_table_const;
|
||||
TableDevice* d_table;
|
||||
TableHost* h_table;
|
||||
|
||||
int **tabindex;
|
||||
F_FLOAT m_cutsq[MAX_TYPES_STACKPARAMS+1][MAX_TYPES_STACKPARAMS+1];
|
||||
|
||||
typename ArrayTypes<DeviceType>::t_ffloat_2d d_cutsq;
|
||||
|
||||
void allocate();
|
||||
|
||||
typename ArrayTypes<DeviceType>::t_x_array_randomread x;
|
||||
typename ArrayTypes<DeviceType>::t_x_array_const c_x;
|
||||
typename ArrayTypes<DeviceType>::t_f_array f;
|
||||
typename ArrayTypes<DeviceType>::t_int_1d_randomread type;
|
||||
typename ArrayTypes<DeviceType>::t_efloat_1d d_eatom;
|
||||
typename ArrayTypes<DeviceType>::t_virial_array d_vatom;
|
||||
|
||||
protected:
|
||||
int nlocal,nall,eflag,vflag,neighflag,newton_pair;
|
||||
|
||||
int update_table;
|
||||
void create_kokkos_tables();
|
||||
void cleanup_copy();
|
||||
|
||||
template<bool STACKPARAMS, class Specialisation>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
F_FLOAT compute_fpair(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const;
|
||||
|
||||
template<bool STACKPARAMS, class Specialisation>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
F_FLOAT compute_evdwl(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const;
|
||||
|
||||
template<bool STACKPARAMS, class Specialisation>
|
||||
KOKKOS_INLINE_FUNCTION
|
||||
F_FLOAT compute_ecoul(const F_FLOAT& rsq, const int& i, const int&j, const int& itype, const int& jtype) const {
|
||||
return 0;
|
||||
}
|
||||
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,true,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,true,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,true,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,true,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,true,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,false,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,false,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,false,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,false,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,false,S_TableRXCompute<DeviceType,LOOKUP> >;
|
||||
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,true,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,true,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,true,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,true,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,true,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,false,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,false,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,false,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,false,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,false,S_TableRXCompute<DeviceType,LINEAR> >;
|
||||
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,true,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,true,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,true,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,true,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,true,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,false,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,false,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,false,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,false,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,false,S_TableRXCompute<DeviceType,SPLINE> >;
|
||||
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,true,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,true,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,true,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,true,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,true,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULL,false,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALF,false,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,HALFTHREAD,false,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,N2,false,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
friend class PairComputeFunctor<PairTableRXKokkos,FULLCLUSTER,false,S_TableRXCompute<DeviceType,BITMAP> >;
|
||||
|
||||
friend void pair_virial_fdotr_compute<PairTableRXKokkos>(PairTableRXKokkos*);
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
|
||||
/* ERROR/WARNING messages:
|
||||
|
||||
E: Pair distance < table inner cutoff
|
||||
|
||||
Two atoms are closer together than the pairwise table allows.
|
||||
|
||||
E: Pair distance > table outer cutoff
|
||||
|
||||
Two atoms are further apart than the pairwise table allows.
|
||||
|
||||
E: Illegal ... command
|
||||
|
||||
Self-explanatory. Check the input script syntax and compare to the
|
||||
documentation for the command. You can use -echo screen as a
|
||||
command-line option when running LAMMPS to see the offending line.
|
||||
|
||||
E: Unknown table style in pair_style command
|
||||
|
||||
Style of table is invalid for use with pair_style table command.
|
||||
|
||||
E: Illegal number of pair table entries
|
||||
|
||||
There must be at least 2 table entries.
|
||||
|
||||
E: Invalid pair table length
|
||||
|
||||
Length of read-in pair table is invalid
|
||||
|
||||
E: Invalid pair table cutoff
|
||||
|
||||
Cutoffs in pair_coeff command are not valid with read-in pair table.
|
||||
|
||||
E: Bitmapped table in file does not match requested table
|
||||
|
||||
Setting for bitmapped table in pair_coeff command must match table
|
||||
in file exactly.
|
||||
|
||||
E: All pair coeffs are not set
|
||||
|
||||
All pair coefficients must be set in the data file or by the
|
||||
pair_coeff command before running a simulation.
|
||||
|
||||
E: Cannot open file %s
|
||||
|
||||
The specified file cannot be opened. Check that the path and name are
|
||||
correct. If the file is a compressed file, also check that the gzip
|
||||
executable can be found and run.
|
||||
|
||||
E: Did not find keyword in table file
|
||||
|
||||
Keyword used in pair_coeff command was not found in table file.
|
||||
|
||||
E: Bitmapped table is incorrect length in table file
|
||||
|
||||
Number of table entries is not a correct power of 2.
|
||||
|
||||
E: Invalid keyword in pair table parameters
|
||||
|
||||
Keyword used in list of table parameters is not recognized.
|
||||
|
||||
E: Pair table parameters did not set N
|
||||
|
||||
List of pair table parameters must include N setting.
|
||||
|
||||
E: Pair table cutoffs must all be equal to use with KSpace
|
||||
|
||||
When using pair style table with a long-range KSpace solver, the
|
||||
cutoffs for all atom type pairs must all be the same, since the
|
||||
long-range solver starts at that cutoff.
|
||||
|
||||
E: Cannot use chosen neighbor list style with lj/cut/kk
|
||||
|
||||
That style is not supported by Kokkos.
|
||||
|
||||
|
||||
|
||||
|
||||
*/
|
||||
@ -18,7 +18,7 @@ PairStyle(multi/lucy,PairMultiLucy)
|
||||
#else
|
||||
|
||||
#ifndef LMP_PAIR_MULTI_LUCY_H
|
||||
#define LMP_PAIR_MUTLI_LUCY_H
|
||||
#define LMP_PAIR_MULTI_LUCY_H
|
||||
|
||||
#include "pair.h"
|
||||
|
||||
|
||||
@ -59,8 +59,7 @@ static const char cite_pair_multi_lucy_rx[] =
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
PairMultiLucyRX::PairMultiLucyRX(LAMMPS *lmp) : Pair(lmp),
|
||||
ntables(0), tables(NULL), tabindex(NULL), site1(NULL), site2(NULL)
|
||||
PairMultiLucyRX::PairMultiLucyRX(LAMMPS *lmp) : Pair(lmp)
|
||||
{
|
||||
if (lmp->citeme) lmp->citeme->add(cite_pair_multi_lucy_rx);
|
||||
|
||||
@ -69,6 +68,9 @@ PairMultiLucyRX::PairMultiLucyRX(LAMMPS *lmp) : Pair(lmp),
|
||||
ntables = 0;
|
||||
tables = NULL;
|
||||
|
||||
tabindex = NULL;
|
||||
site1 = site2 = NULL;
|
||||
|
||||
comm_forward = 1;
|
||||
comm_reverse = 1;
|
||||
|
||||
|
||||
@ -18,7 +18,7 @@ PairStyle(multi/lucy/rx,PairMultiLucyRX)
|
||||
#else
|
||||
|
||||
#ifndef LMP_PAIR_MULTI_LUCY_RX_H
|
||||
#define LMP_PAIR_MUTLI_LUCY_RX_H
|
||||
#define LMP_PAIR_MULTI_LUCY_RX_H
|
||||
|
||||
#include "pair.h"
|
||||
|
||||
|
||||
@ -30,9 +30,9 @@ class PairTable : public Pair {
|
||||
virtual ~PairTable();
|
||||
|
||||
virtual void compute(int, int);
|
||||
void settings(int, char **);
|
||||
virtual void settings(int, char **);
|
||||
void coeff(int, char **);
|
||||
double init_one(int, int);
|
||||
virtual double init_one(int, int);
|
||||
void write_restart(FILE *);
|
||||
void read_restart(FILE *);
|
||||
void write_restart_settings(FILE *);
|
||||
@ -58,7 +58,7 @@ class PairTable : public Pair {
|
||||
|
||||
int **tabindex;
|
||||
|
||||
void allocate();
|
||||
virtual void allocate();
|
||||
void read_table(Table *, char *, char *);
|
||||
void param_extract(Table *, char *);
|
||||
void bcast_table(Table *);
|
||||
|
||||
Reference in New Issue
Block a user