Added some more tests

Author: lisawim

pySDC/projects/PinTSimE/battery_2condensators_model.py

@@ -8,6 +8,7 @@
 from pySDC.implementations.sweeper_classes.imex_1st_order import imex_1st_order
 from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI
 from pySDC.implementations.transfer_classes.TransferMesh import mesh_to_mesh
+from pySDC.projects.PinTSimE.battery_model import get_recomputed
 from pySDC.projects.PinTSimE.piline_model import setup_mpl
 import pySDC.helpers.plot_helper as plt_helper
 from pySDC.core.Hooks import hooks
@@ -167,12 +168,14 @@ def main(use_switch_estimator=True):
     assert np.mean(niters) <= 4, "Mean number of iterations is too high, got %s" % np.mean(niters)
     f.close()
 
-    plot_voltages(description, use_switch_estimator)
+    recomputed = False
+
+    plot_voltages(description, recomputed, use_switch_estimator)
 
     return np.mean(niters)
 
 
-def plot_voltages(description, use_switch_estimator, cwd='./'):
+def plot_voltages(description, recomputed, use_switch_estimator, cwd='./'):
     """
     Routine to plot the numerical solution of the model
     """
@@ -182,9 +185,9 @@ def plot_voltages(description, use_switch_estimator, cwd='./'):
     f.close()
 
     # convert filtered statistics to list of iterations count, sorted by process
-    cL = get_sorted(stats, type='current L', sortby='time')
-    vC1 = get_sorted(stats, type='voltage C1', sortby='time')
-    vC2 = get_sorted(stats, type='voltage C2', sortby='time')
+    cL = get_sorted(stats, type='current L', recomputed=recomputed, sortby='time')
+    vC1 = get_sorted(stats, type='voltage C1', recomputed=recomputed, sortby='time')
+    vC2 = get_sorted(stats, type='voltage C2', recomputed=recomputed, sortby='time')
 
     times = [v[0] for v in cL]
 
@@ -195,7 +198,10 @@ def plot_voltages(description, use_switch_estimator, cwd='./'):
     ax.plot(times, [v[1] for v in vC2], label='$v_{C_2}$')
 
     if use_switch_estimator:
-        switches = get_sorted(stats, type='switch', sortby='time')
+        switches = get_recomputed(stats, type='switch', sortby='time')
+
+        if recomputed is not None:
+            assert len(switches) >= 1 and len(switches) >= 2, "No switches found"
         t_switches = [v[1] for v in switches]
 
         for i in range(len(t_switches)):
@@ -222,11 +228,11 @@ def proof_assertions_description(description, use_switch_estimator):
         description['problem_params']['alpha'] > description['problem_params']['V_ref'][1]
     ), 'Please set "alpha" greater than "V_ref2"'
 
+    assert type(description['problem_params']['V_ref']) == np.ndarray, '"V_ref" needs to be an array (of type float)'
+
     assert description['problem_params']['V_ref'][0] > 0, 'Please set "V_ref1" greater than 0'
     assert description['problem_params']['V_ref'][1] > 0, 'Please set "V_ref2" greater than 0'
 
-    assert type(description['problem_params']['V_ref']) == np.ndarray, '"V_ref" needs to be an array (of type float)'
-
     assert 'errtol' not in description['step_params'].keys(), 'No exact solution known to compute error'
     assert 'alpha' in description['problem_params'].keys(), 'Please supply "alpha" in the problem parameters'
     assert 'V_ref' in description['problem_params'].keys(), 'Please supply "V_ref" in the problem parameters'

pySDC/projects/PinTSimE/battery_model.py

@@ -132,6 +132,8 @@ def main(dt, problem, sweeper, use_switch_estimator, use_adaptivity):
     t0 = 0.0
     Tend = 0.3
 
+    assert dt < Tend, "Time step is too large for the time domain!"
+
     # instantiate controller
     controller = controller_nonMPI(num_procs=1, controller_params=controller_params, description=description)
 
@@ -167,7 +169,7 @@ def main(dt, problem, sweeper, use_switch_estimator, use_adaptivity):
         min_iter = min(min_iter, item[1])
         max_iter = max(max_iter, item[1])
 
-    assert np.mean(niters) <= 5, "Mean number of iterations is too high, got %s" % np.mean(niters)
+    assert np.mean(niters) <= 4, "Mean number of iterations is too high, got %s" % np.mean(niters)
     f.close()
 
     return description
@@ -182,6 +184,7 @@ def run():
     dt = 1e-3
     problem_classes = [battery, battery_implicit]
     sweeper_classes = [imex_1st_order, generic_implicit]
+    recomputed = False
     use_switch_estimator = [True]
     use_adaptivity = [True]
 
@@ -196,10 +199,10 @@ def run():
                     use_adaptivity=use_A,
                 )
 
-            plot_voltages(description, problem.__name__, sweeper.__name__, use_SE, use_A)
+            plot_voltages(description, problem.__name__, sweeper.__name__, recomputed, use_SE, use_A)
 
 
-def plot_voltages(description, problem, sweeper, use_switch_estimator, use_adaptivity, cwd='./'):
+def plot_voltages(description, problem, sweeper, recomputed, use_switch_estimator, use_adaptivity, cwd='./'):
     """
     Routine to plot the numerical solution of the model
     """
@@ -221,8 +224,10 @@ def plot_voltages(description, problem, sweeper, use_switch_estimator, use_adapt
     ax.plot(times, [v[1] for v in vC], label=r'$v_C$')
 
     if use_switch_estimator:
-        val_switch = get_sorted(stats, type='switch', sortby='time')
-        t_switch = [v[1] for v in val_switch]
+        switches = get_recomputed(stats, type='switch', sortby='time')
+
+        assert len(switches) >= 1, 'No switches found!'
+        t_switch = [v[1] for v in switches]
         ax.axvline(x=t_switch[-1], linestyle='--', linewidth=0.8, color='r', label='Switch')
 
     if use_adaptivity:
@@ -241,7 +246,7 @@ def plot_voltages(description, problem, sweeper, use_switch_estimator, use_adapt
 
     fig.savefig('data/{}_model_solution_{}.png'.format(problem, sweeper), dpi=300, bbox_inches='tight')
     plt_helper.plt.close(fig)
-    
+
 def get_recomputed(stats, type, sortby):
     """
     Function that filters statistics after a recomputation

pySDC/projects/PinTSimE/estimation_check.py

@@ -96,6 +96,8 @@ def run(dt, problem, sweeper, use_switch_estimator, use_adaptivity, V_ref):
     t0 = 0.0
     Tend = 0.3
 
+    assert dt < Tend, "Time step is too large for the time domain!"
+
     # instantiate controller
     controller = controller_nonMPI(num_procs=1, controller_params=controller_params, description=description)
 
@@ -156,6 +158,8 @@ def check(cwd='./'):
                         V_ref=V_ref,
                     )
 
+                    assert len(get_recomputed(stats, type='switch', sortby='time')) >= 1, 'No switches found!'
+
                     fname = 'data/battery_dt{}_USE{}_USA{}_{}.dat'.format(dt_item, use_SE, use_A, sweeper.__name__)
                     f = open(fname, 'wb')
                     dill.dump(stats, f)
@@ -225,8 +229,8 @@ def accuracy_check(dt_list, problem, sweeper, V_ref, cwd='./'):
         stats_adapt = dill.load(f4)
         f4.close()
 
-        val_switch_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
-        t_switch_SE_adapt = [v[1] for v in val_switch_SE_adapt]
+        switches_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
+        t_switch_SE_adapt = [v[1] for v in switches_SE_adapt]
         t_switch_SE_adapt = t_switch_SE_adapt[-1]
 
         dt_SE_adapt_val = get_sorted(stats_SE_adapt, type='dt', recomputed=False)
@@ -349,12 +353,12 @@ def differences_around_switch(
         stats_adapt = dill.load(f4)
         f4.close()
 
-        val_switch_SE = get_recomputed(stats_SE, type='switch', sortby='time')
-        t_switch = [v[1] for v in val_switch_SE]
+        switches_SE = get_recomputed(stats_SE, type='switch', sortby='time')
+        t_switch = [v[1] for v in switches_SE]
         t_switch = t_switch[-1]  # battery has only one single switch
 
-        val_switch_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
-        t_switch_SE_adapt = [v[1] for v in val_switch_SE_adapt]
+        switches_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
+        t_switch_SE_adapt = [v[1] for v in switches_SE_adapt]
         t_switch_SE_adapt = t_switch_SE_adapt[-1]
 
         vC_SE = get_sorted(stats_SE, type='voltage C', recomputed=False, sortby='time')
@@ -493,12 +497,12 @@ def differences_over_time(dt_list, problem, sweeper, V_ref, cwd='./'):
         stats_adapt = dill.load(f4)
         f4.close()
 
-        val_switch_SE = get_recomputed(stats_SE, type='switch', sortby='time')
-        t_switch_SE = [v[1] for v in val_switch_SE]
+        switches_SE = get_recomputed(stats_SE, type='switch', sortby='time')
+        t_switch_SE = [v[1] for v in switches_SE]
         t_switch_SE = t_switch_SE[-1]  # battery has only one single switch
 
-        val_switch_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
-        t_switch_SE_adapt = [v[1] for v in val_switch_SE_adapt]
+        switches_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
+        t_switch_SE_adapt = [v[1] for v in switches_SE_adapt]
         t_switch_SE_adapt = t_switch_SE_adapt[-1]
 
         dt_adapt = np.array(get_sorted(stats_adapt, type='dt', recomputed=False, sortby='time'))
@@ -654,12 +658,12 @@ def iterations_over_time(dt_list, maxiter, problem, sweeper, cwd='./'):
         times_adapt.append([v[0] for v in iter_counts_adapt_val])
         times.append([v[0] for v in iter_counts_val])
 
-        val_switch_SE = get_recomputed(stats_SE, type='switch', sortby='time')
-        t_switch_SE = [v[1] for v in val_switch_SE]
+        switches_SE = get_recomputed(stats_SE, type='switch', sortby='time')
+        t_switch_SE = [v[1] for v in switches_SE]
         t_switches_SE.append(t_switch_SE[-1])
 
-        val_switch_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
-        t_switch_SE_adapt = [v[1] for v in val_switch_SE_adapt]
+        switches_SE_adapt = get_recomputed(stats_SE_adapt, type='switch', sortby='time')
+        t_switch_SE_adapt = [v[1] for v in switches_SE_adapt]
         t_switches_SE_adapt.append(t_switch_SE_adapt[-1])
 
     if len(dt_list) > 1:

pySDC/projects/PinTSimE/estimation_check_extended.py

@@ -8,6 +8,7 @@
 from pySDC.implementations.sweeper_classes.imex_1st_order import imex_1st_order
 from pySDC.implementations.controller_classes.controller_nonMPI import controller_nonMPI
 from pySDC.implementations.transfer_classes.TransferMesh import mesh_to_mesh
+from pySDC.projects.PinTSimE.battery_model import get_recomputed
 from pySDC.projects.PinTSimE.piline_model import setup_mpl
 from pySDC.projects.PinTSimE.battery_2condensators_model import log_data, proof_assertions_description
 import pySDC.helpers.plot_helper as plt_helper
@@ -128,6 +129,8 @@ def check(cwd='./'):
         for item in use_switch_estimator:
             stats, description = run(dt=dt_item, use_switch_estimator=item)
 
+            assert len(get_recomputed(stats, type='switch', sortby='time')) >= 1, 'No switches found!'
+
             fname = 'data/battery_2condensators_dt{}_USE{}.dat'.format(dt_item, item)
             f = open(fname, 'wb')
             dill.dump(stats, f)
@@ -159,13 +162,8 @@ def check(cwd='./'):
         stats_false = dill.load(f2)
         f2.close()
 
-        val_switch = get_sorted(stats_true, type='switch', sortby='time')
-        #val_switch2 = get_sorted(stats_true, type='switch2', sortby='time')
-        t_switch = [v[1] for v in val_switch]
-        #t_switch2 = [v[1] for v in val_switch2]
-
-        #t_switch1 = t_switch1[-1]
-        #t_switch2 = t_switch2[-1]
+        switches = get_recomputed(stats_true, type='switch', sortby='time')
+        t_switch = [v[1] for v in switches]
 
         val_switch_all.append([t_switch[0], t_switch[1]])
 

pySDC/projects/PinTSimE/switch_estimator.py

@@ -35,11 +35,10 @@ def setup(self, controller, params, description):
         defaults = {
             'control_order': 100,
             'tol': description['level_params']['dt'],
-            'coll_nodes_local': coll.nodes,
+            'coll_nodes': coll.nodes,
             'switch_detected': False,
             'switch_detected_step': False,
             't_switch': None,
-            'count_switches': 0,
             'dt_initial': description['level_params']['dt'],
         }
         return {**defaults, **params}
@@ -66,7 +65,7 @@ def get_new_step_size(self, controller, S):
 
             if self.params.switch_detected:
                 t_interp = [
-                    L.time + L.dt * self.params.coll_nodes_local[m] for m in range(len(self.params.coll_nodes_local))
+                    L.time + L.dt * self.params.coll_nodes[m] for m in range(len(self.params.coll_nodes))
                 ]
 
                 # only find root if vc_switch[0], vC_switch[-1] have opposite signs (intermediate value theorem)
@@ -158,7 +157,7 @@ def post_step_processing(self, controller, S):
                 self.params.t_switch = None
                 self.params.switch_detected_step = False
 
-                dt_planned = L.status.dt_new if L.status.dt_new is not None else self.params.dt_initial
+                dt_planned = L.status.dt_new if L.status.dt_new is not None else L.params.dt
                 L.status.dt_new = dt_planned
 
         super(SwitchEstimator, self).post_step_processing(controller, S)