From b974ba58841b65c5e8e4ce22bf5ad27f6ad75bd9 Mon Sep 17 00:00:00 2001
From: Sven Marcus <s.marcus@outlook.de>
Date: Mon, 12 Apr 2021 17:32:40 +0200
Subject: [PATCH] L2 Test as separate scripts for SLURM integration

---
 Python/SlurmTests/poiseuille/evaluator.py     | 20 +++++
 .../SlurmTests/poiseuille/result_collector.py | 73 ++++++++++++++++
 Python/SlurmTests/poiseuille/settings.py      | 26 ++++++
 .../poiseuille/simulation_runner.py           | 19 +++++
 Python/SlurmTests/poiseuille/slurm.job        | 25 ++++++
 Python/acousticscaling.py                     | 85 +++++++++++++++++++
 Python/poiseuille/simulation.py               | 12 +--
 Python/test_boundaryconditions.py             |  2 +-
 Python/test_kernel.py                         |  7 ++
 Python/vtk_utilities.py                       |  4 +-
 .../pythonbindings/src/submodules/kernel.cpp  |  2 +-
 11 files changed, 266 insertions(+), 9 deletions(-)
 create mode 100644 Python/SlurmTests/poiseuille/evaluator.py
 create mode 100644 Python/SlurmTests/poiseuille/result_collector.py
 create mode 100644 Python/SlurmTests/poiseuille/settings.py
 create mode 100644 Python/SlurmTests/poiseuille/simulation_runner.py
 create mode 100644 Python/SlurmTests/poiseuille/slurm.job
 create mode 100644 Python/acousticscaling.py

diff --git a/Python/SlurmTests/poiseuille/evaluator.py b/Python/SlurmTests/poiseuille/evaluator.py
new file mode 100644
index 000000000..74602846b
--- /dev/null
+++ b/Python/SlurmTests/poiseuille/evaluator.py
@@ -0,0 +1,20 @@
+import numpy as np
+import scipy.stats as stats
+import errors
+from SlurmTests.poiseuille.result_collector import collect_results
+from SlurmTests.poiseuille.settings import Scaling
+
+analytical_results, numerical_results = collect_results()
+normalized_l2_errors = [errors.normalized_l2_error(analytical, numerical)
+                        for analytical, numerical in zip(analytical_results, numerical_results)]
+
+nodes_in_x3_per_run = []
+for simulation_run in range(0, 3):
+    grid_params, _, _, _ = Scaling.configuration_for_scale_level(simulation_run)
+    nodes_in_x3_per_run.append(grid_params.number_of_nodes_per_direction[2])
+
+nodes_as_log = [np.log10(node) for node in nodes_in_x3_per_run]
+l2_norms_as_log = [np.log10(l2) for l2 in normalized_l2_errors]
+res = stats.linregress(nodes_as_log, l2_norms_as_log)
+
+assert res.slope <= -2, f"Expected slope of l2 error to be <= -2, but was {res.slope}"
diff --git a/Python/SlurmTests/poiseuille/result_collector.py b/Python/SlurmTests/poiseuille/result_collector.py
new file mode 100644
index 000000000..06efa481c
--- /dev/null
+++ b/Python/SlurmTests/poiseuille/result_collector.py
@@ -0,0 +1,73 @@
+from typing import Collection, List
+
+import pyvista as pv
+from poiseuille.analytical import PoiseuilleSettings, poiseuille_at_heights
+from vtk_utilities import vertical_column_from_mesh, get_values_from_indices
+from SlurmTests.poiseuille.settings import Scaling
+
+
+def get_output_file_name(output_folder, runtime_params):
+    timesteps = runtime_params.number_of_timesteps
+    file_name = f"{output_folder}/mq/mq{timesteps}/mq0_{timesteps}.bin.vtu"
+
+    return file_name
+
+
+def get_mesh_for_last_timestep(output_folder, runtime_params):
+    file_name_of_last_timestep = get_output_file_name(output_folder, runtime_params)
+    mesh_of_last_timestep = pv.read(file_name_of_last_timestep)
+    return mesh_of_last_timestep
+
+
+def get_heights_from_indices(mesh, indices):
+    return [mesh.points[index][2] for index in indices]
+
+
+def get_heights(output_folder, runtime_params):
+    mesh_of_last_timestep = get_mesh_for_last_timestep(output_folder, runtime_params)
+    column_indices = vertical_column_from_mesh(mesh_of_last_timestep)
+    heights = get_heights_from_indices(mesh_of_last_timestep, column_indices)
+    return heights
+
+
+def get_numerical_results(runtime_params, output_folder):
+    mesh_of_last_timestep = get_mesh_for_last_timestep(output_folder, runtime_params)
+    velocities_in_x_direction = mesh_of_last_timestep.get_array("Vx")
+    column_indices = vertical_column_from_mesh(mesh_of_last_timestep)
+    numerical_results = get_values_from_indices(velocities_in_x_direction, column_indices)
+
+    return numerical_results
+
+
+def get_analytical_results(grid_params, physical_params, kernel, height_values):
+    channel_height = grid_params.number_of_nodes_per_direction[2]
+    settings = get_analytical_poiseuille_settings(channel_height, physical_params, kernel)
+    max_grid_height = channel_height * grid_params.node_distance
+    adjusted_height_values = [value / max_grid_height * channel_height for value in height_values]
+    analytical_results = poiseuille_at_heights(settings, adjusted_height_values)
+    return analytical_results
+
+
+def get_analytical_poiseuille_settings(height, physical_params, kernel):
+    settings = PoiseuilleSettings()
+    settings.height = height
+    settings.viscosity = physical_params.lattice_viscosity
+    settings.density = 1
+    settings.force = kernel.forcing_in_x1
+
+    return settings
+
+
+def collect_results() -> (List[List[float]], List[List[float]]):
+    analytical_results = []
+    numerical_results = []
+
+    for simulation_run in range(0, 3):
+        output_folder = f"output-{simulation_run}"
+        grid_params, physical_params, runtime_params, kernel = Scaling.configuration_for_scale_level(simulation_run)
+        heights = get_heights(output_folder, runtime_params)
+        analytical_results.append(
+            get_analytical_results(grid_params, physical_params, kernel, heights))
+        numerical_results.append(get_numerical_results(runtime_params, output_folder))
+
+    return analytical_results, numerical_results
diff --git a/Python/SlurmTests/poiseuille/settings.py b/Python/SlurmTests/poiseuille/settings.py
new file mode 100644
index 000000000..7b2f13b83
--- /dev/null
+++ b/Python/SlurmTests/poiseuille/settings.py
@@ -0,0 +1,26 @@
+import os
+from acousticscaling import AcousticScaling
+from pyfluids.kernel import LBMKernel, KernelType
+from pyfluids.parameters import RuntimeParameters, GridParameters, PhysicalParameters
+
+
+grid_params = GridParameters()
+grid_params.node_distance = 1
+grid_params.number_of_nodes_per_direction = [1, 1, 16]
+grid_params.blocks_per_direction = [1, 1, 8]
+grid_params.periodic_boundary_in_x1 = True
+grid_params.periodic_boundary_in_x2 = True
+
+physical_params = PhysicalParameters()
+physical_params.lattice_viscosity = 1e-3
+
+runtime_params = RuntimeParameters()
+runtime_params.number_of_threads = int(os.environ["SLURM_NTASKS"])
+runtime_params.number_of_timesteps = 2_500_000
+runtime_params.timestep_log_interval = 500_000
+
+kernel = LBMKernel(KernelType.CompressibleCumulantFourthOrderViscosity)
+kernel.use_forcing = True
+kernel.forcing_in_x1 = 1e-9
+
+Scaling = AcousticScaling(grid_params, physical_params, runtime_params, kernel)
diff --git a/Python/SlurmTests/poiseuille/simulation_runner.py b/Python/SlurmTests/poiseuille/simulation_runner.py
new file mode 100644
index 000000000..0b75de40b
--- /dev/null
+++ b/Python/SlurmTests/poiseuille/simulation_runner.py
@@ -0,0 +1,19 @@
+import os
+
+from SlurmTests.poiseuille.settings import Scaling
+from poiseuille.simulation import run_simulation
+from pyfluids.writer import Writer, OutputFormat
+
+
+scale_level = int(os.environ["PYFLUIDS_SCALE_LEVEL"])
+grid_params, physical_params, runtime_params, kernel = Scaling.configuration_for_scale_level(scale_level)
+
+writer = Writer()
+writer.output_format = OutputFormat.BINARY
+writer.output_path = "./output-" + str(scale_level)
+
+run_simulation(grid_params=grid_params,
+               physical_params=physical_params,
+               runtime_params=runtime_params,
+               kernel=kernel,
+               writer=writer)
diff --git a/Python/SlurmTests/poiseuille/slurm.job b/Python/SlurmTests/poiseuille/slurm.job
new file mode 100644
index 000000000..1d7072123
--- /dev/null
+++ b/Python/SlurmTests/poiseuille/slurm.job
@@ -0,0 +1,25 @@
+#!/bin/bash
+#SBATCH -J PyFluidsTest
+#SBATCH --nodes=2
+#SBATCH --ntasks-per-node=1
+#SBATCH --cpus-per-task=8
+
+#SBATCH --mem-per-cpu=30000
+#SBATCH --time=02:00:00
+#SBATCH --partition=shortrun_small
+
+source .bashrc
+
+echo "PyFluids Poiseuille Test Case"
+echo "Number of tasks: ${SLURM_NTASKS}"
+
+export PYFLUIDS_SCALE_LEVEL=0
+mpiexec -np "$SLURM_NTASKS" singularity exec SingularityImages/VirtualFluidsMPICH.sif python3 /VirtualFluids_dev/Python/poiseuille/simulation_runner.py
+
+export PYFLUIDS_SCALE_LEVEL=1
+mpiexec -np "$SLURM_NTASKS" singularity exec SingularityImages/VirtualFluidsMPICH.sif python3 /VirtualFluids_dev/Python/poiseuille/simulation_runner.py
+
+export PYFLUIDS_SCALE_LEVEL=2
+mpiexec -np "$SLURM_NTASKS" singularity exec SingularityImages/VirtualFluidsMPICH.sif python3 /VirtualFluids_dev/Python/poiseuille/simulation_runner.py
+
+singularity exec SingularityImages/VirtualFluidsMPICH.sif python3 /VirtualFluids_dev/Python/poiseuille/SlurmTests/poiseuille/evaluator.py
diff --git a/Python/acousticscaling.py b/Python/acousticscaling.py
new file mode 100644
index 000000000..ff10cbd49
--- /dev/null
+++ b/Python/acousticscaling.py
@@ -0,0 +1,85 @@
+from pyfluids.kernel import LBMKernel
+from pyfluids.parameters import GridParameters, PhysicalParameters, RuntimeParameters
+
+
+class AcousticScaling:
+
+    def __init__(self, grid_parameters: GridParameters,
+                 physical_parameters: PhysicalParameters,
+                 runtime_parameters: RuntimeParameters,
+                 kernel: LBMKernel):
+        self._grid_params = grid_parameters
+        self._physical_params = physical_parameters
+        self._runtime_params = runtime_parameters
+        self._kernel = kernel
+
+    def configuration_for_scale_level(self, level: int = 1) -> (GridParameters,
+                                                                PhysicalParameters,
+                                                                RuntimeParameters,
+                                                                LBMKernel):
+        if level < 0:
+            raise ValueError("level must be >= 0")
+
+        grid_params = self.clone_grid_params_for_level(level)
+        physical_params = self.clone_physical_parameters(level)
+        runtime_params = self.clone_runtime_params_for_level(level)
+        kernel = self.clone_kernel_for_level(level)
+
+        return grid_params, physical_params, runtime_params, kernel
+
+    def clone_grid_params_for_level(self, level) -> GridParameters:
+        grid_params = GridParameters()
+        grid_params.reference_direction_index = self._grid_params.reference_direction_index
+        grid_params.periodic_boundary_in_x1 = self._grid_params.periodic_boundary_in_x1
+        grid_params.periodic_boundary_in_x2 = self._grid_params.periodic_boundary_in_x2
+        grid_params.periodic_boundary_in_x3 = self._grid_params.periodic_boundary_in_x3
+
+        grid_params.number_of_nodes_per_direction = list(self._grid_params.number_of_nodes_per_direction)
+        grid_params.blocks_per_direction = list(self._grid_params.blocks_per_direction)
+        grid_params.node_distance = self._grid_params.node_distance
+
+        if level > 0:
+            grid_params.node_distance /= (level * 2)
+            grid_params.number_of_nodes_per_direction = [grid_params.number_of_nodes_per_direction[0],
+                                                         grid_params.number_of_nodes_per_direction[1],
+                                                         grid_params.number_of_nodes_per_direction[2] * (level * 2)]
+
+            grid_params.blocks_per_direction = [grid_params.blocks_per_direction[0],
+                                                grid_params.blocks_per_direction[1],
+                                                grid_params.blocks_per_direction[2] * (level * 2)]
+
+        return grid_params
+
+    def clone_physical_parameters(self, level):
+        physical_params = PhysicalParameters()
+        physical_params.lattice_viscosity = self._physical_params.lattice_viscosity
+
+        if level > 0:
+            physical_params.lattice_viscosity /= (level * 2)
+
+        return physical_params
+
+    def clone_runtime_params_for_level(self, level):
+        runtime_params = RuntimeParameters()
+        runtime_params.number_of_timesteps = self._runtime_params.number_of_timesteps
+        runtime_params.number_of_threads = self._runtime_params.number_of_threads
+        runtime_params.timestep_log_interval = self._runtime_params.timestep_log_interval
+
+        if level > 0:
+            runtime_params.number_of_timesteps *= (level * 2)
+
+        return runtime_params
+
+    def clone_kernel_for_level(self, level):
+        kernel = LBMKernel(self._kernel.type)
+        kernel.use_forcing = self._kernel.use_forcing
+        kernel.forcing_in_x1 = self._kernel.forcing_in_x1
+        kernel.forcing_in_x2 = self._kernel.forcing_in_x2
+        kernel.forcing_in_x3 = self._kernel.forcing_in_x3
+
+        if level > 0:
+            kernel.forcing_in_x1 *= (level * 2)
+            kernel.forcing_in_x2 *= (level * 2)
+            kernel.forcing_in_x3 *= (level * 2)
+
+        return kernel
diff --git a/Python/poiseuille/simulation.py b/Python/poiseuille/simulation.py
index 898809fbb..31ceb1ab9 100644
--- a/Python/poiseuille/simulation.py
+++ b/Python/poiseuille/simulation.py
@@ -20,11 +20,14 @@ default_runtime_params.number_of_threads = 4
 default_runtime_params.number_of_timesteps = 10000
 default_runtime_params.timestep_log_interval = 1000
 
-
 default_kernel = LBMKernel(KernelType.CompressibleCumulantFourthOrderViscosity)
 default_kernel.use_forcing = True
 default_kernel.forcing_in_x1 = 1e-8
 
+default_writer = Writer()
+default_writer.output_path = "./output"
+default_writer.output_format = OutputFormat.BINARY
+
 
 default_kernel = LBMKernel(KernelType.CompressibleCumulantFourthOrderViscosity)
 default_kernel.use_forcing = True
@@ -34,13 +37,10 @@ default_kernel.forcing_in_x1 = 1e-8
 def run_simulation(physical_params=default_physical_params,
                    grid_params=default_grid_params,
                    runtime_params=default_runtime_params,
-                   kernel=default_kernel):
+                   kernel=default_kernel,
+                   writer=default_writer):
     simulation = Simulation()
 
-    writer = Writer()
-    writer.output_path = "./output"
-    writer.output_format = OutputFormat.BINARY
-
     simulation.set_kernel_config(kernel)
     simulation.set_physical_parameters(physical_params)
     simulation.set_grid_parameters(grid_params)
diff --git a/Python/test_boundaryconditions.py b/Python/test_boundaryconditions.py
index 4fd51680e..5a7d61f36 100644
--- a/Python/test_boundaryconditions.py
+++ b/Python/test_boundaryconditions.py
@@ -58,4 +58,4 @@ class BoundaryConditionsTest(unittest.TestCase):
         Should be able to create NonReflectingOutflow
         """
 
-        sut = NonReflectingOutflow()
\ No newline at end of file
+        sut = NonReflectingOutflow()
diff --git a/Python/test_kernel.py b/Python/test_kernel.py
index 5e5487903..b1016f153 100644
--- a/Python/test_kernel.py
+++ b/Python/test_kernel.py
@@ -51,3 +51,10 @@ class TestLBMKernel(unittest.TestCase):
         self.assertEqual(self.sut.forcing_in_x2, 8)
         self.assertEqual(self.sut.forcing_in_x3, 5)
 
+    def test_lbm_kernel__when_getting_type__should_equal_kernel_type_enum_value(self) -> None:
+        """
+        WHEN getting the kernel type IT should equal the corresponding KernelType enum value
+        """
+
+        actual = self.sut.type
+        self.assertEqual(KernelType.BGK, actual)
diff --git a/Python/vtk_utilities.py b/Python/vtk_utilities.py
index 55b83aed1..14d13a40d 100644
--- a/Python/vtk_utilities.py
+++ b/Python/vtk_utilities.py
@@ -1,4 +1,6 @@
 import math
+from typing import List
+
 import pyvista as pv
 
 
@@ -21,7 +23,7 @@ def vertical_column_from_mesh(mesh):
     return relevant_indices
 
 
-def get_values_from_indices(array, indices):
+def get_values_from_indices(array, indices) -> List[float]:
     return [array[index] for index in indices]
 
 
diff --git a/src/cpu/pythonbindings/src/submodules/kernel.cpp b/src/cpu/pythonbindings/src/submodules/kernel.cpp
index 2418aa1d0..fb2917906 100644
--- a/src/cpu/pythonbindings/src/submodules/kernel.cpp
+++ b/src/cpu/pythonbindings/src/submodules/kernel.cpp
@@ -16,9 +16,9 @@ namespace kernel
                 .value("CompressibleCumulantFourthOrderViscosity",
                        KernelFactory::COMPRESSIBLE_CUMULANT_4TH_ORDER_VISCOSITY);
 
-
         py::class_<LBMKernelConfiguration, std::shared_ptr<LBMKernelConfiguration>>(kernelModule, "LBMKernel")
                 .def(py::init<KernelFactory::KernelType>())
+                .def_readwrite("type", &LBMKernelConfiguration::kernelType)
                 .def_readwrite("use_forcing", &LBMKernelConfiguration::useForcing)
                 .def_readwrite("forcing_in_x1", &LBMKernelConfiguration::forcingX1)
                 .def_readwrite("forcing_in_x2", &LBMKernelConfiguration::forcingX2)
-- 
GitLab