diff --git a/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp b/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp
index adb3c98232..b480b644a1 100644
--- a/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp
+++ b/src/KOKKOS/fix_acks2_reaxff_kokkos.cpp
@@ -205,7 +205,6 @@ void FixACKS2ReaxFFKokkos<DeviceType>::pre_force(int vflag)
   type = atomKK->k_type.view<DeviceType>();
   mask = atomKK->k_mask.view<DeviceType>();
   nlocal = atomKK->nlocal;
-  nall = atom->nlocal + atom->nghost;
   newton_pair = force->newton_pair;
 
   k_params.template sync<DeviceType>();
diff --git a/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp b/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp
index 76e1f1e94d..43fd832700 100644
--- a/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp
+++ b/src/KOKKOS/fix_qeq_reaxff_kokkos.cpp
@@ -871,7 +871,8 @@ void FixQEqReaxFFKokkos<DeviceType>::sparse_matvec_kokkos(typename AT::t_ffloat2
   }
 
   if (neighflag != FULL) {
-    Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType,TagQEqZeroQGhosts>(atom->nlocal,atom->nghost),*this);
+    int nall = nlocal + atomKK->nghost;
+    Kokkos::parallel_for(Kokkos::RangePolicy<DeviceType,TagQEqZeroQGhosts>(atom->nlocal,nall),*this);
 
     if (need_dup)
       dup_o.reset_except(d_o);
diff --git a/src/OPENMP/fix_qeq_reaxff_omp.cpp b/src/OPENMP/fix_qeq_reaxff_omp.cpp
index e0714fc636..ee9748d354 100644
--- a/src/OPENMP/fix_qeq_reaxff_omp.cpp
+++ b/src/OPENMP/fix_qeq_reaxff_omp.cpp
@@ -254,8 +254,6 @@ void FixQEqReaxFFOMP::pre_force(int /* vflag */)
 {
   if (update->ntimestep % nevery) return;
 
-  int n = atom->nlocal;
-
   if (reaxff) {
     nn = reaxff->list->inum;
     ilist = reaxff->list->ilist;
@@ -272,7 +270,7 @@ void FixQEqReaxFFOMP::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();
@@ -486,6 +484,9 @@ void FixQEqReaxFFOMP::sparse_matvec(sparse_matrix *A, double *x, double *b)
     int i, j, itr_j;
     int ii;
     int nthreads = comm->nthreads;
+    int nlocal = atom->nlocal;
+    int nall = atom->nlocal + atom->nghost;
+
 #if defined(_OPENMP)
     int tid = omp_get_thread_num();
 #else
@@ -503,14 +504,14 @@ void FixQEqReaxFFOMP::sparse_matvec(sparse_matrix *A, double *x, double *b)
 #if defined(_OPENMP)
 #pragma omp for schedule(dynamic,50)
 #endif
-    for (i = atom->nlocal; i < atom->nghost; ++i) {
+    for (i = nlocal; i < nall; ++i) {
       if (atom->mask[i] & groupbit) b[i] = 0;
     }
 
 #if defined(_OPENMP)
 #pragma omp for schedule(dynamic,50)
 #endif
-    for (i = 0; i < atom->nghost; ++i)
+    for (i = 0; i < nall; ++i)
       for (int t=0; t<nthreads; t++) b_temp[t][i] = 0.0;
 
     // Wait for b accumulated and b_temp zeroed.
@@ -535,7 +536,7 @@ void FixQEqReaxFFOMP::sparse_matvec(sparse_matrix *A, double *x, double *b)
 #pragma omp barrier
 #pragma omp for schedule(dynamic,50)
 #endif
-    for (i = 0; i < atom->nghost; ++i)
+    for (i = 0; i < nall; ++i)
       for (int t = 0; t < nthreads; ++t) b[i] += b_temp[t][i];
 
   } //end omp parallel
@@ -835,7 +836,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 = atom->nlocal; i < atom->nghost; ++i) {
+    for (i = nlocal; i < nall; ++i) {
       if (atom->mask[i] & groupbit) {
         indxI = 2 * i;
         b[indxI]   = 0;
@@ -846,7 +847,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 < atom->nghost; ++i) {
+    for (i = 0; i < nall; ++i) {
       indxI = 2 * i;
       for (int t=0; t<nthreads; t++) {
         b_temp[t][indxI] = 0.0;
@@ -880,7 +881,7 @@ void FixQEqReaxFFOMP::dual_sparse_matvec(sparse_matrix *A, double *x1, double *x
 #pragma omp barrier
 #pragma omp for schedule(dynamic,50)
 #endif
-    for (i = 0; i < atom->nghost; ++i) {
+    for (i = 0; i < nall; ++i) {
       indxI = 2 * i;
       for (int t = 0; t < nthreads; ++t) {
         b[indxI] += b_temp[t][indxI];
@@ -925,7 +926,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 = atom->nlocal; i < atom->nghost; ++i) {
+    for (i = nlocal; i < nall; ++i) {
       if (atom->mask[i] & groupbit) {
         indxI = 2 * i;
         b[indxI]   = 0;
@@ -936,7 +937,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 < atom->nghost; ++i) {
+    for (i = 0; i < nall; ++i) {
       indxI = 2 * i;
       for (int t=0; t<nthreads; t++) {
         b_temp[t][indxI] = 0.0;
@@ -970,7 +971,7 @@ void FixQEqReaxFFOMP::dual_sparse_matvec(sparse_matrix *A, double *x, double *b)
 #pragma omp barrier
 #pragma omp for schedule(dynamic,50)
 #endif
-    for (i = 0; i < atom->nghost; ++i) {
+    for (i = 0; i < nall; ++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 ad026102e0..d9cf4c13b8 100644
--- a/src/REAXFF/fix_acks2_reaxff.cpp
+++ b/src/REAXFF/fix_acks2_reaxff.cpp
@@ -624,7 +624,7 @@ void FixACKS2ReaxFF::sparse_matvec_acks2(sparse_matrix *H, sparse_matrix *X, dou
     }
   }
 
-  for (i = atom->nlocal; i < atom->nghost; ++i) {
+  for (i = atom->nlocal; i < NN; ++i) {
     if (atom->mask[i] & groupbit) {
       b[i] = 0;
       b[NN + i] = 0;
diff --git a/src/REAXFF/fix_qeq_reaxff.cpp b/src/REAXFF/fix_qeq_reaxff.cpp
index 12300d7054..b6c53aa977 100644
--- a/src/REAXFF/fix_qeq_reaxff.cpp
+++ b/src/REAXFF/fix_qeq_reaxff.cpp
@@ -786,7 +786,8 @@ void FixQEqReaxFF::sparse_matvec(sparse_matrix *A, double *x, double *b)
       b[i] = eta[atom->type[i]] * x[i];
   }
 
-  for (i = atom->nlocal; i < atom->nghost; ++i)
+  int nall = atom->nlocal + atom->nghost;
+  for (i = atom->nlocal; i < nall; ++i)
       b[i] = 0;
 
   for (ii = 0; ii < nn; ++ii) {