Add a helper file to simplify PDI usage

See merge request gysela-developpers/gyselalibxx!704

--------------------------------------------

Co-authored-by: Virginie Grandgirard <virginie.grandgirard@cea.fr>

simulations/geometry5D/testcollisions/testcollisions.cpp

@@ -18,6 +18,7 @@
 #include "geometry.hpp"
 #include "input.hpp"
 #include "paraconfpp.hpp"
+#include "pdi_helper.hpp"
 #include "pdi_out.yml.hpp"
 #include "simpson_quadrature.hpp"
 #include "testcollisions.yaml.hpp"
@@ -65,76 +66,33 @@ int main(int argc, char** argv)
             PDI_OUT,
             NULL);
 
-    std::vector<size_t> grid_tor1_extents(1);
-    std::vector<size_t> grid_tor2_extents(1);
-    std::vector<size_t> grid_tor3_extents(1);
-    std::vector<size_t> grid_vpar_extents(1);
-    std::vector<size_t> grid_mu_extents(1);
-    std::vector<size_t> species_extents(1);
-    std::vector<size_t> charges_extents(1);
-    std::vector<size_t> masses_extents(1);
-    PDI_multi_expose(
-            "read_grid_extents",
-            "grid_tor1_extents",
-            grid_tor1_extents.data(),
-            PDI_INOUT,
-            "grid_tor2_extents",
-            grid_tor2_extents.data(),
-            PDI_INOUT,
-            "grid_tor3_extents",
-            grid_tor3_extents.data(),
-            PDI_INOUT,
-            "grid_vpar_extents",
-            grid_vpar_extents.data(),
-            PDI_INOUT,
-            "grid_mu_extents",
-            grid_mu_extents.data(),
-            PDI_INOUT,
-            "species_extents",
-            species_extents.data(),
-            PDI_INOUT,
-            "charges_extents",
-            charges_extents.data(),
-            PDI_INOUT,
-            "masses_extents",
-            masses_extents.data(),
-            PDI_INOUT,
-            NULL);
-    std::vector<double> grid_tor1(grid_tor1_extents[0]);
-    std::vector<double> grid_tor2(grid_tor2_extents[0]);
-    std::vector<double> grid_tor3(grid_tor3_extents[0]);
-    std::vector<double> grid_vpar(grid_vpar_extents[0]);
-    std::vector<double> grid_mu(grid_mu_extents[0]);
-    std::vector<int> species(species_extents[0]);
-    std::vector<double> charges(charges_extents[0]);
-    std::vector<double> masses(masses_extents[0]);
-    PDI_multi_expose(
+    std::vector<double> grid_tor1;
+    std::vector<double> grid_tor2;
+    std::vector<double> grid_tor3;
+    std::vector<double> grid_vpar;
+    std::vector<double> grid_mu;
+    std::vector<int> species;
+    std::vector<double> charges;
+    std::vector<double> masses;
+    PDI_get_arrays(
             "read_grid",
             "grid_tor1",
-            grid_tor1.data(),
-            PDI_INOUT,
+            grid_tor1,
             "grid_tor2",
-            grid_tor2.data(),
-            PDI_INOUT,
+            grid_tor2,
             "grid_tor3",
-            grid_tor3.data(),
-            PDI_INOUT,
+            grid_tor3,
             "grid_vpar",
-            grid_vpar.data(),
-            PDI_INOUT,
+            grid_vpar,
             "grid_mu",
-            grid_mu.data(),
-            PDI_INOUT,
+            grid_mu,
             "species",
-            species.data(),
-            PDI_INOUT,
+            species,
             "charges",
-            charges.data(),
-            PDI_INOUT,
+            charges,
             "masses",
-            masses.data(),
-            PDI_INOUT,
-            NULL);
+            masses);
+
     ddc::init_discrete_space<GridTor1>(grid_tor1);
     IdxRangeTor1 idx_range_tor1(IdxTor1(0), IdxStepTor1(grid_tor1.size()));
     ddc::init_discrete_space<GridTor2>(grid_tor2);

simulations/geometry5D/testcollisions/testcollisions.yaml.hpp

@@ -4,8 +4,8 @@
 
 constexpr char const* const params_yaml = R"PDI_CFG(
 InputFileNames:
-  read_restart : GyselaX_restart_00000.h5
-  write_restart : GyselaX_restart_00001.h5
+  read_restart : GyselaX_restart5D_00000.h5
+  write_restart : GyselaX_restart5D_00001.h5
  
 Algorithm:
   deltat: 0.125

src/io/pdi_helper.hpp

@@ -0,0 +1,62 @@
+// SPDX-License-Identifier: MIT
+#pragma once
+
+#include <pdi.h>
+
+namespace detail {
+
+template <class NameTuple, class OutputTuple, size_t... I>
+void PDI_get_array(
+        std::string event_name,
+        NameTuple names,
+        OutputTuple args,
+        std::integer_sequence<size_t, I...>)
+{
+    std::array<size_t, sizeof...(I)> sizes;
+    // Put sizes into PDI
+    ((PDI_share((std::get<I>(names) + "_extents").c_str(), &sizes[I], PDI_INOUT)), ...);
+    // Collect sizes from file
+    PDI_event((event_name + "_extents").c_str());
+    // Collect sizes from PDI
+    ((PDI_reclaim((std::get<I>(names) + "_extents").c_str())), ...);
+
+    // Update the size of the vectors
+    ((std::get<I>(args).resize(sizes[I])), ...);
+    // Put vector into PDI
+    ((PDI_share(std::get<I>(names).c_str(), std::get<I>(args).data(), PDI_INOUT)), ...);
+    // Collect vector from file
+    PDI_event(event_name.c_str());
+    // Collect vector from PDI
+    ((PDI_reclaim(std::get<I>(names).c_str())), ...);
+}
+
+template <class TupleType, size_t... I>
+auto get_name_tuple(TupleType input_args, std::integer_sequence<size_t, I...>)
+{
+    return std::make_tuple(std::string(std::get<I * 2>(input_args))...);
+}
+
+template <class TupleType, size_t... I>
+auto get_vector_tuple(TupleType input_args, std::integer_sequence<size_t, I...>)
+{
+    return std::tie(std::get<I * 2 + 1>(input_args)...);
+}
+
+} // namespace detail
+
+/**
+ * A helper function to read an unknown number of arrays from a file using PDI.
+ */
+template <class T, class... Args>
+void PDI_get_arrays(
+        std::string event_name,
+        std::string name,
+        std::vector<T>& out_vector,
+        Args&... input_args)
+{
+    auto idx_sequence = std::make_index_sequence<sizeof...(Args) / 2 + 1> {};
+    auto arg_tuple = std::tie(name, out_vector, input_args...);
+    auto names = detail::get_name_tuple(arg_tuple, idx_sequence);
+    auto out_vectors = detail::get_vector_tuple(arg_tuple, idx_sequence);
+    detail::PDI_get_array(event_name, names, out_vectors, idx_sequence);
+}

tests/geometry5D/collisions/input_params_twospecies_1x1x1x128x64.yaml → tests/geometry5D/collisions/input_params_twospecies_geom5D.yaml

@@ -1,9 +1,9 @@
 Mesh:
-  ncell_tor1: 1       # number of cells in the radial direction
-  ncell_tor2: 1       # number of cells in the poloidal direction
-  ncell_tor3: 1       # number of cells in the toroidal direction
-  ncell_vpar: 128     # number of cells in the parallel velocity direction
-  ncell_mu: 64        # number of cells in the adiabatic invariant direction
+  ncell_bkpoints_tor1: 8       # number of cells in the radial direction
+  ncell_bkpoints_tor2: 4       # number of cells in the poloidal direction
+  ncell_bkpoints_tor3: 16      # number of cells in the toroidal direction
+  ncell_bkpoints_vpar: 128     # number of cells in the parallel velocity direction
+  ncell_bkpoints_mu: 64        # number of cells in the adiabatic invariant direction
 
   min_tor1: 0.0
   max_tor1: 1.0

tests/geometry5D/collisions/test_collisions.sh

@@ -1,15 +1,30 @@
 #!/bin/bash
 set -xe
 
+# Collisions test
+#
+#Script to automise :
+# - the creation of the initial restart file with the python script `init_distribution.py`
+# - the creation of the input YAML file required as input of the C++ simulation `testcollision`
+#
+#For instance, the following command
+#`./testcollisions.sh ${GYSELALIBXX_SRC}/build/simulations/geometry5D/testcollisions input_params_twospecies_geom5D.yaml`
+#
+#creates the folder `D_INPUT_PARAMS_TWOSPECIES_GEOM5D.YAML` containing:
+#  - `GysX_rst_00000.h5` : output of the python script `init_distribution.py`
+#  - `GysX_rst_00001.h5` : output of the C++ collision executable
+#  - `coll_ref.yml` : input for C++ collision executable automatically created by the bash script `testcollision.sh`
+#  - `diff_f_vpar_mu_itor1eq0_itor2eq0_itor3eq0_ispeq0.png` : output figure to compare the results between `GysX_rst_00000.h5` and `GysX_rst_00001.h5`
+
+
 # module load cray-python
 
-if [ $# -ne 1 ]; then
-    echo "Usage: $0 <INPUT_YAML_FILE>"
+if [ $# -ne 2 ]; then
+    echo "Usage: $0 <TESTCOLLISION_EXE> <INPUT_YAML_FILE>"
     exit 1
 fi
 
-# TESTCOLL_EXE="$(pwd)/../../../build.mi250.hipcc.adastra.spack/simulations/geometry5D/testcollisions/testcollisions"
-TESTCOLL_EXE="$(pwd)/../../../build-v100/simulations/geometry5D/testcollisions/testcollisions"
+TESTCOLL_EXE=$(realpath $1)
 
 GYSELA_IO_PATH="$(pwd)/gysela_io"
 
@@ -23,7 +38,7 @@ else
 fi
 export PYTHONPATH="${GYSELA_IO_PATH}:${PYTHONPATH}"
 
-INPUT_YAML_FILE="${1}"
+INPUT_YAML_FILE="$2"
 
 CASENAME=$(basename "${INPUT_YAML_FILE}" | tr '[:lower:]' '[:upper:]')
 RESDIR="D_${CASENAME}"
@@ -51,7 +66,7 @@ EOF
 cp -- "${INPUT_YAML_FILE}" "${RESDIR}"
 mv -- "${INPUT_CXX_YAML_FILE}" "${RESDIR}"
 cd -- "${RESDIR}"
-python3 "${GYSELA_IO_PATH}/initialisation/init_distribution.py" -i "${INPUT_YAML_FILE}" -o "${INIT_RST_FILE}"
+python3 "${GYSELA_IO_PATH}/initialisation/init_distribution.py"  -i "${INPUT_YAML_FILE}" -o "${INIT_RST_FILE}"
 if [ -f ${SAVE_RST_FILE} ]; then
   rm -- "${SAVE_RST_FILE}"
 fi