From 67574601ed8bfadb5e4a4139ae52b89399e080b7 Mon Sep 17 00:00:00 2001 From: Trung Nguyen Date: Sun, 15 Jan 2023 15:41:54 -0600 Subject: [PATCH] Cleaned up commented-out and debugging stuffs, removed irrelevant changes to lj/cut/dipole/cut, reverted unwanted changes in the PPPMGPU destructor, fixed unresolved conflicts in tinker.py, updated the userbinsize==0 case in atom.cpp and using Force::pair_match() as suggested. Internal timing stuffs need work. --- cmake/CMakeLists.txt | 1 - examples/amoeba/in.ubiquitin | 2 -- src/AMOEBA/amoeba_induce.cpp | 40 ++---------------------- src/AMOEBA/amoeba_kspace.cpp | 8 ----- src/AMOEBA/fix_amoeba_bitorsion.cpp | 6 ++-- src/AMOEBA/improper_amoeba.cpp | 7 ++--- src/AMOEBA/pair_amoeba.cpp | 7 ----- src/DIPOLE/pair_lj_cut_dipole_cut.cpp | 19 ++---------- src/Depend.sh | 4 +++ src/GPU/pppm_gpu.cpp | 2 ++ src/atom.cpp | 11 +++---- tools/tinker/tinker2lmp.py | 44 --------------------------- 12 files changed, 19 insertions(+), 132 deletions(-) diff --git a/cmake/CMakeLists.txt b/cmake/CMakeLists.txt index d7137c3672..0223750ace 100644 --- a/cmake/CMakeLists.txt +++ b/cmake/CMakeLists.txt @@ -395,7 +395,6 @@ endif() pkg_depends(ML-IAP ML-SNAP) pkg_depends(MPIIO MPI) pkg_depends(ATC MANYBODY) -pkg_depends(AMOEBA KSPACE) pkg_depends(LATBOLTZ MPI) pkg_depends(SCAFACOS MPI) pkg_depends(AMOEBA KSPACE) diff --git a/examples/amoeba/in.ubiquitin b/examples/amoeba/in.ubiquitin index cacb7d3571..cb789a19f8 100644 --- a/examples/amoeba/in.ubiquitin +++ b/examples/amoeba/in.ubiquitin @@ -34,8 +34,6 @@ pair_coeff * * amoeba_ubiquitin.prm amoeba_ubiquitin.key special_bonds lj/coul 0.5 0.5 0.5 one/five yes -# setup force components this way so can dump them (AMOEBA or HIPPO also needs them for now) - fix fhal all store/state 0 fx fy fz fix frepulse all store/state 0 fx fy fz fix fdisp all store/state 0 fx fy fz diff --git a/src/AMOEBA/amoeba_induce.cpp b/src/AMOEBA/amoeba_induce.cpp index 17c4df326d..031173060c 100644 --- a/src/AMOEBA/amoeba_induce.cpp +++ b/src/AMOEBA/amoeba_induce.cpp @@ -86,17 +86,6 @@ void PairAmoeba::induce() crstyle = FIELD; comm->reverse_comm(this); - // DEBUG statements - - /* - for (i = 0; i < nlocal; i++) - if (atom->tag[i] == 1) - printf("AAA FIELD atom %d: field %g %g %g: fieldp %g %g %g\n", - atom->tag[i], - field[i][0],field[i][1],field[i][2], - fieldp[i][0],fieldp[i][1],fieldp[i][2]); - */ - // set induced dipoles to polarizability times direct field for (i = 0; i < nlocal; i++) { @@ -213,16 +202,7 @@ void PairAmoeba::induce() cfstyle = INDUCE; comm->forward_comm(this); -/* - if (comm->me == 0) { - printf("CPU: cutghost = %f\n", comm->cutghost[0]); - for (i = 0; i < 20; i++) { - printf("i = %d: uind = %f %f %f; udirp = %f %f %f\n", - i, uind[i][0], uind[i][1], uind[i][2], - uinp[i][0], uinp[i][1], uinp[i][2]); - } - } -*/ + ufield0c(field,fieldp); crstyle = FIELD; @@ -284,18 +264,6 @@ void PairAmoeba::induce() crstyle = FIELD; comm->reverse_comm(this); - - //error->all(FLERR,"STOP"); -/* - if (comm->me == 0) { - printf("CPU: iter = %d\n", iter); - for (i = 0; i < 10; i++) { - printf("i = %d: field = %f %f %f; fieldp = %f %f %f\n", - i, field[i][0], field[i][1], field[i][2], - fieldp[i][0], fieldp[i][1], fieldp[i][2]); - } - } -*/ for (i = 0; i < nlocal; i++) { for (j = 0; j < 3; j++) { @@ -413,8 +381,6 @@ void PairAmoeba::induce() } } - // if (comm->me == 0) printf("CG iteration count = %d\n",iter); - // terminate the calculation if dipoles failed to converge // NOTE: could make this an error @@ -1033,9 +999,7 @@ void PairAmoeba::umutual2b(double **field, double **fieldp) j = jlist[jj]; uindj = uind[j]; uinpj = uinp[j]; - //if (i==0 && j == 10) - // printf("i = %d: j = %d: tdipdip %f %f %f %f %f %f\n", - // i, j,tdipdip[0],tdipdip[1],tdipdip[2],tdipdip[3],tdipdip[4],tdipdip[5]); + fid[0] = tdipdip[0]*uindj[0] + tdipdip[1]*uindj[1] + tdipdip[2]*uindj[2]; fid[1] = tdipdip[1]*uindj[0] + tdipdip[3]*uindj[1] + tdipdip[4]*uindj[2]; fid[2] = tdipdip[2]*uindj[0] + tdipdip[4]*uindj[1] + tdipdip[5]*uindj[2]; diff --git a/src/AMOEBA/amoeba_kspace.cpp b/src/AMOEBA/amoeba_kspace.cpp index 76d13da780..9213b96042 100644 --- a/src/AMOEBA/amoeba_kspace.cpp +++ b/src/AMOEBA/amoeba_kspace.cpp @@ -68,8 +68,6 @@ void PairAmoeba::moduli() int maxfft = MAX(nfft1,nfft2); maxfft = MAX(maxfft,nfft3); - //double *array = new double[bsorder]; - //double *bsarray = new double[maxfft]; if (maxfft > _nfft_max) { memory->destroy(_moduli_bsarray); _nfft_max = maxfft; @@ -79,7 +77,6 @@ void PairAmoeba::moduli() // compute and load the moduli values double x = 0.0; - //bspline(x,bsorder,array); bspline(x,bsorder,_moduli_array); for (i = 0; i < maxfft; i++) _moduli_bsarray[i] = 0.0; @@ -88,11 +85,6 @@ void PairAmoeba::moduli() dftmod(bsmod1,_moduli_bsarray,nfft1,bsorder); dftmod(bsmod2,_moduli_bsarray,nfft2,bsorder); dftmod(bsmod3,_moduli_bsarray,nfft3,bsorder); - - // perform deallocation of local arrays - - //delete[] array; - //delete[] bsarray; } /* ---------------------------------------------------------------------- diff --git a/src/AMOEBA/fix_amoeba_bitorsion.cpp b/src/AMOEBA/fix_amoeba_bitorsion.cpp index 6c3c31eec8..cb8c62819d 100644 --- a/src/AMOEBA/fix_amoeba_bitorsion.cpp +++ b/src/AMOEBA/fix_amoeba_bitorsion.cpp @@ -194,10 +194,8 @@ void FixAmoebaBiTorsion::init() // error check that PairAmoeba or PairHiippo exist pair = nullptr; - pair = force->pair_match("amoeba",1,0); - if (!pair) pair = force->pair_match("amoeba/gpu",1,0); - if (!pair) pair = force->pair_match("hippo",1,0); - if (!pair) pair = force->pair_match("hippo/gpu",1,0); + pair = force->pair_match("^amoeba",0,0); + if (!pair) pair = force->pair_match("^hippo",0,0); if (!pair) error->all(FLERR,"Cannot use fix amoeba/bitorsion w/out pair amoeba/hippo"); diff --git a/src/AMOEBA/improper_amoeba.cpp b/src/AMOEBA/improper_amoeba.cpp index 136857e74b..cb9db01b59 100644 --- a/src/AMOEBA/improper_amoeba.cpp +++ b/src/AMOEBA/improper_amoeba.cpp @@ -285,10 +285,9 @@ void ImproperAmoeba::init_style() // check if PairAmoeba disabled improper terms Pair *pair = nullptr; - pair = force->pair_match("amoeba",1,0); - if (!pair) pair = force->pair_match("amoeba/gpu",1,0); - if (!pair) pair = force->pair_match("hippo",1,0); - if (!pair) pair = force->pair_match("hippo/gpu",1,0); + pair = force->pair_match("^amoeba",0,0); + if (!pair) pair = force->pair_match("^hippo",0,0); + if (!pair) error->all(FLERR,"Improper amoeba could not find pair amoeba/hippo"); int tmp; diff --git a/src/AMOEBA/pair_amoeba.cpp b/src/AMOEBA/pair_amoeba.cpp index bb06ecb4a4..2a1a10075c 100644 --- a/src/AMOEBA/pair_amoeba.cpp +++ b/src/AMOEBA/pair_amoeba.cpp @@ -1055,13 +1055,6 @@ void PairAmoeba::init_style() // request standard neighbor list - -// int irequest = neighbor->request(this,instance_me); - - // for DEBUGGING with GPU - //neighbor->requests[irequest]->half = 0; - //neighbor->requests[irequest]->full = 1; - neighbor->add_request(this); } diff --git a/src/DIPOLE/pair_lj_cut_dipole_cut.cpp b/src/DIPOLE/pair_lj_cut_dipole_cut.cpp index 2047eb8b9c..a7e5674a88 100644 --- a/src/DIPOLE/pair_lj_cut_dipole_cut.cpp +++ b/src/DIPOLE/pair_lj_cut_dipole_cut.cpp @@ -90,8 +90,6 @@ void PairLJCutDipoleCut::compute(int eflag, int vflag) numneigh = list->numneigh; firstneigh = list->firstneigh; - int maxsize = 10; - // loop over neighbors of my atoms for (ii = 0; ii < inum; ii++) { @@ -104,13 +102,6 @@ void PairLJCutDipoleCut::compute(int eflag, int vflag) jlist = firstneigh[i]; jnum = numneigh[i]; - double scale_dipole = 1.0; - if (jnum > maxsize) { - scale_dipole = maxsize; //1.0/(double)maxsize; - } else { - scale_dipole = jnum; //1.0/(double)jnum; - } - for (jj = 0; jj < jnum; jj++) { j = jlist[jj]; factor_lj = special_lj[sbmask(j)]; @@ -216,7 +207,7 @@ void PairLJCutDipoleCut::compute(int eflag, int vflag) // total force - fq = scale_dipole*factor_coul*qqrd2e; + fq = factor_coul*qqrd2e; fx = fq*forcecoulx + delx*forcelj; fy = fq*forcecouly + dely*forcelj; fz = fq*forcecoulz + delz*forcelj; @@ -230,7 +221,7 @@ void PairLJCutDipoleCut::compute(int eflag, int vflag) torque[i][1] += fq*tiycoul; torque[i][2] += fq*tizcoul; - if (newton_pair) { + if (newton_pair || j < nlocal) { f[j][0] -= fx; f[j][1] -= fy; f[j][2] -= fz; @@ -371,13 +362,7 @@ void PairLJCutDipoleCut::init_style() if (!atom->q_flag || !atom->mu_flag || !atom->torque_flag) error->all(FLERR,"Pair dipole/cut requires atom attributes q, mu, torque"); -<<<<<<< HEAD - int irequest = neighbor->request(this,instance_me); - neighbor->requests[irequest]->half = 0; - neighbor->requests[irequest]->full = 1; -======= neighbor->add_request(this); ->>>>>>> amoeba } /* ---------------------------------------------------------------------- diff --git a/src/Depend.sh b/src/Depend.sh index 6cf613cde7..28ac78d9af 100755 --- a/src/Depend.sh +++ b/src/Depend.sh @@ -45,6 +45,10 @@ depend () { # add one if statement per parent package # add one depend() call per child package that depends on that parent +if (test $1 = "AMOEBA") then + depend GPU +fi + if (test $1 = "ASPHERE") then depend GPU depend OPENMP diff --git a/src/GPU/pppm_gpu.cpp b/src/GPU/pppm_gpu.cpp index 4019eb467d..a2a2b0eed8 100644 --- a/src/GPU/pppm_gpu.cpp +++ b/src/GPU/pppm_gpu.cpp @@ -102,6 +102,8 @@ PPPMGPU::PPPMGPU(LAMMPS *lmp) : PPPM(lmp) PPPMGPU::~PPPMGPU() { PPPM_GPU_API(clear)(poisson_time); + destroy_3d_offset(density_brick_gpu,nzlo_out,nylo_out); + destroy_3d_offset(vd_brick,nzlo_out,nylo_out); } /* ---------------------------------------------------------------------- diff --git a/src/atom.cpp b/src/atom.cpp index 8b78b4f8f7..0de44e50ca 100644 --- a/src/atom.cpp +++ b/src/atom.cpp @@ -2358,16 +2358,13 @@ void Atom::setup_sort_bins() } #ifdef LMP_GPU - if (userbinsize == 0.0) { - int ifix = modify->find_fix("package_gpu"); - if (ifix >= 0) { + if (userbinsize == 0.0) { + auto ifix = dynamic_cast(modify->get_fix_by_id("package_gpu")); + if (ifix) { const double subx = domain->subhi[0] - domain->sublo[0]; const double suby = domain->subhi[1] - domain->sublo[1]; const double subz = domain->subhi[2] - domain->sublo[2]; - - FixGPU *fix = static_cast(modify->fix[ifix]); - binsize = fix->binsize(subx, suby, subz, atom->nlocal, - 0.5 * neighbor->cutneighmax); + binsize = ifix->binsize(subx, suby, subz, atom->nlocal, 0.5 * neighbor->cutneighmax); } } #endif diff --git a/tools/tinker/tinker2lmp.py b/tools/tinker/tinker2lmp.py index d376593ea3..e3ae59748c 100644 --- a/tools/tinker/tinker2lmp.py +++ b/tools/tinker/tinker2lmp.py @@ -227,11 +227,7 @@ class XYZfile(object): print(i+1,label[i],x[i],y[i],z[i],type[i], end=' ', file=fp) for j in bonds[i]: print(j, end=' ', file=fp) print(file=fp) -<<<<<<< HEAD - -======= ->>>>>>> develop fp.close() # triplet of atoms in an angle = atom 1,2,3 @@ -1098,16 +1094,6 @@ for i,one in enumerate(alist): elif len(params[3]) == 2: nbonds,hcount = xyz.angle_hbond_count(atom1,atom2,atom3,lmptype,lmpmass) -<<<<<<< HEAD - - if nbonds != 3: - print("Center angle atom has wrong bond count") - print(" angle atom IDs:",atom1,atom2,atom3) - print(" angle atom classes:",c1,c2,c3) - print(" Tinker FF file param options:",len(params[3])) - print(" Nbonds and hydrogen count:",nbonds,hcount) - #sys.exit() NOTE: allow this for now -======= #if nbonds != 3: #print("Center angle atom has wrong bond count") @@ -1117,33 +1103,12 @@ for i,one in enumerate(alist): #print(" Nbonds and hydrogen count:",nbonds,hcount) # NOTE: allow this for now #sys.exit() ->>>>>>> develop if hcount == 0: which = 1 elif hcount == 1: which = 2 m += 1 -<<<<<<< HEAD - print("3-bond angle") - print(" angle atom IDs:",atom1,atom2,atom3) - print(" angle atom classes:",c1,c2,c3) - print(" Tinker FF file param options:",len(params[3])) - print(" Nbonds and hydrogen count:",nbonds,hcount) - print(" which:",which,m) - - elif len(params[3]) == 3: - nbonds,hcount = xyz.angle_hbond_count(atom1,atom2,atom3,lmptype,lmpmass) - - if nbonds != 4: - print("Center angle atom has wrong bond count") - print(" angle atom IDs:",atom1,atom2,atom3) - print(" angle atom classes:",c1,c2,c3) - print(" Tinker FF file param options:",len(params[3])) - print(" Nbonds and hydrogen count:",nbonds,hcount) - #sys.exit() NOTE: allow this for now - -======= #print("3-bond angle") #print(" angle atom IDs:",atom1,atom2,atom3) #print(" angle atom classes:",c1,c2,c3) @@ -1163,7 +1128,6 @@ for i,one in enumerate(alist): # NOTE: allow this for now #sys.exit() ->>>>>>> develop if hcount == 0: which = 1 elif hcount == 1: which = 2 @@ -1207,12 +1171,8 @@ for itype in range(len(aparams)): elif (c3,c2,c1) in badict: n1,n2,r1,r2 = badict[(c3,c2,c1)] else: -<<<<<<< HEAD - print("Bond-stretch angle triplet not found: %d %d %d" % (c1,c2,c3)) -======= # NOTE: just for debugging #print("Bond-stretch angle triplet not found: %d %d %d" % (c1,c2,c3)) ->>>>>>> develop n1,n2,r1,r2 = 4*[0.0] baparams.append((n1,n2,r1,r2)) @@ -1670,11 +1630,7 @@ print("Natoms =",natoms) print("Ntypes =",ntypes) print("Tinker XYZ types =",len(tink2lmp)) print("Tinker PRM types =",prm.ntypes) -<<<<<<< HEAD -#print "Tinker groups =",ngroups -======= #print("Tinker groups =",ngroups) ->>>>>>> develop print("Nmol =",nmol) print("Nbonds =",nbonds) print("Nangles =",nangles)