From ee00c0f654115629bb9de70254d73d86e888a9a9 Mon Sep 17 00:00:00 2001 From: sjplimp Date: Wed, 4 Mar 2015 15:24:51 +0000 Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@13169 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- src/Makefile | 2 +- src/QEQ/fix_qeq_slater.cpp | 4 +-- src/QEQ/fix_qeq_slater.h | 2 +- src/compute_chunk_atom.cpp | 4 +-- src/compute_chunk_atom.h | 2 +- src/pair_coul_streitz.cpp | 34 ++++++++++------------- src/velocity.cpp | 56 +++++++++++++++++++++++--------------- 7 files changed, 55 insertions(+), 49 deletions(-) diff --git a/src/Makefile b/src/Makefile index 6fafae1fd8..82d1182897 100755 --- a/src/Makefile +++ b/src/Makefile @@ -13,7 +13,7 @@ OBJ = $(SRC:.cpp=.o) # Package variables -PACKAGE = asphere body class2 colloid dipole fld gpu granular kim \ +PACKAGE = asphere body class2 colloid coreshell dipole fld gpu granular kim \ kokkos kspace manybody mc meam misc molecule mpiio opt peri poems \ qeq reax replica rigid shock snap srd voronoi xtc diff --git a/src/QEQ/fix_qeq_slater.cpp b/src/QEQ/fix_qeq_slater.cpp index 95d0079fe0..36b3b2c65f 100644 --- a/src/QEQ/fix_qeq_slater.cpp +++ b/src/QEQ/fix_qeq_slater.cpp @@ -57,7 +57,7 @@ FixQEqSlater::FixQEqSlater(LAMMPS *lmp, int narg, char **arg) : } else error->all(FLERR,"Illegal fix qeq/slater command"); } - if (streitz_flag) extract(); + if (streitz_flag) extract_streitz(); } /* ---------------------------------------------------------------------- */ @@ -88,7 +88,7 @@ void FixQEqSlater::init() /* ---------------------------------------------------------------------- */ -void FixQEqSlater::extract() +void FixQEqSlater::extract_streitz() { int ntypes = atom->ntypes; diff --git a/src/QEQ/fix_qeq_slater.h b/src/QEQ/fix_qeq_slater.h index c0174ee7c1..5bb31cfde4 100644 --- a/src/QEQ/fix_qeq_slater.h +++ b/src/QEQ/fix_qeq_slater.h @@ -37,7 +37,7 @@ class FixQEqSlater : public FixQEq { void compute_H(); double calculate_H(double, double, double, double, double &); double calculate_H_wolf(double, double, double, double, double &); - void extract(); + void extract_streitz(); class PairCoulStreitz *streitz; }; diff --git a/src/compute_chunk_atom.cpp b/src/compute_chunk_atom.cpp index e2bfc1cc27..d83643b047 100644 --- a/src/compute_chunk_atom.cpp +++ b/src/compute_chunk_atom.cpp @@ -373,7 +373,7 @@ ComputeChunkAtom::ComputeChunkAtom(LAMMPS *lmp, int narg, char **arg) : id_fix = NULL; fixstore = NULL; - if (compress) hash = new std::map(); + if (compress) hash = new std::map(); else hash = NULL; maxvar = 0; @@ -1010,7 +1010,7 @@ void ComputeChunkAtom::compress_chunk_ids() void ComputeChunkAtom::idring(int n, char *cbuf) { tagint *list = (tagint *) cbuf; - std::map *hash = cptr->hash; + std::map *hash = cptr->hash; for (int i = 0; i < n; i++) (*hash)[list[i]] = 0; } diff --git a/src/compute_chunk_atom.h b/src/compute_chunk_atom.h index 8f69735d92..aa42a1f0c4 100644 --- a/src/compute_chunk_atom.h +++ b/src/compute_chunk_atom.h @@ -86,7 +86,7 @@ class ComputeChunkAtom : public Compute { int molcheck; // one-time check if all molecule atoms in chunk int *exclude; // 1 if atom is not assigned to any chunk - std::map *hash; // store original chunks IDs before compression + std::map *hash; // store original chunks IDs before compression // static variable for ring communication callback to access class data // callback functions for ring communication diff --git a/src/pair_coul_streitz.cpp b/src/pair_coul_streitz.cpp index a715c1ebc8..f501be61cf 100644 --- a/src/pair_coul_streitz.cpp +++ b/src/pair_coul_streitz.cpp @@ -286,7 +286,7 @@ void PairCoulStreitz::read_file(char *file) // strip comment, skip line if blank - if (ptr = strchr(line,'#')) *ptr = '\0'; + if ((ptr = strchr(line,'#'))) *ptr = '\0'; nwords = atom->count_words(line); if (nwords == 0) continue; @@ -305,7 +305,7 @@ void PairCoulStreitz::read_file(char *file) if (eof) break; MPI_Bcast(&n,1,MPI_INT,0,world); MPI_Bcast(line,n,MPI_CHAR,0,world); - if (ptr = strchr(line,'#')) *ptr = '\0'; + if ((ptr = strchr(line,'#'))) *ptr = '\0'; nwords = atom->count_words(line); } @@ -316,7 +316,7 @@ void PairCoulStreitz::read_file(char *file) nwords = 0; words[nwords++] = strtok(line," \t\n\r\f"); - while (words[nwords++] = strtok(NULL," \t\n\r\f")) continue; + while ((words[nwords++] = strtok(NULL," \t\n\r\f"))) continue; // ielement = 1st args @@ -355,7 +355,7 @@ void PairCoulStreitz::read_file(char *file) void PairCoulStreitz::setup() { - int i,j,k,m,n; + int i,m,n; // set elem2param @@ -391,29 +391,27 @@ void PairCoulStreitz::setup() void PairCoulStreitz::compute(int eflag, int vflag) { int i, j, ii, jj, inum, jnum; - int itag, jtag, itype, jtype, iparam_i,iparam_j; + int itype, jtype, iparam_i,iparam_j; int *ilist, *jlist, *numneigh, **firstneigh; - int *tag = atom->tag; int *type = atom->type; int nlocal = atom->nlocal; int newton_pair = force->newton_pair; - double xtmp, ytmp, ztmp, evdwl, ecoul, fpair; + double xtmp, ytmp, ztmp, ecoul, fpair; double qi, qj, selfion, r, rsq, delr[3]; - double zei, zej, zi, zj, ci_jfi, dci_jfi, ci_fifj, dci_fifj; - double prefactor, forcecoul, factor_coul; + double zei, zej, zj, ci_jfi, dci_jfi, ci_fifj, dci_fifj; + double forcecoul, factor_coul; double **x = atom->x; double **f = atom->f; double *q = atom->q; - double qqrd2e = force->qqrd2e; double *special_coul = force->special_coul; - evdwl = ecoul = 0.0; + ecoul = 0.0; selfion = fpair = 0.0; ci_jfi = ci_fifj = dci_jfi = dci_fifj = 0.0; - prefactor = forcecoul = 0.0; + forcecoul = 0.0; if (eflag || vflag) ev_setup(eflag,vflag); else evflag = vflag_fdotr = vflag_atom = 0; @@ -429,7 +427,6 @@ void PairCoulStreitz::compute(int eflag, int vflag) for (ii = 0; ii < inum; ii++) { i = ilist[ii]; - itag = tag[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -437,7 +434,6 @@ void PairCoulStreitz::compute(int eflag, int vflag) iparam_i = elem2param[itype]; qi = q[i]; zei = params[iparam_i].zeta; - zi = params[iparam_i].zcore; // self energy: ionization + wolf sum @@ -454,7 +450,6 @@ void PairCoulStreitz::compute(int eflag, int vflag) j = jlist[jj]; j &= NEIGHMASK; - jtag = tag[j]; jtype = map[type[j]]; iparam_j = elem2param[jtype]; qj = q[j]; @@ -505,7 +500,6 @@ void PairCoulStreitz::compute(int eflag, int vflag) for (ii = 0; ii < inum; ii++) { i = ilist[ii]; - itag = tag[i]; xtmp = x[i][0]; ytmp = x[i][1]; ztmp = x[i][2]; @@ -513,7 +507,6 @@ void PairCoulStreitz::compute(int eflag, int vflag) iparam_i = elem2param[itype]; qi = q[i]; zei = params[iparam_i].zeta; - zi = params[iparam_i].zcore; // self ionizition energy, only on i atom @@ -529,7 +522,6 @@ void PairCoulStreitz::compute(int eflag, int vflag) for (jj = 0; jj < jnum; jj++) { j = jlist[jj]; j &= NEIGHMASK; - jtag = tag[j]; jtype = map[type[j]]; iparam_j = elem2param[jtype]; qj = q[j]; @@ -587,6 +579,8 @@ double PairCoulStreitz::self(Param *param, double qi) if (kspacetype == 1) return 1.0*qi*(s1+qi*(0.50*s2 - qqrd2e*woself)); if (kspacetype == 2) return 1.0*qi*(s1+qi*(0.50*s2)); + + return 0.0; } /* ---------------------------------------------------------------------- */ @@ -673,8 +667,8 @@ void PairCoulStreitz::wolf_sum(double qi, double qj, double zj, double r, double derfcr = exp(-a*a*r*r); double erfcrc = erfc(a*rc); - double etmp1, etmp2, etmp3, etmp4; - double ftmp1, ftmp2, ftmp3, ftmp4; + double etmp1, etmp2, etmp3; + double ftmp1, ftmp2, ftmp3; etmp = etmp1 = etmp2 = etmp3 = 0.0; ftmp = ftmp1 = ftmp2 = ftmp3 = 0.0; diff --git a/src/velocity.cpp b/src/velocity.cpp index 0ea6b01673..ba1bc06bd0 100644 --- a/src/velocity.cpp +++ b/src/velocity.cpp @@ -98,6 +98,28 @@ void Velocity::command(int narg, char **arg) else if (style == RAMP) options(narg-8,&arg[8]); else if (style == ZERO) options(narg-3,&arg[3]); + // special cases where full init and border communication must be done first + // for CREATE/SET if compute temp/cs is used + // for ZERO if fix rigid/small is used + // b/c methods invoked in the compute/fix perform forward/reverse comm + + int initcomm = 0; + if (style == ZERO && rfix >= 0 && + strcmp(modify->fix[rfix]->style,"rigid/small") == 0) initcomm = 1; + if ((style == CREATE || style == SET) && temperature && + strcmp(temperature->style,"temp/cs2") == 0) initcomm = 1; + + if (initcomm) { + lmp->init(); + if (domain->triclinic) domain->x2lamda(atom->nlocal); + domain->pbc(); + domain->reset_box(); + comm->setup(); + comm->exchange(); + comm->borders(); + if (domain->triclinic) domain->lamda2x(atom->nlocal+atom->nghost); + } + // initialize velocities based on style // create() invoked differently, so can be called externally @@ -670,37 +692,27 @@ void Velocity::ramp(int narg, char **arg) /* ---------------------------------------------------------------------- zero linear or angular momentum of a group - if using rigid/small requires init of entire system since - its methods perform forward/reverse comm, - comm::init needs neighbor::init needs pair::init needs kspace::init, etc - also requires setup_pre_neighbor call to setup bodies ------------------------------------------------------------------------- */ void Velocity::zero(int narg, char **arg) { if (strcmp(arg[0],"linear") == 0) { if (rfix < 0) zero_momentum(); - else { - if (strcmp(modify->fix[rfix]->style,"rigid/small") == 0) { - lmp->init(); - modify->fix[rfix]->setup_pre_neighbor(); - modify->fix[rfix]->zero_momentum(); - } else if (strstr(modify->fix[rfix]->style,"rigid")) { - modify->fix[rfix]->zero_momentum(); - } else error->all(FLERR,"Velocity rigid used with non-rigid fix-ID"); - } + else if (strcmp(modify->fix[rfix]->style,"rigid/small") == 0) { + modify->fix[rfix]->setup_pre_neighbor(); + modify->fix[rfix]->zero_momentum(); + } else if (strstr(modify->fix[rfix]->style,"rigid")) { + modify->fix[rfix]->zero_momentum(); + } else error->all(FLERR,"Velocity rigid used with non-rigid fix-ID"); } else if (strcmp(arg[0],"angular") == 0) { if (rfix < 0) zero_rotation(); - else { - if (strcmp(modify->fix[rfix]->style,"rigid/small") == 0) { - lmp->init(); - modify->fix[rfix]->setup_pre_neighbor(); - modify->fix[rfix]->zero_rotation(); - } else if (strstr(modify->fix[rfix]->style,"rigid")) { - modify->fix[rfix]->zero_rotation(); - } else error->all(FLERR,"Velocity rigid used with non-rigid fix-ID"); - } + else if (strcmp(modify->fix[rfix]->style,"rigid/small") == 0) { + modify->fix[rfix]->setup_pre_neighbor(); + modify->fix[rfix]->zero_rotation(); + } else if (strstr(modify->fix[rfix]->style,"rigid")) { + modify->fix[rfix]->zero_rotation(); + } else error->all(FLERR,"Velocity rigid used with non-rigid fix-ID"); } else error->all(FLERR,"Illegal velocity command"); }