cpp20 concepts

CMakeLists.txt

@@ -387,8 +387,7 @@ set(BOUT_SOURCES
     ./src/sys/range.cxx
     ./src/sys/slepclib.cxx
     ./src/sys/timer.cxx
-    ./src/sys/type_name.cxx
-    ./src/sys/utils.cxx
+        ./src/sys/utils.cxx
     ${CMAKE_CURRENT_BINARY_DIR}/include/bout/revision.hxx
     ${CMAKE_CURRENT_BINARY_DIR}/include/bout/version.hxx
 )

include/bout/field.hxx

@@ -258,8 +258,11 @@ inline T filledFrom(const T& f, Function func, std::string region_string = "RGN_
   return result;
 }
 
+template<typename T>
+concept IsField = std::is_base_of_v<Field, T>;
+
 /// Unary + operator. This doesn't do anything
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 T operator+(const T& f) {
   return f;
 }
@@ -341,7 +344,7 @@ inline T fromFieldAligned(const T& f, const std::string& region = "RGN_ALL") {
 /// @param[in] f      Input field
 /// @param[in] allpe  Minimum over all processors?
 /// @param[in] rgn    The region to calculate the result over
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline BoutReal min(const T& f, bool allpe = false,
                     const std::string& rgn = "RGN_NOBNDRY") {
 
@@ -390,7 +393,7 @@ inline BoutReal min(const BinaryExpr<ResT, L, R, Func>& f, bool allpe = false,
 /// @param[in] f       The field to check
 /// @param[in] allpe   Check over all processors
 /// @param[in] region  The region to check for uniformity over
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline bool isUniform(const T& f, bool allpe = false,
                       const std::string& region = "RGN_ALL") {
   bool result = true;
@@ -418,7 +421,7 @@ inline bool isUniform(const T& f, bool allpe = false,
 /// @param[in] f       The field to check
 /// @param[in] allpe   Check over all processors
 /// @param[in] region  The region to assume is uniform
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline BoutReal getUniform(const T& f, [[maybe_unused]] bool allpe = false,
                            const std::string& region = "RGN_ALL") {
 #if CHECK > 1
@@ -443,7 +446,7 @@ inline BoutReal getUniform(const T& f, [[maybe_unused]] bool allpe = false,
 /// @param[in] f      Input field
 /// @param[in] allpe  Maximum over all processors?
 /// @param[in] rgn    The region to calculate the result over
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline BoutReal max(const T& f, bool allpe = false,
                     const std::string& rgn = "RGN_NOBNDRY") {
 
@@ -494,7 +497,7 @@ inline BoutReal max(const BinaryExpr<ResT, L, R, Func>& f, bool allpe = false,
 /// @param[in] f      Input field
 /// @param[in] allpe  Mean over all processors?
 /// @param[in] rgn    The region to calculate the result over
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline BoutReal mean(const T& f, bool allpe = false,
                      const std::string& rgn = "RGN_NOBNDRY") {
 
@@ -544,7 +547,7 @@ inline BoutReal mean(const BinaryExpr<ResT, L, R, Func>& f, bool allpe = false,
 /// This loops over the entire domain, including guard/boundary cells by
 /// default (can be changed using the \p rgn argument)
 /// If CHECK >= 3 then the result will be checked for non-finite numbers
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 T pow(const T& lhs, const T& rhs, const std::string& rgn = "RGN_ALL") {
 
   ASSERT1(areFieldsCompatible(lhs, rhs));
@@ -557,7 +560,7 @@ T pow(const T& lhs, const T& rhs, const std::string& rgn = "RGN_ALL") {
   return result;
 }
 
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 T pow(const T& lhs, BoutReal rhs, const std::string& rgn = "RGN_ALL") {
 
   // Check if the inputs are allocated
@@ -572,7 +575,7 @@ T pow(const T& lhs, BoutReal rhs, const std::string& rgn = "RGN_ALL") {
   return result;
 }
 
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 T pow(BoutReal lhs, const T& rhs, const std::string& rgn = "RGN_ALL") {
 
   // Check if the inputs are allocated
@@ -639,7 +642,7 @@ std::optional<int> getPerpYIndex(const BinaryExpr<ResT, L, R, Func>& expr) {
     }                                                                                   \
   };                                                                                    \
   };                                                                                    \
-  template <typename T, typename = bout::utils::EnableIfField<T>>                       \
+  template <IsField T>                       \
   inline auto name(const T& f, const std::string& rgn = "RGN_ALL") {                    \
     if constexpr (std::is_same_v<T, Field3DParallel>) {                                 \
       /* Check if the input is allocated */                                             \
@@ -696,7 +699,7 @@ struct Square {
 };
 }; // namespace bout::op
 
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline auto SQ(const T& f, const std::string& rgn = "RGN_ALL") {
   if constexpr (std::is_same_v<T, Field3DParallel>) {
     checkData(f);
@@ -840,7 +843,7 @@ FIELD_FUNC(tanh, ::tanh)
 /// Check if all values of a field \p var are finite.
 /// Loops over all points including the boundaries by
 /// default (can be changed using the \p rgn argument
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline bool finite(const T& f, const std::string& rgn = "RGN_ALL") {
 
   if (!f.isAllocated()) {
@@ -858,7 +861,7 @@ inline bool finite(const T& f, const std::string& rgn = "RGN_ALL") {
 
 /// Makes a copy of a field \p f, ensuring that the underlying data is
 /// not shared.
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 T copy(const T& f) {
   T result = f;
   result.allocate();
@@ -873,7 +876,7 @@ class Field3DParallel;
 /// @param[in] var  Variable to apply floor to
 /// @param[in] f    The floor value
 /// @param[in] rgn  The region to calculate the result over
-template <typename T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 inline T floor(const T& var, BoutReal f, const std::string& rgn = "RGN_ALL") {
   checkData(var);
   T result = copy(var);

include/bout/field3d.hxx

@@ -669,7 +669,7 @@ protected:
       _track(change, operation);
     }
   }
-  template <typename T, typename = bout::utils::EnableIfField<T>>
+  template <IsField T>
   void _track(const T& change, std::string operation);
   void _track(const BoutReal& change, std::string operation);
 

include/bout/fieldgroup.hxx

@@ -1,7 +1,6 @@
 #ifndef BOUT_FIELDGROUP_H
 #define BOUT_FIELDGROUP_H
 
-#include <bout/traits.hxx>
 #include <bout/vector2d.hxx>
 #include <bout/vector3d.hxx>
 
@@ -62,7 +61,7 @@ public:
   /// to FieldGroup, leading to an infinite loop.
   template <typename... Ts>
   explicit FieldGroup(Ts&... ts) {
-    add(ts...);
+    (add(ts), ...);
   }
 
   /// Copy contents of another FieldGroup \p other into this group.
@@ -121,27 +120,8 @@ public:
 
   /// Add multiple fields to this group
   template <typename... Ts>
-  void add(Field& t, Ts&... ts) {
-    add(t);     // Add the first using functions above
-    add(ts...); // Add the rest
-  }
-
-  template <typename... Ts>
-  void add(Field3D& t, Ts&... ts) {
-    add(t);     // Add the first using functions above
-    add(ts...); // Add the rest
-  }
-
-  template <typename... Ts>
-  void add(Vector3D& t, Ts&... ts) {
-    add(t);     // Add the first using functions above
-    add(ts...); // Add the rest
-  }
-
-  template <typename... Ts>
-  void add(Vector2D& t, Ts&... ts) {
-    add(t);     // Add the first using functions above
-    add(ts...); // Add the rest
+  requires(sizeof...(Ts) > 1) void add(Ts&... ts) {
+    (add(ts), ...);
   }
 
   /// Return number of fields

include/bout/sys/type_name.hxx

@@ -8,8 +8,10 @@
 #include "bout/bout_types.hxx"
 
 #include <string>
+#include <type_traits> // Required for std::is_same_v
 #include <typeinfo>
 
+// Forward declarations
 class Field2D;
 class Field3D;
 class FieldPerp;
@@ -23,44 +25,39 @@ namespace utils {
 
 template <typename T>
 std::string typeName() {
-  return typeid(T).name();
+  if constexpr (std::is_same_v<T, bool>)
+    return "bool";
+  else if constexpr (std::is_same_v<T, int>)
+    return "int";
+  else if constexpr (std::is_same_v<T, std::string>)
+    return "string";
+  // Specialised for BOUT++ types to ensure that the result is human-readable
+  else if constexpr (std::is_same_v<T, BoutReal>)
+    return "BoutReal";
+  else if constexpr (std::is_same_v<T, Field2D>)
+    return "Field2D";
+  else if constexpr (std::is_same_v<T, Field3D>)
+    return "Field3D";
+  else if constexpr (std::is_same_v<T, FieldPerp>)
+    return "FieldPerp";
+  else if constexpr (std::is_same_v<T, Array<int>>)
+    return "Array<int>";
+  else if constexpr (std::is_same_v<T, Array<BoutReal>>)
+    return "Array<BoutReal>";
+  else if constexpr (std::is_same_v<T, Matrix<int>>)
+    return "Matrix<int>";
+  else if constexpr (std::is_same_v<T, Matrix<BoutReal>>)
+    return "Matrix<BoutReal>";
+  else if constexpr (std::is_same_v<T, Tensor<int>>)
+    return "Tensor<int>";
+  else if constexpr (std::is_same_v<T, Tensor<BoutReal>>)
+    return "Tensor<BoutReal>";
+  else {
+    return typeid(T).name();
+  }
 }
 
-template <>
-std::string typeName<bool>();
-
-template <>
-std::string typeName<int>();
-
-template <>
-std::string typeName<std::string>();
-
-// Specialised for BOUT++ types to ensure that the result is human-readable
-template <>
-std::string typeName<BoutReal>();
-
-template <>
-std::string typeName<Field2D>();
-
-template <>
-std::string typeName<Field3D>();
-
-template <>
-std::string typeName<FieldPerp>();
-
-template <>
-std::string typeName<Array<int>>();
-template <>
-std::string typeName<Array<BoutReal>>();
-template <>
-std::string typeName<Matrix<int>>();
-template <>
-std::string typeName<Matrix<BoutReal>>();
-template <>
-std::string typeName<Tensor<int>>();
-template <>
-std::string typeName<Tensor<BoutReal>>();
 } // namespace utils
 } // namespace bout
 
-#endif //TYPE_NAME_HXX
+#endif // TYPE_NAME_HXX

include/bout/where.hxx

@@ -31,14 +31,18 @@
 #include "bout/field.hxx"
 #include "bout/field2d.hxx"
 #include "bout/field3d.hxx"
+#include <type_traits>
 
 /// For each point, choose between two inputs based on a third input
 ///
 /// @param[in] test   The value which determines which input to use
 /// @param[in] gt0    Uses this value if test > 0.0
 /// @param[in] le0    Uses this value if test <= 0.0
-template <class T, class U, class V,
-          class ResultType = typename bout::utils::EnableIfField<T, U, V>>
+
+// Overload 1: Three fields
+template <typename T, typename U, typename V,
+          typename ResultType = std::common_type_t<T, U, V>>
+requires IsField<T>&& IsField<U>&& IsField<V>
 auto where(const T& test, const U& gt0, const V& le0) -> ResultType {
   ASSERT1_FIELDS_COMPATIBLE(test, gt0);
   ASSERT1_FIELDS_COMPATIBLE(test, le0);
@@ -51,7 +55,9 @@ auto where(const T& test, const U& gt0, const V& le0) -> ResultType {
   return result;
 }
 
-template <class T, class U, class ResultType = typename bout::utils::EnableIfField<T, U>>
+// Overload 2: Two fields, one BoutReal (le0)
+template <typename T, typename U, typename ResultType = std::common_type_t<T, U>>
+requires IsField<T>&& IsField<U>
 auto where(const T& test, const U& gt0, BoutReal le0) -> ResultType {
   ASSERT1_FIELDS_COMPATIBLE(test, gt0);
 
@@ -63,7 +69,9 @@ auto where(const T& test, const U& gt0, BoutReal le0) -> ResultType {
   return result;
 }
 
-template <class T, class V, class ResultType = typename bout::utils::EnableIfField<T, V>>
+// Overload 3: Two fields, one BoutReal (gt0)
+template <typename T, typename V, typename ResultType = std::common_type_t<T, V>>
+requires IsField<T>&& IsField<V>
 auto where(const T& test, BoutReal gt0, const V& le0) -> ResultType {
   ASSERT1_FIELDS_COMPATIBLE(test, le0);
 
@@ -75,7 +83,9 @@ auto where(const T& test, BoutReal gt0, const V& le0) -> ResultType {
   return result;
 }
 
-template <class T, class ResultType = T>
+// Overload 4: One field, two BoutReals
+template <typename T, typename ResultType = T>
+requires IsField<T>
 auto where(const T& test, BoutReal gt0, BoutReal le0) -> ResultType {
   ResultType result{emptyFrom(test)};
 

src/field/field3d.cxx

@@ -982,7 +982,7 @@ Field3D& Field3D::enableTracking(const std::string& name,
   return *this;
 }
 
-template <typename T, typename>
+template <IsField T>
 void Field3D::_track(const T& change, std::string operation) {
   if (tracking_state == 0) {
     return;
@@ -1010,12 +1010,10 @@ void Field3D::_track(const T& change, std::string operation) {
   });
 }
 
-template void
-Field3D::_track<Field3D, bout::utils::EnableIfField<Field3D>>(const Field3D&,
-                                                              std::string);
+template void Field3D::_track<Field3D>(const Field3D&, std::string);
 template void Field3D::_track<Field3DParallel>(const Field3DParallel&, std::string);
 template void Field3D::_track<Field2D>(const Field2D&, std::string);
-template void Field3D::_track<>(const FieldPerp&, std::string);
+template void Field3D::_track<FieldPerp>(const FieldPerp&, std::string);
 
 void Field3D::_track(const BoutReal& change, std::string operation) {
   if (tracking_state == 0) {

src/mesh/coordinates.cxx

@@ -56,7 +56,7 @@ namespace {
 // initializing yet, leading to an infinite recursion.
 // Also, here we interpolate for the boundary points at xstart/ystart and
 // (xend+1)/(yend+1) instead of extrapolating.
-template <class T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 void fillGuards_impl(T& result, CELL_LOC location, const T& f, bool extrapolate_x,
                      bool extrapolate_y, bool no_extra_interpolate = false) {
   const auto* localmesh = result.getMesh();
@@ -172,7 +172,7 @@ void fillGuards_impl(T& result, CELL_LOC location, const T& f, bool extrapolate_
 /// Boundary guard cells are set by extrapolating from the grid, like
 /// 'free_o3' boundary conditions
 /// Corner guard cells are set to BoutNaN
-template <class T, typename = bout::utils::EnableIfField<T>>
+template <IsField T>
 auto interpolateAndExtrapolate(const T& f_, CELL_LOC location, bool extrapolate_x,
                                bool extrapolate_y, bool no_extra_interpolate,
                                ParallelTransform* pt_ = nullptr) -> T {