Stuff

#!!!!!!!!!! WARNING: RUFF FAILED !!!!!!!!!!: 

#pySDC/projects/GPU/analysis_scripts/compare_RBC3D.py:298:5: F841 Local variable `data` is assigned to but never used
#    |
Parallel-in-Time#296 |     res = 32
Parallel-in-Time#297 |
Parallel-in-Time#298 |     data = []
#    |     ^^^^ F841
Parallel-in-Time#299 |     labels = []
Parallel-in-Time#300 |     linestyles = []
#    |
#    = help: Remove assignment to unused variable `data`
#
#pySDC/projects/GPU/analysis_scripts/compare_RBC3D.py:299:5: F841 Local variable `labels` is assigned to but never used
#    |
Parallel-in-Time#298 |     data = []
Parallel-in-Time#299 |     labels = []
#    |     ^^^^^^ F841
Parallel-in-Time#300 |     linestyles = []
#    |
#    = help: Remove assignment to unused variable `labels`
#
#pySDC/projects/GPU/analysis_scripts/compare_RBC3D.py:300:5: F841 Local variable `linestyles` is assigned to but never used
#    |
Parallel-in-Time#298 |     data = []
Parallel-in-Time#299 |     labels = []
Parallel-in-Time#300 |     linestyles = []
#    |     ^^^^^^^^^^ F841
Parallel-in-Time#301 |
Parallel-in-Time#302 |     ref_data = get_pySDC_data(Ra, res=res, dt=0.01, config_name='RBC3DG4R4')
#    |
#    = help: Remove assignment to unused variable `linestyles`
#
#pySDC/projects/GPU/analysis_scripts/compare_RBC3D.py:312:9: F841 Local variable `t` is assigned to but never used
#    |
Parallel-in-Time#310 |         ax.set_title(title)
Parallel-in-Time#311 |         Nu_ref = np.array(ref['Nu']['V'])
Parallel-in-Time#312 |         t = ref['t']
#    |         ^ F841
Parallel-in-Time#313 |
Parallel-in-Time#314 |         for dt in dts:
#    |
#    = help: Remove assignment to unused variable `t`
#
#pySDC/projects/GPU/analysis_scripts/compare_RBC3D.py:587:9: E722 Do not use bare `except`
#    |
Parallel-in-Time#585 |                 k[_s > 1e-16], _s[_s > 1e-16], color=last_line.get_color(), ls=last_line.get_linestyle(), label=label
Parallel-in-Time#586 |             )
Parallel-in-Time#587 |         except:
#    |         ^^^^^^ E722
Parallel-in-Time#588 |             pass
#    |
#
#pySDC/projects/GPU/analysis_scripts/compare_RBC3D.py:686:5: F841 Local variable `Delta_Nu` is assigned to but never used
#    |
#684 |     t = data['t']
#685 |     avg_Nu = np.array([np.mean(Nu[40 : 40 + i + 1]) for i in range(len(Nu[40:]))])
#686 |     Delta_Nu = np.array([abs(avg_Nu[i + 1] - avg_Nu[i]) for i in range(len(avg_Nu) - 1)])
#    |     ^^^^^^^^ F841
#687 |     # ax.plot(data['t'][40:-1], Delta_Nu / avg_Nu[:-1])
#688 |     # ax.plot(data['t'], np.abs(avg_Nu - avg_Nu[-1]) / avg_Nu[-1])
#    |
#    = help: Remove assignment to unused variable `Delta_Nu`:

Author: brownbaerchen

pySDC/helpers/plot_helper.py

@@ -45,13 +45,15 @@ def figsize_by_journal(journal, scale, ratio):  # pragma: no cover
         'Springer_proceedings': 347.12354,
         'JSC_thesis': 434.26027,
         'TUHH_thesis': 426.79135,
+        'Nature_CS': 372.0,
     }
     # store text height in points here, get this from LaTeX using \the\textheight
     textheights = {
         'JSC_beamer': 214.43411,
         'JSC_thesis': 635.5,
         'TUHH_thesis': 631.65118,
         'Springer_proceedings': 549.13828,
+        'Nature_CS': 552.69478,
     }
     assert (
         journal in textwidths.keys()

pySDC/projects/GPU/analysis_scripts/RBC3D_order.py

@@ -73,12 +73,13 @@ def compare_order(Ra):
     fig, ax = plt.subplots()
     ls = {'SDC': '-', 'RK': '--', 'Euler': '-.'}
     if Ra == 1e5:
-        paths = [f'./data/RBC3DG4R4{me}Ra1e5-res-1-order.pickle' for me in ['', 'RK', 'Euler']]
-        labels = ['SDC', 'RK', 'Euler']
+        paths = [f'./data/RBC3DG4R4{me}Ra1e5-res-1-order.pickle' for me in ['', 'RK', 'Euler', 'O3']]
+        labels = ['SDC', 'RK', 'Euler', 'SDC']
+
     else:
         raise NotImplementedError
 
-    for path, label in zip(paths, labels):
+    for path, label in zip(paths, labels, strict=True):
         with open(path, 'rb') as file:
             errors = pickle.load(file)
 

pySDC/projects/GPU/analysis_scripts/compare_RBC3D.py

@@ -2,6 +2,9 @@
 import matplotlib.pyplot as plt
 import numpy as np
 from scipy import integrate
+from pySDC.helpers.plot_helper import figsize_by_journal, setup_mpl
+
+setup_mpl()
 
 
 def get_pySDC_data(Ra, RK=False, res=-1, dt=-1, config_name='RBC3DG4'):
@@ -282,11 +285,12 @@ def interpolate_NuV_to_reference_times(data, reference_data, order=12):
 
 
 def compare_Nusselt_over_time1e5():
-    fig, Nu_ax = plt.subplots()
-    _, axs = plt.subplots(1, 3, figsize=(10, 3))
-    prof_ax = axs[0]
-    rms_ax = axs[1]
-    spectrum_ax = axs[2]
+    Nu_fig, Nu_axs = plt.subplots(4, 1, sharex=True, sharey=True, figsize=figsize_by_journal('Nature_CS', 1, 1.4))
+    # _, axs = plt.subplots(1, 3, figsize=(10, 3))
+    # prof_ax = axs[0]
+    # rms_ax = axs[1]
+    # spectrum_ax = axs[2]
+    # deviation_fig, deviation_ax = plt.subplots(figsize=figsize_by_journal('Nature_CS', 0.7, 0.7))
 
     Ra = '1e5'
     res = 32
@@ -296,89 +300,151 @@ def compare_Nusselt_over_time1e5():
     linestyles = []
 
     ref_data = get_pySDC_data(Ra, res=res, dt=0.01, config_name='RBC3DG4R4')
+    # deviations = {'SDC 3': [], 'SDC': [], 'RK': [], 'Euler': []}
+    # Nu_all = {'SDC 3': [], 'SDC': [], 'RK': [], 'Euler': []}
+
+    _Nu_axs = {'SDC 3': Nu_axs[1], 'SDC': Nu_axs[0], 'RK': Nu_axs[2], 'Euler': Nu_axs[3]}
+
+    def _plot_Nu_(Ra, res, dts, config_name, ref, ax, title):
+        ax.plot(ref['t'], ref['Nu']['V'], color='black', ls='--')
+        ax.set_title(title)
+        Nu_ref = np.array(ref['Nu']['V'])
+        t = ref['t']
+
+        for dt in dts:
+            data = get_pySDC_data(Ra=Ra, res=res, dt=dt, config_name=config_name)
+            t_i, Nu_i = interpolate_NuV_to_reference_times(data, ref)
+            ax.plot(t_i, Nu_i, label=rf'$\Delta t={{{dt}}}$')
+
+            # error = np.maximum.accumulate(np.abs(Nu_ref[: Nu_i.shape[0]] - Nu_i) / np.abs(Nu_ref[: Nu_i.shape[0]]))
+            error = np.abs(Nu_ref[: Nu_i.shape[0]] - Nu_i) / np.abs(Nu_ref[: Nu_i.shape[0]])
+
+            # compute mean Nu
+            mask = np.logical_and(t_i >= 100, t_i <= 200)
+            Nu_mean = np.mean(Nu_i[mask])
+            Nu_std = np.std(Nu_i[mask])
+
+            last_line = ax.get_lines()[-1]
+            if any(error > 1e-2):
+                deviates = min(t_i[error > 1e-2])
+                ax.axvline(deviates, color=last_line.get_color(), ls=':')
+                print(f'{title} dt={dt} Nu={Nu_mean:.3f}+={Nu_std:.3f}, deviates more than 1% from t={deviates:.2f}')
+            else:
+                print(f'{title} dt={dt} Nu={Nu_mean:.3f}+={Nu_std:.3f}')
+            ax.legend(frameon=True, loc='upper left')
+
+    _plot_Nu_('1e5', 32, [0.06, 0.04, 0.02], 'RBC3DG4R4', ref_data, Nu_axs[0], 'SDC44')
+    _plot_Nu_(
+        '1e5',
+        32,
+        [
+            0.06,
+            0.04,
+            0.02,
+        ],
+        'RBC3DG4R4O4',
+        ref_data,
+        Nu_axs[0],
+        'SDC34',
+    )
+    _plot_Nu_('1e5', 32, [0.06, 0.05, 0.02, 0.01], 'RBC3DG4R4O3', ref_data, Nu_axs[1], 'SDC23')
+    _plot_Nu_('1e5', 32, [0.05, 0.04, 0.02, 0.01, 0.005], 'RBC3DG4R4RK', ref_data, Nu_axs[2], 'RK443')
+    _plot_Nu_('1e5', 32, [0.02, 0.01, 0.005], 'RBC3DG4R4Euler', ref_data, Nu_axs[3], 'RK111')
+
+    # for dt in [0.06, 0.02, 0.01]:
+    #     data.append(get_pySDC_data(Ra, res=res, dt=dt, config_name='RBC3DG4R4'))
+    #     labels.append(f'SDC, dt={dt:.4f}')
+    #     linestyles.append('-')
+
+    # for dt in [0.06, 0.05, 0.02, 0.01]:
+    #     data.append(get_pySDC_data(Ra, res=res, dt=dt, config_name='RBC3DG4R4O3'))
+    #     labels.append(f'SDC 3, dt={dt:.4f}')
+    #     linestyles.append('-')
 
-    _, ting_ax = plt.subplots()
-
-    for dt in [0.06, 0.02, 0.01]:
-        data.append(get_pySDC_data(Ra, res=res, dt=dt, config_name='RBC3DG4R4'))
-        labels.append(f'SDC, dt={dt:.4f}')
-        linestyles.append('-')
-
-    # ----------------- RK ------------------------
-
-    for dt in [0.05, 0.04, 0.02, 0.01, 0.005]:
-        data.append(get_pySDC_data(Ra, res=res, dt=dt, RK=True, config_name='RBC3DG4R4'))
-        labels.append(f'RK, dt={dt:.3f}')
-        linestyles.append('--')
-
-    # ----------------- Euler ---------------------
-
-    for dt in [0.02, 0.005]:
-        data.append(get_pySDC_data(Ra, res=res, dt=dt, config_name='RBC3DG4R4Euler'))
-        labels.append(f'Euler, dt={dt:.3f}')
-        linestyles.append('-.')
-
-    for dat, label, linestyle in zip(data, labels, linestyles):
-        Nu = np.array(dat['Nu']['V'])
-        t = dat['t']
-        Nu_ref = np.array(ref_data['Nu']['V'])
-        t_i, Nu_i = interpolate_NuV_to_reference_times(dat, ref_data)
-        Nu_ax.plot(t, Nu, label=label, ls=linestyle)
-
-        error = np.maximum.accumulate(np.abs(Nu_ref[: Nu_i.shape[0]] - Nu_i) / np.abs(Nu_ref[: Nu_i.shape[0]]))
-
-        # compute mean Nu
-        mask = np.logical_and(t >= 20, t <= 200)
-        Nu_mean = np.mean(Nu[mask])
-        Nu_std = np.std(Nu[mask])
-
-        last_line = Nu_ax.get_lines()[-1]
-        if any(error > 1e-2):
-            deviates = min(t_i[error > 1e-2])
-            Nu_ax.axvline(deviates, color=last_line.get_color(), ls=last_line.get_linestyle())
-            print(f'{label} Nu={Nu_mean:.3f}+={Nu_std:.3f}, deviates more than 1% from t={deviates:.2f}')
-        else:
-            print(f'{label} Nu={Nu_mean:.3f}+={Nu_std:.3f}')
-
-        k = dat['k']
-        spectrum = np.array(dat['spectrum'])
-        u_spectrum = np.mean(spectrum, axis=0)[1]
-        idx = dat['res_in_boundary_layer']
-        _s = u_spectrum[idx]
-        spectrum_ax.loglog(
-            k[_s > 1e-16], _s[_s > 1e-16], color=last_line.get_color(), ls=last_line.get_linestyle(), label=label
-        )
-
-        prof_ax.plot(dat['profile_T'], dat['z'], color=last_line.get_color(), ls=last_line.get_linestyle(), label=label)
-        rms_ax.plot(
-            dat['rms_profile_T'], dat['z'], color=last_line.get_color(), ls=last_line.get_linestyle(), label=label
-        )
-
-        print(dat['avg_Nu'])
-        try:
-            ting_ax.scatter(abs(dat['avg_Nu']['V'] - dat['avg_Nu']['thermal']), dat['avg_Nu']['b'], label=label)
-        except:
-            pass
-        ting_ax.set_xscale('log')
-        ting_ax.legend()
-
-        # if 'dissipation' in dat.keys():
-
-        #     mean_ke = np.mean([me[1] for me in dat['dissipation']])
-        #     print(
-        #         f'Energy error: {(dat["dissipation"][-1][1] - np.sum([me[0] for me in dat["dissipation"]])/2)/mean_ke:.2e}'
-        #     )
-        #     ting_ax.plot(t, np.cumsum([me[0] for me in dat['dissipation']]) / 2, label=label)
-        #     ting_ax.plot(t, [me[1] for me in dat['dissipation']], label=label, ls='--')
-        #     ting_ax.legend()
-
-    Nu_ax.legend(frameon=True)
-    Nu_ax.set_xlabel('t')
-    Nu_ax.set_ylabel('Nu')
+    # # ----------------- RK ------------------------
 
-    spectrum_ax.legend(frameon=False)
-    spectrum_ax.set_xlabel('$k$')
-    spectrum_ax.set_ylabel(r'$\|\hat{u}_x\|$')
+    # for dt in [0.05, 0.04, 0.02, 0.01, 0.005]:
+    #     data.append(get_pySDC_data(Ra, res=res, dt=dt, RK=True, config_name='RBC3DG4R4'))
+    #     labels.append(f'RK, dt={dt:.3f}')
+    #     linestyles.append('--')
+
+    # # ----------------- Euler ---------------------
+
+    # for dt in [0.02, 0.01, 0.005]:
+    #     data.append(get_pySDC_data(Ra, res=res, dt=dt, config_name='RBC3DG4R4Euler'))
+    #     labels.append(f'Euler, dt={dt:.3f}')
+    #     linestyles.append('-.')
+
+    # for dat, label, linestyle in zip(data, labels, linestyles):
+    #     Nu_ax = _Nu_axs[label[:label.index('dt')-2]]
+    #     Nu = np.array(dat['Nu']['V'])
+    #     t = dat['t']
+    #     Nu_ref = np.array(ref_data['Nu']['V'])
+    #     t_i, Nu_i = interpolate_NuV_to_reference_times(dat, ref_data)
+    #     Nu_ax.plot(t, Nu, label=label, ls=linestyle)
+    #     dt = float(label[label.index('dt')+3:])
+
+    #     error = np.maximum.accumulate(np.abs(Nu_ref[: Nu_i.shape[0]] - Nu_i) / np.abs(Nu_ref[: Nu_i.shape[0]]))
+
+    #     # compute mean Nu
+    #     mask = np.logical_and(t >= 20, t <= 200)
+    #     Nu_mean = np.mean(Nu[mask])
+    #     Nu_std = np.std(Nu[mask])
+
+    #     for key in Nu_all.keys():
+    #         if str(key) in label:
+    #             Nu_all[key].append((Nu_mean, Nu_std, dt))
+    #             break
+
+    #     last_line = Nu_ax.get_lines()[-1]
+    #     if any(error > 1e-2):
+    #         deviates = min(t_i[error > 1e-2])
+    #         for key in deviations.keys():
+    #             if str(key) in label:
+    #                 deviations[key].append((dt, deviates))
+    #                 break
+    #         Nu_ax.axvline(deviates, color=last_line.get_color(), ls=last_line.get_linestyle())
+    #         print(f'{label} Nu={Nu_mean:.3f}+={Nu_std:.3f}, deviates more than 1% from t={deviates:.2f}')
+    #     else:
+    #         print(f'{label} Nu={Nu_mean:.3f}+={Nu_std:.3f}')
+
+    #     k = dat['k']
+    #     spectrum = np.array(dat['spectrum'])
+    #     u_spectrum = np.mean(spectrum, axis=0)[1]
+    #     idx = dat['res_in_boundary_layer']
+    #     _s = u_spectrum[idx]
+    #     spectrum_ax.loglog(
+    #         k[_s > 1e-16], _s[_s > 1e-16], color=last_line.get_color(), ls=last_line.get_linestyle(), label=label
+    #     )
+
+    #     prof_ax.plot(dat['profile_T'], dat['z'], color=last_line.get_color(), ls=last_line.get_linestyle(), label=label)
+    #     rms_ax.plot(
+    #         dat['rms_profile_T'], dat['z'], color=last_line.get_color(), ls=last_line.get_linestyle(), label=label
+    #     )
+
+    # for key, value in deviations.items():
+    #     deviation_ax.plot([me[0] for me in value], [me[1] for me in value], label=key, marker='.')
+    # deviation_ax.set_xscale('log')
+    # deviation_ax.legend(frameon=False)
+    # deviation_ax.set_xlabel(r'$\Delta t$')
+    # deviation_ax.set_ylabel('Nu deviates from reference')
+
+    # _, Nu_summary_ax = plt.subplots()
+    # for key in list(Nu_all.keys())[::-1]:
+    #     value = Nu_all[key]
+    #     Nu_summary_ax.errorbar([me[2] for me in value], [me[0] for me in value], yerr=[me[1] for me in value], label=key, marker='.')
+    # Nu_summary_ax.legend(frameon=False)
+
+    # spectrum_ax.legend(frameon=False)
+    # spectrum_ax.set_xlabel('$k$')
+    # spectrum_ax.set_ylabel(r'$\|\hat{u}_x\|$')
+
+    Nu_axs[-1].set_xlabel('$t$')
+    Nu_axs[-1].set_ylabel('$Nu$')
+
+    Nu_fig.tight_layout()
+    Nu_fig.savefig('./plots/Nu_over_time_Ra1e5.pdf', bbox_inches='tight')
+    # deviation_fig.savefig('./plots/Nu_deviation_Ra1e5.pdf', bbox_inches='tight')
 
 
 def compare_Nusselt_over_time1e6():