Hardcoded references for battery_2condensators

Author: lisawim

pySDC/projects/PinTSimE/battery_2condensators_model.py

@@ -7,7 +7,6 @@
 from pySDC.implementations.problem_classes.Battery_2Condensators import battery_2condensators
 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
@@ -80,6 +79,8 @@ def main(use_switch_estimator=True):
     level_params['restol'] = -1
     level_params['dt'] = 1e-2
 
+    assert level_params['dt'] == 1e-2, 'Error! Do not use the time step dt != 1e-2!'
+
     # initialize sweeper parameters
     sweeper_params = dict()
     sweeper_params['quad_type'] = 'LOBATTO'
@@ -122,7 +123,6 @@ def main(use_switch_estimator=True):
     description['sweeper_params'] = sweeper_params  # pass sweeper parameters
     description['level_params'] = level_params  # pass level parameters
     description['step_params'] = step_params
-    description['space_transfer_class'] = mesh_to_mesh  # pass spatial transfer class
 
     if use_switch_estimator:
         description['convergence_controllers'] = convergence_controllers
@@ -176,6 +176,8 @@ def main(use_switch_estimator=True):
 
     recomputed = False
 
+    check_solution(stats, use_switch_estimator)
+
     plot_voltages(description, recomputed, use_switch_estimator)
 
     return description
@@ -228,6 +230,68 @@ def plot_voltages(description, recomputed, use_switch_estimator, cwd='./'):
     plt_helper.plt.close(fig)
 
 
+def check_solution(stats, use_switch_estimator):
+    """
+    Function that checks the solution based on a hardcoded reference solution. Based on check_solution function from @brownbaerchen.
+
+    Args:
+        stats (dict): Raw statistics from a controller run
+        use_switch_estimator (bool): flag if the switch estimator wants to be used or not
+    """
+
+    data = get_data_dict(stats, use_switch_estimator)
+
+    if use_switch_estimator:
+        msg = 'Error when using the switch estimator for battery_2condensators:'
+        expected = {
+            'cL': 1.1597046304825833,
+            'vC1': 1.000472118416924,
+            'vC2': 1.000226101799117,
+            'switch1': 1.6094379124373626,
+            'switch2': 3.2184040405613974,
+            'restarts': 2.0,
+        }
+
+    got = {
+        'cL': data['cL'][-1],
+        'vC1': data['vC1'][-1],
+        'vC2': data['vC2'][-1],
+        'switch1': data['switch1'],
+        'switch2': data['switch2'],
+        'restarts': data['restarts'],
+    }
+
+    for key in expected.keys():
+        assert np.isclose(
+            expected[key], got[key], rtol=1e-4
+        ), f'{msg} Expected {key}={expected[key]:.4e}, got {key}={got[key]:.4e}'
+
+
+def get_data_dict(stats, use_switch_estimator, recomputed=False):
+    """
+    Converts the statistics in a useful data dictionary so that it can be easily checked in the check_solution function.
+    Based on @brownbaerchen's get_data function.
+
+    Args:
+        stats (dict): Raw statistics from a controller run
+        use_switch_estimator (bool): flag if the switch estimator wants to be used or not
+        recomputed (bool): flag if the values after a restart are used or before
+
+    Return:
+        data (dict): contains all information as the statistics dict
+    """
+
+    data = dict()
+    data['cL'] = np.array(get_sorted(stats, type='current L', recomputed=recomputed, sortby='time'))[:, 1]
+    data['vC1'] = np.array(get_sorted(stats, type='voltage C1', recomputed=recomputed, sortby='time'))[:, 1]
+    data['vC2'] = np.array(get_sorted(stats, type='voltage C2', recomputed=recomputed, sortby='time'))[:, 1]
+    data['switch1'] = np.array(get_recomputed(stats, type='switch', sortby='time'))[0, 1]
+    data['switch2'] = np.array(get_recomputed(stats, type='switch', sortby='time'))[-1, 1]
+    data['restarts'] = np.sum(np.array(get_sorted(stats, type='restart', recomputed=None, sortby='time'))[:, 1])
+
+    return data
+
+
 def proof_assertions_description(description, use_switch_estimator):
     """
     Function to proof the assertions (function to get cleaner code)

pySDC/projects/PinTSimE/battery_model.py

@@ -336,9 +336,10 @@ def check_solution(stats, problem, use_adaptivity, use_switch_estimator):
         ), f'{msg} Expected {key}={expected[key]:.4e}, got {key}={got[key]:.4e}'
 
 
-def get_data_dict(stats, use_adaptivity, use_switch_estimator, recomputed=False):
+def get_data_dict(stats, use_adaptivity=True, use_switch_estimator=True, recomputed=False):
     """
     Converts the statistics in a useful data dictionary so that it can be easily checked in the check_solution function.
+    Based on @brownbaerchen's get_data function.
 
     Args:
         stats (dict): Raw statistics from a controller run

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.projects.PinTSimE.battery_model import get_recomputed
+from pySDC.projects.PinTSimE.battery_2condensators_model import get_data_dict
 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
@@ -33,6 +34,10 @@ def run(dt, use_switch_estimator=True):
     level_params['restol'] = -1
     level_params['dt'] = dt
 
+    assert (
+        dt == 4e-1 or dt == 4e-2 or dt == 4e-3
+    ), "Error! Do not use other time steps dt != 4e-1 or dt != 4e-2 or dt != 4e-3 due to hardcoded references!"
+
     # initialize sweeper parameters
     sweeper_params = dict()
     sweeper_params['quad_type'] = 'LOBATTO'
@@ -143,9 +148,10 @@ def check(cwd='./'):
             stats, description = run(dt=dt_item, use_switch_estimator=use_SE)
 
             if use_SE:
-                assert (
-                    len(get_recomputed(stats, type='switch', sortby='time')) >= 1
-                ), 'No switches found for dt={}!'.format(dt_item)
+                switches = get_recomputed(stats, type='switch', sortby='time')
+                assert (len(switches) >= 2), f"Expected at least 2 switches for dt: {dt_item}, got {len(switches)}!"
+
+                check_solution(stats, dt_item, use_SE)
 
             fname = 'data/battery_2condensators_dt{}_USE{}.dat'.format(dt_item, use_SE)
             f = open(fname, 'wb')
@@ -299,5 +305,62 @@ def check(cwd='./'):
     plt_helper.plt.close(fig2)
 
 
+def check_solution(stats, dt, use_switch_estimator):
+    """
+    Function that checks the solution based on a hardcoded reference solution. Based on check_solution function from @brownbaerchen.
+
+    Args:
+        stats (dict): Raw statistics from a controller run
+        dt (float): initial time step
+        use_switch_estimator (bool): flag if the switch estimator wants to be used or not
+    """
+
+    data = get_data_dict(stats, use_switch_estimator)
+
+    if use_switch_estimator:
+        msg = f'Error when using the switch estimator for battery_2condensators for dt={dt:.1e}:'
+        if dt == 4e-1:
+            expected = {
+                'cL': 1.1556732037544801,
+                'vC1': 1.002239522400514,
+                'vC2': 1.0000329223874842,
+                'switch1': 1.607586793484041,
+                'switch2': 3.216645438176962,
+                'restarts': 2.0,
+            }
+        elif dt == 4e-2:
+            expected = {
+                'cL': 1.1492603893091364,
+                'vC1': 1.0005011122241925,
+                'vC2': 1.000015039670507,
+                'switch1': 1.6094074085596919,
+                'switch2': 3.2183750611596893,
+                'restarts': 2.0,
+            }
+        elif dt == 4e-3:
+            expected = {
+                'cL': 1.1476283937778273,
+                'vC1': 1.0001336962511904,
+                'vC2': 1.0000182217245925,
+                'switch1': 1.6093728710270498,
+                'switch2': 3.218742142249058,
+                'restarts': 2.0,
+            }
+
+    got = {
+        'cL': data['cL'][-1],
+        'vC1': data['vC1'][-1],
+        'vC2': data['vC2'][-1],
+        'switch1': data['switch1'],
+        'switch2': data['switch2'],
+        'restarts': data['restarts'],
+    }
+
+    for key in expected.keys():
+        assert np.isclose(
+            expected[key], got[key], rtol=1e-4
+        ), f'{msg} Expected {key}={expected[key]:.4e}, got {key}={got[key]:.4e}'
+
+
 if __name__ == "__main__":
     check()