mdarray is a mess, reverting to std::array

CMakeLists.txt

@@ -4,7 +4,7 @@ cmake_minimum_required(VERSION 3.10)
 project(OpenQCD++)
 
 # Specify the C++ standard
-set(CMAKE_CXX_STANDARD 17)
+set(CMAKE_CXX_STANDARD 20)
 set(CMAKE_CXX_STANDARD_REQUIRED True)
 
 # Add the executables

app/init_print_lattice/main.cpp

@@ -2,21 +2,19 @@
 #include <complex>
 
 int main() {
-    LatticeField<float> latticeFloat(4, 4, 4, 4);
-    latticeFloat.initialize();
-    latticeFloat.print();
+    LatticeField<float> latticeFloat(4, 4, 4, 4, 0.0);
+    latticeFloat.print(1,2,2,3);
 
-    LatticeField<double> latticeDouble(4, 4, 4, 4);
-    latticeDouble.initialize();
-    latticeDouble.print();
+    LatticeField<double> latticeDouble(4, 4, 4, 4, 1.0);
+    latticeFloat.print(1,2,2,3);
 
-    LatticeField<std::complex<float>> latticeComplexFloat(4, 4, 4, 4);
-    latticeComplexFloat.initialize();
-    latticeComplexFloat.print();
+    LatticeField<std::complex<float>> latticeComplexFloat(4, 4, 4, 4, std::complex<float>(0.0, 0.0));
+    latticeComplexFloat.print(0,0,0,0);
+    latticeComplexFloat.print(0,0,3,0);
 
-    LatticeField<std::complex<double>> latticeComplexDouble(4, 4, 4, 4);
-    latticeComplexDouble.initialize();
-    latticeComplexDouble.print();
+    LatticeField<std::complex<double>> latticeComplexDouble(4, 4, 4, 4, std::complex<double>(1.0, 1.0));
+    latticeComplexDouble.print(0,0,0,0);
+    latticeComplexDouble.print(0,0,3,0);
 
     return 0;
-}
+}

include/lattice/LatticeField.h

@@ -10,20 +10,18 @@ template <typename T>
 class LatticeField {
 public:
     LatticeField(
-        const int dim1, const int dim2, const int dim3, const int dim4,
+        const int dim0, const int dim1, const int dim2, const int dim3,
         const T default_value
     );
     ~LatticeField() = default;
 
-    void print() const;
+    void print(size_t i0, size_t i1, size_t i2, size_t i3) const;
     void randomize();
 
 private:
-    int dim1, dim2, dim3, dim4;
-    int volume;
-    // FIXME: Better to use std::vector instead?
-    mdarray<T, dim0, dim1, dim2, dim3> field_values;
-
+    const int dim0, dim1, dim2, dim3;
+    const int volume;
+    std::array<T, volume> field_values;
 };
 
 #endif // LATTICEFIELD_H

include/mdarray/Makefile

@@ -5,18 +5,20 @@ CXX = g++
 CXXFLAGS = -std=c++20 -Wall -Wextra -O2
 
 # Source files
-SRCS = main.cpp
+SRCS = $(wildcard *.cpp)
 
 # Header files
 HDRS = mdarray.h
 
-# Target executable
-TARGET = main
+# Target executables
+TARGETS = $(SRCS:.cpp=)
 
 # Build rule
-$(TARGET): $(SRCS) $(HDRS)
-	$(CXX) $(CXXFLAGS) -o $(TARGET) $(SRCS)
+all: $(TARGETS)
+
+%: %.cpp $(HDRS)
+	$(CXX) $(CXXFLAGS) -o $@ $<
 
 # Clean rule
 clean:
-	rm -f $(TARGET)
+	rm -f $(TARGETS)

include/mdarray/mdarray.h

@@ -1,28 +1,112 @@
 // TODO: FIXME: This header file templates a multidimensional array for use in lattice fields.
 // This should be replaced with a fuller std::mdarray in a hopefully not so distant C++ standard.
 // Cf. https://www.open-std.org/jtc1/sc22/wg21/docs/papers/2023/p1684r5.html
-
 #include <array>
-#include <cstddef>
+#include <algorithm>
+
+
+// Type alias for multi-dimensional array
+template <typename T, std::size_t... sizes>
+struct mdarray_impl;
+
 
 // This recursively instantiates the template arguments and creates a nested std::array type.
-template <class T, std::size_t size, std::size_t... sizes>
-struct ArrayHelper {
-    using type = std::array<typename ArrayHelper<T, sizes...>::type, size>;
+template <typename T, std::size_t dim_i, std::size_t... dim_ixx>
+struct mdarray_impl<T, dim_i, dim_ixx...> {
+    using type = std::array<typename mdarray_impl<T, dim_ixx...>::type, dim_i>;
 };
 
+
 // Base case.
-template <class T, std::size_t size>
-struct ArrayHelper<T, size> {
-    using type = std::array<T, size>;
+template <class T, std::size_t dim_i>
+struct mdarray_impl<T, dim_i> {
+    using type = std::array<T, dim_i>;
 };
 
+
 // TODO: This would be changed to the std::mdarray type in a future C++ standard.
-template <class T, std::size_t... sizes>
-using mdarray = typename ArrayHelper<T, sizes...>::type;
+template <typename T, std::size_t... sizes>
+using mdarray = typename mdarray_impl<T, sizes...>::type;
+
+
+// Function template to fill the mdarray with a stream of numbers
+// template <typename T, std::size_t size, std::size_t... sizes, typename InputIt>
+// void fill_mdarray_with_iterator(mdarray<T, size, sizes...>& arr, InputIt& it) {
+//     for (auto& sub_arr : arr) {
+//         fill_mdarray_with_iterator(sub_arr, it);
+//     }
+// }
+
+// // Base case for the innermost array with iterator
+// template <typename T, std::size_t size, typename InputIt>
+// void fill_mdarray_with_iterator(mdarray<T, size>& arr, InputIt& it) {
+//     for (auto& elem : arr) {
+//         elem = *it;
+//         ++it;
+//     }
+// }
+
+// Function template to fill the mdarray with a default value
+// FIXME: This is buggy.
+// template <typename T, std::size_t dim_i, std::size_t... dim_ixx>
+// void fill_mdarray_with_value(mdarray_impl<T, dim_i, dim_ixx...>& arr, const T& default_value) {
+//     for (auto& sub_arr : arr) {
+//         fill_mdarray_with_value<T, dim_ixx...>(sub_arr, default_value); // Recursive call
+//     }
+// }
 
+// Function to fill with default value on 4D array
+template <typename T, std::size_t dim_0, std::size_t dim_1, std::size_t dim_2, std::size_t dim_3>
+void fill_mdarray_with_value_4d(mdarray<T, dim_0, dim_1, dim_2, dim_3>& arr, const T& default_value) {
+    for (auto& sub_arr_3d : arr) {
+        for (auto& sub_arr_2d : sub_arr_3d) {
+            for (auto& sub_arr_1d : sub_arr_2d) {
+                std::fill(sub_arr_1d.begin(), sub_arr_1d.end(), default_value); // Fill innermost array
+            }
+        }
+    }
+}
 
-// Ex. use case
-// int main() {
-//     mdarray<int, 5, 3, 4, 3> arr;
+// Base case for the innermost array
+template <typename T, std::size_t dim_i>
+void fill_mdarray_with_value(std::array<T, dim_i>& arr, const T& default_value) {
+    std::fill(arr.begin(), arr.end(), default_value); // Fill innermost array
+}
+
+
+
+// Function to fill a 4D array with a stream of numbers (e.g. for randomly generated data)
+template <typename T, std::size_t dim_0, std::size_t dim_1, std::size_t dim_2, std::size_t dim_3, typename InputIt>
+void fill_mdarray_with_iterator_4d(mdarray<T, dim_0, dim_1, dim_2, dim_3>& arr, InputIt& it) {
+    for (auto& sub_arr_3d : arr) {
+        for (auto& sub_arr_2d : sub_arr_3d) {
+            for (auto& sub_arr_1d : sub_arr_2d) {
+                for (auto& elem : sub_arr_1d) {
+                    elem = *it;
+                    ++it;
+                }
+            }
+        }
+    }
+}
+
+
+
+
+
+
+
+// // Function template to fill the mdarray with a default value
+// template <typename T, std::size_t size, std::size_t... sizes>
+// void fill_mdarray_with_value(mdarray<T, size, sizes...>& arr, const T& default_value) {
+//     for (auto& sub_arr : arr) {
+//         fill_mdarray_with_value(sub_arr, default_value); // Recursive call
+//     }
 // }
+
+// // Base case for the innermost array
+// template <typename T, std::size_t size>
+// void fill_mdarray_with_value(std::array<T, size>& arr, const T& default_value) {
+//     std::fill(arr.begin(), arr.end(), default_value); // Fill innermost array
+// }
+

include/mdarray/test2.cpp

@@ -0,0 +1,35 @@
+// Ex. 2 use case
+#include <vector>
+#include <iostream>
+#include "mdarray.h"
+#include <cstddef>
+#include <iterator>
+#include <algorithm>
+#include <numeric>
+
+
+int main() {
+
+    // Define the dimensions
+    constexpr std::size_t dim_0 = 4,
+                          dim_1 = 4,
+                          dim_2 = 3,
+                          dim_3 = 2;
+
+    mdarray<int, dim_0, dim_1, dim_2, dim_3> arr;
+
+    // Fill the mdarray with a default value
+    fill_mdarray_with_value_4d<int, dim_0, dim_1, dim_2, dim_3>(arr, 42);
+
+    // Print the first element to verify
+    std::cout << arr[0][0][0][0] << std::endl; // Should print 42
+
+    // Now with a 1D array
+    constexpr std::size_t dim_4 = 5;
+
+    mdarray<int, dim_4> arr_1d;
+    int default_value = 33;
+    fill_mdarray_with_value(arr_1d, default_value);
+    std::cout << arr_1d[0] << std::endl; // Should print 33
+    return 0;
+}

include/mdarray/test3.cpp

@@ -0,0 +1,40 @@
+// This test checks that we can fill an mdarray with a stream of randomly generated numbers.
+#include <iostream>
+#include <random>
+#include <vector>
+#include "mdarray.h"
+
+int main() {
+    // Define the dimensions
+    constexpr std::size_t dim_0 = 4,
+                          dim_1 = 4,
+                          dim_2 = 3,
+                          dim_3 = 25;
+
+    mdarray<int, dim_0, dim_1, dim_2, dim_3> arr;
+
+    // Create a random number generator
+    std::random_device rd;
+    std::mt19937 gen(rd());
+    std::uniform_int_distribution<> dis(1, 100);
+
+    // Generate a stream of random numbers
+    std::vector<int> random_numbers;
+    for (std::size_t i = 0; i < dim_0 * dim_1 * dim_2 * dim_3; ++i) {
+        random_numbers.push_back(dis(gen));
+    }
+
+    // Fill the mdarray with the random numbers
+    auto it = random_numbers.begin();
+    fill_mdarray_with_iterator_4d<int, dim_0, dim_1, dim_2, dim_3>(arr, it);
+
+    // Print the first couple of element to verify
+    for (size_t i3 = 0; i3 < dim_3; i3++)
+    {
+        std::cout << arr[0][0][0][i3] << std::endl;
+    }
+
+    return 0;
+}
+
+

include/mdarray/test4.cpp

@@ -0,0 +1,36 @@
+#include <iostream>
+#include <random>
+#include <vector>
+#include "mdarray.h"
+
+int main() {
+    // Define the dimensions
+    constexpr std::size_t dim_0 = 4,
+                          dim_1 = 4,
+                          dim_2 = 3,
+                          dim_3 = 25;
+
+    mdarray<double, dim_0, dim_1, dim_2, dim_3> arr;
+
+    // Create a random number generator
+    std::random_device rd;
+    std::mt19937 gen(rd());
+    std::normal_distribution<> dis(0.0, 1.0); // mu = 0, stdev = 1
+
+    // Generate a stream of random numbers
+    std::vector<double> random_numbers;
+    for (std::size_t i = 0; i < dim_0 * dim_1 * dim_2 * dim_3; ++i) {
+        random_numbers.push_back(dis(gen));
+    }
+
+    // Fill the mdarray with the random numbers
+    auto it = random_numbers.begin();
+    fill_mdarray_with_iterator_4d<double, dim_0, dim_1, dim_2, dim_3>(arr, it);
+
+    // Print the first couple of elements to verify
+    for (size_t i3 = 0; i3 < dim_3; i3++) {
+        std::cout << arr[0][0][0][i3] << std::endl;
+    }
+
+    return 0;
+}