From 76724addf4d2e6f63625cce88e8fd83b5815eb9e Mon Sep 17 00:00:00 2001
From: Thomas Baumann <39156931+brownbaerchen@users.noreply.github.com>
Date: Mon, 10 Nov 2025 16:56:48 +0100
Subject: [PATCH 1/8] Added order script

---
 .../GPU/analysis_scripts/RBC3D_order.py       | 134 ++++++++++++++++++
 1 file changed, 134 insertions(+)
 create mode 100644 pySDC/projects/GPU/analysis_scripts/RBC3D_order.py

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
new file mode 100644
index 0000000000..124505b8a6
--- /dev/null
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
@@ -0,0 +1,134 @@
+import os
+import pickle
+import numpy as np
+from pySDC.helpers.fieldsIO import FieldsIO
+from pySDC.projects.GPU.configs.base_config import get_config
+from pySDC.implementations.problem_classes.RayleighBenard3D import RayleighBenard3D
+from mpi4py import MPI
+import matplotlib.pyplot as plt
+
+step_sizes = {
+    'RBC3DG4R4Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
+    'RBC3DG4R4SDC23Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3, 2.5e-3],
+    'RBC3DG4R4SDC34Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
+    'RBC3DG4R4RKRa1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3, 2.5e-3],
+    # 'RBC3DG4R4EulerRa1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
+    'RBC3DG4R4EulerRa1e5': [1e-3 * 2**i for i in range(5)],
+}
+n_freefall_times = {}
+
+
+def no_logging_hook(*args, **kwargs):
+    return None
+
+
+def get_path(args):
+    config = get_config(args)
+    fname = config.get_file_name()
+    return f'{fname[:fname.index('dt')]}order.pickle'
+
+
+def compute_errors(args, dts, Tend):
+    errors = {'u': [], 'v': [], 'w': [], 'T': [], 'p': [], 'dt': []}
+    prob = RayleighBenard3D(nx=4, ny=4, nz=4)
+
+    dts = np.sort(dts)[::-1]
+    ref = run(args, dts[-1], Tend)
+    for dt in dts[:-1]:
+        u = run(args, dt, Tend)
+        e = u - ref
+        for comp in ['u', 'v', 'w', 'T', 'p']:
+            i = prob.index(comp)
+            e_comp = np.max(np.abs(e[i]))
+            e_comp = MPI.COMM_WORLD.allreduce(e_comp, op=MPI.MAX)
+            errors[comp].append(e_comp)
+        errors['dt'].append(dt)
+
+    path = get_path(args)
+    if MPI.COMM_WORLD.rank == 0:
+        with open(path, 'wb') as file:
+            pickle.dump(errors, file)
+            print(f'Saved errors to {path}', flush=True)
+
+
+def plot_error_all_components(args):
+    fig, ax = plt.subplots()
+    with open(get_path(args), 'rb') as file:
+        errors = pickle.load(file)
+
+    for comp in ['u', 'v', 'w', 'T', 'p']:
+        e = np.array(errors[comp])
+        dt = np.array(errors['dt'])
+        order = np.log(e[1:] / e[:-1]) / np.log(dt[1:] / dt[:-1])
+        ax.loglog(errors['dt'], errors[comp], label=f'{comp} order {np.mean(order):.1f}')
+
+    ax.loglog(errors['dt'], np.array(errors['dt']) ** 4, label='Order 4', ls='--')
+    ax.loglog(errors['dt'], np.array(errors['dt']) ** 2, label='Order 2', ls='--')
+    ax.legend()
+    ax.set_xlabel(r'$\Delta t$')
+    ax.set_ylabel(r'$e$')
+
+
+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', 'SDC23']]
+        labels = ['SDC', 'RK', 'Euler', 'SDC']
+
+    else:
+        raise NotImplementedError
+
+    for path, label in zip(paths, labels, strict=True):
+        with open(path, 'rb') as file:
+            errors = pickle.load(file)
+
+        e = np.array(errors['T'])
+        dt = np.array(errors['dt'])
+        order = np.log(e[1:] / e[:-1]) / np.log(dt[1:] / dt[:-1])
+        ax.loglog(dt, e, label=f'{label} order {np.mean(order):.1f}', ls=ls[label])
+
+    ax.legend(frameon=False)
+    ax.set_xlabel(r'$\Delta t$')
+    ax.set_ylabel(r'$e$')
+
+
+def run(args, dt, Tend):
+    from pySDC.projects.GPU.run_experiment import run_experiment
+    from pySDC.core.errors import ConvergenceError
+
+    args['mode'] = 'run'
+    args['dt'] = dt
+
+    config = get_config(args)
+    config.Tend = n_freefall_times.get(type(config).__name__, 3)
+
+    ic_config_name = type(config).__name__
+    for name in ['RK', 'Euler', 'O3', 'O4', 'SDC23', 'SDC34']:
+        ic_config_name = ic_config_name.replace(name, '')
+    ic_config = get_config({**args, 'config': ic_config_name})
+    config.ic_config['config'] = type(ic_config)
+    config.ic_config['res'] = ic_config.res
+    config.ic_config['dt'] = ic_config.dt
+
+    config.get_LogToFile = no_logging_hook
+    config.Tend = Tend
+
+    u = run_experiment(args, config)
+    return u
+
+
+if __name__ == '__main__':
+    from pySDC.projects.GPU.run_experiment import parse_args
+
+    args = parse_args()
+    config = get_config(args)
+
+    dts = step_sizes[type(config).__name__]
+    if args['mode'] == 'run':
+        compute_errors(args, dts, max(dts))
+
+    plot_error_all_components(args)
+    compare_order(1e5)
+    if MPI.COMM_WORLD.rank == 0:
+        plt.show()

From b54285dd814c21eab6347f8f7f525400c9906466 Mon Sep 17 00:00:00 2001
From: Thomas Baumann <t.baumann@fz-juelich.de>
Date: Tue, 11 Nov 2025 16:56:01 +0100
Subject: [PATCH 2/8] Fixes

---
 .../projects/GPU/analysis_scripts/RBC3D_order.py | 12 ++++++++----
 pySDC/projects/GPU/configs/RBC3D_configs.py      | 16 ++++++++++++++++
 2 files changed, 24 insertions(+), 4 deletions(-)

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
index 124505b8a6..b5a0a0ded1 100644
--- a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
@@ -9,8 +9,12 @@
 
 step_sizes = {
     'RBC3DG4R4Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
-    'RBC3DG4R4SDC23Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3, 2.5e-3],
-    'RBC3DG4R4SDC34Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
+    # 'RBC3DG4R4SDC23Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3, 2.5e-3],
+    # 'RBC3DG4R4SDC34Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
+    'RBC3DG4R4SDC23Ra1e5': [5e-3 * 2**i for i in range(5)],
+    'RBC3DG4R4SDC44Ra1e5': [5e-3 * 2**i for i in range(5)],
+    # 'RBC3DG4R4SDC44Ra1e5': [1e-3 * 2**i for i in range(8)],
+    'RBC3DG4R4SDC34Ra1e5': [1e-3 * 2**i for i in range(10)],
     'RBC3DG4R4RKRa1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3, 2.5e-3],
     # 'RBC3DG4R4EulerRa1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
     'RBC3DG4R4EulerRa1e5': [1e-3 * 2**i for i in range(5)],
@@ -30,7 +34,7 @@ def get_path(args):
 
 def compute_errors(args, dts, Tend):
     errors = {'u': [], 'v': [], 'w': [], 'T': [], 'p': [], 'dt': []}
-    prob = RayleighBenard3D(nx=4, ny=4, nz=4)
+    prob = RayleighBenard3D(nx=4, ny=4, nz=4, comm=MPI.COMM_SELF)
 
     dts = np.sort(dts)[::-1]
     ref = run(args, dts[-1], Tend)
@@ -104,7 +108,7 @@ def run(args, dt, Tend):
     config.Tend = n_freefall_times.get(type(config).__name__, 3)
 
     ic_config_name = type(config).__name__
-    for name in ['RK', 'Euler', 'O3', 'O4', 'SDC23', 'SDC34']:
+    for name in ['RK', 'Euler', 'O3', 'O4', 'SDC23', 'SDC34', 'SDC44']:
         ic_config_name = ic_config_name.replace(name, '')
     ic_config = get_config({**args, 'config': ic_config_name})
     config.ic_config['config'] = type(ic_config)
diff --git a/pySDC/projects/GPU/configs/RBC3D_configs.py b/pySDC/projects/GPU/configs/RBC3D_configs.py
index e12eb2e192..064ed60ea9 100644
--- a/pySDC/projects/GPU/configs/RBC3D_configs.py
+++ b/pySDC/projects/GPU/configs/RBC3D_configs.py
@@ -221,6 +221,15 @@ def get_description(self, *args, **kwargs):
         return desc
 
 
+class RBC3DM4K4(RBC3Dverification):
+
+    def get_description(self, *args, **kwargs):
+        desc = super().get_description(*args, **kwargs)
+        desc['level_params']['nsweeps'] = 4
+        desc['sweeper_params']['num_nodes'] = 4
+        return desc
+
+
 class RBC3DverificationRK(RBC3Dverification):
 
     def get_description(self, *args, res=-1, dt=-1, **kwargs):
@@ -268,6 +277,13 @@ class RBC3DG4R4SDC34Ra1e5(RBC3DM3K4):
     converged = 50
 
 
+class RBC3DG4R4SDC44Ra1e5(RBC3DM4K4):
+    Tend = 200
+    dt = 6e-2
+    res = 32
+    converged = 50
+
+
 class RBC3DG4R4RKRa1e5(RBC3DverificationRK):
     Tend = 200
     dt = 8e-2

From 5a0d1a7f2b88c071d595ec56be68774c236671af Mon Sep 17 00:00:00 2001
From: Thomas Baumann <39156931+brownbaerchen@users.noreply.github.com>
Date: Wed, 12 Nov 2025 16:34:56 +0100
Subject: [PATCH 3/8] Fixed plotting

---
 .../GPU/analysis_scripts/RBC3D_order.py       | 30 ++++++++++++++-----
 1 file changed, 23 insertions(+), 7 deletions(-)

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
index b5a0a0ded1..8b177606ec 100644
--- a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
@@ -6,6 +6,7 @@
 from pySDC.implementations.problem_classes.RayleighBenard3D import RayleighBenard3D
 from mpi4py import MPI
 import matplotlib.pyplot as plt
+from pySDC.helpers.plot_helper import setup_mpl, figsize_by_journal
 
 step_sizes = {
     'RBC3DG4R4Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
@@ -56,6 +57,7 @@ def compute_errors(args, dts, Tend):
 
 
 def plot_error_all_components(args):
+    setup_mpl()
     fig, ax = plt.subplots()
     with open(get_path(args), 'rb') as file:
         errors = pickle.load(file)
@@ -74,27 +76,37 @@ def plot_error_all_components(args):
 
 
 def compare_order(Ra):
-    fig, ax = plt.subplots()
-    ls = {'SDC': '-', 'RK': '--', 'Euler': '-.'}
+    setup_mpl()
+    fig, ax = plt.subplots(figsize=figsize_by_journal('Nature_CS', 1, 0.6))
+    ls = {'SD': '-', 'RK': '--', 'Eu': '-.'}
     if Ra == 1e5:
-        paths = [f'./data/RBC3DG4R4{me}Ra1e5-res-1-order.pickle' for me in ['', 'RK', 'Euler', 'SDC23']]
-        labels = ['SDC', 'RK', 'Euler', 'SDC']
+        labels = ['RK443', 'Euler', 'SDC23', 'SDC34', 'SDC44'][::-1]
+        names = ['RK', 'Euler', 'SDC23', 'SDC34', 'SDC44'][::-1]
+        markers = ['>', '.', 'o', '<', 'x'][::-1]
+        paths = [f'./data/RBC3DG4R4{me}Ra1e5-res-1-order.pickle' for me in names]
 
     else:
         raise NotImplementedError
 
-    for path, label in zip(paths, labels, strict=True):
+    for path, label, marker in zip(paths, labels, markers, strict=True):
         with open(path, 'rb') as file:
             errors = pickle.load(file)
 
         e = np.array(errors['T'])
         dt = np.array(errors['dt'])
         order = np.log(e[1:] / e[:-1]) / np.log(dt[1:] / dt[:-1])
-        ax.loglog(dt, e, label=f'{label} order {np.mean(order):.1f}', ls=ls[label])
+        print(f'{label}: order: mean={np.mean(order):.1f} / median={np.median(order):.1f}')
+        ax.loglog(dt, e, label=f'{label}', ls=ls[label[:2]], marker=marker, markersize=6)
+
+    for _dt in dt:
+        for i in [1, 3, 4]:
+            ax.text(_dt, _dt**i, i, fontweight='bold', fontsize=14, ha='center', va='center')
+            ax.loglog(dt, dt**i, ls=':', color='black')
 
     ax.legend(frameon=False)
     ax.set_xlabel(r'$\Delta t$')
     ax.set_ylabel(r'$e$')
+    fig.savefig('plots/SDC_order_Ra1e5.pdf', bbox_inches='tight')
 
 
 def run(args, dt, Tend):
@@ -107,6 +119,9 @@ def run(args, dt, Tend):
     config = get_config(args)
     config.Tend = n_freefall_times.get(type(config).__name__, 3)
 
+    desc = config.get_description()
+    prob = desc['problem_class'](**desc['problem_params'])
+
     ic_config_name = type(config).__name__
     for name in ['RK', 'Euler', 'O3', 'O4', 'SDC23', 'SDC34', 'SDC44']:
         ic_config_name = ic_config_name.replace(name, '')
@@ -118,7 +133,8 @@ def run(args, dt, Tend):
     config.get_LogToFile = no_logging_hook
     config.Tend = Tend
 
-    u = run_experiment(args, config)
+    u_hat = run_experiment(args, config)
+    u = prob.itransform(u_hat)
     return u
 
 

From 8773225828408b4e615ec62772ac2abf53214909 Mon Sep 17 00:00:00 2001
From: Thomas Baumann <t.baumann@fz-juelich.de>
Date: Wed, 12 Nov 2025 16:36:39 +0100
Subject: [PATCH 4/8] Fixed recording

---
 .../GPU/analysis_scripts/RBC3D_order.py         | 17 +++++++----------
 1 file changed, 7 insertions(+), 10 deletions(-)

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
index 8b177606ec..5723b31526 100644
--- a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
@@ -10,15 +10,11 @@
 
 step_sizes = {
     'RBC3DG4R4Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
-    # 'RBC3DG4R4SDC23Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3, 2.5e-3],
-    # 'RBC3DG4R4SDC34Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
-    'RBC3DG4R4SDC23Ra1e5': [5e-3 * 2**i for i in range(5)],
-    'RBC3DG4R4SDC44Ra1e5': [5e-3 * 2**i for i in range(5)],
-    # 'RBC3DG4R4SDC44Ra1e5': [1e-3 * 2**i for i in range(8)],
-    'RBC3DG4R4SDC34Ra1e5': [1e-3 * 2**i for i in range(10)],
-    'RBC3DG4R4RKRa1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3, 2.5e-3],
-    # 'RBC3DG4R4EulerRa1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
-    'RBC3DG4R4EulerRa1e5': [1e-3 * 2**i for i in range(5)],
+    'RBC3DG4R4SDC23Ra1e5': [5e-3 * 2**i for i in range(8)],
+    'RBC3DG4R4SDC34Ra1e5': [5e-3 * 2**i for i in range(8)],
+    'RBC3DG4R4SDC44Ra1e5': [5e-3 * 2**i for i in range(8)],
+    'RBC3DG4R4RKRa1e5': [5e-3 * 2**i for i in range(8)],
+    'RBC3DG4R4EulerRa1e5': [5e-3 * 2**i for i in range(8)],
 }
 n_freefall_times = {}
 
@@ -44,7 +40,7 @@ def compute_errors(args, dts, Tend):
         e = u - ref
         for comp in ['u', 'v', 'w', 'T', 'p']:
             i = prob.index(comp)
-            e_comp = np.max(np.abs(e[i]))
+            e_comp = np.max(np.abs(e[i])) / np.max(np.abs(ref[i]))
             e_comp = MPI.COMM_WORLD.allreduce(e_comp, op=MPI.MAX)
             errors[comp].append(e_comp)
         errors['dt'].append(dt)
@@ -135,6 +131,7 @@ def run(args, dt, Tend):
 
     u_hat = run_experiment(args, config)
     u = prob.itransform(u_hat)
+
     return u
 
 

From d19a4ee8963b201ef9a59f7ab1ac6e42f9059dee Mon Sep 17 00:00:00 2001
From: Thomas Baumann <39156931+brownbaerchen@users.noreply.github.com>
Date: Wed, 12 Nov 2025 18:02:03 +0100
Subject: [PATCH 5/8] Added order plot for 3D RBC

---
 .../GPU/analysis_scripts/RBC3D_order.py       | 11 ++++----
 .../GPU/tests/test_RBC_3D_analysis.py         | 25 +++++++++++++++++--
 2 files changed, 29 insertions(+), 7 deletions(-)

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
index 5723b31526..582eb09b96 100644
--- a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
@@ -42,7 +42,7 @@ def compute_errors(args, dts, Tend):
             i = prob.index(comp)
             e_comp = np.max(np.abs(e[i])) / np.max(np.abs(ref[i]))
             e_comp = MPI.COMM_WORLD.allreduce(e_comp, op=MPI.MAX)
-            errors[comp].append(e_comp)
+            errors[comp].append(float(e_comp))
         errors['dt'].append(dt)
 
     path = get_path(args)
@@ -52,7 +52,7 @@ def compute_errors(args, dts, Tend):
             print(f'Saved errors to {path}', flush=True)
 
 
-def plot_error_all_components(args):
+def plot_error_all_components(args):  # pragma: no cover
     setup_mpl()
     fig, ax = plt.subplots()
     with open(get_path(args), 'rb') as file:
@@ -71,7 +71,7 @@ def plot_error_all_components(args):
     ax.set_ylabel(r'$e$')
 
 
-def compare_order(Ra):
+def compare_order(Ra):  # pragma: no cover
     setup_mpl()
     fig, ax = plt.subplots(figsize=figsize_by_journal('Nature_CS', 1, 0.6))
     ls = {'SD': '-', 'RK': '--', 'Eu': '-.'}
@@ -115,12 +115,13 @@ def run(args, dt, Tend):
     config = get_config(args)
     config.Tend = n_freefall_times.get(type(config).__name__, 3)
 
-    desc = config.get_description()
+    desc = config.get_description(res=args['res'])
     prob = desc['problem_class'](**desc['problem_params'])
 
     ic_config_name = type(config).__name__
     for name in ['RK', 'Euler', 'O3', 'O4', 'SDC23', 'SDC34', 'SDC44']:
-        ic_config_name = ic_config_name.replace(name, '')
+        ic_config_name = ic_config_name.replace(name, 'SDC34')
+
     ic_config = get_config({**args, 'config': ic_config_name})
     config.ic_config['config'] = type(ic_config)
     config.ic_config['res'] = ic_config.res
diff --git a/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py b/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
index 696f5840e1..847f3b0dbb 100644
--- a/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
+++ b/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
@@ -26,10 +26,10 @@ def get_config(args):
     return config
 
 
-def generate_simulation_file(path):
+def generate_simulation_file(path, args=None):
     from pySDC.projects.GPU.run_experiment import run_experiment
 
-    args = get_args(path)
+    args = {**get_args(path), **args} if args is not None else get_args(path)
     config = get_config(args)
 
     run_experiment(args, config)
@@ -112,3 +112,24 @@ def _get_Nu(_t):
     tI, NuI = interpolate_NuV_to_reference_times(data, ref_data, order=4)
     assert not np.allclose(NuI, ref_data['Nu']['V'])
     assert np.allclose(tI, ref_data['t'])
+
+
+def test_order_computation(tmp_sim_data, tmp_path):
+    from pySDC.projects.GPU.analysis_scripts.RBC3D_order import compute_errors, get_path
+    from pySDC.projects.GPU.configs.RBC3D_configs import RBC3DG4R4SDC34Ra1e5
+    import numpy as np
+    import pickle
+
+    args = get_args(tmp_path)
+
+    generate_simulation_file(tmp_path, args={'config': args['config'].replace('SDC23', 'SDC34')})
+    RBC3DG4R4SDC34Ra1e5.res = args['res']
+    RBC3DG4R4SDC34Ra1e5.dt = args['dt']
+
+    dts = [1e-2, 5e-4]
+    compute_errors(args, dts, np.max(dts))
+
+    with open(get_path(args), 'rb') as file:
+        errors = pickle.load(file)
+    for comp in ['T', 'p']:
+        assert max(errors[comp]) < 1e-12

From 6cbb58beb6ecce6a80fee7e9ffd7d869a8bfc550 Mon Sep 17 00:00:00 2001
From: Thomas Baumann <39156931+brownbaerchen@users.noreply.github.com>
Date: Thu, 13 Nov 2025 15:18:45 +0100
Subject: [PATCH 6/8] Added spectrum plot

---
 .../GPU/analysis_scripts/RBC3D_order.py       | 25 +++++------
 .../analysis_scripts/RBC3D_plotting_utils.py  | 45 +++++++++++++++++++
 .../GPU/analysis_scripts/RBC3D_spectrum.py    | 33 ++++++++++++++
 .../projects/GPU/analysis_scripts/plot_Nu.py  |  9 +---
 .../analysis_scripts/process_RBC3D_data.py    |  8 ++++
 .../GPU/tests/test_RBC_3D_analysis.py         |  2 +-
 6 files changed, 100 insertions(+), 22 deletions(-)
 create mode 100644 pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py
 create mode 100644 pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
index 582eb09b96..9bae3c1121 100644
--- a/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_order.py
@@ -6,7 +6,8 @@
 from pySDC.implementations.problem_classes.RayleighBenard3D import RayleighBenard3D
 from mpi4py import MPI
 import matplotlib.pyplot as plt
-from pySDC.helpers.plot_helper import setup_mpl, figsize_by_journal
+from pySDC.helpers.plot_helper import figsize_by_journal
+from pySDC.projects.GPU.analysis_scripts.RBC3D_plotting_utils import get_plotting_style, savefig
 
 step_sizes = {
     'RBC3DG4R4Ra1e5': [8e-2, 4e-2, 2e-2, 1e-2, 5e-3],
@@ -53,7 +54,6 @@ def compute_errors(args, dts, Tend):
 
 
 def plot_error_all_components(args):  # pragma: no cover
-    setup_mpl()
     fig, ax = plt.subplots()
     with open(get_path(args), 'rb') as file:
         errors = pickle.load(file)
@@ -72,27 +72,24 @@ def plot_error_all_components(args):  # pragma: no cover
 
 
 def compare_order(Ra):  # pragma: no cover
-    setup_mpl()
     fig, ax = plt.subplots(figsize=figsize_by_journal('Nature_CS', 1, 0.6))
-    ls = {'SD': '-', 'RK': '--', 'Eu': '-.'}
     if Ra == 1e5:
-        labels = ['RK443', 'Euler', 'SDC23', 'SDC34', 'SDC44'][::-1]
         names = ['RK', 'Euler', 'SDC23', 'SDC34', 'SDC44'][::-1]
-        markers = ['>', '.', 'o', '<', 'x'][::-1]
+        configs = [f'RBC3DG4R4{me}Ra1e5' for me in names]
         paths = [f'./data/RBC3DG4R4{me}Ra1e5-res-1-order.pickle' for me in names]
 
     else:
         raise NotImplementedError
 
-    for path, label, marker in zip(paths, labels, markers, strict=True):
+    for path, config in zip(paths, configs, strict=True):
         with open(path, 'rb') as file:
             errors = pickle.load(file)
 
         e = np.array(errors['T'])
         dt = np.array(errors['dt'])
         order = np.log(e[1:] / e[:-1]) / np.log(dt[1:] / dt[:-1])
-        print(f'{label}: order: mean={np.mean(order):.1f} / median={np.median(order):.1f}')
-        ax.loglog(dt, e, label=f'{label}', ls=ls[label[:2]], marker=marker, markersize=6)
+        print(f'{config}: order: mean={np.mean(order):.1f} / median={np.median(order):.1f}')
+        ax.loglog(dt, e, **get_plotting_style(config))
 
     for _dt in dt:
         for i in [1, 3, 4]:
@@ -102,7 +99,7 @@ def compare_order(Ra):  # pragma: no cover
     ax.legend(frameon=False)
     ax.set_xlabel(r'$\Delta t$')
     ax.set_ylabel(r'$e$')
-    fig.savefig('plots/SDC_order_Ra1e5.pdf', bbox_inches='tight')
+    savefig(fig, 'RBC3D_order_Ra1e5')
 
 
 def run(args, dt, Tend):
@@ -140,13 +137,15 @@ def run(args, dt, Tend):
     from pySDC.projects.GPU.run_experiment import parse_args
 
     args = parse_args()
-    config = get_config(args)
 
-    dts = step_sizes[type(config).__name__]
     if args['mode'] == 'run':
+        config = get_config(args)
+        dts = step_sizes[type(config).__name__]
         compute_errors(args, dts, max(dts))
 
-    plot_error_all_components(args)
+    if args['config'] is not None:
+        plot_error_all_components(args)
+
     compare_order(1e5)
     if MPI.COMM_WORLD.rank == 0:
         plt.show()
diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py
new file mode 100644
index 0000000000..6b7886f894
--- /dev/null
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py
@@ -0,0 +1,45 @@
+from pySDC.helpers.plot_helper import figsize_by_journal, setup_mpl
+import warnings
+
+setup_mpl()
+
+
+def get_plotting_style(config):
+
+    args = {'color': None, 'ls': None, 'marker': None, 'markersize': 6, 'label': None}
+
+    if config == 'RBC3DG4R4SDC23Ra1e5':
+        args['color'] = 'tab:blue'
+        args['ls'] = '-'
+        args['marker'] = 'o'
+        args['label'] = 'SDC23'
+    elif config == 'RBC3DG4R4SDC34Ra1e5':
+        args['color'] = 'tab:orange'
+        args['ls'] = '-'
+        args['marker'] = '<'
+        args['label'] = 'SDC34'
+    elif config in ['RBC3DG4R4SDC44Ra1e5', 'RBC3DG4R4Ra1e5']:
+        args['color'] = 'tab:green'
+        args['ls'] = '-'
+        args['marker'] = 'x'
+        args['label'] = 'SDC44'
+    elif config == 'RBC3DG4R4EulerRa1e5':
+        args['color'] = 'tab:purple'
+        args['ls'] = '--'
+        args['marker'] = '.'
+        args['label'] = 'Euler'
+    elif config == 'RBC3DG4R4RKRa1e5':
+        args['color'] = 'tab:red'
+        args['ls'] = '--'
+        args['marker'] = '>'
+        args['label'] = 'RK444'
+    else:
+        warnings.warn(f'No plotting style for {config=!r}')
+
+    return args
+
+
+def savefig(fig, name, format='pdf', base_path='./plots', **kwargs):
+    path = f'{base_path}/{name}.{format}'
+    fig.savefig(path, bbox_inches='tight', **kwargs)
+    print(f'Saved figure {path!r}')
diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py
new file mode 100644
index 0000000000..9faaee634b
--- /dev/null
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py
@@ -0,0 +1,33 @@
+from pySDC.projects.GPU.analysis_scripts.process_RBC3D_data import get_pySDC_data
+from pySDC.projects.GPU.analysis_scripts.RBC3D_plotting_utils import figsize_by_journal, get_plotting_style, savefig
+import matplotlib.pyplot as plt
+
+
+def plot_spectrum(res, dt, config_name, ax):  # pragma: no cover
+    data = get_pySDC_data(res=res, dt=dt, config_name=config_name)
+
+    spectrum = data['avg_spectrum']
+    k = data['k']
+    ax.loglog(k[spectrum > 1e-16], spectrum[spectrum > 1e-16], **get_plotting_style(config_name), markevery=5)
+    ax.set_xlabel('$k$')
+    ax.set_ylabel(r'$\|\hat{u}_x\|$')
+
+
+def plot_spectra_Ra1e5():  # pragma: no cover
+    fig, ax = plt.subplots(figsize=figsize_by_journal('Nature_CS', 1, 0.6))
+
+    configs = [f'RBC3DG4R4{name}Ra1e5' for name in ['SDC34', 'SDC23', '', 'Euler', 'RK']]
+    dts = [0.06, 0.06, 0.06, 0.02, 0.04]
+    res = 32
+
+    for config, dt in zip(configs, dts):
+        plot_spectrum(res, dt, config, ax)
+
+    ax.legend(frameon=False)
+    savefig(fig, 'RBC3D_spectrum_Ra1e5')
+
+
+if __name__ == '__main__':
+    plot_spectra_Ra1e5()
+
+    plt.show()
diff --git a/pySDC/projects/GPU/analysis_scripts/plot_Nu.py b/pySDC/projects/GPU/analysis_scripts/plot_Nu.py
index f0b9ba0e62..b63cd8762e 100644
--- a/pySDC/projects/GPU/analysis_scripts/plot_Nu.py
+++ b/pySDC/projects/GPU/analysis_scripts/plot_Nu.py
@@ -3,18 +3,11 @@
 import numpy as np
 from scipy import integrate
 from pySDC.helpers.plot_helper import figsize_by_journal, setup_mpl
+from pySDC.projects.GPU.analysis_scripts.process_RBC3D_data import get_pySDC_data
 
 setup_mpl()
 
 
-def get_pySDC_data(res=-1, dt=-1, config_name='RBC3DG4', base_path='data/RBC_time_averaged'):
-    path = f'{base_path}/{config_name}-res{res}-dt{dt:.0e}.pickle'
-    with open(path, 'rb') as file:
-        data = pickle.load(file)
-
-    return data
-
-
 def interpolate_NuV_to_reference_times(data, reference_data, order=12):
     from qmat.lagrange import getSparseInterpolationMatrix
 
diff --git a/pySDC/projects/GPU/analysis_scripts/process_RBC3D_data.py b/pySDC/projects/GPU/analysis_scripts/process_RBC3D_data.py
index 5fb4e3a501..25e153207e 100644
--- a/pySDC/projects/GPU/analysis_scripts/process_RBC3D_data.py
+++ b/pySDC/projects/GPU/analysis_scripts/process_RBC3D_data.py
@@ -204,5 +204,13 @@ def process_RBC3D_data(base_path='./data/RBC_time_averaged', plot=True, args=Non
     return path
 
 
+def get_pySDC_data(res=-1, dt=-1, config_name='RBC3DG4', base_path='data/RBC_time_averaged'):
+    path = f'{base_path}/{config_name}-res{res}-dt{dt:.0e}.pickle'
+    with open(path, 'rb') as file:
+        data = pickle.load(file)
+
+    return data
+
+
 if __name__ == '__main__':
     process_RBC3D_data()
diff --git a/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py b/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
index 847f3b0dbb..739e2206d5 100644
--- a/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
+++ b/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
@@ -80,7 +80,7 @@ def test_processing(tmp_processed_data):
 
 
 def test_get_pySDC_data(tmp_processed_data, tmp_path):
-    from pySDC.projects.GPU.analysis_scripts.plot_Nu import get_pySDC_data
+    from pySDC.projects.GPU.analysis_scripts.process_RBC3D_data import get_pySDC_data
 
     args = get_args(tmp_path)
     data = get_pySDC_data(res=args['res'], dt=args['dt'], config_name=args['config'], base_path=tmp_path)

From 10efb2320b19c9177b6964f6e9628214c71dc156 Mon Sep 17 00:00:00 2001
From: Thomas Baumann <39156931+brownbaerchen@users.noreply.github.com>
Date: Thu, 13 Nov 2025 15:21:24 +0100
Subject: [PATCH 7/8] Refactor

---
 pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py | 4 ++--
 pySDC/projects/GPU/tests/test_RBC_3D_analysis.py            | 2 +-
 2 files changed, 3 insertions(+), 3 deletions(-)

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py
index 6b7886f894..6403c0524f 100644
--- a/pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_plotting_utils.py
@@ -4,7 +4,7 @@
 setup_mpl()
 
 
-def get_plotting_style(config):
+def get_plotting_style(config):  # pragma: no cover
 
     args = {'color': None, 'ls': None, 'marker': None, 'markersize': 6, 'label': None}
 
@@ -39,7 +39,7 @@ def get_plotting_style(config):
     return args
 
 
-def savefig(fig, name, format='pdf', base_path='./plots', **kwargs):
+def savefig(fig, name, format='pdf', base_path='./plots', **kwargs):  # pragma: no cover
     path = f'{base_path}/{name}.{format}'
     fig.savefig(path, bbox_inches='tight', **kwargs)
     print(f'Saved figure {path!r}')
diff --git a/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py b/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
index 739e2206d5..a005851f8a 100644
--- a/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
+++ b/pySDC/projects/GPU/tests/test_RBC_3D_analysis.py
@@ -114,7 +114,7 @@ def _get_Nu(_t):
     assert np.allclose(tI, ref_data['t'])
 
 
-def test_order_computation(tmp_sim_data, tmp_path):
+def test_error_computation(tmp_sim_data, tmp_path):
     from pySDC.projects.GPU.analysis_scripts.RBC3D_order import compute_errors, get_path
     from pySDC.projects.GPU.configs.RBC3D_configs import RBC3DG4R4SDC34Ra1e5
     import numpy as np

From 32d081e4351640d87f66881f4494a94250e54911 Mon Sep 17 00:00:00 2001
From: Thomas Baumann <39156931+brownbaerchen@users.noreply.github.com>
Date: Thu, 13 Nov 2025 15:22:21 +0100
Subject: [PATCH 8/8] Refactor

---
 pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py b/pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py
index 9faaee634b..921bb756c6 100644
--- a/pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py
+++ b/pySDC/projects/GPU/analysis_scripts/RBC3D_spectrum.py
@@ -20,7 +20,7 @@ def plot_spectra_Ra1e5():  # pragma: no cover
     dts = [0.06, 0.06, 0.06, 0.02, 0.04]
     res = 32
 
-    for config, dt in zip(configs, dts):
+    for config, dt in zip(configs, dts, strict=True):
         plot_spectrum(res, dt, config, ax)
 
     ax.legend(frameon=False)