diff --git a/kharma/b_ct/b_ct_boundaries.cpp b/kharma/b_ct/b_ct_boundaries.cpp
index 355605b6..242f657c 100644
--- a/kharma/b_ct/b_ct_boundaries.cpp
+++ b/kharma/b_ct/b_ct_boundaries.cpp
@@ -1,25 +1,25 @@
-/* 
+/*
  *  File: b_ct_boundaries.cpp
- *  
+ *
  *  BSD 3-Clause License
- *  
+ *
  *  Copyright (c) 2020, AFD Group at UIUC
  *  All rights reserved.
- *  
+ *
  *  Redistribution and use in source and binary forms, with or without
  *  modification, are permitted provided that the following conditions are met:
- *  
+ *
  *  1. Redistributions of source code must retain the above copyright notice, this
  *     list of conditions and the following disclaimer.
- *  
+ *
  *  2. Redistributions in binary form must reproduce the above copyright notice,
  *     this list of conditions and the following disclaimer in the documentation
  *     and/or other materials provided with the distribution.
- *  
+ *
  *  3. Neither the name of the copyright holder nor the names of its
  *     contributors may be used to endorse or promote products derived from
  *     this software without specific prior written permission.
- *  
+ *
  *  THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
  *  AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  *  IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
@@ -121,7 +121,7 @@ void B_CT::AverageBoundaryEMF(MeshBlockData<Real> *rc, IndexDomain domain, const
                     inner_dir = X1DIR;
                 }
             }
-            parthenon::par_for_outer(DEFAULT_OUTER_LOOP_PATTERN, "reduce_EMF1_" + bname, pmb->exec_space,
+            parthenon::par_for_outer(DEFAULT_OUTER_LOOP_PATTERN, "reduce_EMF_" + bname, pmb->exec_space,
                 0, 1, outer.s, outer.e,
                 KOKKOS_LAMBDA(parthenon::team_mbr_t member, const int& o) {
                     double emf_sum = 0.;
@@ -157,6 +157,7 @@ void B_CT::AverageBoundaryEMF(MeshBlockData<Real> *rc, IndexDomain domain, const
                             }
                         , sum_reducer);
                     }
+                    member.team_barrier();
 
                     // Calculate the average
                     const double emf_av = emf_sum / len;
diff --git a/kharma/grmhd/grmhd.cpp b/kharma/grmhd/grmhd.cpp
index 5d119ece..831070ec 100644
--- a/kharma/grmhd/grmhd.cpp
+++ b/kharma/grmhd/grmhd.cpp
@@ -519,6 +519,7 @@ void CancelBoundaryU3(MeshBlockData<Real> *rc, IndexDomain domain, bool coarse)
                     }
                 );
             }
+            member.team_barrier();
 
             // Sum the first rank of U3
             Real U3_sum = 0.;
@@ -528,6 +529,7 @@ void CancelBoundaryU3(MeshBlockData<Real> *rc, IndexDomain domain, bool coarse)
                     local_result += isnan(P(m_p.U3, k, jf, i)) ? 0. : P(m_p.U3, k, jf, i);
                 }
             , sum_reducer);
+            member.team_barrier();
 
             // Subtract the average, floor, restore conserved vars, update ctop
             const Real U3_avg = U3_sum / (bi.ke - bi.ks + 1);
@@ -591,6 +593,7 @@ void CancelBoundaryT3(MeshBlockData<Real> *rc, IndexDomain domain, bool coarse)
                     }
                 );
             }
+            member.team_barrier();
 
             // Sum the first rank of the angular momentum T3
             Real T3_sum = 0.;
@@ -600,6 +603,7 @@ void CancelBoundaryT3(MeshBlockData<Real> *rc, IndexDomain domain, bool coarse)
                     local_result += isnan(U(m_u.U3, k, jf, i)) ? 0. : U(m_u.U3, k, jf, i);
                 }
             , sum_reducer);
+            member.team_barrier();
 
             // Calculate the average and subtract it
             const Real T3_avg = T3_sum / (bi.ke - bi.ks + 1);
diff --git a/pars/bondi/bondi_b_vertical.par b/pars/bondi/bondi_b_vertical.par
index 31380300..6e165083 100644
--- a/pars/bondi/bondi_b_vertical.par
+++ b/pars/bondi/bondi_b_vertical.par
@@ -55,7 +55,9 @@ on = true
 # Transmitting boundaries, since there will be material at the pole
 inner_x2 = transmitting
 outer_x2 = transmitting
-excise_polar_flux = true
+# This might additionally help, but the option is unstable
+#excise_polar_flux = true
+
 # Allow in material from boundary
 check_inflow_outer_x1 = false
 
diff --git a/pars/tori_3d/mad.par b/pars/tori_3d/mad.par
index 47579314..809217c9 100644
--- a/pars/tori_3d/mad.par
+++ b/pars/tori_3d/mad.par
@@ -66,8 +66,9 @@ two_sync = false
 <boundaries>
 inner_x2 = transmitting
 outer_x2 = transmitting
-# This reduces the polar "wake," but may shrink timestep
-excise_polar_flux = true
+# This reduces the polar "wake," but shrinks the timestep,
+# and has caused instabilities. YMMV.
+#excise_polar_flux = true
 
 <floors>
 frame        = drift
diff --git a/pars/tori_3d/sane.par b/pars/tori_3d/sane.par
index b1f05cbb..eccdf7a4 100644
--- a/pars/tori_3d/sane.par
+++ b/pars/tori_3d/sane.par
@@ -64,8 +64,9 @@ two_sync = false
 <boundaries>
 inner_x2 = transmitting
 outer_x2 = transmitting
-# This reduces the polar "wake," but may shrink timestep
-excise_polar_flux = true
+# This reduces the polar "wake," but shrinks the timestep,
+# and has caused instabilities. YMMV.
+#excise_polar_flux = true
 
 <floors>
 # Now with FOFC it might work to set "frame = normal",