test verifies that do_coll_update=False leads to simple copy of last stage; can now switch between collocation update and copy be setting do_coll_update as sweeper parameter

Author: Daniel Ruprecht

pySDC/Sweeper.py

@@ -45,7 +45,7 @@ def __init__(self,params):
         coll = params['collocation_class'](params['num_nodes'],0,1)
         assert isinstance(coll, CollBase)
         if not coll.right_is_node:
-          assert params['do_coll_update'], "For nodes where the right end point is not a node, do_coll_update has to be set to True"
+          assert self.params.do_coll_update, "For nodes where the right end point is not a node, do_coll_update has to be set to True"
 
         # This will be set as soon as the sweeper is instantiated at the level
         self.__level = None

pySDC/sweeper_classes/generic_LU.py

@@ -135,10 +135,10 @@ def compute_end_point(self):
         L = self.level
         P = L.prob
 
--        # check if Mth node is equal to right point and do_coll_update is false, perform a simple copy
--        if (self.coll.right_is_node and not self.params['do_coll_update']):
--            # a copy is sufficient
--            L.uend = P.dtype_u(L.u[-1])
+         # check if Mth node is equal to right point and do_coll_update is false, perform a simple copy
+        if (self.coll.right_is_node and not self.params.do_coll_update):
+          # a copy is sufficient
+          L.uend = P.dtype_u(L.u[-1])
         else:
           # start with u0 and add integral over the full interval (using coll.weights)
           L.uend = P.dtype_u(L.u[0])

pySDC/sweeper_classes/imex_1st_order.py

@@ -134,19 +134,24 @@ def compute_end_point(self):
         """
         Compute u at the right point of the interval
 
-        The value uend computed here is a full evaluation of the Picard formulation (always!)
+        The value uend computed here is a full evaluation of the Picard formulation unless do_full_update==False
         """
 
         # get current level and problem description
         L = self.level
         P = L.prob
 
-        # start with u0 and add integral over the full interval (using coll.weights)
-        L.uend = P.dtype_u(L.u[0])
-        for m in range(self.coll.num_nodes):
-            L.uend += L.dt*self.coll.weights[m]*(L.f[m+1].impl + L.f[m+1].expl)
-        # add up tau correction of the full interval (last entry)
-        if L.tau is not None:
-            L.uend += L.tau[-1]
+        # check if Mth node is equal to right point and do_coll_update is false, perform a simple copy
+        if (self.coll.right_is_node and not self.params.do_coll_update):
+           # a copy is sufficient
+           L.uend = P.dtype_u(L.u[-1])
+        else:
+          # start with u0 and add integral over the full interval (using coll.weights)
+          L.uend = P.dtype_u(L.u[0])
+          for m in range(self.coll.num_nodes):
+              L.uend += L.dt*self.coll.weights[m]*(L.f[m+1].impl + L.f[m+1].expl)
+          # add up tau correction of the full interval (last entry)
+          if L.tau is not None:
+              L.uend += L.tau[-1]
 
         return None

tests/test_imexsweeper.py

@@ -239,3 +239,41 @@ def test_norightnode_collupdate_fails(self):
     # Has to throw an exception
     with self.assertRaises(AssertionError):
       step, level, problem, nnodes = self.setupLevelStepProblem()   
+
+  #
+  # Make sure the update with do_coll_update=False reproduces last stage
+  # 
+  def test_update_nocollupdate_laststage(self):
+    self.swparams['do_coll_update'] = False
+    step, level, problem, nnodes = self.setupLevelStepProblem()
+    level.sweep.predict()
+    ulaststage = np.random.rand()
+    level.u[nnodes].values = ulaststage
+    level.sweep.compute_end_point()
+    uend = level.uend.values
+    assert abs(uend-ulaststage)<1e-14, "compute_end_point with do_coll_update=False did not reproduce last stage value"
+
+  #
+  # Make sure that update with do_coll_update=False is identical to update formula with q=(0,...,0,1)
+  #
+  def test_updateformula_no_coll_update(self):
+    self.swparams['do_coll_update'] = False
+    step, level, problem, nnodes = self.setupLevelStepProblem()
+    level.sweep.predict()
+    u0full = np.array([ level.u[l].values.flatten() for l in range(1,nnodes+1) ])
+
+    # Perform update step in sweeper
+    level.sweep.update_nodes()
+    ustages = np.array([ level.u[l].values.flatten() for l in range(1,nnodes+1) ])
+
+    # Compute end value through provided function
+    level.sweep.compute_end_point()
+    uend_sweep = level.uend.values
+    # Compute end value from matrix formulation
+    q = np.zeros(nnodes)
+    q[nnodes-1] = 1.0
+    print q
+    uend_mat   = q.dot(ustages)
+    print uend_sweep
+    print uend_mat
+    assert np.linalg.norm(uend_sweep - uend_mat, np.infty)<1e-14, "For do_coll_update=False, update formula in sweeper gives different result than matrix update formula with q=(0,..,0,1)"