jianmu
/
lammps
forked from lijiext/lammps
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@7644 f3b2605a-c512-4ea7-a41b-209d697bcdaa

This commit is contained in:
sjplimp 2012-02-01 16:07:44 +00:00
parent 7e11d0e7ba
commit 2c42373bc8
28 changed files with 300 additions and 206 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

20
lib/gpu/lal_table.cpp
View File

@ -122,7 +122,7 @@ int TableT::init(const int ntypes,
host_write[ix*lj_types+iy].w = (numtyp)0.0;
}
ucl_copy(coeff2,host_write,false);
// Allocate tablength arrays
UCL_H_Vec<numtyp4> host_write2(_ntables*_tablength,*(this->ucl_device),
UCL_WRITE_OPTIMIZED);
@ -307,29 +307,29 @@ void TableT::loop(const bool _eflag, const bool _vflag) {
&this->ans->dev_engv.begin(), &eflag, &vflag,
&ainum, &nbor_pitch, &this->_threads_per_atom,
&_tablength);
}
}
} else {
if (_tabstyle == LOOKUP) {
this->k_pair.set_size(GX,BX);
this->k_pair.run(&this->atom->dev_x.begin(), &tabindex.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &cutsq.begin(),
&_lj_types, &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &_lj_types,
&cutsq.begin(), &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&this->_nbor_data->begin(), &this->ans->dev_ans.begin(),
&this->ans->dev_engv.begin(), &eflag, &vflag, &ainum,
&nbor_pitch, &this->_threads_per_atom, &_tablength);
} else if (_tabstyle == LINEAR) {
this->k_pair_linear.set_size(GX,BX);
this->k_pair_linear.run(&this->atom->dev_x.begin(), &tabindex.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &cutsq.begin(),
&_lj_types, &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &_lj_types,
&cutsq.begin(), &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&this->_nbor_data->begin(), &this->ans->dev_ans.begin(),
&this->ans->dev_engv.begin(), &eflag, &vflag, &ainum,
&nbor_pitch, &this->_threads_per_atom, &_tablength);
} else if (_tabstyle == SPLINE) {
this->k_pair_spline.set_size(GX,BX);
this->k_pair_spline.run(&this->atom->dev_x.begin(), &tabindex.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &cutsq.begin(),
&_lj_types, &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &_lj_types,
&cutsq.begin(), &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&this->_nbor_data->begin(), &this->ans->dev_ans.begin(),
&this->ans->dev_engv.begin(), &eflag, &vflag, &ainum,
&nbor_pitch, &this->_threads_per_atom, &_tablength);
@ -337,8 +337,8 @@ void TableT::loop(const bool _eflag, const bool _vflag) {
this->k_pair_bitmap.set_size(GX,BX);
this->k_pair_bitmap.run(&this->atom->dev_x.begin(), &tabindex.begin(),
&nshiftbits.begin(), &nmask.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &cutsq.begin(),
&_lj_types, &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&coeff2.begin(), &coeff3.begin(), &coeff4.begin(), &_lj_types,
&cutsq.begin(), &sp_lj.begin(), &this->nbor->dev_nbor.begin(),
&this->_nbor_data->begin(), &this->ans->dev_ans.begin(),
&this->ans->dev_engv.begin(), &eflag, &vflag, &ainum,
&nbor_pitch, &this->_threads_per_atom, &_tablength);

2
src/GPU/pair_buck_coul_long_gpu.cpp
View File

@ -248,7 +248,7 @@ void PairBuckCoulLongGPU::cpu_compute(int start, int inum, int eflag,
if (rsq < cutsq[itype][jtype]) {
r2inv = 1.0/rsq;
r = sqrt(rsq);
if (rsq < cut_coulsq) {
grij = g_ewald * r;
expm2 = exp(-grij*grij);

9
src/STUBS/mpi.c
View File

@ -125,6 +125,15 @@ int MPI_Send(void *buf, int count, MPI_Datatype datatype,
/* ---------------------------------------------------------------------- */
int MPI_Isend(void *buf, int count, MPI_Datatype datatype,
int source, int tag, MPI_Comm comm, MPI_Request *request)
{
printf("MPI Stub WARNING: Should not send message to self\n");
return 0;
}
/* ---------------------------------------------------------------------- */
int MPI_Rsend(void *buf, int count, MPI_Datatype datatype,
int dest, int tag, MPI_Comm comm)
{

2
src/STUBS/mpi.h
View File

@ -73,6 +73,8 @@ int MPI_Type_size(int, int *);
int MPI_Send(void *buf, int count, MPI_Datatype datatype,
int dest, int tag, MPI_Comm comm);
int MPI_Isend(void *buf, int count, MPI_Datatype datatype,
int source, int tag, MPI_Comm comm, MPI_Request *request);
int MPI_Rsend(void *buf, int count, MPI_Datatype datatype,
int dest, int tag, MPI_Comm comm);
int MPI_Recv(void *buf, int count, MPI_Datatype datatype,

6
src/USER-MISC/dihedral_table.h
View File

@ -36,15 +36,15 @@ namespace LAMMPS_NS {
class DihedralTable : public Dihedral {
public:
DihedralTable(class LAMMPS *);
~DihedralTable();
void compute(int, int);
virtual ~DihedralTable();
virtual void compute(int, int);
void settings(int, char **);
void coeff(int, char **);
void write_restart(FILE *);
void read_restart(FILE *);
double single(int type, int i1, int i2, int i3, int i4);
private:
protected:
int tabstyle,tablength;
// double *phi0; <- equilibrium angles not supported

7
src/USER-OMP/bond_fene_expand_omp.cpp
View File

@ -81,6 +81,7 @@ void BondFENEExpandOMP::eval(int nfrom, int nto, ThrData * const thr)
double * const * const f = thr->get_f();
const int * const * const bondlist = neighbor->bondlist;
const int nlocal = atom->nlocal;
const int tid = thr->get_tid();
for (n = nfrom; n < nto; n++) {
i1 = bondlist[n][0];
@ -105,10 +106,14 @@ void BondFENEExpandOMP::eval(int nfrom, int nto, ThrData * const thr)
if (rlogarg < 0.1) {
char str[128];
sprintf(str,"FENE bond too long: " BIGINT_FORMAT " %d %d %g",
update->ntimestep,atom->tag[i1],atom->tag[i2],sqrt(rsq));
error->warning(FLERR,str,0);
if (rlogarg <= -3.0) error->one(FLERR,"Bad FENE bond");
if (check_error_thr((rlogarg <= -3.0),tid,FLERR,"Bad FENE bond"))
return;
rlogarg = 0.1;
}

7
src/USER-OMP/bond_fene_omp.cpp
View File

@ -80,6 +80,7 @@ void BondFENEOMP::eval(int nfrom, int nto, ThrData * const thr)
double * const * const f = thr->get_f();
const int * const * const bondlist = neighbor->bondlist;
const int nlocal = atom->nlocal;
const int tid = thr->get_tid();
for (n = nfrom; n < nto; n++) {
i1 = bondlist[n][0];
@ -101,10 +102,14 @@ void BondFENEOMP::eval(int nfrom, int nto, ThrData * const thr)
if (rlogarg < 0.1) {
char str[128];
sprintf(str,"FENE bond too long: " BIGINT_FORMAT " %d %d %g",
update->ntimestep,atom->tag[i1],atom->tag[i2],sqrt(rsq));
error->warning(FLERR,str,0);
if (rlogarg <= -3.0) error->one(FLERR,"Bad FENE bond");
if (check_error_thr((rlogarg <= -3.0),tid,FLERR,"Bad FENE bond"))
return;
rlogarg = 0.1;
}

8
src/USER-OMP/neigh_derive_omp.cpp
View File

@ -86,8 +86,8 @@ void Neighbor::half_from_full_no_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = inum_full;
@ -175,8 +175,8 @@ void Neighbor::half_from_full_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = inum_full;

22
src/USER-OMP/neigh_full_omp.cpp
View File

@ -105,8 +105,8 @@ void Neighbor::full_nsq_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -213,8 +213,8 @@ void Neighbor::full_nsq_ghost_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -314,8 +314,8 @@ void Neighbor::full_bin_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -324,7 +324,7 @@ void Neighbor::full_bin_omp(NeighList *list)
/* ----------------------------------------------------------------------
binned neighbor list construction for all neighbors
include neighbors of ghost atoms (no "special neighbors" for ghosts)
include neighbors of ghost atoms, but no "special neighbors" for ghosts
every neighbor pair appears in list of both atoms i and j
------------------------------------------------------------------------- */
@ -446,8 +446,8 @@ void Neighbor::full_bin_ghost_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -555,8 +555,8 @@ void Neighbor::full_multi_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;

20
src/USER-OMP/neigh_gran_omp.cpp
View File

@ -163,8 +163,8 @@ void Neighbor::granular_nsq_no_newton_omp(NeighList *list)
firstshear[i] = shearptr;
}
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -266,8 +266,8 @@ void Neighbor::granular_nsq_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -427,8 +427,8 @@ void Neighbor::granular_bin_no_newton_omp(NeighList *list)
firstshear[i] = shearptr;
}
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -544,8 +544,8 @@ void Neighbor::granular_bin_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -649,8 +649,8 @@ void Neighbor::granular_bin_newton_tri_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;

12
src/USER-OMP/neigh_half_bin_omp.cpp
View File

@ -118,8 +118,8 @@ void Neighbor::half_bin_no_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -245,8 +245,8 @@ void Neighbor::half_bin_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -356,8 +356,8 @@ void Neighbor::half_bin_newton_tri_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;

12
src/USER-OMP/neigh_half_multi_omp.cpp
View File

@ -126,8 +126,8 @@ void Neighbor::half_multi_no_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -263,8 +263,8 @@ void Neighbor::half_multi_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -384,8 +384,8 @@ void Neighbor::half_multi_newton_tri_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;

8
src/USER-OMP/neigh_half_nsq_omp.cpp
View File

@ -105,8 +105,8 @@ void Neighbor::half_nsq_no_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;
@ -214,8 +214,8 @@ void Neighbor::half_nsq_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
NEIGH_OMP_CLOSE;
list->inum = nlocal;

80
src/USER-OMP/neigh_respa_omp.cpp
View File

@ -52,7 +52,7 @@ void Neighbor::respa_nsq_no_newton_omp(NeighList *list)
#endif
NEIGH_OMP_SETUP(nlocal);
int i,j,n,itype,jtype,which,n_inner,n_middle;
int i,j,n,itype,jtype,n_inner,n_middle;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
int *neighptr,*neighptr_inner,*neighptr_middle;
@ -91,6 +91,8 @@ void Neighbor::respa_nsq_no_newton_omp(NeighList *list)
int npage_middle = tid;
int npnt_middle = 0;
int which = 0;
for (i = ifrom; i < ito; i++) {
#if defined(_OPENMP)
@ -168,23 +170,23 @@ void Neighbor::respa_nsq_no_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
ilist_inner[i] = i;
firstneigh_inner[i] = neighptr_inner;
numneigh_inner[i] = n_inner;
npnt_inner += n_inner;
if (npnt_inner >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_inner > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
if (respamiddle) {
ilist_middle[i] = i;
firstneigh_middle[i] = neighptr_middle;
numneigh_middle[i] = n_middle;
npnt_middle += n_middle;
if (npnt_middle >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_middle > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
}
NEIGH_OMP_CLOSE;
@ -225,7 +227,7 @@ void Neighbor::respa_nsq_newton_omp(NeighList *list)
#endif
NEIGH_OMP_SETUP(nlocal);
int i,j,n,itype,jtype,itag,jtag,which,n_inner,n_middle;
int i,j,n,itype,jtype,itag,jtag,n_inner,n_middle;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
int *neighptr,*neighptr_inner,*neighptr_middle;
@ -264,6 +266,8 @@ void Neighbor::respa_nsq_newton_omp(NeighList *list)
int npage_middle = tid;
int npnt_middle = 0;
int which = 0;
for (i = ifrom; i < ito; i++) {
#if defined(_OPENMP)
@ -359,23 +363,23 @@ void Neighbor::respa_nsq_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
ilist_inner[i] = i;
firstneigh_inner[i] = neighptr_inner;
numneigh_inner[i] = n_inner;
npnt_inner += n_inner;
if (npnt_inner >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_inner > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
if (respamiddle) {
ilist_middle[i] = i;
firstneigh_middle[i] = neighptr_middle;
numneigh_middle[i] = n_middle;
npnt_middle += n_middle;
if (npnt_middle >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_middle > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
}
NEIGH_OMP_CLOSE;
@ -419,7 +423,7 @@ void Neighbor::respa_bin_no_newton_omp(NeighList *list)
#endif
NEIGH_OMP_SETUP(nlocal);
int i,j,k,n,itype,jtype,ibin,which,n_inner,n_middle;
int i,j,k,n,itype,jtype,ibin,n_inner,n_middle;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
int *neighptr,*neighptr_inner,*neighptr_middle;
@ -459,6 +463,8 @@ void Neighbor::respa_bin_no_newton_omp(NeighList *list)
int npage_middle = tid;
int npnt_middle = 0;
int which = 0;
for (i = ifrom; i < ito; i++) {
#if defined(_OPENMP)
@ -545,23 +551,23 @@ void Neighbor::respa_bin_no_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
ilist_inner[i] = i;
firstneigh_inner[i] = neighptr_inner;
numneigh_inner[i] = n_inner;
npnt_inner += n_inner;
if (npnt_inner >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_inner > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
if (respamiddle) {
ilist_middle[i] = i;
firstneigh_middle[i] = neighptr_middle;
numneigh_middle[i] = n_middle;
npnt_middle += n_middle;
if (npnt_middle >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_middle > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
}
NEIGH_OMP_CLOSE;
@ -604,7 +610,7 @@ void Neighbor::respa_bin_newton_omp(NeighList *list)
#endif
NEIGH_OMP_SETUP(nlocal);
int i,j,k,n,itype,jtype,ibin,which,n_inner,n_middle;
int i,j,k,n,itype,jtype,ibin,n_inner,n_middle;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
int *neighptr,*neighptr_inner,*neighptr_middle;
@ -644,6 +650,8 @@ void Neighbor::respa_bin_newton_omp(NeighList *list)
int npage_middle = tid;
int npnt_middle = 0;
int which = 0;
for (i = ifrom; i < ito; i++) {
#if defined(_OPENMP)
@ -766,23 +774,23 @@ void Neighbor::respa_bin_newton_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
ilist_inner[i] = i;
firstneigh_inner[i] = neighptr_inner;
numneigh_inner[i] = n_inner;
npnt_inner += n_inner;
if (npnt_inner >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_inner > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
if (respamiddle) {
ilist_middle[i] = i;
firstneigh_middle[i] = neighptr_middle;
numneigh_middle[i] = n_middle;
npnt_middle += n_middle;
if (npnt_middle >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_middle > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
}
NEIGH_OMP_CLOSE;
@ -825,7 +833,7 @@ void Neighbor::respa_bin_newton_tri_omp(NeighList *list)
#endif
NEIGH_OMP_SETUP(nlocal);
int i,j,k,n,itype,jtype,ibin,which,n_inner,n_middle;
int i,j,k,n,itype,jtype,ibin,n_inner,n_middle;
double xtmp,ytmp,ztmp,delx,dely,delz,rsq;
int *neighptr,*neighptr_inner,*neighptr_middle;
@ -865,6 +873,8 @@ void Neighbor::respa_bin_newton_tri_omp(NeighList *list)
int npage_middle = tid;
int npnt_middle = 0;
int which = 0;
for (i = ifrom; i < ito; i++) {
#if defined(_OPENMP)
@ -959,23 +969,23 @@ void Neighbor::respa_bin_newton_tri_omp(NeighList *list)
firstneigh[i] = neighptr;
numneigh[i] = n;
npnt += n;
if (n > oneatom || npnt >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
ilist_inner[i] = i;
firstneigh_inner[i] = neighptr_inner;
numneigh_inner[i] = n_inner;
npnt_inner += n_inner;
if (npnt_inner >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_inner > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
if (respamiddle) {
ilist_middle[i] = i;
firstneigh_middle[i] = neighptr_middle;
numneigh_middle[i] = n_middle;
npnt_middle += n_middle;
if (npnt_middle >= pgsize)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one or page");
if (n_middle > oneatom)
error->one(FLERR,"Neighbor list overflow, boost neigh_modify one");
}
}
NEIGH_OMP_CLOSE;

7
src/USER-OMP/neighbor_omp.h
View File

@ -1,4 +1,4 @@
/* ----------------------------------------------------------------------
/* -*- c++ -*- ----------------------------------------------------------
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
http://lammps.sandia.gov, Sandia National Laboratories
Steve Plimpton, sjplimp@sandia.gov
@ -35,9 +35,8 @@ namespace LAMMPS_NS {
const int tid = omp_get_thread_num(); \
const int idelta = 1 + num/nthreads; \
const int ifrom = tid*idelta; \
int ito = ifrom + idelta; \
if (ito > num) \
ito = num
const int ito = ((ifrom + idelta) > num) \
? num : (ifrom+idelta); \
#define NEIGH_OMP_CLOSE }

6
src/USER-OMP/pair_airebo_omp.cpp
View File

@ -2653,7 +2653,6 @@ void PairAIREBOOMP::TORSION_thr(int ifrom, int ito,
void PairAIREBOOMP::REBO_neigh_thr()
{
const int nlocal = atom->nlocal;
const int nall = nlocal + atom->nghost;
const int nthreads = comm->nthreads;
if (atom->nmax > maxlocal) {
@ -2697,9 +2696,7 @@ void PairAIREBOOMP::REBO_neigh_thr()
const int iidelta = 1 + allnum/nthreads;
const int iifrom = tid*iidelta;
int iito = iifrom + iidelta;
if (iito > allnum)
iito = allnum;
const int iito = ((iifrom+iidelta)>allnum) ? allnum : (iifrom+iidelta);
// store all REBO neighs of owned and ghost atoms
// scan full neighbor list of I
@ -2757,6 +2754,7 @@ void PairAIREBOOMP::REBO_neigh_thr()
}
}
/* ---------------------------------------------------------------------- */
double PairAIREBOOMP::memory_usage()

7
src/USER-OMP/pair_colloid_omp.cpp
View File

@ -85,6 +85,7 @@ void PairColloidOMP::eval(int iifrom, int iito, ThrData * const thr)
double * const * const f = thr->get_f();
const int * const type = atom->type;
const int nlocal = atom->nlocal;
const int tid = thr->get_tid();
const double * const special_lj = force->special_lj;
double fxtmp,fytmp,fztmp;
@ -146,7 +147,11 @@ void PairColloidOMP::eval(int iifrom, int iito, ThrData * const thr)
evdwl = 2.0/9.0*fR *
(1.0-(K[1]*(K[1]*(K[1]/3.0+3.0*K[2])+4.2*K[4])+K[2]*K[4]) *
sigma6[itype][jtype]/K[6]) - offset[itype][jtype];
if (rsq <= K[1]) error->one(FLERR,"Overlapping small/large in pair colloid");
if (check_error_thr((rsq <= K[1]),tid,FLERR,
"Overlapping small/large in pair colloid"))
return;
break;
case LARGE_LARGE:

101
src/USER-OMP/pair_comb_omp.cpp
View File

@ -22,12 +22,10 @@
#include "neighbor.h"
#include "neigh_list.h"
#if defined(_OPENMP)
#include <omp.h>
#endif
using namespace LAMMPS_NS;
#define MAXNEIGH 24
/* ---------------------------------------------------------------------- */
PairCombOMP::PairCombOMP(LAMMPS *lmp) :
@ -48,17 +46,8 @@ void PairCombOMP::compute(int eflag, int vflag)
const int nthreads = comm->nthreads;
const int inum = list->inum;
// grow coordination array if necessary
if (atom->nmax > nmax) {
memory->destroy(NCo);
memory->destroy(bbij);
nmax = atom->nmax;
memory->create(NCo,nmax,"pair:NCo");
memory->create(bbij,nmax,nmax,"pair:bbij");
}
// Build short range neighbor list
Short_neigh_thr();
#if defined(_OPENMP)
@ -102,7 +91,7 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
double yaself;
double potal,fac11,fac11e;
double vionij,fvionij,sr1,sr2,sr3,Eov,Fov;
int sht_jnum, *sht_jlist;
int sht_jnum, *sht_jlist, nj;
evdwl = ecoul = 0.0;
@ -139,7 +128,8 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
fxtmp = fytmp = fztmp = 0.0;
iq = q[i];
NCo[i] = 0;
NCo[i] = 0;
nj = 0;
iparam_i = elem2param[itype][itype][itype];
// self energy, only on i atom
@ -242,7 +232,6 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
int numcoor = 0;
for (jj = 0; jj < sht_jnum; jj++) {
j = sht_jlist[jj];
j &= NEIGHMASK;
jtype = map[type[j]];
iparam_ij = elem2param[itype][jtype][jtype];
@ -264,7 +253,6 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
for (jj = 0; jj < sht_jnum; jj++) {
j = sht_jlist[jj];
j &= NEIGHMASK;
jtype = map[type[j]];
iparam_ij = elem2param[itype][jtype][jtype];
@ -279,6 +267,7 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
rsq1 = vec3_dot(delr1,delr1);
if (rsq1 > params[iparam_ij].cutsq) continue;
nj ++;
// accumulate bondorder zeta for each i-j interaction via loop over k
@ -289,7 +278,6 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
for (kk = 0; kk < sht_jnum; kk++) {
k = sht_jlist[kk];
if (j == k) continue;
k &= NEIGHMASK;
ktype = map[type[k]];
iparam_ijk = elem2param[itype][jtype][ktype];
@ -309,7 +297,7 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
if (cuo_flag1 && cuo_flag2) cuo_flag = 1;
else cuo_flag = 0;
force_zeta(&params[iparam_ij],EFLAG,i,j,rsq1,zeta_ij,
force_zeta(&params[iparam_ij],EFLAG,i,nj,rsq1,zeta_ij,
iq,jq,fpair,prefactor,evdwl);
// over-coordination correction for HfO2
@ -334,7 +322,6 @@ void PairCombOMP::eval(int iifrom, int iito, ThrData * const thr)
for (kk = 0; kk < sht_jnum; kk++) {
k = sht_jlist[kk];
if (j == k) continue;
k &= NEIGHMASK;
ktype = map[type[k]];
iparam_ijk = elem2param[itype][jtype][ktype];
@ -429,6 +416,7 @@ double PairCombOMP::yasu_char(double *qf_fix, int &igroup)
const int i = ilist[ii];
const int itag = tag[i];
int nj = 0;
if (mask[i] & groupbit) {
fqi = fqj = fqij = fqji = fqjj = 0.0; // should not be needed.
@ -511,10 +499,11 @@ double PairCombOMP::yasu_char(double *qf_fix, int &igroup)
const int iparam_ij = elem2param[itype][jtype][jtype];
if (rsq1 > params[iparam_ij].cutsq) continue;
nj ++;
// charge force in Aij and Bij
qfo_short(&params[iparam_ij],i,j,rsq1,iq,jq,fqij,fqjj);
qfo_short(&params[iparam_ij],i,nj,rsq1,iq,jq,fqij,fqjj);
fqi += fqij;
#if defined(_OPENMP)
#pragma omp atomic
@ -522,7 +511,7 @@ double PairCombOMP::yasu_char(double *qf_fix, int &igroup)
qf[j] += fqjj;
}
#if defined(_OPENMP) && 0
#if defined(_OPENMP)
#pragma omp atomic
#endif
qf[i] += fqi;
@ -559,12 +548,13 @@ void PairCombOMP::Short_neigh_thr()
{
if (atom->nmax > nmax) {
nmax = int(1.0 * atom->nmax);
memory->sfree(sht_num);
nmax = atom->nmax;
memory->sfree(sht_first);
memory->create(sht_num,nmax,"pair:sht_num");
sht_first = (int **) memory->smalloc(nmax*sizeof(int *),
"pair:sht_first");
memory->grow(sht_num,nmax,"pair:sht_num");
memory->grow(NCo,nmax,"pair:NCo");
memory->grow(bbij,nmax,nmax,"pair:bbij");
}
const int nthreads = comm->nthreads;
@ -602,46 +592,47 @@ void PairCombOMP::Short_neigh_thr()
int npage = tid;
npntj = 0;
for (ii = iifrom; ii < iito; ii++) {
i = ilist[ii];
itype = type[i];
if (iifrom < inum) {
for (ii = iifrom; ii < iito; ii++) {
i = ilist[ii];
itype = type[i];
#if defined(_OPENMP)
#pragma omp critical
#endif
if(pgsize - npntj < oneatom) {
npntj = 0;
npage++;
if (npage == maxpage) add_pages(nthreads);
}
if(pgsize - npntj < oneatom) {
npntj = 0;
npage += nthreads;
if (npage >= maxpage) add_pages(nthreads);
}
neighptrj = &pages[npage][npntj];
nj = 0;
neighptrj = &pages[npage][npntj];
nj = 0;
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
xtmp = x[i][0];
ytmp = x[i][1];
ztmp = x[i][2];
jlist = firstneigh[i];
jnum = numneigh[i];
jlist = firstneigh[i];
jnum = numneigh[i];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
iparam_ij = elem2param[itype][jtype][jtype];
for (jj = 0; jj < jnum; jj++) {
j = jlist[jj];
j &= NEIGHMASK;
jtype = type[j];
delrj[0] = xtmp - x[j][0];
delrj[1] = ytmp - x[j][1];
delrj[2] = ztmp - x[j][2];
rsq = vec3_dot(delrj,delrj);
delrj[0] = xtmp - x[j][0];
delrj[1] = ytmp - x[j][1];
delrj[2] = ztmp - x[j][2];
rsq = vec3_dot(delrj,delrj);
if (rsq > cutmin) continue;
neighptrj[nj++] = j;
if (rsq > cutmin) continue;
neighptrj[nj++] = j;
}
sht_first[i] = neighptrj;
sht_num[i] = nj;
npntj += nj;
}
sht_first[i] = neighptrj;
sht_num[i] = nj;
npntj += nj;
}
}
}

8
src/USER-OMP/pair_lj_cut_coul_long_tip4p_omp.cpp
View File

@ -151,6 +151,7 @@ void PairLJCutCoulLongTIP4POMP::eval(int iifrom, int iito, ThrData * const thr)
const double * const q = atom->q;
const int * const type = atom->type;
const int nlocal = atom->nlocal;
const int tid = thr->get_tid();
const double * const special_coul = force->special_coul;
const double * const special_lj = force->special_lj;
const double qqrd2e = force->qqrd2e;
@ -221,8 +222,11 @@ void PairLJCutCoulLongTIP4POMP::eval(int iifrom, int iito, ThrData * const thr)
x2 = mpos[j];
jH1 = h1idx[j];
jH2 = h2idx[j];
if (jtype == typeO && ( jH1 < 0 || jH2 < 0 ))
error->one(FLERR,"TIP4P hydrogen is missing");
if (check_error_thr((jtype == typeO && ( jH1 < 0 || jH2 < 0 )),
tid, FLERR,"TIP4P hydrogen is missing"))
return;
delx = x1[0] - x2[0];
dely = x1[1] - x2[1];
delz = x1[2] - x2[2];

8
src/USER-OMP/pair_lj_cut_coul_pppm_tip4p_omp.cpp
View File

@ -179,6 +179,7 @@ void PairLJCutCoulPPPMTIP4POMP::eval(int iifrom, int iito, ThrData * const thr)
const double * const q = atom->q;
const int * const type = atom->type;
const int nlocal = atom->nlocal;
const int tid = thr->get_tid();
const double * const special_coul = force->special_coul;
const double * const special_lj = force->special_lj;
const double qqrd2e = force->qqrd2e;
@ -249,8 +250,11 @@ void PairLJCutCoulPPPMTIP4POMP::eval(int iifrom, int iito, ThrData * const thr)
x2 = mpos[j];
jH1 = h1idx[j];
jH2 = h2idx[j];
if (jtype == typeO && ( jH1 < 0 || jH2 < 0 ))
error->one(FLERR,"TIP4P hydrogen is missing");
if (check_error_thr((jtype == typeO && ( jH1 < 0 || jH2 < 0 )),
tid, FLERR,"TIP4P hydrogen is missing"))
return;
delx = x1[0] - x2[0];
dely = x1[1] - x2[1];
delz = x1[2] - x2[2];

10
src/USER-OMP/pair_peri_lps_omp.cpp
View File

@ -213,16 +213,15 @@ void PairPeriLPSOMP::eval(int iifrom, int iito, ThrData * const thr)
// each thread works on a fixed chunk of atoms.
const int idelta = 1 + nlocal/comm->nthreads;
iifrom = thr->get_tid()*idelta;
iito = iifrom + idelta;
if (iito > nlocal)
iito = nlocal;
iito = ((iifrom + idelta) > nlocal) ? nlocal : (iifrom + idelta);
#else
iifrom = 0;
iito = nlocal;
#endif
// Compute the dilatation on each particle
compute_dilatation_thr(iifrom, iito);
if (iifrom < nlocal)
compute_dilatation_thr(iifrom, iito);
// wait until all threads are done before communication
sync_threads();
@ -350,7 +349,8 @@ void PairPeriLPSOMP::eval(int iifrom, int iito, ThrData * const thr)
sync_threads();
// store new s0 (in parallel)
for (i = iifrom; i < iito; i++) s0[i] = s0_new[i];
if (iifrom < nlocal)
for (i = iifrom; i < iito; i++) s0[i] = s0_new[i];
}
/* ---------------------------------------------------------------------- */

17
src/USER-OMP/pair_peri_pmb_omp.cpp
View File

@ -203,9 +203,7 @@ void PairPeriPMBOMP::eval(int iifrom, int iito, ThrData * const thr)
// each thread works on a fixed chunk of atoms.
const int idelta = 1 + nlocal/comm->nthreads;
iifrom = thr->get_tid()*idelta;
iito = iifrom + idelta;
if (iito > nlocal)
iito = nlocal;
iito = ((iifrom + idelta) > nlocal) ? nlocal : (iifrom + idelta);
#else
iifrom = 0;
iito = nlocal;
@ -234,8 +232,8 @@ void PairPeriPMBOMP::eval(int iifrom, int iito, ThrData * const thr)
// check if lost a partner without first breaking bond
if (j < 0) {
partner[i][jj] = 0;
continue;
partner[i][jj] = 0;
continue;
}
// compute force density, add to PD equation of motion
@ -257,7 +255,7 @@ void PairPeriPMBOMP::eval(int iifrom, int iito, ThrData * const thr)
// scale vfrac[j] if particle j near the horizon
if ((fabs(r0[i][jj] - delta)) <= half_lc)
vfrac_scale = (-1.0/(2*half_lc))*(r0[i][jj]) +
vfrac_scale = (-1.0/(2*half_lc))*(r0[i][jj]) +
(1.0 + ((delta - half_lc)/(2*half_lc) ) );
else vfrac_scale = 1.0;
@ -285,9 +283,9 @@ void PairPeriPMBOMP::eval(int iifrom, int iito, ThrData * const thr)
// update s0 for next timestep
if (first)
s0_new[i] = s00[itype][jtype] - (alpha[itype][jtype] * stretch);
s0_new[i] = s00[itype][jtype] - (alpha[itype][jtype] * stretch);
else
s0_new[i] = MAX(s0_new[i],s00[itype][jtype] - (alpha[itype][jtype] * stretch));
s0_new[i] = MAX(s0_new[i],s00[itype][jtype] - (alpha[itype][jtype] * stretch));
first = false;
}
@ -296,7 +294,8 @@ void PairPeriPMBOMP::eval(int iifrom, int iito, ThrData * const thr)
sync_threads();
// store new s0 (in parallel)
for (i = iifrom; i < iito; i++) s0[i] = s0_new[i];
if (iifrom < nlocal)
for (i = iifrom; i < iito; i++) s0[i] = s0_new[i];
}
/* ---------------------------------------------------------------------- */

28
src/USER-OMP/pair_table_omp.cpp
View File

@ -88,6 +88,7 @@ void PairTableOMP::eval(int iifrom, int iito, ThrData * const thr)
double * const * const f = thr->get_f();
const int * const type = atom->type;
const int nlocal = atom->nlocal;
const int tid = thr->get_tid();
const double * const special_lj = force->special_lj;
double fxtmp,fytmp,fztmp;
@ -121,25 +122,36 @@ void PairTableOMP::eval(int iifrom, int iito, ThrData * const thr)
if (rsq < cutsq[itype][jtype]) {
tb = &tables[tabindex[itype][jtype]];
if (rsq < tb->innersq)
error->one(FLERR,"Pair distance < table inner cutoff");
if (check_error_thr((rsq < tb->innersq),tid,
FLERR,"Pair distance < table inner cutoff"))
return;
if (tabstyle == LOOKUP) {
itable = static_cast<int> ((rsq - tb->innersq) * tb->invdelta);
if (itable >= tlm1)
error->one(FLERR,"Pair distance > table outer cutoff");
if (check_error_thr((itable >= tlm1),tid,
FLERR,"Pair distance > table outer cutoff"))
return;
fpair = 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");
if (check_error_thr((itable >= tlm1),tid,
FLERR,"Pair distance > table outer cutoff"))
return;
fraction = (rsq - tb->rsq[itable]) * tb->invdelta;
value = tb->f[itable] + fraction*tb->df[itable];
fpair = 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");
if (check_error_thr((itable >= tlm1),tid,
FLERR,"Pair distance > table outer cutoff"))
return;
b = (rsq - tb->rsq[itable]) * tb->invdelta;
a = 1.0 - b;
value = a * tb->f[itable] + b * tb->f[itable+1] +

8
src/USER-OMP/pair_tersoff_table_omp.cpp
View File

@ -16,6 +16,7 @@
#include "pair_tersoff_table_omp.h"
#include "atom.h"
#include "comm.h"
#include "error.h"
#include "force.h"
#include "neighbor.h"
#include "neigh_list.h"
@ -89,6 +90,7 @@ void PairTersoffTableOMP::eval(int iifrom, int iito, ThrData * const thr)
double * const * const f = thr->get_f();
const int * const type = atom->type;
const int nlocal = atom->nlocal;
const int tid = thr->get_tid();
inum = list->inum;
ilist = list->ilist;
@ -108,6 +110,12 @@ void PairTersoffTableOMP::eval(int iifrom, int iito, ThrData * const thr)
jlist = firstneigh[i];
jnum = numneigh[i];
if (check_error_thr((jnum > leadingDimensionInteractionList), tid,
FLERR,"Too many neighbors for interaction list.\n"
"Check your system or increase 'leadingDimension"
"InteractionList'"))
return;
// Pre-calculate gteta and cutoff function
for (int neighbor_j = 0; neighbor_j < jnum; neighbor_j++) {

24
src/USER-OMP/thr_omp.cpp
View File

@ -42,7 +42,8 @@ using namespace MathConst;
/* ---------------------------------------------------------------------- */
ThrOMP::ThrOMP(LAMMPS *ptr, int style) : lmp(ptr), fix(NULL), thr_style(style)
ThrOMP::ThrOMP(LAMMPS *ptr, int style) : lmp(ptr), fix(NULL),
thr_style(style), thr_error(0)
{
// register fix omp with this class
int ifix = lmp->modify->find_fix("package_omp");
@ -66,6 +67,7 @@ void ThrOMP::ev_setup_thr(int eflag, int vflag, int nall, double *eatom,
double **vatom, ThrData *thr)
{
const int tid = thr->get_tid();
if (tid == 0) thr_error = 0;
if (thr_style & THR_PAIR) {
if (eflag & 2) {
@ -171,6 +173,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
double **f = lmp->atom->f;
double **x = lmp->atom->x;
if (evflag)
sync_threads();
switch (thr_style) {
case THR_PAIR: {
@ -201,11 +206,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
}
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(pair->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(pair->vatom[0][0]), nall, nthreads, 6, tid);
}
}
@ -248,11 +251,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
}
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(pair->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(pair->vatom[0][0]), nall, nthreads, 6, tid);
}
}
@ -272,11 +273,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
bond->virial[i] += thr->virial_bond[i];
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(bond->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(bond->vatom[0][0]), nall, nthreads, 6, tid);
}
}
@ -295,11 +294,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
angle->virial[i] += thr->virial_angle[i];
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(angle->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(angle->vatom[0][0]), nall, nthreads, 6, tid);
}
}
@ -318,11 +315,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
dihedral->virial[i] += thr->virial_dihed[i];
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(dihedral->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(dihedral->vatom[0][0]), nall, nthreads, 6, tid);
}
}
@ -351,12 +346,10 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
}
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(dihedral->eatom[0]), nall, nthreads, 1, tid);
data_reduce_thr(&(pair->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(dihedral->vatom[0][0]), nall, nthreads, 6, tid);
data_reduce_thr(&(pair->vatom[0][0]), nall, nthreads, 6, tid);
}
@ -376,11 +369,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
improper->virial[i] += thr->virial_imprp[i];
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(improper->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(improper->vatom[0][0]), nall, nthreads, 6, tid);
}
}
@ -401,11 +392,9 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
kspace->virial[i] += thr->virial_kspce[i];
}
if (eflag & 2) {
sync_threads();
data_reduce_thr(&(kspace->eatom[0]), nall, nthreads, 1, tid);
}
if (vflag & 4) {
sync_threads();
data_reduce_thr(&(kspace->vatom[0][0]), nall, nthreads, 6, tid);
}
}
@ -418,7 +407,6 @@ void ThrOMP::reduce_thr(void *style, const int eflag, const int vflag,
}
if (style == fix->last_omp_style) {
sync_threads();
data_reduce_thr(&(f[0][0]), nall, nthreads, 3, tid);
if (lmp->atom->torque)
data_reduce_thr(&(lmp->atom->torque[0][0]), nall, nthreads, 3, tid);

31
src/USER-OMP/thr_omp.h
View File

@ -19,6 +19,7 @@
#define LMP_THR_OMP_H
#include "pointers.h"
#include "error.h"
#include "fix_omp.h"
#include "thr_data.h"
@ -40,6 +41,7 @@ class ThrOMP {
FixOMP *fix; // pointer to fix_omp;
const int thr_style;
int thr_error;
public:
ThrOMP(LAMMPS *, int);
@ -68,7 +70,34 @@ class ThrOMP {
const int, const int, const int);
// reduce per thread data as needed
void reduce_thr(void * const style, const int eflag, const int vflag, ThrData * const thr, const int nproxy=0);
void reduce_thr(void * const style, const int eflag, const int vflag,
ThrData * const thr, const int nproxy=0);
// thread safe variant error abort support.
// signals an error condition in any thread by making
// thr_error > 0, if condition "cond" is true.
// will abort from thread 0 if thr_error is > 0
// otherwise return true.
// returns false if no error found on any thread.
// use return value to jump/return to end of threaded region.
bool check_error_thr(const bool cond, const int tid, const char *fname,
const int line, const char *errmsg) {
if (cond) {
#if defined(_OPENMP)
#pragma omp atomic
++thr_error;
#endif
if (tid > 0) return true;
else lmp->error->one(fname,line,errmsg);
} else {
if (thr_error > 0) {
if (tid == 0) lmp->error->one(fname,line,errmsg);
else return true;
} else return false;
}
return false;
};
protected:

3
src/memory.cpp
View File

@ -34,7 +34,8 @@ void *Memory::smalloc(bigint nbytes, const char *name)
#if defined(LAMMPS_MEMALIGN)
void *ptr;
posix_memalign(&ptr, LAMMPS_MEMALIGN, nbytes);
int retval = posix_memalign(&ptr, LAMMPS_MEMALIGN, nbytes);
if (retval) ptr = NULL;
#else
void *ptr = malloc(nbytes);
#endif

33
tools/restart2data.cpp
View File

@ -475,7 +475,7 @@ int main (int narg, char **arg)
void header(FILE *fp, Data &data)
{
const char *version = "5 Jan 2012";
const char *version = "31 Jan 2012";
data.triclinic = 0;
@ -1506,6 +1506,8 @@ void pair(FILE *fp, Data &data, char *style, int flag)
}
}
} else if (strcmp(style,"comb") == 0) {
} else if (strcmp(style,"coul/diel") == 0) {
m = 1;
double cut_coul = read_double(fp);
@ -2534,6 +2536,16 @@ void angle(FILE *fp, Data &data)
} else if ((strcmp(data.angle_style,"cg/cmm") == 0) ||
(strcmp(data.angle_style,"sdk") == 0)) {
data.angle_harmonic_k = new double[data.nangletypes+1];
data.angle_harmonic_theta0 = new double[data.nangletypes+1];
data.angle_cg_cmm_epsilon = new double[data.nangletypes+1];
fread(&data.angle_harmonic_k[1],sizeof(double),data.nangletypes,fp);
fread(&data.angle_harmonic_theta0[1],sizeof(double),data.nangletypes,fp);
fread(&data.angle_cg_cmm_epsilon[1],sizeof(double),data.nangletypes,fp);
} else if (strcmp(data.angle_style,"cg/cmm/old") == 0) {
data.angle_harmonic_k = new double[data.nangletypes+1];
data.angle_harmonic_theta0 = new double[data.nangletypes+1];
data.angle_cg_cmm_epsilon = new double[data.nangletypes+1];
@ -3046,6 +3058,7 @@ void Data::write(FILE *fp, FILE *fp2)
(strcmp(pair_style,"adp") != 0) &&
(strcmp(pair_style,"airebo") != 0) &&
(strcmp(pair_style,"brownian") != 0) &&
(strcmp(pair_style,"comb") != 0) &&
(strcmp(pair_style,"coul/cut") != 0) &&
(strcmp(pair_style,"coul/debye") != 0) &&
(strcmp(pair_style,"coul/diel") != 0) &&
@ -3181,7 +3194,7 @@ void Data::write(FILE *fp, FILE *fp2)
(strcmp(pair_style,"cg/cmm/coul/long") == 0) ||
(strcmp(pair_style,"lj/sdk") == 0) ||
(strcmp(pair_style,"lj/sdk/coul/long") == 0)) {
printf("ERROR: Cannot write pair_style %s to data file\n",
printf("ERROR: Cannot write pair_style %s to data file alone. please provide an input file, too.\n",
pair_style);
exit(1);
}
@ -3351,7 +3364,13 @@ void Data::write(FILE *fp, FILE *fp2)
fprintf(fp,"%d %g %g\n",i,
angle_harmonic_k[i],angle_harmonic_theta0[i]/PI*180.0);
} else if (strcmp(angle_style,"cg/cmm") == 0) {
} else if ((strcmp(angle_style,"cg/cmm") == 0) ||
(strcmp(angle_style,"sdk") == 0)) {
for (int i = 1; i <= nangletypes; i++)
fprintf(fp,"%d %g %g %g\n",i,angle_harmonic_k[i],
angle_harmonic_theta0[i]/PI*180.0,angle_cg_cmm_epsilon[i]);
} else if (strcmp(angle_style,"cg/cmm/old") == 0) {
for (int i = 1; i <= nangletypes; i++)
fprintf(fp,"%d %g %g %s %g %g\n",i,
angle_harmonic_k[i],angle_harmonic_theta0[i]/PI*180.0,
@ -3376,7 +3395,13 @@ void Data::write(FILE *fp, FILE *fp2)
fprintf(fp2,"angle_coeffs %d %g %g\n",i,
angle_harmonic_k[i],angle_harmonic_theta0[i]/PI*180.0);
} else if (strcmp(angle_style,"cg/cmm") == 0) {
} else if ((strcmp(angle_style,"cg/cmm") == 0) ||
(strcmp(angle_style,"sdk") == 0)) {
for (int i = 1; i <= nangletypes; i++)
fprintf(fp2,"angle_coeffs %d %g %g %g\n",i,angle_harmonic_k[i],
angle_harmonic_theta0[i]/PI*180.0,angle_cg_cmm_epsilon[i]);
} else if (strcmp(angle_style,"cg/cmm/old") == 0) {
for (int i = 1; i <= nangletypes; i++)
fprintf(fp2,"angle_coeffs %d %g %g %s %g %g\n",i,
angle_harmonic_k[i],angle_harmonic_theta0[i]/PI*180.0,