diff --git a/src/KOKKOS/atom_kokkos.h b/src/KOKKOS/atom_kokkos.h
index 6eebbad661..b554442b37 100644
--- a/src/KOKKOS/atom_kokkos.h
+++ b/src/KOKKOS/atom_kokkos.h
@@ -100,18 +100,18 @@ struct SortFunctor {
     dest(i) = source(index(i));
   }
   void operator()(const typename std::enable_if<ViewType::rank==2, int>::type& i) {
-    for(int j=0; j < (int)source.extent(1); j++)
+    for (int j=0; j < (int)source.extent(1); j++)
       dest(i,j) = source(index(i),j);
   }
   void operator()(const typename std::enable_if<ViewType::rank==3, int>::type& i) {
-    for(int j=0; j < (int)source.extent(1); j++)
-      for(int k=0; k < (int)source.extent(2); k++)
+    for (int j=0; j < (int)source.extent(1); j++)
+      for (int k=0; k < (int)source.extent(2); k++)
         dest(i,j,k) = source(index(i),j,k);
   }
   void operator()(const typename std::enable_if<ViewType::rank==4, int>::type& i) {
-    for(int j=0; j < (int)source.extent(1); j++)
-      for(int k=0; k < (int)source.extent(2); k++)
-        for(int l=0; l < (int)source.extent(3); l++)
+    for (int j=0; j < (int)source.extent(1); j++)
+      for (int k=0; k < (int)source.extent(2); k++)
+        for (int l=0; l < (int)source.extent(3); l++)
           dest(i,j,k,l) = source(index(i),j,k,l);
   }
 };
diff --git a/src/KOKKOS/atom_vec_kokkos.h b/src/KOKKOS/atom_vec_kokkos.h
index 9fbf172535..f81b7715ad 100644
--- a/src/KOKKOS/atom_vec_kokkos.h
+++ b/src/KOKKOS/atom_vec_kokkos.h
@@ -178,7 +178,7 @@ class AtomVecKokkos : public AtomVec {
     }
     mirror_type tmp_view((typename ViewType::value_type*)buffer, src.d_view.layout());
 
-    if(space == Device) {
+    if (space == Device) {
       Kokkos::deep_copy(LMPHostType(),tmp_view,src.h_view),
       Kokkos::deep_copy(LMPHostType(),src.d_view,tmp_view);
       src.clear_sync_state();
@@ -191,7 +191,7 @@ class AtomVecKokkos : public AtomVec {
   #else
   template<class ViewType>
   void perform_async_copy(ViewType& src, unsigned int space) {
-    if(space == Device)
+    if (space == Device)
       src.template sync<LMPDeviceType>();
     else
       src.template sync<LMPHostType>();
diff --git a/src/KOKKOS/memory_kokkos.h b/src/KOKKOS/memory_kokkos.h
index 2f9e0cc375..445ecd87f2 100644
--- a/src/KOKKOS/memory_kokkos.h
+++ b/src/KOKKOS/memory_kokkos.h
@@ -92,7 +92,7 @@ void destroy_kokkos(TYPE data, typename TYPE::value_type* &array)
 template <typename TYPE>
 TYPE destroy_kokkos(TYPE &data)
 {
-  /*if(data.data()!=nullptr)
+  /*if (data.data()!=nullptr)
     free(data.data());*/
   data = TYPE();
   return data;
@@ -173,7 +173,7 @@ TYPE create_kokkos(TYPE &data, typename TYPE::value_type **&array,
 
   bigint n = 0;
   for (int i = 0; i < n1; i++) {
-    if(n2==0)
+    if (n2==0)
       array[i] = nullptr;
     else
       array[i] = &data.h_view(i,0);
@@ -194,7 +194,7 @@ template <typename TYPE, typename HTYPE>
 
   bigint n = 0;
   for (int i = 0; i < n1; i++) {
-    if(n2==0)
+    if (n2==0)
       array[i] = nullptr;
     else
       array[i] = &h_data(i,0);
@@ -218,7 +218,7 @@ TYPE grow_kokkos(TYPE &data, typename TYPE::value_type **&array,
   array = (typename TYPE::value_type**) srealloc(array,nbytes,name);
 
   for (int i = 0; i < n1; i++)
-    if(n2==0)
+    if (n2==0)
       array[i] = nullptr;
     else
       array[i] = &data.h_view(i,0);
@@ -235,7 +235,7 @@ TYPE create_kokkos(TYPE &data, typename TYPE::value_type **&array,
   array = (typename TYPE::value_type **) smalloc(nbytes,name);
 
   for (int i = 0; i < n1; i++)
-    if(data.h_view.extent(1)==0)
+    if (data.h_view.extent(1)==0)
       array[i] = nullptr;
     else
       array[i] = &data.h_view(i,0);
@@ -255,7 +255,7 @@ TYPE grow_kokkos(TYPE &data, typename TYPE::value_type **&array,
   array = (typename TYPE::value_type **) smalloc(nbytes,name);
 
   for (int i = 0; i < n1; i++)
-    if(data.h_view.extent(1)==0)
+    if (data.h_view.extent(1)==0)
       array[i] = nullptr;
     else
       array[i] = &data.h_view(i,0);
diff --git a/src/KOKKOS/nbin_ssa_kokkos.h b/src/KOKKOS/nbin_ssa_kokkos.h
index cc98859913..a99e0d448c 100644
--- a/src/KOKKOS/nbin_ssa_kokkos.h
+++ b/src/KOKKOS/nbin_ssa_kokkos.h
@@ -108,20 +108,20 @@ class NBinSSAKokkos : public NBinStandard {
     if (y >= subhi_[1]) iy = 1;
     if (x < sublo_[0]) ix = -1;
     if (x >= subhi_[0]) ix = 1;
-    if(iz < 0){
+    if (iz < 0){
       return -1;
-    } else if(iz == 0){
-      if( iy<0 ) return -1; // bottom left/middle/right
-      if( (iy==0) && (ix<0)  ) return -1; // left atoms
-      if( (iy==0) && (ix==0) ) return 0; // Locally owned atoms
-      if( (iy==0) && (ix>0)  ) return 2; // Right atoms
-      if( (iy>0)  && (ix==0) ) return 1; // Top-middle atoms
-      if( (iy>0)  && (ix!=0) ) return 3; // Top-right and top-left atoms
+    } else if (iz == 0){
+      if (iy<0) return -1; // bottom left/middle/right
+      if ((iy==0) && (ix<0) ) return -1; // left atoms
+      if ((iy==0) && (ix==0)) return 0; // Locally owned atoms
+      if ((iy==0) && (ix>0) ) return 2; // Right atoms
+      if ((iy>0)  && (ix==0)) return 1; // Top-middle atoms
+      if ((iy>0)  && (ix!=0)) return 3; // Top-right and top-left atoms
     } else { // iz > 0
-      if((ix==0) && (iy==0)) return 4; // Back atoms
-      if((ix==0) && (iy!=0)) return 5; // Top-back and bottom-back atoms
-      if((ix!=0) && (iy==0)) return 6; // Left-back and right-back atoms
-      if((ix!=0) && (iy!=0)) return 7; // Back corner atoms
+      if ((ix==0) && (iy==0)) return 4; // Back atoms
+      if ((ix==0) && (iy!=0)) return 5; // Top-back and bottom-back atoms
+      if ((ix!=0) && (iy==0)) return 6; // Left-back and right-back atoms
+      if ((ix!=0) && (iy!=0)) return 7; // Back corner atoms
     }
     return -2;
   }
diff --git a/src/KOKKOS/pair_kokkos.h b/src/KOKKOS/pair_kokkos.h
index cb55dd3141..7dc3e33f74 100644
--- a/src/KOKKOS/pair_kokkos.h
+++ b/src/KOKKOS/pair_kokkos.h
@@ -144,7 +144,7 @@ struct PairComputeFunctor  {
       const int jtype = c.type(j);
       const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
 
-      if(rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
+      if (rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
 
         const F_FLOAT fpair = factor_lj*c.template compute_fpair<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
 
@@ -213,13 +213,13 @@ struct PairComputeFunctor  {
       const int jtype = c.type(j);
       const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
 
-      if(rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
+      if (rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
 
         F_FLOAT fpair = F_FLOAT();
 
-        if(rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype)))
+        if (rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype)))
           fpair+=factor_lj*c.template compute_fpair<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
-        if(rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype)))
+        if (rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype)))
           fpair+=c.template compute_fcoul<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype,factor_coul,qtmp);
 
         fxtmp += delx*fpair;
@@ -236,11 +236,11 @@ struct PairComputeFunctor  {
           F_FLOAT evdwl = 0.0;
           F_FLOAT ecoul = 0.0;
           if (c.eflag) {
-            if(rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype))) {
+            if (rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype))) {
               evdwl = factor_lj * c.template compute_evdwl<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
               ev.evdwl += (((NEIGHFLAG==HALF || NEIGHFLAG==HALFTHREAD)&&(NEWTON_PAIR||(j<c.nlocal)))?1.0:0.5)*evdwl;
             }
-            if(rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype))) {
+            if (rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype))) {
               ecoul = c.template compute_ecoul<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype,factor_coul,qtmp);
               ev.ecoul += (((NEIGHFLAG==HALF || NEIGHFLAG==HALFTHREAD)&&(NEWTON_PAIR||(j<c.nlocal)))?1.0:0.5)*ecoul;
             }
@@ -294,7 +294,7 @@ struct PairComputeFunctor  {
         const int jtype = c.type(j);
         const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
 
-        if(rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
+        if (rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
 
           const F_FLOAT fpair = factor_lj*c.template compute_fpair<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
 
@@ -351,13 +351,13 @@ struct PairComputeFunctor  {
         const int jtype = c.type(j);
         const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
 
-        if(rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
+        if (rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
 
           F_FLOAT fpair = F_FLOAT();
 
-          if(rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype)))
+          if (rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype)))
             fpair+=factor_lj*c.template compute_fpair<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
-          if(rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype)))
+          if (rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype)))
             fpair+=c.template compute_fcoul<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype,factor_coul,qtmp);
 
           ftmp.x += delx*fpair;
@@ -413,7 +413,7 @@ struct PairComputeFunctor  {
         const int jtype = c.type(j);
         const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
 
-        if(rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
+        if (rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
 
           const F_FLOAT fpair = factor_lj*c.template compute_fpair<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
 
@@ -510,13 +510,13 @@ struct PairComputeFunctor  {
         const int jtype = c.type(j);
         const F_FLOAT rsq = delx*delx + dely*dely + delz*delz;
 
-        if(rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
+        if (rsq < (STACKPARAMS?c.m_cutsq[itype][jtype]:c.d_cutsq(itype,jtype))) {
 
           F_FLOAT fpair = F_FLOAT();
 
-          if(rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype)))
+          if (rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype)))
             fpair+=factor_lj*c.template compute_fpair<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
-          if(rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype)))
+          if (rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype)))
             fpair+=c.template compute_fcoul<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype,factor_coul,qtmp);
 
           fev_tmp.f[0] += delx*fpair;
@@ -526,11 +526,11 @@ struct PairComputeFunctor  {
           F_FLOAT evdwl = 0.0;
           F_FLOAT ecoul = 0.0;
           if (c.eflag) {
-            if(rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype))) {
+            if (rsq < (STACKPARAMS?c.m_cut_ljsq[itype][jtype]:c.d_cut_ljsq(itype,jtype))) {
               evdwl = factor_lj * c.template compute_evdwl<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype);
               fev_tmp.evdwl += 0.5*evdwl;
             }
-            if(rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype))) {
+            if (rsq < (STACKPARAMS?c.m_cut_coulsq[itype][jtype]:c.d_cut_coulsq(itype,jtype))) {
               ecoul = c.template compute_ecoul<STACKPARAMS,Specialisation>(rsq,i,j,itype,jtype,factor_coul,qtmp);
               fev_tmp.ecoul += 0.5*ecoul;
             }
@@ -722,7 +722,7 @@ int GetTeamSize(FunctorStyle& KOKKOS_GPU_ARG(functor), int KOKKOS_GPU_ARG(inum),
     else
       team_size_max = Kokkos::TeamPolicy<>(inum,Kokkos::AUTO).team_size_max(functor,Kokkos::ParallelForTag());
 
-    if(team_size*vector_length > team_size_max)
+    if (team_size*vector_length > team_size_max)
       team_size = team_size_max/vector_length;
 #else
     team_size = 1;
@@ -743,7 +743,7 @@ EV_FLOAT pair_compute_neighlist (PairStyle* fpair, typename std::enable_if<(NEIG
     int vector_length = 8;
     int atoms_per_team = 32;
 
-    if(fpair->atom->ntypes > MAX_TYPES_STACKPARAMS) {
+    if (fpair->atom->ntypes > MAX_TYPES_STACKPARAMS) {
       PairComputeFunctor<PairStyle,NEIGHFLAG,false,Specialisation > ff(fpair,list);
       atoms_per_team = GetTeamSize(ff, list->inum, (fpair->eflag || fpair->vflag), atoms_per_team, vector_length);
       Kokkos::TeamPolicy<Kokkos::IndexType<int> > policy(list->inum,atoms_per_team,vector_length);
@@ -757,7 +757,7 @@ EV_FLOAT pair_compute_neighlist (PairStyle* fpair, typename std::enable_if<(NEIG
       else                              Kokkos::parallel_for(policy,ff);
     }
   } else {
-    if(fpair->atom->ntypes > MAX_TYPES_STACKPARAMS) {
+    if (fpair->atom->ntypes > MAX_TYPES_STACKPARAMS) {
       PairComputeFunctor<PairStyle,NEIGHFLAG,false,Specialisation > ff(fpair,list);
       if (fpair->eflag || fpair->vflag) Kokkos::parallel_reduce(list->inum,ff,ev);
       else                              Kokkos::parallel_for(list->inum,ff);
diff --git a/src/KOKKOS/pair_snap_kokkos_impl.h b/src/KOKKOS/pair_snap_kokkos_impl.h
index da99cbd52a..cf4d5851b6 100644
--- a/src/KOKKOS/pair_snap_kokkos_impl.h
+++ b/src/KOKKOS/pair_snap_kokkos_impl.h
@@ -1323,7 +1323,7 @@ void PairSNAPKokkos<DeviceType, real_type, vector_length>::check_team_size_for(i
 
   team_size_max = Kokkos::TeamPolicy<DeviceType,TagStyle>(inum,Kokkos::AUTO).team_size_max(*this,Kokkos::ParallelForTag());
 
-  if(team_size*vector_length > team_size_max)
+  if (team_size*vector_length > team_size_max)
     team_size = team_size_max/vector_length;
 }
 
@@ -1334,7 +1334,7 @@ void PairSNAPKokkos<DeviceType, real_type, vector_length>::check_team_size_reduc
 
   team_size_max = Kokkos::TeamPolicy<DeviceType,TagStyle>(inum,Kokkos::AUTO).team_size_max(*this,Kokkos::ParallelReduceTag());
 
-  if(team_size*vector_length > team_size_max)
+  if (team_size*vector_length > team_size_max)
     team_size = team_size_max/vector_length;
 }
 
diff --git a/src/KOKKOS/sna_kokkos_impl.h b/src/KOKKOS/sna_kokkos_impl.h
index 23c1670bd8..d7aabb22c8 100644
--- a/src/KOKKOS/sna_kokkos_impl.h
+++ b/src/KOKKOS/sna_kokkos_impl.h
@@ -67,7 +67,7 @@ SNAKokkos<DeviceType, real_type, vector_length>::SNAKokkos(real_type rfac0_in,
     auto h_bzero = Kokkos::create_mirror_view(bzero);
 
     double www = wself*wself*wself;
-    for(int j = 0; j <= twojmax; j++)
+    for (int j = 0; j <= twojmax; j++)
       if (bnorm_flag)
         h_bzero[j] = www;
       else
@@ -95,9 +95,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   auto h_idxcg_block = Kokkos::create_mirror_view(idxcg_block);
 
   int idxcg_count = 0;
-  for(int j1 = 0; j1 <= twojmax; j1++)
-    for(int j2 = 0; j2 <= j1; j2++)
-      for(int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
+  for (int j1 = 0; j1 <= twojmax; j1++)
+    for (int j2 = 0; j2 <= j1; j2++)
+      for (int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
         h_idxcg_block(j1,j2,j) = idxcg_count;
         for (int m1 = 0; m1 <= j1; m1++)
           for (int m2 = 0; m2 <= j2; m2++)
@@ -114,10 +114,10 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
 
   int idxu_count = 0;
 
-  for(int j = 0; j <= twojmax; j++) {
+  for (int j = 0; j <= twojmax; j++) {
     h_idxu_block[j] = idxu_count;
-    for(int mb = 0; mb <= j; mb++)
-      for(int ma = 0; ma <= j; ma++)
+    for (int mb = 0; mb <= j; mb++)
+      for (int ma = 0; ma <= j; ma++)
         idxu_count++;
   }
   idxu_max = idxu_count;
@@ -128,10 +128,10 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   auto h_idxu_half_block = Kokkos::create_mirror_view(idxu_half_block);
 
   int idxu_half_count = 0;
-  for(int j = 0; j <= twojmax; j++) {
+  for (int j = 0; j <= twojmax; j++) {
     h_idxu_half_block[j] = idxu_half_count;
-    for(int mb = 0; 2*mb <= j; mb++)
-      for(int ma = 0; ma <= j; ma++)
+    for (int mb = 0; 2*mb <= j; mb++)
+      for (int ma = 0; ma <= j; ma++)
         idxu_half_count++;
   }
   idxu_half_max = idxu_half_count;
@@ -142,10 +142,10 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   auto h_idxu_full_half = Kokkos::create_mirror_view(idxu_full_half);
 
   idxu_count = 0;
-  for(int j = 0; j <= twojmax; j++) {
+  for (int j = 0; j <= twojmax; j++) {
     int jju_half = h_idxu_half_block[j];
-    for(int mb = 0; mb <= j; mb++) {
-      for(int ma = 0; ma <= j; ma++) {
+    for (int mb = 0; mb <= j; mb++) {
+      for (int ma = 0; ma <= j; ma++) {
         FullHalfMapper mapper;
         if (2*mb <= j) {
           mapper.idxu_half = jju_half + mb * (j + 1) + ma;
@@ -169,9 +169,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   auto h_idxu_cache_block = Kokkos::create_mirror_view(idxu_cache_block);
 
   int idxu_cache_count = 0;
-  for(int j = 0; j <= twojmax; j++) {
+  for (int j = 0; j <= twojmax; j++) {
     h_idxu_cache_block[j] = idxu_cache_count;
-    for(int mb = 0; mb < ((j+3)/2); mb++)
+    for (int mb = 0; mb < ((j+3)/2); mb++)
       for (int ma = 0; ma <= j; ma++)
         idxu_cache_count++;
   }
@@ -181,9 +181,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   // index list for beta and B
 
   int idxb_count = 0;
-  for(int j1 = 0; j1 <= twojmax; j1++)
-    for(int j2 = 0; j2 <= j1; j2++)
-      for(int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2)
+  for (int j1 = 0; j1 <= twojmax; j1++)
+    for (int j2 = 0; j2 <= j1; j2++)
+      for (int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2)
         if (j >= j1) idxb_count++;
 
   idxb_max = idxb_count;
@@ -191,9 +191,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   auto h_idxb = Kokkos::create_mirror_view(idxb);
 
   idxb_count = 0;
-  for(int j1 = 0; j1 <= twojmax; j1++)
-    for(int j2 = 0; j2 <= j1; j2++)
-      for(int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2)
+  for (int j1 = 0; j1 <= twojmax; j1++)
+    for (int j2 = 0; j2 <= j1; j2++)
+      for (int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2)
         if (j >= j1) {
           h_idxb(idxb_count,0) = j1;
           h_idxb(idxb_count,1) = j2;
@@ -208,9 +208,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   auto h_idxb_block = Kokkos::create_mirror_view(idxb_block);
 
   idxb_count = 0;
-  for(int j1 = 0; j1 <= twojmax; j1++)
-    for(int j2 = 0; j2 <= j1; j2++)
-      for(int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
+  for (int j1 = 0; j1 <= twojmax; j1++)
+    for (int j2 = 0; j2 <= j1; j2++)
+      for (int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
         if (j >= j1) {
           h_idxb_block(j1,j2,j) = idxb_count;
           idxb_count++;
@@ -222,9 +222,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
 
   int idxz_count = 0;
 
-  for(int j1 = 0; j1 <= twojmax; j1++)
-    for(int j2 = 0; j2 <= j1; j2++)
-      for(int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2)
+  for (int j1 = 0; j1 <= twojmax; j1++)
+    for (int j2 = 0; j2 <= j1; j2++)
+      for (int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2)
         for (int mb = 0; 2*mb <= j; mb++)
           for (int ma = 0; ma <= j; ma++)
             idxz_count++;
@@ -237,9 +237,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::build_indexlist()
   auto h_idxz_block = Kokkos::create_mirror_view(idxz_block);
 
   idxz_count = 0;
-  for(int j1 = 0; j1 <= twojmax; j1++)
-    for(int j2 = 0; j2 <= j1; j2++)
-      for(int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
+  for (int j1 = 0; j1 <= twojmax; j1++)
+    for (int j2 = 0; j2 <= j1; j2++)
+      for (int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
         h_idxz_block(j1,j2,j) = idxz_count;
 
         // find right beta(ii,jjb) entry
@@ -286,7 +286,7 @@ template<class DeviceType, typename real_type, int vector_length>
 inline
 void SNAKokkos<DeviceType, real_type, vector_length>::grow_rij(int newnatom, int newnmax)
 {
-  if(newnatom <= natom && newnmax <= nmax) return;
+  if (newnatom <= natom && newnmax <= nmax) return;
   natom = newnatom;
   nmax = newnmax;
 
@@ -644,7 +644,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_zi(const int& iato
       #ifdef LMP_KK_DEVICE_COMPILE
       #pragma unroll
       #endif
-      for(int ib = 0; ib < nb; ib++) {
+      for (int ib = 0; ib < nb; ib++) {
 
         int ma1 = ma1min;
         int ma2 = ma2max;
@@ -653,7 +653,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_zi(const int& iato
         #ifdef LMP_KK_DEVICE_COMPILE
         #pragma unroll
         #endif
-        for(int ia = 0; ia < na; ia++) {
+        for (int ia = 0; ia < na; ia++) {
           const auto utot1 = ulisttot_pack(iatom_mod, jju1+ma1, elem1, iatom_div);
           const auto utot2 = ulisttot_pack(iatom_mod, jju2+ma2, elem2, iatom_div);
           const auto cgcoeff_a = cgblock[icga];
@@ -717,8 +717,8 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_bi(const int& iato
         double sumzu = 0.0;
         double sumzu_temp = 0.0;
 
-        for(int mb = 0; 2*mb < j; mb++) {
-          for(int ma = 0; ma <= j; ma++) {
+        for (int mb = 0; 2*mb < j; mb++) {
+          for (int ma = 0; ma <= j; ma++) {
             const int jju_index = jju+mb*(j+1)+ma;
             const int jjz_index = jjz+mb*(j+1)+ma;
             if (2*mb == j) return; // I think we can remove this?
@@ -734,7 +734,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_bi(const int& iato
           sumzu_temp = 0.;
 
           const int mb = j/2;
-          for(int ma = 0; ma < mb; ma++) {
+          for (int ma = 0; ma < mb; ma++) {
             const int jju_index = jju+(mb-1)*(j+1)+(j+1)+ma;
             const int jjz_index = jjz+(mb-1)*(j+1)+(j+1)+ma;
 
@@ -1161,7 +1161,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_zi_cpu(const int&
       int jju1 = idxu_block[j1] + (j1+1)*mb1min;
       int jju2 = idxu_block[j2] + (j2+1)*mb2max;
       int icgb = mb1min*(j2+1) + mb2max;
-      for(int ib = 0; ib < nb; ib++) {
+      for (int ib = 0; ib < nb; ib++) {
 
         real_type suma1_r = static_cast<real_type>(0.0);
         real_type suma1_i = static_cast<real_type>(0.0);
@@ -1169,7 +1169,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_zi_cpu(const int&
         int ma1 = ma1min;
         int ma2 = ma2max;
         int icga = ma1min * (j2 + 1) + ma2max;
-        for(int ia = 0; ia < na; ia++) {
+        for (int ia = 0; ia < na; ia++) {
           suma1_r += cgblock[icga] * (ulisttot_full(jju1+ma1, elem1, iatom).re * ulisttot_full(jju2+ma2, elem2, iatom).re -
                                       ulisttot_full(jju1+ma1, elem1, iatom).im * ulisttot_full(jju2+ma2, elem2, iatom).im);
           suma1_i += cgblock[icga] * (ulisttot_full(jju1+ma1, elem1, iatom).re * ulisttot_full(jju2+ma2, elem2, iatom).im +
@@ -1222,7 +1222,6 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_bi_cpu(const typen
       for (int elem3 = 0; elem3 < nelements; elem3++) {
         Kokkos::parallel_for(Kokkos::TeamThreadRange(team,idxb_max),
           [&] (const int& jjb) {
-        //for(int jjb = 0; jjb < idxb_max; jjb++) {
           const int j1 = idxb(jjb, 0);
           const int j2 = idxb(jjb, 1);
           const int j = idxb(jjb, 2);
@@ -1234,8 +1233,6 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_bi_cpu(const typen
           const int bound = (j+2)/2;
           Kokkos::parallel_reduce(Kokkos::ThreadVectorRange(team,(j+1)*bound),
               [&] (const int mbma, real_type& sum) {
-              //for(int mb = 0; 2*mb < j; mb++)
-                //for(int ma = 0; ma <= j; ma++) {
               const int ma = mbma % (j + 1);
               const int mb = mbma / (j + 1);
               const int jju_index = jju + mb * (j + 1) + ma;
@@ -1253,7 +1250,6 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_bi_cpu(const typen
             const int mb = j/2;
             Kokkos::parallel_reduce(Kokkos::ThreadVectorRange(team, mb),
                 [&] (const int ma, real_type& sum) {
-            //for(int ma = 0; ma < mb; ma++) {
               const int jju_index = jju+(mb-1)*(j+1)+(j+1)+ma;
               const int jjz_index = jjz+(mb-1)*(j+1)+(j+1)+ma;
               sum +=
@@ -1449,14 +1445,13 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_deidrj_cpu(const t
   t_scalar3<real_type> final_sum;
   const int jelem = element(iatom, jnbor);
 
-  //for(int j = 0; j <= twojmax; j++) {
   Kokkos::parallel_reduce(Kokkos::ThreadVectorRange(team,twojmax+1),
       [&] (const int& j, t_scalar3<real_type>& sum_tmp) {
     int jju_half = idxu_half_block[j];
     int jju_cache = idxu_cache_block[j];
 
-    for(int mb = 0; 2*mb < j; mb++)
-      for(int ma = 0; ma <= j; ma++) {
+    for (int mb = 0; 2*mb < j; mb++)
+      for (int ma = 0; ma <= j; ma++) {
         sum_tmp.x += dulist(jju_cache,iatom,jnbor,0).re * ylist(jju_half,jelem,iatom).re +
                      dulist(jju_cache,iatom,jnbor,0).im * ylist(jju_half,jelem,iatom).im;
         sum_tmp.y += dulist(jju_cache,iatom,jnbor,1).re * ylist(jju_half,jelem,iatom).re +
@@ -1471,7 +1466,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::compute_deidrj_cpu(const t
     if (j%2 == 0) {
 
       int mb = j/2;
-      for(int ma = 0; ma < mb; ma++) {
+      for (int ma = 0; ma < mb; ma++) {
         sum_tmp.x += dulist(jju_cache,iatom,jnbor,0).re * ylist(jju_half,jelem,iatom).re +
                      dulist(jju_cache,iatom,jnbor,0).im * ylist(jju_half,jelem,iatom).im;
         sum_tmp.y += dulist(jju_cache,iatom,jnbor,1).re * ylist(jju_half,jelem,iatom).re +
@@ -2015,9 +2010,9 @@ void SNAKokkos<DeviceType, real_type, vector_length>::init_clebsch_gordan()
   int ifac;
 
   int idxcg_count = 0;
-  for(int j1 = 0; j1 <= twojmax; j1++)
-    for(int j2 = 0; j2 <= j1; j2++)
-      for(int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
+  for (int j1 = 0; j1 <= twojmax; j1++)
+    for (int j2 = 0; j2 <= j1; j2++)
+      for (int j = j1 - j2; j <= MIN(twojmax, j1 + j2); j += 2) {
         for (int m1 = 0; m1 <= j1; m1++) {
           aa2 = 2 * m1 - j1;
 
@@ -2028,7 +2023,7 @@ void SNAKokkos<DeviceType, real_type, vector_length>::init_clebsch_gordan()
             bb2 = 2 * m2 - j2;
             m = (aa2 + bb2 + j) / 2;
 
-            if(m < 0 || m > j) {
+            if (m < 0 || m > j) {
               h_cglist[idxcg_count] = 0.0;
               idxcg_count++;
               continue;
@@ -2120,8 +2115,8 @@ real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_sfac(real_typ
   constexpr real_type onehalf = static_cast<real_type>(0.5);
   if (switch_flag == 0) return one;
   if (switch_flag == 1) {
-    if(r <= rmin0) return one;
-    else if(r > rcut) return zero;
+    if (r <= rmin0) return one;
+    else if (r > rcut) return zero;
     else {
       auto rcutfac = static_cast<real_type>(MY_PI) / (rcut - rmin0);
       return onehalf * (cos((r - rmin0) * rcutfac) + one);
@@ -2140,8 +2135,8 @@ real_type SNAKokkos<DeviceType, real_type, vector_length>::compute_dsfac(real_ty
   constexpr real_type onehalf = static_cast<real_type>(0.5);
   if (switch_flag == 0) return zero;
   if (switch_flag == 1) {
-    if(r <= rmin0) return zero;
-    else if(r > rcut) return zero;
+    if (r <= rmin0) return zero;
+    else if (r > rcut) return zero;
     else {
       auto rcutfac = static_cast<real_type>(MY_PI) / (rcut - rmin0);
       return -onehalf * sin((r - rmin0) * rcutfac) * rcutfac;
diff --git a/src/MANYBODY/pair_eam_cd.h b/src/MANYBODY/pair_eam_cd.h
index 8f306d1019..4d938aad10 100644
--- a/src/MANYBODY/pair_eam_cd.h
+++ b/src/MANYBODY/pair_eam_cd.h
@@ -60,7 +60,7 @@ public:
   // Evaluates the h(x) polynomial for a given local concentration x.
   inline double evalH(double x) const {
     double v = 0.0;
-    for(int i = nhcoeff-1; i >= 1; i--) {
+    for (int i = nhcoeff-1; i >= 1; i--) {
       v = (v + hcoeff[i]) * x;
     }
     return v + hcoeff[0];
@@ -69,7 +69,7 @@ public:
   // Calculates the derivative of the h(x) polynomial.
   inline double evalHprime(double x) const {
     double v = 0.0;
-    for(int i = nhcoeff-1; i >= 2; i--) {
+    for (int i = nhcoeff-1; i >= 2; i--) {
       v = (v + (double)i * hcoeff[i]) * x;
     }
     return v + hcoeff[1];
diff --git a/src/USER-CGDNA/mf_oxdna.h b/src/USER-CGDNA/mf_oxdna.h
index e4ed1bbd03..00df53234d 100644
--- a/src/USER-CGDNA/mf_oxdna.h
+++ b/src/USER-CGDNA/mf_oxdna.h
@@ -96,13 +96,13 @@ inline double MFOxdna::F2(double r, double k, double cut_0, double cut_lc,
                           double b_lo, double b_hi, double cut_c)
 {
 
-  if(r < cut_lc || r > cut_hc){
+  if (r < cut_lc || r > cut_hc){
     return 0;
   }
-  else if(r < cut_lo){
+  else if (r < cut_lo){
     return k * b_lo * (cut_lc - r)*(cut_lc-r);
   }
-  else if(r < cut_hi){
+  else if (r < cut_hi){
     return k * 0.5 * ((r - cut_0)*(r-cut_0) - (cut_0 - cut_c)*(cut_0 - cut_c));
   }
   else{
@@ -118,13 +118,13 @@ inline double MFOxdna::DF2(double r, double k, double cut_0, double cut_lc,
                            double cut_hc, double cut_lo, double cut_hi,
                            double b_lo, double b_hi)
 {
-  if(r < cut_lc || r > cut_hc){
+  if (r < cut_lc || r > cut_hc){
     return 0;
   }
-  else if(r < cut_lo){
+  else if (r < cut_lo){
     return 2*k * b_lo * (r - cut_lc);
   }
-  else if(r < cut_hi){
+  else if (r < cut_hi){
     return k * (r - cut_0);
   }
   else{
@@ -171,7 +171,7 @@ inline double MFOxdna::F4(double theta, double a, double theta_0,
   else if (dtheta > dtheta_ast) {
     return b * (dtheta-dtheta_c)*(dtheta-dtheta_c);
   }
-  else if(dtheta > -dtheta_ast) {
+  else if (dtheta > -dtheta_ast) {
     return 1 - a * dtheta*dtheta;
   }
   else {
@@ -235,13 +235,13 @@ inline double MFOxdna::F5(double x, double a, double x_ast,
 inline double MFOxdna::DF5(double x, double a, double x_ast,
                            double b, double x_c)
 {
-  if(x >= 0) {
+  if (x >= 0) {
     return 0.0;
   }
   else if (x > x_ast) {
     return -2 * a * x;
   }
-  else if(x > x_c) {
+  else if (x > x_c) {
     return 2 * b * (x-x_c);
   }
   else {
diff --git a/src/USER-DPD/random_external_state.h b/src/USER-DPD/random_external_state.h
index de017ee8bd..e7d11a5926 100644
--- a/src/USER-DPD/random_external_state.h
+++ b/src/USER-DPD/random_external_state.h
@@ -141,9 +141,9 @@ namespace random_external_state {
   // RNGs with k calls to genNextParallelState()
   LAMMPS_INLINE
   void es_init(es_RNG_t &serial_state, uint64_t seed) {
-    if(seed==0) seed = uint64_t(1318319);
+    if (seed==0) seed = uint64_t(1318319);
     serial_state = seed;
-    for(int i = 0; i < 17; i++) es_rand(serial_state);
+    for (int i = 0; i < 17; i++) es_rand(serial_state);
   }
 
   // Call genNextParallelState() once for each RNG to generate
diff --git a/src/USER-MANIFOLD/manifold.h b/src/USER-MANIFOLD/manifold.h
index 9adf7055c7..7c3197e161 100644
--- a/src/USER-MANIFOLD/manifold.h
+++ b/src/USER-MANIFOLD/manifold.h
@@ -77,17 +77,17 @@ namespace user_manifold {
   // Some utility functions that are templated, so I implement them
   // here in the header.
   template< unsigned int size > inline
-  double infnorm( double *vect )
+  double infnorm(double *vect)
   {
     double largest = fabs( vect[0] );
-    for( unsigned int i = 1; i < size; ++i ){
+    for (unsigned int i = 1; i < size; ++i){
       double c = fabs( vect[i] );
       largest = ( c > largest ) ? c : largest;
     }
     return largest;
   }
 
-  inline double dot( double *a, double *b ){
+  inline double dot(double *a, double *b){
     return a[0]*b[0] + a[1]*b[1] + a[2]*b[2];
   }
 
diff --git a/src/USER-MANIFOLD/manifold_factory.h b/src/USER-MANIFOLD/manifold_factory.h
index 1eac8bf644..cc23874eba 100644
--- a/src/USER-MANIFOLD/manifold_factory.h
+++ b/src/USER-MANIFOLD/manifold_factory.h
@@ -93,8 +93,8 @@ class manifold;
   void make_manifold_if( manifold **man_ptr, const char *name,
                          LAMMPS *lmp, int narg, char **arg )
   {
-    if( strcmp( m_type::type(), name ) == 0 ){
-      if( *man_ptr == nullptr ){
+    if ( strcmp( m_type::type(), name ) == 0 ){
+      if ( *man_ptr == nullptr ){
         *man_ptr = new m_type(lmp, narg, arg);
       }
     }
diff --git a/src/USER-MISC/pair_ilp_graphene_hbn.h b/src/USER-MISC/pair_ilp_graphene_hbn.h
index 5ca8eb64a7..aadb792cad 100644
--- a/src/USER-MISC/pair_ilp_graphene_hbn.h
+++ b/src/USER-MISC/pair_ilp_graphene_hbn.h
@@ -82,7 +82,7 @@ class PairILPGrapheneHBN : public Pair {
     double Tap_coeff[8] = {1.0,0.0,0.0,0.0,-35.0,84.0,-70.0,20.0};
 
     r = r_ij/Rcut;
-    if(r >= 1.0) {Tap = 0.0;}
+    if (r >= 1.0) {Tap = 0.0;}
     else {
       Tap = Tap_coeff[7] * r + Tap_coeff[6];
       Tap = Tap * r  + Tap_coeff[5];
@@ -101,7 +101,7 @@ class PairILPGrapheneHBN : public Pair {
     double Tap_coeff[8] = {1.0,0.0,0.0,0.0,-35.0,84.0,-70.0,20.0};
 
     r = r_ij/Rcut;
-    if(r >= 1.0) {dTap = 0.0;}
+    if (r >= 1.0) {dTap = 0.0;}
     else {
       dTap = 7.0*Tap_coeff[7] * r + 6.0*Tap_coeff[6];
       dTap = dTap * r  + 5.0*Tap_coeff[5];
diff --git a/src/USER-MISC/pair_kolmogorov_crespi_full.h b/src/USER-MISC/pair_kolmogorov_crespi_full.h
index c579788110..f7d9237be4 100644
--- a/src/USER-MISC/pair_kolmogorov_crespi_full.h
+++ b/src/USER-MISC/pair_kolmogorov_crespi_full.h
@@ -84,7 +84,7 @@ class PairKolmogorovCrespiFull : public Pair {
     double Tap_coeff[8] = {1.0,0.0,0.0,0.0,-35.0,84.0,-70.0,20.0};
 
     r = r_ij/Rcut;
-    if(r >= 1.0) {Tap = 0.0;}
+    if (r >= 1.0) {Tap = 0.0;}
     else{
       Tap = Tap_coeff[7] * r + Tap_coeff[6];
       Tap = Tap * r  + Tap_coeff[5];
@@ -104,7 +104,7 @@ class PairKolmogorovCrespiFull : public Pair {
     double Tap_coeff[8] = {1.0,0.0,0.0,0.0,-35.0,84.0,-70.0,20.0};
 
     r = r_ij/Rcut;
-    if(r >= 1.0) {dTap = 0.0;}
+    if (r >= 1.0) {dTap = 0.0;}
     else {
       dTap = 7.0*Tap_coeff[7] * r + 6.0*Tap_coeff[6];
       dTap = dTap * r  + 5.0*Tap_coeff[5];
diff --git a/src/USER-MISC/pair_meam_spline.h b/src/USER-MISC/pair_meam_spline.h
index 31be826d0d..d1dc5064f5 100644
--- a/src/USER-MISC/pair_meam_spline.h
+++ b/src/USER-MISC/pair_meam_spline.h
@@ -118,10 +118,10 @@ protected:
     inline double eval(double x) const
     {
       x -= xmin;
-      if(x <= 0.0) {  // Left extrapolation.
+      if (x <= 0.0) {  // Left extrapolation.
         return Y[0] + deriv0 * x;
       }
-      else if(x >= xmax_shifted) {  // Right extrapolation.
+      else if (x >= xmax_shifted) {  // Right extrapolation.
         return Y[N-1] + derivN * (x - xmax_shifted);
       }
       else {
@@ -131,7 +131,7 @@ protected:
         int khi = N-1;
         while(khi - klo > 1) {
           int k = (khi + klo) / 2;
-          if(Xs[k] > x) khi = k;
+          if (Xs[k] > x) khi = k;
           else klo = k;
         }
         double h = Xs[khi] - Xs[klo];
@@ -156,11 +156,11 @@ protected:
     inline double eval(double x, double& deriv) const
     {
       x -= xmin;
-      if(x <= 0.0) {  // Left extrapolation.
+      if (x <= 0.0) {  // Left extrapolation.
         deriv = deriv0;
         return Y[0] + deriv0 * x;
       }
-      else if(x >= xmax_shifted) {  // Right extrapolation.
+      else if (x >= xmax_shifted) {  // Right extrapolation.
         deriv = derivN;
         return Y[N-1] + derivN * (x - xmax_shifted);
       }
@@ -171,7 +171,7 @@ protected:
         int khi = N-1;
         while(khi - klo > 1) {
           int k = (khi + klo) / 2;
-          if(Xs[k] > x) khi = k;
+          if (Xs[k] > x) khi = k;
           else klo = k;
         }
         double h = Xs[khi] - Xs[klo];
diff --git a/src/USER-MISC/pair_meam_sw_spline.h b/src/USER-MISC/pair_meam_sw_spline.h
index 1fd59b34a5..7b40d123ed 100644
--- a/src/USER-MISC/pair_meam_sw_spline.h
+++ b/src/USER-MISC/pair_meam_sw_spline.h
@@ -106,10 +106,10 @@ protected:
                 inline double eval(double x) const
                 {
                         x -= xmin;
-                        if(x <= 0.0) {  // Left extrapolation.
+                        if (x <= 0.0) {  // Left extrapolation.
                                 return Y[0] + deriv0 * x;
                         }
-                        else if(x >= xmax_shifted) {  // Right extrapolation.
+                        else if (x >= xmax_shifted) {  // Right extrapolation.
                                 return Y[N-1] + derivN * (x - xmax_shifted);
                         }
                         else {
@@ -119,7 +119,7 @@ protected:
                                 int khi = N-1;
                                 while(khi - klo > 1) {
                                         int k = (khi + klo) / 2;
-                                        if(Xs[k] > x) khi = k;
+                                        if (Xs[k] > x) khi = k;
                                         else klo = k;
                                 }
                                 double h = Xs[khi] - Xs[klo];
@@ -144,11 +144,11 @@ protected:
                 inline double eval(double x, double& deriv) const
                 {
                         x -= xmin;
-                        if(x <= 0.0) {  // Left extrapolation.
+                        if (x <= 0.0) {  // Left extrapolation.
                                 deriv = deriv0;
                                 return Y[0] + deriv0 * x;
                         }
-                        else if(x >= xmax_shifted) {  // Right extrapolation.
+                        else if (x >= xmax_shifted) {  // Right extrapolation.
                                 deriv = derivN;
                                 return Y[N-1] + derivN * (x - xmax_shifted);
                         }
@@ -159,7 +159,7 @@ protected:
                                 int khi = N-1;
                                 while(khi - klo > 1) {
                                         int k = (khi + klo) / 2;
-                                        if(Xs[k] > x) khi = k;
+                                        if (Xs[k] > x) khi = k;
                                         else klo = k;
                                 }
                                 double h = Xs[khi] - Xs[klo];
diff --git a/src/USER-PTM/ptm_voronoi_cell.h b/src/USER-PTM/ptm_voronoi_cell.h
index 69e3b5bdc4..138f290c42 100644
--- a/src/USER-PTM/ptm_voronoi_cell.h
+++ b/src/USER-PTM/ptm_voronoi_cell.h
@@ -328,12 +328,12 @@ class voronoicell_neighbor : public voronoicell_base {
                 inline void n_allocate(int i,int m) {mne[i]=new int[m*i];}
                 inline void n_add_memory_vertices(int i) {
                         int **pp=new int*[i];
-                        for(int j=0;j<current_vertices;j++) pp[j]=ne[j];
+                        for (int j=0;j<current_vertices;j++) pp[j]=ne[j];
                         delete [] ne;ne=pp;
                 }
                 inline void n_add_memory_vorder(int i) {
                         int **p2=new int*[i];
-                        for(int j=0;j<current_vertex_order;j++) p2[j]=mne[j];
+                        for (int j=0;j<current_vertex_order;j++) p2[j]=mne[j];
                         delete [] mne;mne=p2;
                 }
                 inline void n_set_pointer(int p,int n) {
@@ -346,7 +346,7 @@ class voronoicell_neighbor : public voronoicell_base {
                 inline void n_copy_aux1_shift(int a,int b) {paux1[b]=ne[a][b+1];}
                 inline void n_set_aux2_copy(int a,int b) {
                         paux2=mne[b]+b*mec[b];
-                        for(int i=0;i<b;i++) ne[a][i]=paux2[i];
+                        for (int i=0;i<b;i++) ne[a][i]=paux2[i];
                 }
                 inline void n_copy_pointer(int a,int b) {ne[a]=ne[b];}
                 inline void n_set_to_aux1(int j) {ne[j]=paux1;}
diff --git a/src/math_eigen_impl.h b/src/math_eigen_impl.h
index 7d6f0385ae..98027812dd 100644
--- a/src/math_eigen_impl.h
+++ b/src/math_eigen_impl.h
@@ -399,7 +399,7 @@ void Alloc2D(size_t nrows,          // size of the array (number of rows)
   //assert(paaX);
   *paaX = new Entry* [nrows];  //conventional 2D C array (pointer-to-pointer)
   (*paaX)[0] = new Entry [nrows * ncols];  // 1D C array (contiguous memor)
-  for(size_t iy=0; iy<nrows; iy++)
+  for (size_t iy=0; iy<nrows; iy++)
     (*paaX)[iy] = (*paaX)[0] + iy*ncols;
   // The caller can access the contents using (*paaX)[i][j]
 }
@@ -453,7 +453,7 @@ inline real_t<T> l2_norm(const std::vector<T>& vec) {
 template <typename T1, typename T2>
 inline void scalar_mul(T1 a, std::vector<T2>& vec) {
   int n = vec.size();
-  for(int i = 0;i < n;i++)
+  for (int i = 0;i < n;i++)
     vec[i] *= a;
 }
 
@@ -465,7 +465,7 @@ inline void normalize(std::vector<T>& vec) {
 template <typename T>
 inline real_t<T> l1_norm(const std::vector<T>& vec) {
   real_t<T> norm = real_t<T>(); // Zero initialization
-  for(const T& element : vec)
+  for (const T& element : vec)
     norm += std::abs(element);
   return norm;
 }
@@ -746,7 +746,7 @@ template<typename Scalar,typename Vector,typename Matrix,typename ConstMatrix>
 int Jacobi<Scalar, Vector, Matrix, ConstMatrix>::
 MaxEntryRow(Scalar const *const *M, int i) const {
   int j_max = i+1;
-  for(int j = i+2; j < n; j++)
+  for (int j = i+2; j < n; j++)
     if (std::abs(M[i][j]) > std::abs(M[i][j_max]))
       j_max = j;
   return j_max;
@@ -963,9 +963,9 @@ run(real_t<T>& eigvalue, std::vector<T>& eigvec) const
   pev = std::numeric_limits<real_t<T>>::max();
 
   int itern = this->max_iteration;
-  for(int k = 1;k <= this->max_iteration;k++) {
+  for (int k = 1;k <= this->max_iteration;k++) {
     // vk = (A + offset*E)uk, here E is the identity matrix
-    for(int i = 0;i < n;i++) {
+    for (int i = 0;i < n;i++) {
       vk[i] = uk[i]*this->eigenvalue_offset;
     }
     this->mv_mul(uk, vk);
@@ -984,7 +984,7 @@ run(real_t<T>& eigvalue, std::vector<T>& eigvec) const
 
     alpha.push_back(alphak);
 
-    for(int i = 0;i < n; i++) {
+    for (int i = 0;i < n; i++) {
       uk[i] = vk[i] - betak*u[k-1][i] - alphak*u[k][i];
     }
 
@@ -993,14 +993,14 @@ run(real_t<T>& eigvalue, std::vector<T>& eigvec) const
     betak = l2_norm(uk);
     beta.push_back(betak);
 
-    if(this->find_maximum) {
+    if (this->find_maximum) {
       ev = find_maximum_eigenvalue(alpha, beta);
     } else {
       ev = find_minimum_eigenvalue(alpha, beta);
     }
 
     const real_t<T> zero_threshold = minimum_effective_decimal<real_t<T>>()*1e-1;
-    if(betak < zero_threshold) {
+    if (betak < zero_threshold) {
       u.push_back(uk);
       // This element will never be accessed,
       // but this "push" guarantees u to always have one more element than
@@ -1012,7 +1012,7 @@ run(real_t<T>& eigvalue, std::vector<T>& eigvec) const
     normalize(uk);
     u.push_back(uk);
 
-    if(abs(ev-pev) < std::min(abs(ev), abs(pev))*this->eps) {
+    if (abs(ev-pev) < std::min(abs(ev), abs(pev))*this->eps) {
       itern = k;
       break;
     } else {
@@ -1030,18 +1030,18 @@ run(real_t<T>& eigvalue, std::vector<T>& eigvec) const
 
   beta[m-1] = 0.0;
 
-  if(eigvec.size() < n) {
+  if (eigvec.size() < n) {
     eigvec.resize(n);
   }
 
-  for(int i = 0;i < n;i++) {
+  for (int i = 0;i < n;i++) {
     eigvec[i] = cv[m-1]*u[m-1][i];
   }
 
-  for(int k = m-2;k >= 1;k--) {
+  for (int k = m-2;k >= 1;k--) {
     cv[k] = ((ev - alpha[k+1])*cv[k+1] - beta[k+1]*cv[k+2])/beta[k];
 
-    for(int i = 0;i < n;i++) {
+    for (int i = 0;i < n;i++) {
       eigvec[i] += cv[k]*u[k][i];
     }
   }
@@ -1062,9 +1062,9 @@ schmidt_orth(std::vector<T>& uorth, const std::vector<std::vector<T>>& u)
 
   int n = uorth.size();
 
-  for(int k = 0;k < u.size();k++) {
+  for (int k = 0;k < u.size();k++) {
     T innprod = inner_prod(uorth, u[k]);
-    for(int i = 0;i < n;i++)
+    for (int i = 0;i < n;i++)
       uorth[i] -= innprod * u[k][i];
   }
 }
@@ -1083,16 +1083,16 @@ find_minimum_eigenvalue(const std::vector<real_t<T>>& alpha,
   real_t<T> mid;
   int nmid; // Number of eigenvalues smaller than the "mid"
 
-  while(upper-lower > std::min(abs(lower), abs(upper))*eps) {
+  while (upper-lower > std::min(abs(lower), abs(upper))*eps) {
     mid = (lower+upper)/2.0;
     nmid = num_of_eigs_smaller_than(mid, alpha, beta);
-    if(nmid >= 1) {
+    if (nmid >= 1) {
       upper = mid;
     } else {
       lower = mid;
     }
 
-    if(mid == pmid) {
+    if (mid == pmid) {
       break; // This avoids an infinite loop due to zero matrix
     }
     pmid = mid;
@@ -1119,17 +1119,17 @@ find_maximum_eigenvalue(const std::vector<real_t<T>>& alpha,
                              // triangular matrix, which equals the rank of it
 
 
-  while(upper-lower > std::min(abs(lower), abs(upper))*eps) {
+  while (upper-lower > std::min(abs(lower), abs(upper))*eps) {
     mid = (lower+upper)/2.0;
     nmid = num_of_eigs_smaller_than(mid, alpha, beta);
 
-    if(nmid < m) {
+    if (nmid < m) {
       lower = mid;
     } else {
       upper = mid;
     }
 
-    if(mid == pmid) {
+    if (mid == pmid) {
       break; // This avoids an infinite loop due to zero matrix
     }
     pmid = mid;
@@ -1173,12 +1173,12 @@ num_of_eigs_smaller_than(real_t<T> c,
   int count = 0;
   int m = alpha.size();
 
-  for(int i = 1;i < m;i++){
+  for (int i = 1;i < m;i++){
     q_i = alpha[i] - c - beta[i-1]*beta[i-1]/q_i;
-    if(q_i < 0){
+    if (q_i < 0){
       count++;
     }
-    if(q_i == 0){
+    if (q_i == 0){
       q_i = minimum_effective_decimal<real_t<T>>();
     }
   }
@@ -1271,7 +1271,7 @@ init(std::vector<T>& v)
   std::uniform_real_distribution<T> rand((T)(-1.0), (T)(1.0));
 
   int n = v.size();
-  for(int i = 0;i < n;i++) {
+  for (int i = 0;i < n;i++) {
     v[i] = rand(mt);
   }
 
@@ -1288,7 +1288,7 @@ init(std::vector<std::complex<T>>& v)
   std::uniform_real_distribution<T> rand((T)(-1.0), (T)(1.0));
 
   int n = v.size();
-  for(int i = 0;i < n;i++) {
+  for (int i = 0;i < n;i++) {
     v[i] = std::complex<T>(rand(mt), rand(mt));
   }
 
@@ -1319,14 +1319,14 @@ PrincipalEigen(ConstMatrix matrix,
 {
   //assert(n > 0);
   auto matmul = [&](const std::vector<Scalar>& in, std::vector<Scalar>& out) {
-    for(int i = 0; i < n; i++) {
-      for(int j = 0; j < n; j++) {
+    for (int i = 0; i < n; i++) {
+      for (int j = 0; j < n; j++) {
         out[i] += matrix[i][j]*in[j];
       }
     }
   };
   auto init_vec = [&](std::vector<Scalar>& vec) {
-    for(int i = 0; i < n; i++)
+    for (int i = 0; i < n; i++)
       vec[i] = 0.0;
     vec[0] = 1.0;
   };