From 69115c472467c3c5762e0b4df9d829296d5f611f Mon Sep 17 00:00:00 2001
From: Rob Falgout <rfalgout@llnl.gov>
Date: Fri, 29 Oct 2021 17:48:25 -0700
Subject: [PATCH] First attempt at a linear MGRIT implementation (doesn't work
 yet)

---
 braid/braid.h     |  2 ++
 braid/hierarchy.c |  5 +++++
 braid/restrict.c  | 24 ++++++++++++++++++++++++
 3 files changed, 31 insertions(+)

diff --git a/braid/braid.h b/braid/braid.h
index 11e87239..4c1beeae 100644
--- a/braid/braid.h
+++ b/braid/braid.h
@@ -29,6 +29,8 @@
 #ifndef braid_HEADER
 #define braid_HEADER
 
+#define braid_LINEAR 1
+
 /* To enable XBraid without MPI, re-compile with
  *   make sequential=yes 
  * This will compile using the mpistubs.h */
diff --git a/braid/hierarchy.c b/braid/hierarchy.c
index b916f393..d09d38ea 100644
--- a/braid/hierarchy.c
+++ b/braid/hierarchy.c
@@ -504,6 +504,11 @@ _braid_CopyFineToCoarse(braid_Core  core)
          _braid_MapCoarseToFine(index, f_cfactor, f_index);
          _braid_UGetVector(core, level-1, f_index, &u);
          _braid_Coarsen(core, level, f_index, index, u, &va[index-ilower]);
+         if (braid_LINEAR)
+         {
+            /* Set initial guess to zero on coarse grids */
+            _braid_BaseSum(core, app, -1.0, va[index-ilower], 1.0, va[index-ilower]);
+         }
          
          _braid_BaseFree(core, app,  u);
          _braid_BaseClone(core, app,  va[index-ilower], &u);
diff --git a/braid/restrict.c b/braid/restrict.c
index d3a4da40..53a9b8ed 100644
--- a/braid/restrict.c
+++ b/braid/restrict.c
@@ -199,12 +199,22 @@ _braid_FRestrict(braid_Core   core,
          _braid_MapFineToCoarse(ci, cfactor, c_index);
          _braid_Coarsen(core, c_level, ci, c_index, u, &c_va[c_index-c_ilower]);
          _braid_Coarsen(core, c_level, ci, c_index, r, &c_fa[c_index-c_ilower]);
+         if (braid_LINEAR)
+         {
+            /* Set initial guess to zero on coarse grids */
+            _braid_BaseSum(core, app, -1.0, c_va[c_index-c_ilower], 1.0, c_va[c_index-c_ilower]);
+         }
       }
       else if (ci == 0)
       {
          /* Restrict initial condition, coarsening in space if needed */
          _braid_UGetVectorRef(core, level, 0, &u);
          _braid_Coarsen(core, c_level, 0, 0, u, &c_va[0]);
+         if (braid_LINEAR)
+         {
+            /* Set initial guess to zero on coarse grids */
+            _braid_BaseSum(core, app, -1.0, c_va[0], 1.0, c_va[0]);
+         }
       }
 
       if ((flo <= fhi) || (ci > _braid_CoreElt(core, initiali)))
@@ -219,6 +229,12 @@ _braid_FRestrict(braid_Core   core,
    /* Set initial guess on coarse level */
    _braid_InitGuess(core, c_level);
 
+   if (braid_LINEAR)
+   {
+      /* Don't need to update the rhs in the linear case (no Richardson either) */
+   }
+   else
+   {
    /* Initialize update of c_va[-1] boundary */
    if (c_ilower <= c_iupper)
    {
@@ -302,6 +318,7 @@ _braid_FRestrict(braid_Core   core,
       }
    }
    _braid_CommWait(core, &send_handle);
+   }
    
    /* Compute global rnorm (only on level 0) */
    if (level == 0)
@@ -332,12 +349,19 @@ _braid_FRestrict(braid_Core   core,
       _braid_PrintSpatialNorms(core, tnorm_a, ncpoints);
    }
 
+   if (braid_LINEAR)
+   {
+      /* Don't need to update the rhs in the linear case (no Richardson either) */
+   }
+   else
+   {
    /* Need to finalize the error estimates at the F-points */ 
    if ( level == 0 && est_error )
    {
       _braid_FinalizeErrorEstimates( core, estimate , c_iupper-c_ilower + 1 );
       _braid_TFree(estimate);
    } 
+   }
    
    return _braid_error_flag;
 }