Index: param.ccl
===================================================================
--- param.ccl	(revision 33)
+++ param.ccl	(working copy)
@@ -32,6 +32,20 @@
 {
 } "no"
 
+CCTK_REAL refluxing_fraction "Fraction of refluxing that should be applied" STEERABLE=always
+{
+  *:* :: "0 = no refluxing, 1 = full refluxing"
+} 1.0
+
+BOOLEAN halftime_refluxing "Apply some refluxing at half time already (THIS DOES NOT WORK, BECAUSE POSTRESTRICT IS NOT CALLED ON THE FINEST LEVEL)" STEERABLE=always
+{
+} "no"
+
+CCTK_REAL halftime_refluxing_fraction "Fraction of refluxing that should be applied at half time" STEERABLE=always
+{
+  *:* :: "0 = no refluxing, 1 = full refluxing, 2 = over-refluxing"
+} 1.0
+
 STRING debug_variables "List of intermediate variables to output for debugging" STEERABLE=always
 {
   "" :: ""
Index: src/correct.cc
===================================================================
--- src/correct.cc	(revision 33)
+++ src/correct.cc	(working copy)
@@ -418,9 +418,79 @@
 
 
 
+// Reset the fine grid flux register (on the next finer level)
+static
+void flux_register_fine_reset (cGH const * restrict const cctkGH)
+{
+  assert (is_level_mode());
+  
+  SWITCH_TO_LEVEL (cctkGH, reflevel+1) {
+    BEGIN_LOCAL_MAP_LOOP (cctkGH, CCTK_GF) {
+      BEGIN_LOCAL_COMPONENT_LOOP (cctkGH, CCTK_GF) {
+        DECLARE_CCTK_ARGUMENTS;
+        
+        ivect const& lsh = ivect::ref(cctk_lsh);
+        int const np = prod(lsh);
+        
+        vector<CCTK_REAL *> flux_register_fine_ptrs =
+          get_varptrs (cctkGH, reflevel, variables::flux_register_fine);
+        int const nvars = flux_register_fine_ptrs.size();
+        
+        for (int n=0; n<nvars; ++n) {
+          for (int dir=0; dir<3; ++dir) {
+            assert (dim == 3);
+#pragma omp parallel
+            CCTK_LOOP3_ALL(GRHydro_Reflux_fine_reset, cctkGH, i,j,k) {
+              int const ind = CCTK_GFINDEX3D (cctkGH, i, j, k);
+              flux_register_fine_ptrs.AT(n)[ind+dir*np] = 0.0;
+            } CCTK_ENDLOOP3_ALL(GRHydro_Reflux_fine_reset);
+          }
+        }
+        
+      } END_LOCAL_COMPONENT_LOOP;
+    } END_LOCAL_MAP_LOOP;
+  } END_SWITCH_TO_LEVEL;
+}  
+
+
+
+// Reset the coarse grid flux register (on this level)
+static
+void flux_register_coarse_reset (cGH const * restrict const cctkGH)
+{
+  assert (is_level_mode());
+  
+  BEGIN_LOCAL_MAP_LOOP (cctkGH, CCTK_GF) {
+    BEGIN_LOCAL_COMPONENT_LOOP (cctkGH, CCTK_GF) {
+      DECLARE_CCTK_ARGUMENTS;
+      
+      ivect const& lsh = ivect::ref(cctk_lsh);
+      int const np = prod(lsh);
+      
+      vector<CCTK_REAL *> flux_register_coarse_ptrs =
+        get_varptrs (cctkGH, reflevel, variables::flux_register_coarse);
+      int const nvars = flux_register_coarse_ptrs.size();
+      
+      for (int n=0; n<nvars; ++n) {
+        for (int dir=0; dir<3; ++dir) {
+          assert (dim == 3);
+#pragma omp parallel
+          CCTK_LOOP3_ALL(GRHydro_Reflux_coarse_reset, cctkGH, i,j,k) {
+            int const ind = CCTK_GFINDEX3D (cctkGH, i, j, k);
+            flux_register_coarse_ptrs.AT(n)[ind+dir*np] = 0.0;
+          } CCTK_ENDLOOP3_ALL(GRHydro_Reflux_coarse_reset);
+        }
+      }
+      
+    } END_LOCAL_COMPONENT_LOOP;
+  } END_LOCAL_MAP_LOOP;
+}
+
+
+
 // Reflux
 static
-void reflux (cGH const * restrict const cctkGH)
+void reflux (cGH const * restrict const cctkGH, bool const halftime = false)
 {
   DECLARE_CCTK_PARAMETERS;
   
@@ -540,6 +610,9 @@
   
   // Scale the restricted fine grid register (in the correction), and
   // subtract the coarse grid register; then update the state
+  CCTK_REAL const halftime_factor = halftime ? 0.5 : 1.0;
+  CCTK_REAL const the_refluxing_fraction =
+    halftime ? halftime_refluxing_fraction : refluxing_fraction;
   BEGIN_LOCAL_MAP_LOOP (cctkGH, CCTK_GF) {
     BEGIN_LOCAL_COMPONENT_LOOP (cctkGH, CCTK_GF) {
       DECLARE_CCTK_ARGUMENTS;
@@ -573,7 +646,7 @@
           // lower) if on the lower face, or the cell to the right
           // (same index) if on the upper face
           int const ioff = index (lsh, (face ? 0 : -1) * idir);
-          CCTK_REAL const factor = face ? +1 : -1;
+          CCTK_REAL const face_factor = face ? +1 : -1;
           ibset const coarse_boundary =
             local_box.coarse_boundary[dir][face].shift(idir, 2);
           LOOP_OVER_BSET(cctkGH, coarse_boundary, box, imin, imax) {
@@ -593,23 +666,24 @@
                          cctk_lsh[0],cctk_lsh[1],cctk_lsh[2])
               {
                 int const ind = CCTK_GFINDEX3D (cctkGH, i, j, k);
+                // Update the state
                 flux_correction_ptrs.AT(n)[ind+dir*np] -=
+                  halftime_factor *
                   flux_register_coarse_ptrs.AT(n)[ind+dir*np];
-#warning "TODO"
-                assert (not isnan(flux_correction_ptrs.AT(n)[ind+dir*np]));
-                // Update the state
                 CCTK_REAL const difference =
-                  factor * flux_correction_ptrs.AT(n)[ind+dir*np] / delta[dir];
+                  the_refluxing_fraction *
+                  face_factor *
+                  flux_correction_ptrs.AT(n)[ind+dir*np] / delta[dir];
                 if (refluxing_debug_variables) {
                   // Keep a total of the refluxing changes
                   correction_total_ptrs.AT(n)[ind+dir*np+ioff] += difference;
-#warning "TODO"
-                  assert (not isnan(correction_total_ptrs.AT(n)[ind+dir*np+ioff]));
                 }
                 if (not suppress_refluxing) {
                   var_ptrs.AT(n)[ind+ioff] += difference;
-#warning "TODO"
-                  assert (not isnan(var_ptrs.AT(n)[ind+ioff]));
+                  // Correct the coarse grid flux register (this is
+                  // irrevelant except for halftime refluxing)
+                  flux_register_coarse_ptrs.AT(n)[ind+dir*np] +=
+                    face_factor * difference * delta[dir];
                 }
               } CCTK_ENDLOOP3(GRHydro_Reflux_correction_calculate);
             }
@@ -620,76 +694,50 @@
       
     } END_LOCAL_COMPONENT_LOOP;
   } END_LOCAL_MAP_LOOP;
+  
+  if (not halftime) {
+    flux_register_fine_reset (cctkGH);
+    flux_register_coarse_reset (cctkGH);
+  }
 }
 
 
 
-// Reset the fine grid flux register (on the next finer level)
 static
-void flux_register_fine_reset (cGH const * restrict const cctkGH)
+void HalfTime_CorrectState (CCTK_ARGUMENTS)
 {
-  assert (is_level_mode());
+  DECLARE_CCTK_PARAMETERS;
   
-  SWITCH_TO_LEVEL (cctkGH, reflevel+1) {
-    BEGIN_LOCAL_MAP_LOOP (cctkGH, CCTK_GF) {
-      BEGIN_LOCAL_COMPONENT_LOOP (cctkGH, CCTK_GF) {
-        DECLARE_CCTK_ARGUMENTS;
-        
-        ivect const& lsh = ivect::ref(cctk_lsh);
-        int const np = prod(lsh);
-        
-        vector<CCTK_REAL *> flux_register_fine_ptrs =
-          get_varptrs (cctkGH, reflevel, variables::flux_register_fine);
-        int const nvars = flux_register_fine_ptrs.size();
-        
-        for (int n=0; n<nvars; ++n) {
-          for (int dir=0; dir<3; ++dir) {
-            assert (dim == 3);
-#pragma omp parallel
-            CCTK_LOOP3_ALL(GRHydro_Reflux_fine_reset, cctkGH, i,j,k) {
-              int const ind = CCTK_GFINDEX3D (cctkGH, i, j, k);
-              flux_register_fine_ptrs.AT(n)[ind+dir*np] = 0.0;
-            } CCTK_ENDLOOP3_ALL(GRHydro_Reflux_fine_reset);
-          }
-        }
-        
-      } END_LOCAL_COMPONENT_LOOP;
-    } END_LOCAL_MAP_LOOP;
+  // Sanity check
+  assert (mglevel>=0 and reflevel>=0 and Carpet::map==-1 and component==-1);
+  
+  // We reflux on level L by taking a correction from level L+1 into
+  // account. (This corresponds to restriction, where level L is
+  // "corrected" from level L+1.)
+  
+  // In this case, we apply (part of) the refluxing at half-time
+  // already. At half-time we are called on the fine level L+1, so we
+  // have to switch to the corresponding coarse level L first.
+  
+  // If there is no coarser level; do nothing
+  if (reflevel == 0) return;
+  
+  // If the next coarser level is aligned; do nothing
+  assert (maxtimereflevelfact % timereffacts.AT(reflevel-1) == 0);
+  int const iter_step = maxtimereflevelfact / timereffacts.AT(reflevel-1);
+  if (cctkGH->cctk_iteration % iter_step == 0) return;
+  
+  SWITCH_TO_LEVEL(cctkGH, reflevel-1) {
+    
+    if (verbose or veryverbose) {
+      CCTK_VInfo (CCTK_THORNSTRING,
+                  "Half-time refluxing at iteration %d on level %d of %d",
+                  cctkGH->cctk_iteration, reflevel, reflevels);
+    }
+    
+    reflux (cctkGH, true);
+    
   } END_SWITCH_TO_LEVEL;
-}  
-
-
-
-// Reset the coarse grid flux register (on this level)
-static
-void flux_register_coarse_reset (cGH const * restrict const cctkGH)
-{
-  assert (is_level_mode());
-  
-  BEGIN_LOCAL_MAP_LOOP (cctkGH, CCTK_GF) {
-    BEGIN_LOCAL_COMPONENT_LOOP (cctkGH, CCTK_GF) {
-      DECLARE_CCTK_ARGUMENTS;
-      
-      ivect const& lsh = ivect::ref(cctk_lsh);
-      int const np = prod(lsh);
-      
-      vector<CCTK_REAL *> flux_register_coarse_ptrs =
-        get_varptrs (cctkGH, reflevel, variables::flux_register_coarse);
-      int const nvars = flux_register_coarse_ptrs.size();
-      
-      for (int n=0; n<nvars; ++n) {
-        for (int dir=0; dir<3; ++dir) {
-          assert (dim == 3);
-#pragma omp parallel
-          CCTK_LOOP3_ALL(GRHydro_Reflux_coarse_reset, cctkGH, i,j,k) {
-            int const ind = CCTK_GFINDEX3D (cctkGH, i, j, k);
-            flux_register_coarse_ptrs.AT(n)[ind+dir*np] = 0.0;
-          } CCTK_ENDLOOP3_ALL(GRHydro_Reflux_coarse_reset);
-        }
-      }
-      
-    } END_LOCAL_COMPONENT_LOOP;
-  } END_LOCAL_MAP_LOOP;
 }
 
 
@@ -702,6 +750,10 @@
   // Sanity check
   assert (mglevel>=0 and reflevel>=0 and Carpet::map==-1 and component==-1);
   
+  if (halftime_refluxing) {
+    HalfTime_CorrectState (CCTK_PASS_CTOC);
+  }
+  
   // We reflux on level L by taking a correction from level L+1 into
   // account. (This corresponds to restriction, where level L is
   // "corrected" from level L+1.)
@@ -733,9 +785,6 @@
   }
   
   reflux (cctkGH);
-  
-  flux_register_fine_reset (cctkGH);
-  flux_register_coarse_reset (cctkGH);
 }