diff --git a/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp b/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp index d8fbe32b9a..adb3c98232 100644 --- a/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp +++ b/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp @@ -219,7 +219,7 @@ void FixACKS2ReaxFFKokkos::pre_force(int vflag) d_ilist = k_list->d_ilist; nn = list->inum; - NN = list->inum + list->gnum; + NN = atom->nlocal + atom->nghost; copymode = 1; @@ -526,7 +526,7 @@ void FixACKS2ReaxFFKokkos::allocate_array() if (efield) get_chi_field(); // init_storage - Kokkos::parallel_for(Kokkos::RangePolicy(0,NN),*this); + Kokkos::parallel_for(Kokkos::RangePolicy(0,nn),*this); } @@ -1378,18 +1378,15 @@ template void FixACKS2ReaxFFKokkos::calculate_Q() { - Kokkos::parallel_for(Kokkos::RangePolicy(0,nn),*this); + Kokkos::parallel_for(Kokkos::RangePolicy(0,nn),*this); - pack_flag = 2; - //comm->forward_comm( this ); //Dist_vector( s ); - k_s.modify(); - k_s.sync(); + pack_flag = 4; + //comm->forward_comm( this ); //Dist_vector( q ); + atomKK->k_q.modify(); + atomKK->k_q.sync(); comm->forward_comm(this); - k_s.modify(); - k_s.sync(); - - Kokkos::parallel_for(Kokkos::RangePolicy(0,NN),*this); - + atomKK->k_q.modify(); + atomKK->k_q.sync(); } /* ---------------------------------------------------------------------- */ @@ -1822,11 +1819,13 @@ void FixACKS2ReaxFFKokkos::operator() (TagACKS2Norm3, const int &ii, template KOKKOS_INLINE_FUNCTION -void FixACKS2ReaxFFKokkos::operator() (TagACKS2CalculateQ1, const int &ii) const +void FixACKS2ReaxFFKokkos::operator() (TagACKS2CalculateQ, const int &ii) const { const int i = d_ilist[ii]; if (mask[i] & groupbit) { + q(i) = d_s(i); + /* backup s */ for (int k = nprev-1; k > 0; --k) { d_s_hist(i,k) = d_s_hist(i,k-1); @@ -1848,17 +1847,6 @@ void FixACKS2ReaxFFKokkos::operator() (TagACKS2CalculateQ1, const in /* ---------------------------------------------------------------------- */ -template -KOKKOS_INLINE_FUNCTION -void FixACKS2ReaxFFKokkos::operator() (TagACKS2CalculateQ2, const int &ii) const -{ - const int i = d_ilist[ii]; - if (mask[i] & groupbit) - q(i) = d_s(i); -} - -/* ---------------------------------------------------------------------- */ - template void FixACKS2ReaxFFKokkos::cleanup_copy() { diff --git a/src/KOKKOS/fix_acks2_reaxff_kokkos.h b/src/KOKKOS/fix_acks2_reaxff_kokkos.h index 735b478f1b..e81479d3ba 100644 --- a/src/KOKKOS/fix_acks2_reaxff_kokkos.h +++ b/src/KOKKOS/fix_acks2_reaxff_kokkos.h @@ -54,8 +54,7 @@ struct TagACKS2Precon1B{}; struct TagACKS2Precon2{}; struct TagACKS2Add{}; struct TagACKS2ZeroQGhosts{}; -struct TagACKS2CalculateQ1{}; -struct TagACKS2CalculateQ2{}; +struct TagACKS2CalculateQ{}; template class FixACKS2ReaxFFKokkos : public FixACKS2ReaxFF { @@ -152,10 +151,7 @@ class FixACKS2ReaxFFKokkos : public FixACKS2ReaxFF { void operator()(TagACKS2ZeroQGhosts, const int&) const; KOKKOS_INLINE_FUNCTION - void operator()(TagACKS2CalculateQ1, const int&) const; - - KOKKOS_INLINE_FUNCTION - void operator()(TagACKS2CalculateQ2, const int&) const; + void operator()(TagACKS2CalculateQ, const int&) const; KOKKOS_INLINE_FUNCTION double calculate_H_k(const F_FLOAT &r, const F_FLOAT &shld) const; diff --git a/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp b/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp index 4aa23b77be..76e1f1e94d 100644 --- a/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp +++ b/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp @@ -208,7 +208,6 @@ void FixQEqReaxFFKokkos::pre_force(int /*vflag*/) d_neighbors = k_list->d_neighbors; d_ilist = k_list->d_ilist; nn = list->inum; - NN = list->inum + list->gnum; copymode = 1; @@ -375,7 +374,7 @@ void FixQEqReaxFFKokkos::allocate_array() if (efield) get_chi_field(); - Kokkos::parallel_for(Kokkos::RangePolicy(0,NN),*this); + Kokkos::parallel_for(Kokkos::RangePolicy(0,nn),*this); } /* ---------------------------------------------------------------------- */ @@ -872,7 +871,7 @@ void FixQEqReaxFFKokkos::sparse_matvec_kokkos(typename AT::t_ffloat2 } if (neighflag != FULL) { - Kokkos::parallel_for(Kokkos::RangePolicy(nn,NN),*this); + Kokkos::parallel_for(Kokkos::RangePolicy(atom->nlocal,atom->nghost),*this); if (need_dup) dup_o.reset_except(d_o); diff --git a/src/KOKKOS/pair_reaxff_kokkos.cpp b/src/KOKKOS/pair_reaxff_kokkos.cpp index 714c0c4914..2408f6649d 100644 --- a/src/KOKKOS/pair_reaxff_kokkos.cpp +++ b/src/KOKKOS/pair_reaxff_kokkos.cpp @@ -715,7 +715,7 @@ void PairReaxFFKokkos::compute(int eflag_in, int vflag_in) newton_pair = force->newton_pair; nn = list->inum; - NN = list->inum + list->gnum; + NN = atom->nlocal + atom->nghost; const int inum = list->inum; const int ignum = inum + list->gnum; diff --git a/src/OPENMP/fix_qeq_reaxff_omp.cpp b/src/OPENMP/fix_qeq_reaxff_omp.cpp index 2bcf639eae..e0714fc636 100644 --- a/src/OPENMP/fix_qeq_reaxff_omp.cpp +++ b/src/OPENMP/fix_qeq_reaxff_omp.cpp @@ -237,7 +237,7 @@ void FixQEqReaxFFOMP::init_storage() #if defined(_OPENMP) #pragma omp parallel for schedule(static) #endif - for (int i = 0; i < NN; i++) { + for (int i = 0; i < nn; i++) { Hdia_inv[i] = 1. / eta[atom->type[i]]; b_s[i] = -chi[atom->type[i]]; if (efield) b_s[i] -= chi_field[i]; @@ -258,13 +258,11 @@ void FixQEqReaxFFOMP::pre_force(int /* vflag */) if (reaxff) { nn = reaxff->list->inum; - NN = reaxff->list->inum + reaxff->list->gnum; ilist = reaxff->list->ilist; numneigh = reaxff->list->numneigh; firstneigh = reaxff->list->firstneigh; } else { nn = list->inum; - NN = list->inum + list->gnum; ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; @@ -505,15 +503,14 @@ void FixQEqReaxFFOMP::sparse_matvec(sparse_matrix *A, double *x, double *b) #if defined(_OPENMP) #pragma omp for schedule(dynamic,50) #endif - for (ii = nn; ii < NN; ++ii) { - i = ilist[ii]; + for (i = atom->nlocal; i < atom->nghost; ++i) { if (atom->mask[i] & groupbit) b[i] = 0; } #if defined(_OPENMP) #pragma omp for schedule(dynamic,50) #endif - for (i = 0; i < NN; ++i) + for (i = 0; i < atom->nghost; ++i) for (int t=0; tnghost; ++i) for (int t = 0; t < nthreads; ++t) b[i] += b_temp[t][i]; } //end omp parallel @@ -838,8 +835,7 @@ void FixQEqReaxFFOMP::dual_sparse_matvec(sparse_matrix *A, double *x1, double *x #if defined(_OPENMP) #pragma omp for schedule(dynamic,50) #endif - for (ii = nn; ii < NN; ++ii) { - i = ilist[ii]; + for (i = atom->nlocal; i < atom->nghost; ++i) { if (atom->mask[i] & groupbit) { indxI = 2 * i; b[indxI] = 0; @@ -850,7 +846,7 @@ void FixQEqReaxFFOMP::dual_sparse_matvec(sparse_matrix *A, double *x1, double *x #if defined(_OPENMP) #pragma omp for schedule(dynamic,50) #endif - for (i = 0; i < NN; ++i) { + for (i = 0; i < atom->nghost; ++i) { indxI = 2 * i; for (int t=0; tnghost; ++i) { indxI = 2 * i; for (int t = 0; t < nthreads; ++t) { b[indxI] += b_temp[t][indxI]; @@ -929,8 +925,7 @@ void FixQEqReaxFFOMP::dual_sparse_matvec(sparse_matrix *A, double *x, double *b) #if defined(_OPENMP) #pragma omp for schedule(dynamic,50) #endif - for (ii = nn; ii < NN; ++ii) { - i = ilist[ii]; + for (i = atom->nlocal; i < atom->nghost; ++i) { if (atom->mask[i] & groupbit) { indxI = 2 * i; b[indxI] = 0; @@ -941,7 +936,7 @@ void FixQEqReaxFFOMP::dual_sparse_matvec(sparse_matrix *A, double *x, double *b) #if defined(_OPENMP) #pragma omp for schedule(dynamic,50) #endif - for (i = 0; i < NN; ++i) { + for (i = 0; i < atom->nghost; ++i) { indxI = 2 * i; for (int t=0; tnghost; ++i) { indxI = 2 * i; for (int t = 0; t < nthreads; ++t) { b[indxI] += b_temp[t][indxI]; diff --git a/src/REAXFF/fix_acks2_reaxff.cpp b/src/REAXFF/fix_acks2_reaxff.cpp index c3457d1712..ad026102e0 100644 --- a/src/REAXFF/fix_acks2_reaxff.cpp +++ b/src/REAXFF/fix_acks2_reaxff.cpp @@ -306,7 +306,7 @@ void FixACKS2ReaxFF::init_storage() { if (efield) get_chi_field(); - for (int ii = 0; ii < NN; ii++) { + for (int ii = 0; ii < nn; ii++) { int i = ilist[ii]; if (atom->mask[i] & groupbit) { b_s[i] = -chi[atom->type[i]]; @@ -329,17 +329,15 @@ void FixACKS2ReaxFF::pre_force(int /*vflag*/) { if (update->ntimestep % nevery) return; - int n = atom->nlocal; + NN = atom->nlocal + atom->nghost; if (reaxff) { nn = reaxff->list->inum; - NN = reaxff->list->inum + reaxff->list->gnum; ilist = reaxff->list->ilist; numneigh = reaxff->list->numneigh; firstneigh = reaxff->list->firstneigh; } else { nn = list->inum; - NN = list->inum + list->gnum; ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; @@ -349,7 +347,7 @@ void FixACKS2ReaxFF::pre_force(int /*vflag*/) // need to be atom->nmax in length if (atom->nmax > nmax) reallocate_storage(); - if (n > n_cap*DANGER_ZONE || m_fill > m_cap*DANGER_ZONE) + if (atom->nlocal > n_cap*DANGER_ZONE || m_fill > m_cap*DANGER_ZONE) reallocate_matrix(); if (efield) get_chi_field(); @@ -626,8 +624,7 @@ void FixACKS2ReaxFF::sparse_matvec_acks2(sparse_matrix *H, sparse_matrix *X, dou } } - for (ii = nn; ii < NN; ++ii) { - i = ilist[ii]; + for (i = atom->nlocal; i < atom->nghost; ++i) { if (atom->mask[i] & groupbit) { b[i] = 0; b[NN + i] = 0; @@ -680,6 +677,8 @@ void FixACKS2ReaxFF::calculate_Q() i = ilist[ii]; if (atom->mask[i] & groupbit) { + atom->q[i] = s[i]; + /* backup s */ for (k = nprev-1; k > 0; --k) { s_hist[i][k] = s_hist[i][k-1]; @@ -698,14 +697,8 @@ void FixACKS2ReaxFF::calculate_Q() } } - pack_flag = 2; - comm->forward_comm(this); //Dist_vector(s); - - for (int ii = 0; ii < NN; ++ii) { - i = ilist[ii]; - if (atom->mask[i] & groupbit) - atom->q[i] = s[i]; - } + pack_flag = 4; + comm->forward_comm(this); //Dist_vector(q); } /* ---------------------------------------------------------------------- */ @@ -733,6 +726,11 @@ int FixACKS2ReaxFF::pack_forward_comm(int n, int *list, double *buf, buf[m++] = q_hat[j]; buf[m++] = q_hat[NN+j]; } + } else if (pack_flag == 4) { + for(int i = 0; i < n; i++) { + int j = list[i]; + buf[m++] = atom->q[j]; + } } return m; } @@ -761,6 +759,10 @@ void FixACKS2ReaxFF::unpack_forward_comm(int n, int first, double *buf) q_hat[i] = buf[m++]; q_hat[NN+i] = buf[m++]; } + } else if (pack_flag == 4) { + for(i = first; i < last; i++) { + atom->q[i] = buf[m++]; + } } } diff --git a/src/REAXFF/fix_acks2_reaxff.h b/src/REAXFF/fix_acks2_reaxff.h index 2d21f80fe0..27ba196bbd 100644 --- a/src/REAXFF/fix_acks2_reaxff.h +++ b/src/REAXFF/fix_acks2_reaxff.h @@ -37,7 +37,7 @@ class FixACKS2ReaxFF : public FixQEqReaxFF { double *get_s() { return s; } protected: - int last_rows_rank, last_rows_flag; + int NN, last_rows_rank, last_rows_flag; double **s_hist_X, **s_hist_last; double *bcut_acks2, bond_softness, **bcut; // acks2 parameters diff --git a/src/REAXFF/fix_qeq_reaxff.cpp b/src/REAXFF/fix_qeq_reaxff.cpp index c2469ee9eb..12300d7054 100644 --- a/src/REAXFF/fix_qeq_reaxff.cpp +++ b/src/REAXFF/fix_qeq_reaxff.cpp @@ -104,7 +104,7 @@ FixQEqReaxFF::FixQEqReaxFF(LAMMPS *lmp, int narg, char **arg) : shld = nullptr; nn = n_cap = 0; - NN = nmax = 0; + nmax = 0; m_fill = m_cap = 0; pack_flag = 0; s = nullptr; @@ -319,7 +319,7 @@ void FixQEqReaxFF::reallocate_storage() void FixQEqReaxFF::allocate_matrix() { - int i,ii,n,m; + int i,ii,m; int mincap; double safezone; @@ -332,8 +332,7 @@ void FixQEqReaxFF::allocate_matrix() safezone = REAX_SAFE_ZONE; } - n = atom->nlocal; - n_cap = MAX((int)(n * safezone), mincap); + n_cap = MAX((int)(atom->nlocal * safezone), mincap); // determine the total space for the H matrix @@ -493,13 +492,11 @@ void FixQEqReaxFF::setup_pre_force(int vflag) { if (reaxff) { nn = reaxff->list->inum; - NN = reaxff->list->inum + reaxff->list->gnum; ilist = reaxff->list->ilist; numneigh = reaxff->list->numneigh; firstneigh = reaxff->list->firstneigh; } else { nn = list->inum; - NN = list->inum + list->gnum; ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; @@ -537,7 +534,7 @@ void FixQEqReaxFF::init_storage() { if (efield) get_chi_field(); - for (int ii = 0; ii < NN; ii++) { + for (int ii = 0; ii < nn; ii++) { int i = ilist[ii]; if (atom->mask[i] & groupbit) { Hdia_inv[i] = 1. / eta[atom->type[i]]; @@ -561,13 +558,11 @@ void FixQEqReaxFF::pre_force(int /*vflag*/) if (reaxff) { nn = reaxff->list->inum; - NN = reaxff->list->inum + reaxff->list->gnum; ilist = reaxff->list->ilist; numneigh = reaxff->list->numneigh; firstneigh = reaxff->list->firstneigh; } else { nn = list->inum; - NN = list->inum + list->gnum; ilist = list->ilist; numneigh = list->numneigh; firstneigh = list->firstneigh; @@ -791,11 +786,8 @@ void FixQEqReaxFF::sparse_matvec(sparse_matrix *A, double *x, double *b) b[i] = eta[atom->type[i]] * x[i]; } - for (ii = nn; ii < NN; ++ii) { - i = ilist[ii]; - if (atom->mask[i] & groupbit) + for (i = atom->nlocal; i < atom->nghost; ++i) b[i] = 0; - } for (ii = 0; ii < nn; ++ii) { i = ilist[ii]; diff --git a/src/REAXFF/fix_qeq_reaxff.h b/src/REAXFF/fix_qeq_reaxff.h index 5ebb410a10..f0f46ed4ce 100644 --- a/src/REAXFF/fix_qeq_reaxff.h +++ b/src/REAXFF/fix_qeq_reaxff.h @@ -59,7 +59,7 @@ class FixQEqReaxFF : public Fix { protected: int nevery, reaxflag; int matvecs; - int nn, NN, m_fill; + int nn, m_fill; int n_cap, nmax, m_cap; int pack_flag; int nlevels_respa;