From a892f5f54e09b6b508fa757fd714cf88b8c5b369 Mon Sep 17 00:00:00 2001 From: Daniel Ruprecht Date: Tue, 22 Dec 2015 09:00:57 +0000 Subject: [PATCH 1/4] introduced flag do_coll_update as parameter for sweeper objects. If True (default), a proper collocation update is performed in compute_end_point. If False, the stage corresponding to the end point is copied. Throws an exception of set to False for collocation nodes where endpoint is not a stage --- pySDC/Sweeper.py | 7 +++++-- pySDC/sweeper_classes/generic_LU.py | 23 ++++++++++++++--------- 2 files changed, 19 insertions(+), 11 deletions(-) diff --git a/pySDC/Sweeper.py b/pySDC/Sweeper.py index fdba6341f1..a19e7398af 100644 --- a/pySDC/Sweeper.py +++ b/pySDC/Sweeper.py @@ -31,7 +31,8 @@ def __init__(self,params): defaults = dict() defaults['do_LU'] = False - + defaults['do_coll_update'] = True + for k,v in defaults.items(): setattr(self,k,v) @@ -43,6 +44,8 @@ def __init__(self,params): coll = params['collocation_class'](params['num_nodes'],0,1) assert isinstance(coll, CollBase) + if not coll.right_is_node: + 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 @@ -151,4 +154,4 @@ def update_nodes(self): """ Abstract interface to node update """ - return None \ No newline at end of file + return None diff --git a/pySDC/sweeper_classes/generic_LU.py b/pySDC/sweeper_classes/generic_LU.py index 00697a7223..77d7b5542c 100644 --- a/pySDC/sweeper_classes/generic_LU.py +++ b/pySDC/sweeper_classes/generic_LU.py @@ -135,12 +135,17 @@ def compute_end_point(self): 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] - # add up tau correction of the full interval (last entry) - if L.tau is not None: - L.uend += L.tau[-1] - - return None \ No newline at end of file +- # 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] + # add up tau correction of the full interval (last entry) + if L.tau is not None: + L.uend += L.tau[-1] + + return None From 5bb103a61983d67cdb58516f77c18bbd4a35c6b3 Mon Sep 17 00:00:00 2001 From: Daniel Ruprecht Date: Tue, 22 Dec 2015 09:11:15 +0000 Subject: [PATCH 2/4] Swepper correctly throws AssertionError if do_coll_update is set to False for collocation object with right_is_node=False --- pySDC/Sweeper.py | 2 +- tests/test_imexsweeper.py | 12 +++++++++++- 2 files changed, 12 insertions(+), 2 deletions(-) diff --git a/pySDC/Sweeper.py b/pySDC/Sweeper.py index a19e7398af..ff2e7850bd 100644 --- a/pySDC/Sweeper.py +++ b/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 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" + 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" # This will be set as soon as the sweeper is instantiated at the level self.__level = None diff --git a/tests/test_imexsweeper.py b/tests/test_imexsweeper.py index 6f0e40166a..2de74557a3 100644 --- a/tests/test_imexsweeper.py +++ b/tests/test_imexsweeper.py @@ -204,7 +204,7 @@ def test_manysweepsequalmatrix(self): # # Make sure that update function for K sweeps computed from K-sweep matrix gives same result as K sweeps in node-to-node form plus compute_end_point # - def test_maysweepupdate(self): + def test_manysweepupdate(self): step, level, problem, nnodes = self.setupLevelStepProblem() step.levels[0].sweep.predict() @@ -229,3 +229,13 @@ def test_maysweepupdate(self): # Multiply u0 by value of update function to get end value directly uend_matrix = update*self.pparams['u0'] assert abs(uend_matrix - uend_sweep)<1e-14, "Node-to-node sweep plus update yields different result than update function computed through K-sweep matrix" + + # + # Make sure that creating a sweeper object with a collocation object with right_is_node=False and do_coll_update=False throws an exception + # + def test_norightnode_collupdate_fails(self): + self.swparams['collocation_class'] = collclass.CollGaussLegendre + self.swparams['do_coll_update'] = False + # Has to throw an exception + with self.assertRaises(AssertionError): + step, level, problem, nnodes = self.setupLevelStepProblem() From 9cf2c4a5afa14a31bfc7ddd3a2fbc2ecf8b27552 Mon Sep 17 00:00:00 2001 From: Daniel Ruprecht Date: Tue, 22 Dec 2015 09:38:21 +0000 Subject: [PATCH 3/4] 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 --- pySDC/Sweeper.py | 2 +- pySDC/sweeper_classes/generic_LU.py | 8 +++--- pySDC/sweeper_classes/imex_1st_order.py | 21 ++++++++------ tests/test_imexsweeper.py | 38 +++++++++++++++++++++++++ 4 files changed, 56 insertions(+), 13 deletions(-) diff --git a/pySDC/Sweeper.py b/pySDC/Sweeper.py index ff2e7850bd..cd90e1873d 100644 --- a/pySDC/Sweeper.py +++ b/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 diff --git a/pySDC/sweeper_classes/generic_LU.py b/pySDC/sweeper_classes/generic_LU.py index 77d7b5542c..9f1b91f27a 100644 --- a/pySDC/sweeper_classes/generic_LU.py +++ b/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]) diff --git a/pySDC/sweeper_classes/imex_1st_order.py b/pySDC/sweeper_classes/imex_1st_order.py index d35edfc270..70abaa320c 100644 --- a/pySDC/sweeper_classes/imex_1st_order.py +++ b/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 diff --git a/tests/test_imexsweeper.py b/tests/test_imexsweeper.py index 2de74557a3..36ef52a985 100644 --- a/tests/test_imexsweeper.py +++ b/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)" From 9bc47ca2e8a211d551c615a1673f10af50658e4c Mon Sep 17 00:00:00 2001 From: Daniel Ruprecht Date: Tue, 22 Dec 2015 09:45:05 +0000 Subject: [PATCH 4/4] removed leftover print statements --- tests/test_imexsweeper.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/tests/test_imexsweeper.py b/tests/test_imexsweeper.py index 36ef52a985..e619ab26d9 100644 --- a/tests/test_imexsweeper.py +++ b/tests/test_imexsweeper.py @@ -272,8 +272,5 @@ def test_updateformula_no_coll_update(self): # 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)"