Merge pull request #2203 from akohlmey/collected-small-changes

Various small fixes and cleanups. Bugfixes for atom style hybrid

src/BODY/pair_body_rounded_polygon.cpp

@@ -236,7 +236,7 @@ void PairBodyRoundedPolygon::compute(int eflag, int vflag)
         edge[jefirst+nj][4] = 0;
       }
 
-      int interact, num_contacts, done;
+      int num_contacts, done;
       double delta_a, j_a;
       Contact contact_list[MAX_CONTACTS];
 
@@ -244,15 +244,13 @@ void PairBodyRoundedPolygon::compute(int eflag, int vflag)
 
       // check interaction between i's vertices and j' edges
 
-      interact = vertex_against_edge(i, j, k_nij, k_naij,
+      vertex_against_edge(i, j, k_nij, k_naij, x, f, torque, tag,
-                                     x, f, torque, tag, contact_list,
+                          contact_list, num_contacts, evdwl, facc);
-                                     num_contacts, evdwl, facc);
 
       // check interaction between j's vertices and i' edges
 
-      interact = vertex_against_edge(j, i, k_nij, k_naij,
+      vertex_against_edge(j, i, k_nij, k_naij, x, f, torque, tag,
-                                     x, f, torque, tag, contact_list,
+                          contact_list, num_contacts, evdwl, facc);
-                                     num_contacts, evdwl, facc);
 
       if (num_contacts >= 2) {
 
@@ -595,7 +593,7 @@ void PairBodyRoundedPolygon::sphere_against_sphere(int i, int j,
 {
   double rradi,rradj;
   double vr1,vr2,vr3,vnnr,vn1,vn2,vn3,vt1,vt2,vt3;
-  double rij,rsqinv,R,fx,fy,fz,fpair,shift,energy;
+  double rij,rsqinv,R,fx,fy,fz,fn[3],ft[3],fpair,shift,energy;
   int nlocal = atom->nlocal;
   int newton_pair = force->newton_pair;
 
@@ -641,6 +639,23 @@ void PairBodyRoundedPolygon::sphere_against_sphere(int i, int j,
     vt1 = vr1 - vn1;
     vt2 = vr2 - vn2;
     vt3 = vr3 - vn3;
+
+    // normal friction term at contact
+
+    fn[0] = -c_n * vn1;
+    fn[1] = -c_n * vn2;
+    fn[2] = -c_n * vn3;
+
+    // tangential friction term at contact,
+    // excluding the tangential deformation term for now
+
+    ft[0] = -c_t * vt1;
+    ft[1] = -c_t * vt2;
+    ft[2] = -c_t * vt3;
+
+    fx += fn[0] + ft[0];
+    fy += fn[1] + ft[1];
+    fz += fn[2] + ft[2];
   }
 
   f[i][0] += fx;
@@ -1349,4 +1364,3 @@ void PairBodyRoundedPolygon::distance(const double* x2, const double* x1,
     + (x2[1] - x1[1]) * (x2[1] - x1[1])
     + (x2[2] - x1[2]) * (x2[2] - x1[2]));
 }
-

src/KOKKOS/atom_vec_kokkos.h

@@ -103,14 +103,14 @@ class AtomVecKokkos : public AtomVec {
                          ExecutionSpace space) = 0;
 
   virtual int
-    pack_border_vel_kokkos(int n, DAT::tdual_int_2d k_sendlist,
+    pack_border_vel_kokkos(int /*n*/, DAT::tdual_int_2d /*k_sendlist*/,
-                           DAT::tdual_xfloat_2d buf,int iswap,
+                           DAT::tdual_xfloat_2d /*buf*/,int /*iswap*/,
-                           int pbc_flag, int *pbc, ExecutionSpace space) { return 0; }
+                           int /*pbc_flag*/, int * /*pbc*/, ExecutionSpace /*space*/) { return 0; }
 
   virtual void
-    unpack_border_vel_kokkos(const int &n, const int &nfirst,
+    unpack_border_vel_kokkos(const int &/*n*/, const int & /*nfirst*/,
-                             const DAT::tdual_xfloat_2d &buf,
+                             const DAT::tdual_xfloat_2d & /*buf*/,
-                             ExecutionSpace space) {}
+                             ExecutionSpace /*space*/) {}
 
   virtual int
     pack_exchange_kokkos(const int &nsend, DAT::tdual_xfloat_2d &buf,

src/KOKKOS/kokkos_type.h

@@ -272,7 +272,7 @@ template<typename T1, typename T2>
 class ScatterViewHelper<Kokkos::Experimental::ScatterDuplicated,T1,T2> {
 public:
   KOKKOS_INLINE_FUNCTION
-  static T1 get(const T1 &dup, const T2 &nondup) {
+  static T1 get(const T1 &dup, const T2 & /*nondup*/) {
     return dup;
   }
 };
@@ -1025,7 +1025,7 @@ struct params_lj_coul {
   KOKKOS_INLINE_FUNCTION
   params_lj_coul(){cut_ljsq=0;cut_coulsq=0;lj1=0;lj2=0;lj3=0;lj4=0;offset=0;};
   KOKKOS_INLINE_FUNCTION
-  params_lj_coul(int i){cut_ljsq=0;cut_coulsq=0;lj1=0;lj2=0;lj3=0;lj4=0;offset=0;};
+  params_lj_coul(int /*i*/){cut_ljsq=0;cut_coulsq=0;lj1=0;lj2=0;lj3=0;lj4=0;offset=0;};
   F_FLOAT cut_ljsq,cut_coulsq,lj1,lj2,lj3,lj4,offset;
 };
 

src/KOKKOS/sna_kokkos_impl.h

@@ -368,7 +368,6 @@ void SNAKokkos<DeviceType>::compute_ui(const typename Kokkos::TeamPolicy<DeviceT
 
   for (int j = 1; j <= twojmax; j++) {
     const int jju = idxu_block[j];
-    const int jjup = idxu_block[j-1];
 
     // fill in left side of matrix layer from previous layer
 
@@ -1080,7 +1079,6 @@ void SNAKokkos<DeviceType>::compute_bi_cpu(const typename Kokkos::TeamPolicy<Dev
         Kokkos::parallel_for(Kokkos::TeamThreadRange(team,idxb_max),
           [&] (const int& jjb) {
         //for(int jjb = 0; jjb < idxb_max; jjb++) {
-          const auto jjballoy = itriple;
           const int j1 = idxb(jjb, 0);
           const int j2 = idxb(jjb, 1);
           const int j = idxb(jjb, 2);
@@ -2024,7 +2022,6 @@ double SNAKokkos<DeviceType>::memory_usage()
 {
   int jdimpq = twojmax + 2;
   int jdim = twojmax + 1;
-  int natom_pad = ((natom + 32 - 1) / 32) * 32; // for AoSoA layouts
   double bytes;
 
   bytes = 0;
@@ -2032,10 +2029,9 @@ double SNAKokkos<DeviceType>::memory_usage()
   bytes += jdimpq*jdimpq * sizeof(double);               // pqarray
   bytes += idxcg_max * sizeof(double);                   // cglist
 
-
-
 #ifdef KOKKOS_ENABLE_CUDA
   if (std::is_same<DeviceType,Kokkos::Cuda>::value) {
+    int natom_pad = ((natom + 32 - 1) / 32) * 32; // for AoSoA layouts
 
     bytes += natom * idxu_max * nelements * sizeof(double);          // ulisttot_re
     bytes += natom * idxu_max * nelements * sizeof(double);          // ulisttot_im

src/KSPACE/pppm_disp.cpp

@@ -316,7 +316,6 @@ void PPPMDisp::init()
   for (int i=0; i<=EWALD_MAXORDER; ++i)                 // transcribe order
     if (ewald_order&(1<<i)) {                           // from pair_style
       int  k=0;
-      char str[128];
       switch (i) {
         case 1:
           k = 0; break;

src/MLIAP/mliap_descriptor_snap.cpp

@@ -149,15 +149,13 @@ void MLIAPDescriptorSNAP::forward(int* map, NeighList* list, double **descriptor
 void MLIAPDescriptorSNAP::backward(PairMLIAP* pairmliap, NeighList* list, double **beta, int vflag)
 {
   int i,j,jnum,ninside;
-  double delx,dely,delz,evdwl,rsq;
+  double delx,dely,delz,rsq;
   double fij[3];
   int *jlist,*numneigh,**firstneigh;
 
   double **x = atom->x;
   double **f = atom->f;
   int *type = atom->type;
-  int nlocal = atom->nlocal;
-  int newton_pair = force->newton_pair;
 
   numneigh = list->numneigh;
   firstneigh = list->firstneigh;

src/MLIAP/mliap_model.cpp

@@ -142,7 +142,6 @@ void MLIAPModel::read_coeffs(char *coefffilename)
   }
 
   if (comm->me == 0) fclose(fpcoeff);
-
 }
 
 /* ----------------------------------------------------------------------
@@ -153,9 +152,7 @@ double MLIAPModel::memory_usage()
 {
   double bytes = 0;
 
-  int n = atom->ntypes+1;
   bytes += nelements*nparams*sizeof(double);  // coeffelem
-
   return bytes;
 }
 

src/RIGID/fix_rigid_small.cpp

@@ -410,9 +410,9 @@ FixRigidSmall::FixRigidSmall(LAMMPS *lmp, int narg, char **arg) :
 
   if (comm->me == 0) {
     if (screen)
-      fprintf(screen,"  create bodies CPU = %g secs\n",time2-time1);
+      fprintf(screen,"  create bodies CPU = %g seconds\n",time2-time1);
     if (logfile)
-      fprintf(logfile,"  create bodies CPU = %g secs\n",time2-time1);
+      fprintf(logfile,"  create bodies CPU = %g seconds\n",time2-time1);
   }
 
   // set nlocal_body and allocate bodies I own

src/RIGID/fix_shake.cpp

@@ -229,9 +229,9 @@ FixShake::FixShake(LAMMPS *lmp, int narg, char **arg) :
 
   if (comm->me == 0) {
     if (screen)
-      fprintf(screen,"  find clusters CPU = %g secs\n",time2-time1);
+      fprintf(screen,"  find clusters CPU = %g seconds\n",time2-time1);
     if (logfile)
-      fprintf(logfile,"  find clusters CPU = %g secs\n",time2-time1);
+      fprintf(logfile,"  find clusters CPU = %g seconds\n",time2-time1);
   }
 
   // initialize list of SHAKE clusters to constrain

src/SNAP/compute_sna_atom.cpp

@@ -259,7 +259,7 @@ void ComputeSNAAtom::compute_peratom()
         int jtype = type[j];
         int jelem = 0;
         if (chemflag)
-          int jelem = map[jtype];
+          jelem = map[jtype];
         if (rsq < cutsq[itype][jtype] && rsq>1e-20) {
           snaptr->rij[ninside][0] = delx;
           snaptr->rij[ninside][1] = dely;

src/USER-CGDNA/bond_oxdna_fene.cpp

@@ -328,29 +328,9 @@ void BondOxdnaFene::coeff(int narg, char **arg)
 
 void BondOxdnaFene::init_style()
 {
-  /* special bonds have to be lj = 0 1 1 and coul = 1 1 1 to exclude
-     the ss excluded volume interaction between nearest neighbors   */
-
-  force->special_lj[1] = 0.0;
-  force->special_lj[2] = 1.0;
-  force->special_lj[3] = 1.0;
-  force->special_coul[1] = 1.0;
-  force->special_coul[2] = 1.0;
-  force->special_coul[3] = 1.0;
-
-  fprintf(screen,"Finding 1-2 1-3 1-4 neighbors ...\n"
-                 " Special bond factors lj:   %-10g %-10g %-10g\n"
-                 " Special bond factors coul: %-10g %-10g %-10g\n",
-                 force->special_lj[1],force->special_lj[2],force->special_lj[3],
-                 force->special_coul[1],force->special_coul[2],force->special_coul[3]);
-
   if (force->special_lj[1] != 0.0 || force->special_lj[2] != 1.0 || force->special_lj[3] != 1.0 ||
       force->special_coul[1] != 1.0 || force->special_coul[2] != 1.0 || force->special_coul[3] != 1.0)
-  {
+    error->all(FLERR,"Must use 'special_bonds lj 0 1 1 coul 1 1 1' with bond style oxdna/fene");
-    if (comm->me == 0)
-      error->warning(FLERR,"Use special bonds lj = 0,1,1 and coul = 1,1,1 with bond style oxdna/fene");
-  }
-
 }
 
 /* ---------------------------------------------------------------------- */

src/atom_vec_hybrid.cpp

@@ -57,8 +57,6 @@ AtomVecHybrid::~AtomVecHybrid()
   delete [] styles;
   for (int k = 0; k < nstyles; k++) delete [] keywords[k];
   delete [] keywords;
-
-  for (int k = 0; k < nstyles_bonus; k++) delete styles_bonus[k];
   delete [] styles_bonus;
 
   if (!fields_allocated) return;
@@ -225,6 +223,7 @@ void AtomVecHybrid::process_args(int narg, char **arg)
 
   // set bonus_flag if any substyle has bonus data
   // set nstyles_bonus & styles_bonus
+  // sum two sizes over contributions from each substyle with bonus data.
 
   nstyles_bonus = 0;
   for (int k = 0; k < nstyles; k++)
@@ -234,9 +233,14 @@ void AtomVecHybrid::process_args(int narg, char **arg)
     bonus_flag = 1;
     styles_bonus = new AtomVec*[nstyles_bonus];
     nstyles_bonus = 0;
+    size_forward_bonus = 0;
+    size_border_bonus = 0;
     for (int k = 0; k < nstyles; k++) {
-      if (styles[k]->bonus_flag)
+      if (styles[k]->bonus_flag) {
         styles_bonus[nstyles_bonus++] = styles[k];
+        size_forward_bonus += styles[k]->size_forward_bonus;
+        size_border_bonus += styles[k]->size_border_bonus;
+      }
     }
   }
 

src/domain.cpp

@@ -1939,12 +1939,13 @@ void Domain::print_box(const std::string &prefix)
   if (comm->me == 0) {
     std::string mesg = prefix;
     if (triclinic == 0) {
-      mesg += fmt::format("orthogonal box = ({} {} {}) to ({} {} {})\n",
+      mesg += fmt::format("orthogonal box = ({:.8g} {:.8g} {:.8g}) to "
-                          boxlo[0],boxlo[1],boxlo[2],
+                          "({:.8g} {:.8g} {:.8g})\n",boxlo[0],boxlo[1],
-                          boxhi[0],boxhi[1],boxhi[2]);
+                          boxlo[2],boxhi[0],boxhi[1],boxhi[2]);
     } else {
-      mesg += fmt::format("triclinic box = ({} {} {}) to ({} {} {}) "
+      mesg += fmt::format("triclinic box = ({:.8g} {:.8g} {:.8g}) to "
-                          "with tilt ({} {} {})\n",boxlo[0],boxlo[1],
+                          "({:.8g} {:.8g} {:.8g}) with tilt "
+                          "({:.8g} {:.8g} {:.8g})\n",boxlo[0],boxlo[1],
                           boxlo[2],boxhi[0],boxhi[1],boxhi[2],xy,xz,yz);
     }
     utils::logmesg(lmp,mesg);

src/fix_adapt.cpp

@@ -588,7 +588,7 @@ void FixAdapt::change_settings()
       // for scaleflag, previous_diam_scale is the scale factor on previous step
 
       if (ad->aparam == DIAMETER) {
-        double density,scale;
+        double scale;
         double *radius = atom->radius;
         double *rmass = atom->rmass;
         int *mask = atom->mask;
@@ -691,7 +691,7 @@ void FixAdapt::restore_settings()
 
     } else if (ad->which == ATOM) {
       if (diamflag) {
-        double density,scale;
+        double scale;
 
         double *vec = fix_diam->vstore;
         double *radius = atom->radius;

src/info.cpp

@@ -793,7 +793,6 @@ bool Info::is_active(const char *category, const char *name)
 {
   if ((category == NULL) || (name == NULL)) return false;
   const char *style = "none";
-  const int len = strlen(name);
 
   if (strcmp(category,"package") == 0) {
     if (strcmp(name,"gpu") == 0) {

src/lattice.cpp

@@ -304,7 +304,8 @@ Lattice::Lattice(LAMMPS *lmp, int narg, char **arg) : Pointers(lmp)
   // print lattice spacings
 
   if (comm->me == 0)
-    utils::logmesg(lmp,fmt::format("Lattice spacing in x,y,z = {} {} {}\n",
+    utils::logmesg(lmp,fmt::format("Lattice spacing in x,y,z = "
+                                   "{:.8g} {:.8g} {:.8g}\n",
                                    xlattice,ylattice,zlattice));
 }
 

src/neighbor.cpp

@@ -1438,7 +1438,7 @@ void Neighbor::init_topology()
 
 void Neighbor::print_pairwise_info()
 {
-  int i,m;
+  int i;
   NeighRequest *rq;
 
   const double cutghost = MAX(cutneighmax,comm->cutghostuser);
@@ -1466,10 +1466,10 @@ void Neighbor::print_pairwise_info()
                      every,delay,dist_check ? "yes" : "no");
   out += fmt::format("  max neighbors/atom: {}, page size: {}\n",
                      oneatom, pgsize);
-  out += fmt::format("  master list distance cutoff = {}\n",cutneighmax);
+  out += fmt::format("  master list distance cutoff = {:.8g}\n",cutneighmax);
-  out += fmt::format("  ghost atom cutoff = {}\n",cutghost);
+  out += fmt::format("  ghost atom cutoff = {:.8g}\n",cutghost);
   if (style != Neighbor::NSQ)
-    out += fmt::format("  binsize = {}, bins = {} {} {}\n",binsize,
+    out += fmt::format("  binsize = {:.8g}, bins = {:g} {:g} {:g}\n",binsize,
                        ceil(bbox[0]/binsize), ceil(bbox[1]/binsize),
                        ceil(bbox[2]/binsize));
 

src/pair_table.cpp

@@ -70,7 +70,6 @@ void PairTable::compute(int eflag, int vflag)
   int i,j,ii,jj,inum,jnum,itype,jtype,itable;
   double xtmp,ytmp,ztmp,delx,dely,delz,evdwl,fpair;
   double rsq,factor_lj,fraction,value,a,b;
-  char estr[128];
   int *ilist,*jlist,*numneigh,**firstneigh;
   Table *tb;
 

src/read_data.cpp

@@ -897,7 +897,7 @@ void ReadData::command(int narg, char **arg)
   MPI_Barrier(world);
 
   if (comm->me == 0)
-    utils::logmesg(lmp,fmt::format("  read_data CPU = {:.3f} secs\n",
+    utils::logmesg(lmp,fmt::format("  read_data CPU = {:.3f} seconds\n",
                                    MPI_Wtime()-time1));
 }
 

src/read_restart.cpp

@@ -511,7 +511,7 @@ void ReadRestart::command(int narg, char **arg)
   MPI_Barrier(world);
 
   if (comm->me == 0)
-    utils::logmesg(lmp,fmt::format("  read_restart CPU = {:.3f} secs\n",
+    utils::logmesg(lmp,fmt::format("  read_restart CPU = {:.3f} seconds\n",
                                    MPI_Wtime()-time1));
 }
 

src/replicate.cpp

@@ -779,6 +779,6 @@ void Replicate::command(int narg, char **arg)
   MPI_Barrier(world);
 
   if (me == 0)
-    utils::logmesg(lmp,fmt::format("  replicate CPU = {:.3f} secs\n",
+    utils::logmesg(lmp,fmt::format("  replicate CPU = {:.3f} seconds\n",
                                    MPI_Wtime()-time1));
 }

src/special.cpp

@@ -63,8 +63,8 @@ void Special::build()
     const double * const special_lj   = force->special_lj;
     const double * const special_coul = force->special_coul;
     fmt::print(screen,"Finding 1-2 1-3 1-4 neighbors ...\n"
-               "  special bond factors lj:   {:<10g} {:<10g} {:<10g}\n"
+               "  special bond factors lj:   {:<8g} {:<8g} {:<8g}\n"
-               "  special bond factors coul: {:<10g} {:<10g} {:<10g}\n",
+               "  special bond factors coul: {:<8g} {:<8g} {:<8g}\n",
                special_lj[1],special_lj[2],special_lj[3],
                special_coul[1],special_coul[2],special_coul[3]);
   }
@@ -1314,6 +1314,6 @@ void Special::fix_alteration()
 void Special::timer_output(double time1)
 {
   if (comm->me == 0)
-    utils::logmesg(lmp,fmt::format("  special bonds CPU = {:.3f} secs\n",
+    utils::logmesg(lmp,fmt::format("  special bonds CPU = {:.3f} seconds\n",
                                    MPI_Wtime()-time1));
 }

src/utils.cpp

@@ -549,9 +549,9 @@ bool utils::is_double(const std::string & str) {
 
 std::string utils::path_basename(const std::string & path) {
 #if defined(_WIN32)
-  size_t start = path.find_last_of('/\\');
+  size_t start = path.find_last_of("/\\");
 #else
-  size_t start = path.find_last_of('/');
+  size_t start = path.find_last_of("/");
 #endif
 
   if (start == std::string::npos) {

unittest/force-styles/pair_style.cpp

@@ -1007,6 +1007,9 @@ TEST(PairStyle, single)
     char **argv = (char **)args;
     int argc    = sizeof(args) / sizeof(char *);
 
+    // need to add this dependency
+    test_config.prerequisites.push_back(std::make_pair("atom", "full"));
+
     // create a LAMMPS instance with standard settings to detect the number of atom types
     if (!verbose) ::testing::internal::CaptureStdout();
     LAMMPS *lmp = init_lammps(argc, argv, test_config);
@@ -1018,8 +1021,10 @@ TEST(PairStyle, single)
         for (auto prerequisite : test_config.prerequisites) {
             std::cerr << prerequisite.first << "_style " << prerequisite.second << "\n";
         }
+        test_config.prerequisites.pop_back();
         GTEST_SKIP();
     }
+    test_config.prerequisites.pop_back();
 
     // gather some information and skip if unsupported
     int ntypes    = lmp->atom->ntypes;

unittest/force-styles/tests/atomic-pair-eam.yaml

@@ -1,7 +1,7 @@
 ---
 lammps_version: 15 Jun 2020
 date_generated: Thu Jun 25 12:25:55 202
-epsilon: 5e-12
+epsilon: 6e-12
 prerequisites: ! |
   pair eam
 pre_commands: ! |