Create test for the predictor corrector in RTheta geometry.

CHANGELOG.md

@@ -16,7 +16,6 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 - Allow a `DerivField` to be stored in a `MultipatchField`.
 - Add `DerivFieldOnPatch` and `IdxRangeSliceOnPatch` aliases for the `MultipatchField` in `types.hpp`.
 - Allow `min` and `max` from `math_tools.hpp` to be called from GPU.
-- Add a `periodic_strips_non_uniform_2d_9patches` geometry.
 
 ### Fixed
 

tests/geometryRTheta/CMakeLists.txt

@@ -3,10 +3,12 @@
 # Allow the possibility of deactivating slow polar poisson tests
 option(POISSON_2D_BUILD_TESTING "Build polar poisson tests" ON)
 
+
+add_subdirectory(advection_field_rtheta)
+add_subdirectory(advection_rtheta)
 if("${POISSON_2D_BUILD_TESTING}")
     add_subdirectory(polar_poisson)
 endif()
-add_subdirectory(spline_interpolator_rtheta)
 add_subdirectory(quadrature)
-add_subdirectory(advection_rtheta)
-add_subdirectory(advection_field_rtheta)
+add_subdirectory(spline_interpolator_rtheta)
+add_subdirectory(time_solver)

tests/geometryRTheta/time_solver/CMakeLists.txt

@@ -0,0 +1,22 @@
+# SPDX-License-Identifier: MIT
+
+find_package(Python3 REQUIRED COMPONENTS Interpreter)
+
+function (diocotron_test_executable PREDCORR_METHOD TIME_METHOD )
+    set(test_name "TestSimulationDiocotron_${PREDCORR_METHOD}_${TIME_METHOD}")
+
+    add_test(NAME "${test_name}"
+        COMMAND bash "${CMAKE_CURRENT_SOURCE_DIR}/test_diocotron.sh"
+            "${PROJECT_SOURCE_DIR}"
+            "$<TARGET_FILE:diocotron_${PREDCORR_METHOD}_${TIME_METHOD}>"
+            "$<TARGET_FILE:Python3::Interpreter>")
+
+    # The test should take less than 15 min.         
+    set_property(TEST "${test_name}" PROPERTY TIMEOUT 900000)
+    set_property(TEST "${test_name}" PROPERTY COST 100)
+endfunction()
+
+# Select the parameters of the simulation:
+diocotron_test_executable(PREDCORR EULER_METHOD)
+diocotron_test_executable(EXPLICIT_PREDCORR EULER_METHOD)
+diocotron_test_executable(IMPLICIT_PREDCORR EULER_METHOD)

tests/geometryRTheta/time_solver/growth_rate_test.py

@@ -0,0 +1,52 @@
+#!/usr/bin/python3
+
+# SPDX-License-Identifier: MIT
+
+"""
+File which tests whether the growth rate of the perturbation of a diocotron instabilities test case
+follows the predicted slope. 
+"""
+
+from argparse import ArgumentParser
+from pathlib import Path
+import numpy as np
+
+from gysdata import DiskStore
+
+
+if __name__ == '__main__':
+    parser = ArgumentParser(description="Check growth rate.")
+    parser.add_argument('data_dir',
+                        action='store',
+                        nargs='?',
+                        default=Path.cwd(),
+                        type=Path,
+                        help='Location of the results (output folder).')
+    args = parser.parse_args()
+
+    path_data_structure = Path('data_structure_RTheta.yaml')
+    ds = DiskStore(args.data_dir, data_structure=path_data_structure)
+
+    # Get initial data
+    rho_eq = np.array(ds['density_eq'])
+    phi_eq = np.array(ds['electrical_potential_eq'])
+
+    jacobian = np.array(ds["jacobian"])
+
+    T = float(ds["final_T"])
+    omega_Im = float(ds["slope"])
+
+    # Get the data at each time step
+    Time = np.array(ds["density"].coords["time"])
+    rho = np.array(ds['density'])
+    phi = np.array(ds['electrical_potential'])
+
+    # Compute norms
+    L2norms_rho = np.linalg.norm((rho - rho_eq[None,:,:])*abs(jacobian[None,:,:]), 2, axis=(1,2))
+    L2norms_phi = np.linalg.norm((phi - phi_eq[None,:,:])*abs(jacobian[None,:,:]), 2, axis=(1,2))
+
+    slope_rho = np.polyfit(Time, L2norms_rho, 1)[0]
+    slope_phi = np.polyfit(Time, L2norms_phi, 1)[0]
+
+    assert abs(slope_rho - omega_Im) / abs(omega_Im) < 1
+    assert abs(slope_phi - omega_Im) / abs(omega_Im) < 1

tests/geometryRTheta/time_solver/test_diocotron.sh

@@ -0,0 +1,48 @@
+#!/bin/bash
+set -xe
+
+if [ $# -ne 3 ]; then
+    echo "Usage: $0 <GYSELALIBXX_SRCDIR> <GYSELALIBXX_EXEC> <PYTHON3_EXE>"
+    exit 1
+fi
+GYSELALIBXX_SRCDIR="$1"
+GYSELALIBXX_EXEC="$2"
+PYTHON3_EXE="$3"
+
+OUTDIR="${PWD}"
+
+TMPDIR="$(mktemp -p "${PWD}" -d run-XXXXXXXXXX)"
+function finish {
+    rm -rf "${TMPDIR}"
+}
+trap finish EXIT QUIT ABRT KILL SEGV TERM STOP
+
+cd "${TMPDIR}"
+
+# Create a parameter file with the default values. 
+"${GYSELALIBXX_EXEC}" "--dump-config" "${PWD}/diocotron_params.yaml"
+
+# Modify the default parameter files for a faster test.
+MODIFY_PARAMETERS=$(cat <<EOF
+import yaml
+
+with open('diocotron_params.yaml') as f:
+    data = yaml.safe_load(f)
+
+data['SplineMesh']['r_ncells'] = 32
+data['SplineMesh']['theta_ncells'] = 64
+data['Time']['final_T'] = 40.
+data['Output']['time_step_diag'] = 40
+
+with open('diocotron_params.yaml', 'w') as f:
+    yaml.dump(data, f)
+EOF
+)
+
+"${PYTHON3_EXE}" -c "$MODIFY_PARAMETERS"
+
+# Launch the test with the modified parameter file. 
+"${GYSELALIBXX_EXEC}" "${PWD}/diocotron_params.yaml"
+
+export PYTHONPATH="${GYSELALIBXX_SRCDIR}/post-process/PythonScripts:${PYTHONPATH}"
+"${PYTHON3_EXE}" -B "${GYSELALIBXX_SRCDIR}/tests/geometryRTheta/time_solver/growth_rate_test.py" "${TMPDIR}/output"