diff --git a/include/bout/field.hxx b/include/bout/field.hxx index 3f7f591528..8372549122 100644 --- a/include/bout/field.hxx +++ b/include/bout/field.hxx @@ -6,7 +6,7 @@ * Copyright 2010 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu * * Contact: Ben Dudson, bd512@york.ac.uk - * + * * This file is part of BOUT++. * * BOUT++ is free software: you can redistribute it and/or modify diff --git a/include/bout/field2d.hxx b/include/bout/field2d.hxx index 906af38a1d..c278fbfd4e 100644 --- a/include/bout/field2d.hxx +++ b/include/bout/field2d.hxx @@ -250,8 +250,7 @@ public: Field2D& operator/=(BoutReal rhs); // FieldData virtual functions - - bool is3D() const override { return false; } + FieldType field_type() const override { return FieldType::field2d; } #if CHECK > 0 void doneComms() override { bndry_xin = bndry_xout = bndry_yup = bndry_ydown = true; } diff --git a/include/bout/field3d.hxx b/include/bout/field3d.hxx index 54d9fb85c9..0cc5f2fa3c 100644 --- a/include/bout/field3d.hxx +++ b/include/bout/field3d.hxx @@ -2,7 +2,7 @@ * Copyright 2010 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu * * Contact: Ben Dudson, bd512@york.ac.uk - * + * * This file is part of BOUT++. * * BOUT++ is free software: you can redistribute it and/or modify @@ -361,7 +361,7 @@ public: * Direct access to the underlying data array * * If CHECK > 2 then bounds checking is performed - * + * * If CHECK <= 2 then no checks are performed, to * allow inlining and optimisation of inner loops */ @@ -473,7 +473,7 @@ public: ///@} // FieldData virtual functions - bool is3D() const override { return true; } + FieldType field_type() const override { return FieldType::field3d; } #if CHECK > 0 void doneComms() override { bndry_xin = bndry_xout = bndry_yup = bndry_ydown = true; } diff --git a/include/bout/field_data.hxx b/include/bout/field_data.hxx index 185dcabf2d..2f0d17a36d 100644 --- a/include/bout/field_data.hxx +++ b/include/bout/field_data.hxx @@ -6,7 +6,7 @@ * Copyright 2010 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu * * Contact: Ben Dudson, bd512@york.ac.uk - * + * * This file is part of BOUT++. * * BOUT++ is free software: you can redistribute it and/or modify @@ -33,6 +33,7 @@ class FieldData; #include "bout/bout_types.hxx" #include "bout/unused.hxx" +#include #include #include #include @@ -71,13 +72,22 @@ public: /// Get variable location virtual CELL_LOC getLocation() const; + /// Enum to distinguish the different kinds of Fields + enum class FieldType : std::uint8_t { field3d, field2d, fieldperp }; + /// Is this an instance of `Field3D`, `Field2D`, or `FieldPerp`? + virtual FieldType field_type() const = 0; + // Defines interface which must be implemented /// True if variable is 3D - virtual bool is3D() const = 0; + [[deprecated("Use `field_type()` instead")]] + bool is3D() const { + return field_type() == FieldType::field3d; + } + /// Number of BoutReals in one element virtual int elementSize() const { return 1; } - virtual void doneComms(){}; // Notifies that communications done + virtual void doneComms() {}; // Notifies that communications done // Boundary conditions void setBoundary(const std::string& name); ///< Set the boundary conditions @@ -86,13 +96,13 @@ public: copyBoundary(const FieldData& f); ///< Copy the boundary conditions from another field virtual void applyBoundary(bool UNUSED(init) = false) {} - virtual void applyTDerivBoundary(){}; + virtual void applyTDerivBoundary() {}; - virtual void applyParallelBoundary(){}; - virtual void applyParallelBoundary(BoutReal UNUSED(t)){}; - virtual void applyParallelBoundary(const std::string& UNUSED(condition)){}; + virtual void applyParallelBoundary() {}; + virtual void applyParallelBoundary(BoutReal UNUSED(t)) {}; + virtual void applyParallelBoundary(const std::string& UNUSED(condition)) {}; virtual void applyParallelBoundary(const std::string& UNUSED(region), - const std::string& UNUSED(condition)){}; + const std::string& UNUSED(condition)) {}; // JMAD void addBndryFunction(FuncPtr userfunc, BndryLoc location); void addBndryGenerator(FieldGeneratorPtr gen, BndryLoc location); diff --git a/include/bout/fieldgroup.hxx b/include/bout/fieldgroup.hxx index 184766c6b8..440655e7ce 100644 --- a/include/bout/fieldgroup.hxx +++ b/include/bout/fieldgroup.hxx @@ -1,22 +1,23 @@ #ifndef BOUT_FIELDGROUP_H #define BOUT_FIELDGROUP_H -#include "bout/field_data.hxx" -#include - +#include #include #include #include -#include +class Field2D; +class Field3D; +class FieldPerp; +class Field; /// Group together fields for easier communication /// -/// Note: The FieldData class is used as a base class, -/// which is inherited by Field2D, Field3D, Vector2D and Vector3D -/// however Vector2D and Vector3D are stored by reference to their -/// components (x,y,z) as Field2D or Field3D objects. +/// Note: The `Field` class is used as a base class, +/// which is inherited by `Field2D`, `Field3D`, `FieldPerp`; +/// however `Vector2D` and `Vector3D` are stored by reference to their +/// components ``(x, y, z)`` as `Field2D` or `Field3D` objects. class FieldGroup { public: FieldGroup() = default; @@ -24,11 +25,9 @@ public: FieldGroup(FieldGroup&& other) = default; FieldGroup& operator=(const FieldGroup& other) = default; FieldGroup& operator=(FieldGroup&& other) = default; + ~FieldGroup() = default; - /// Constructor with a single FieldData \p f - FieldGroup(FieldData& f) { fvec.push_back(&f); } - - /// Constructor with a single Field3D \p f + FieldGroup(Field& f) { fvec.push_back(&f); } FieldGroup(Field3D& f) { fvec.push_back(&f); f3vec.push_back(&f); @@ -56,7 +55,7 @@ public: } /// Variadic constructor. Allows an arbitrary number of - /// FieldData arguments + /// Field arguments /// /// The explicit keyword prevents FieldGroup being constructed with arbitrary /// types. In particular arguments to add() cannot be implicitly converted @@ -78,12 +77,12 @@ public: return *this; } - /// Add a FieldData \p f to the group. + /// Add a Field \p f to the group. /// /// A pointer to this field will be stored internally, /// so the lifetime of this variable should be longer /// than the lifetime of this group. - void add(FieldData& f) { fvec.push_back(&f); } + void add(Field& f) { fvec.push_back(&f); } // Add a 3D field \p f, which goes into both vectors. // @@ -121,12 +120,8 @@ public: } /// Add multiple fields to this group - /// - /// This is a variadic template which allows Field3D objects to be - /// treated as a special case. An arbitrary number of fields can be - /// added. template - void add(FieldData& t, Ts&... ts) { + void add(Field& t, Ts&... ts) { add(t); // Add the first using functions above add(ts...); // Add the rest } @@ -165,16 +160,14 @@ public: } /// Iteration over all fields - using iterator = std::vector::iterator; - iterator begin() { return fvec.begin(); } - iterator end() { return fvec.end(); } + auto begin() { return fvec.begin(); } + auto end() { return fvec.end(); } /// Const iteration over all fields - using const_iterator = std::vector::const_iterator; - const_iterator begin() const { return fvec.begin(); } - const_iterator end() const { return fvec.end(); } + auto begin() const { return fvec.cbegin(); } + auto end() const { return fvec.cend(); } - const std::vector& get() const { return fvec; } + const std::vector& get() const { return fvec; } /// Iteration over 3D fields const std::vector& field3d() const { return f3vec; } @@ -183,8 +176,8 @@ public: void makeUnique(); private: - std::vector fvec; // Vector of fields - std::vector f3vec; // Vector of 3D fields + std::vector fvec; // Vector of fields + std::vector f3vec; // Vector of 3D fields }; /// Combine two FieldGroups diff --git a/include/bout/fieldperp.hxx b/include/bout/fieldperp.hxx index 1f70089990..9c1806b5ce 100644 --- a/include/bout/fieldperp.hxx +++ b/include/bout/fieldperp.hxx @@ -5,7 +5,7 @@ * Copyright 2010 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu * * Contact: Ben Dudson, bd512@york.ac.uk - * + * * This file is part of BOUT++. * * BOUT++ is free software: you can redistribute it and/or modify @@ -45,8 +45,8 @@ class Field3D; // #include "bout/field3d.hxx" /*! * Represents a 2D field perpendicular to the magnetic field - * at a particular index in Y, which only varies in X-Z. - * + * at a particular index in Y, which only varies in X-Z. + * * Primarily used inside field solvers */ class FieldPerp : public Field { @@ -205,7 +205,7 @@ public: /*! * Access to the underlying data array at a given x,z index - * + * * If CHECK > 2 then bounds checking is performed, otherwise * no checks are performed */ @@ -241,9 +241,9 @@ public: } /*! - * Access to the underlying data array. (X,Y,Z) indices for consistency with + * Access to the underlying data array. (X,Y,Z) indices for consistency with * other field types - * + * */ BoutReal& operator()(int jx, int UNUSED(jy), int jz) { return (*this)(jx, jz); } @@ -252,7 +252,7 @@ public: } /*! - * Addition, modifying in-place. + * Addition, modifying in-place. * This loops over the entire domain, including guard/boundary cells */ FieldPerp& operator+=(const FieldPerp& rhs); @@ -261,7 +261,7 @@ public: FieldPerp& operator+=(BoutReal rhs); /*! - * Subtraction, modifying in place. + * Subtraction, modifying in place. * This loops over the entire domain, including guard/boundary cells */ FieldPerp& operator-=(const FieldPerp& rhs); @@ -270,7 +270,7 @@ public: FieldPerp& operator-=(BoutReal rhs); /*! - * Multiplication, modifying in place. + * Multiplication, modifying in place. * This loops over the entire domain, including guard/boundary cells */ FieldPerp& operator*=(const FieldPerp& rhs); @@ -279,7 +279,7 @@ public: FieldPerp& operator*=(BoutReal rhs); /*! - * Division, modifying in place. + * Division, modifying in place. * This loops over the entire domain, including guard/boundary cells */ FieldPerp& operator/=(const FieldPerp& rhs); @@ -300,7 +300,7 @@ public: */ int getNz() const override { return nz; }; - bool is3D() const override { return false; } + FieldType field_type() const override { return FieldType::fieldperp; } friend void swap(FieldPerp& first, FieldPerp& second) noexcept; diff --git a/include/bout/mesh.hxx b/include/bout/mesh.hxx index a1ed6a9011..5797617a3e 100644 --- a/include/bout/mesh.hxx +++ b/include/bout/mesh.hxx @@ -3,7 +3,7 @@ * * Interface for mesh classes. Contains standard variables and useful * routines. - * + * * Changelog * ========= * @@ -11,7 +11,7 @@ * * Removing coordinate system into separate * Coordinates class * * Adding index derivative functions from derivs.cxx - * + * * 2010-06 Ben Dudson, Sean Farley * * Initial version, adapted from GridData class * * Incorporates code from topology.cpp and Communicator @@ -20,7 +20,7 @@ * Copyright 2010-2025 BOUT++ contributors * * Contact: Ben Dudson, dudson2@llnl.gov - * + * * This file is part of BOUT++. * * BOUT++ is free software: you can redistribute it and/or modify @@ -65,6 +65,7 @@ class Mesh; #include #include #include +#include class BoundaryRegion; class BoundaryRegionPar; @@ -301,11 +302,6 @@ public: /// @param g The group of fields to communicate. Guard cells will be modified void communicateYZ(FieldGroup& g); - /*! - * Communicate an X-Z field - */ - virtual void communicate(FieldPerp& f); - /*! * Send a list of FieldData objects * Packs arguments into a FieldGroup and passes @@ -815,8 +811,8 @@ protected: const std::vector readInts(const std::string& name, int n); /// Calculates the size of a message for a given x and y range - int msg_len(const std::vector& var_list, int xge, int xlt, int yge, - int ylt); + int msg_len(const std::vector& var_list, int xge, int xlt, int yge, + int ylt) const; /// Initialise derivatives void derivs_init(Options* options); diff --git a/include/bout/vector2d.hxx b/include/bout/vector2d.hxx index bdc375e698..1e19241a8e 100644 --- a/include/bout/vector2d.hxx +++ b/include/bout/vector2d.hxx @@ -13,7 +13,7 @@ * Copyright 2010 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu * * Contact: Ben Dudson, bd512@york.ac.uk - * + * * This file is part of BOUT++. * * BOUT++ is free software: you can redistribute it and/or modify @@ -142,7 +142,7 @@ public: CELL_LOC getLocation() const override; // FieldData virtual functions - bool is3D() const override { return false; } + FieldType field_type() const override { return FieldType::field2d; } int elementSize() const override { return 3; } /// Apply boundary condition to all fields diff --git a/include/bout/vector3d.hxx b/include/bout/vector3d.hxx index 0c71dcffa5..655a85ca73 100644 --- a/include/bout/vector3d.hxx +++ b/include/bout/vector3d.hxx @@ -9,7 +9,7 @@ * Copyright 2010 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu * * Contact: Ben Dudson, bd512@york.ac.uk - * + * * This file is part of BOUT++. * * BOUT++ is free software: you can redistribute it and/or modify @@ -45,10 +45,10 @@ class Vector2D; * ------- * * Vector3D f; - * + * * a.x; // Error! a.x not allocated * - * + * */ class Vector3D : public FieldData { public: @@ -70,7 +70,7 @@ public: ~Vector3D() override; /*! - * The components of the vector. These can be + * The components of the vector. These can be * either co- or contra-variant, depending on * the boolean flag "covariant" */ @@ -92,8 +92,8 @@ public: bool covariant{true}; /*! - * In-place conversion to covariant form. - * + * In-place conversion to covariant form. + * * If already covariant (covariant = true) then does nothing * If contravariant, multiplies by metric tensor g_{ij} */ @@ -104,7 +104,7 @@ public: * * If already contravariant (covariant = false) then does nothing * If covariant, multiplies by metric tensor g^{ij} - * + * */ void toContravariant(); @@ -114,14 +114,14 @@ public: * The first time this is called, a new Vector3D object is created. * Subsequent calls return a pointer to this same object * - * For convenience, a standalone function "ddt" exists, so that + * For convenience, a standalone function "ddt" exists, so that * * ddt(v) is equivalent to *(v.timeDeriv()) - * + * * This does some book-keeping to ensure that the time derivative * of the components is the same as the components of the time derivative * - * ddt(v).x == ddt(v.x) + * ddt(v).x == ddt(v.x) */ Vector3D* timeDeriv(); @@ -172,7 +172,7 @@ public: CELL_LOC getLocation() const override; // FieldData virtual functions - bool is3D() const override { return true; } + FieldType field_type() const override { return FieldType::field3d; } int elementSize() const override { return 3; } void applyBoundary(bool init = false) override; @@ -203,7 +203,7 @@ const Vector3D cross(const Vector3D& lhs, const Vector2D& rhs); /*! * Absolute magnitude (modulus) of a vector |v| - * + * * sqrt( v.x^2 + v.y^2 + v.z^2 ) */ const Field3D abs(const Vector3D& v, const std::string& region = "RGN_ALL"); diff --git a/src/mesh/impls/bout/boutmesh.cxx b/src/mesh/impls/bout/boutmesh.cxx index ea1b6a41b2..6c0a559972 100644 --- a/src/mesh/impls/bout/boutmesh.cxx +++ b/src/mesh/impls/bout/boutmesh.cxx @@ -2218,9 +2218,9 @@ void BoutMesh::topology() { } for (int i = 0; i < limiter_count; ++i) { - int const yind = limiter_yinds[i]; - int const xstart = limiter_xstarts[i]; - int const xend = limiter_xends[i]; + const int yind = limiter_yinds[i]; + const int xstart = limiter_xstarts[i]; + const int xend = limiter_xends[i]; output_info.write("Adding a limiter between y={} and {}. X indices {} to {}\n", yind, yind + 1, xstart, xend); add_target(yind, xstart, xend); @@ -2397,72 +2397,112 @@ void BoutMesh::overlapHandleMemory(BoutMesh* yup, BoutMesh* ydown, BoutMesh* xin * Communication utilities ****************************************************************/ -int BoutMesh::pack_data(const std::vector& var_list, int xge, int xlt, - int yge, int ylt, BoutReal* buffer) { +int BoutMesh::pack_data(const std::vector& var_list, int xge, int xlt, int yge, + int ylt, BoutReal* buffer) const { + using enum Field::FieldType; int len = 0; + const int zge = 0; + const int zlt = LocalNz; - /// Loop over variables for (const auto& var : var_list) { - if (var->is3D()) { - // 3D variable - auto* var3d_ref_ptr = dynamic_cast(var); + switch (var->field_type()) { + case field3d: { + const auto* var3d_ref_ptr = dynamic_cast(var); ASSERT0(var3d_ref_ptr != nullptr); - auto& var3d_ref = *var3d_ref_ptr; + const auto& var3d_ref = *var3d_ref_ptr; ASSERT2(var3d_ref.isAllocated()); - for (int jx = xge; jx != xlt; jx++) { + for (int jx = xge; jx < xlt; jx++) { for (int jy = yge; jy < ylt; jy++) { - for (int jz = 0; jz < LocalNz; jz++, len++) { + for (int jz = zge; jz < zlt; jz++, len++) { buffer[len] = var3d_ref(jx, jy, jz); } } } - } else { - // 2D variable - auto* var2d_ref_ptr = dynamic_cast(var); + break; + } + case field2d: { + const auto* var2d_ref_ptr = dynamic_cast(var); ASSERT0(var2d_ref_ptr != nullptr); - auto& var2d_ref = *var2d_ref_ptr; + const auto& var2d_ref = *var2d_ref_ptr; ASSERT2(var2d_ref.isAllocated()); - for (int jx = xge; jx != xlt; jx++) { + for (int jx = xge; jx < xlt; jx++) { for (int jy = yge; jy < ylt; jy++, len++) { buffer[len] = var2d_ref(jx, jy); } } + break; + } + case fieldperp: { + const auto* varperp_ref_ptr = dynamic_cast(var); + ASSERT0(varperp_ref_ptr != nullptr); + const auto& varperp_ref = *varperp_ref_ptr; + ASSERT2(varperp_ref.isAllocated()); + for (int jx = xge; jx < xlt; jx++) { + for (int jz = zge; jz < zlt; jz++, len++) { + buffer[len] = varperp_ref(jx, jz); + } + } + break; + } } } - return (len); + return len; } -int BoutMesh::unpack_data(const std::vector& var_list, int xge, int xlt, - int yge, int ylt, BoutReal* buffer) { +int BoutMesh::unpack_data(const std::vector& var_list, int xge, int xlt, int yge, + int ylt, const BoutReal* buffer) const { + using enum Field::FieldType; int len = 0; + const int zge = 0; + const int zlt = LocalNz; - /// Loop over variables for (const auto& var : var_list) { - if (var->is3D()) { - // 3D variable - auto& var3d_ref = *dynamic_cast(var); - for (int jx = xge; jx != xlt; jx++) { + switch (var->field_type()) { + case field3d: { + auto* var3d_ref_ptr = dynamic_cast(var); + ASSERT0(var3d_ref_ptr != nullptr); + auto& var3d_ref = *var3d_ref_ptr; + ASSERT2(var3d_ref.isAllocated()); + for (int jx = xge; jx < xlt; jx++) { for (int jy = yge; jy < ylt; jy++) { - for (int jz = 0; jz < LocalNz; jz++, len++) { + for (int jz = zge; jz < zlt; jz++, len++) { var3d_ref(jx, jy, jz) = buffer[len]; } } } - } else { - // 2D variable - auto& var2d_ref = *dynamic_cast(var); - for (int jx = xge; jx != xlt; jx++) { + break; + } + case field2d: { + auto* var2d_ref_ptr = dynamic_cast(var); + ASSERT0(var2d_ref_ptr != nullptr); + auto& var2d_ref = *var2d_ref_ptr; + ASSERT2(var2d_ref.isAllocated()); + for (int jx = xge; jx < xlt; jx++) { for (int jy = yge; jy < ylt; jy++, len++) { var2d_ref(jx, jy) = buffer[len]; } } + break; + } + case fieldperp: { + auto* varperp_ref_ptr = dynamic_cast(var); + ASSERT0(varperp_ref_ptr != nullptr); + auto& varperp_ref = *varperp_ref_ptr; + ASSERT2(varperp_ref.isAllocated()); + for (int jx = xge; jx < xlt; jx++) { + for (int jz = zge; jz < zlt; jz++, len++) { + varperp_ref(jx, jz) = buffer[len]; + } + } + break; + } } } - return (len); + return len; } /**************************************************************** diff --git a/src/mesh/impls/bout/boutmesh.hxx b/src/mesh/impls/bout/boutmesh.hxx index b42bf325b5..cfb91c64da 100644 --- a/src/mesh/impls/bout/boutmesh.hxx +++ b/src/mesh/impls/bout/boutmesh.hxx @@ -4,6 +4,7 @@ #include "mpi.h" +#include "bout/bout_types.hxx" #include "bout/unused.hxx" #include @@ -12,6 +13,8 @@ #include #include +class Field; + /// Implementation of Mesh (mostly) compatible with BOUT /// /// Topology and communications compatible with BOUT @@ -476,12 +479,11 @@ private: void post_receiveY(CommHandle& ch); /// Take data from objects and put into a buffer - int pack_data(const std::vector& var_list, int xge, int xlt, int yge, - int ylt, BoutReal* buffer); + int pack_data(const std::vector& var_list, int xge, int xlt, int yge, int ylt, + BoutReal* buffer) const; /// Copy data from a buffer back into the fields - - int unpack_data(const std::vector& var_list, int xge, int xlt, int yge, - int ylt, BoutReal* buffer); + int unpack_data(const std::vector& var_list, int xge, int xlt, int yge, int ylt, + const BoutReal* buffer) const; }; namespace { diff --git a/src/mesh/mesh.cxx b/src/mesh/mesh.cxx index 829c150a15..585c0f7257 100644 --- a/src/mesh/mesh.cxx +++ b/src/mesh/mesh.cxx @@ -376,38 +376,28 @@ void Mesh::communicate(FieldGroup& g) { } } -/// This is a bit of a hack for now to get FieldPerp communications -/// The FieldData class needs to be changed to accomodate FieldPerp objects -void Mesh::communicate(FieldPerp& f) { - comm_handle recv[2]; - - int nin = xstart; // Number of x points in inner guard cell - int nout = LocalNx - xend - 1; // Number of x points in outer guard cell - - // Post receives for guard cell regions - - recv[0] = irecvXIn(f[0], nin * LocalNz, 0); - recv[1] = irecvXOut(f[xend + 1], nout * LocalNz, 1); +int Mesh::msg_len(const std::vector& var_list, int xge, int xlt, int yge, + int ylt) const { + int len = 0; - // Send data - sendXIn(f[xstart], nin * LocalNz, 1); - sendXOut(f[xend - nout + 1], nout * LocalNz, 0); + using enum Field::FieldType; - // Wait for receive - wait(recv[0]); - wait(recv[1]); -} - -int Mesh::msg_len(const std::vector& var_list, int xge, int xlt, int yge, - int ylt) { - int len = 0; + const auto x_length = xlt - xge; + const auto y_length = ylt - yge; + const auto z_length = LocalNz; /// Loop over variables for (const auto& var : var_list) { - if (var->is3D()) { - len += (xlt - xge) * (ylt - yge) * LocalNz * var->elementSize(); - } else { - len += (xlt - xge) * (ylt - yge) * var->elementSize(); + switch (var->field_type()) { + case field3d: + len += x_length * y_length * z_length * var->elementSize(); + break; + case field2d: + len += x_length * y_length * var->elementSize(); + break; + case fieldperp: + len += x_length * z_length * var->elementSize(); + break; } } diff --git a/tests/unit/fake_parallel_mesh.hxx b/tests/unit/fake_parallel_mesh.hxx index 80406fb676..b3e57c4968 100644 --- a/tests/unit/fake_parallel_mesh.hxx +++ b/tests/unit/fake_parallel_mesh.hxx @@ -6,8 +6,10 @@ #include #include #include +#include #include "../../src/mesh/impls/bout/boutmesh.hxx" +#include "bout/assert.hxx" #include "bout/boundary_op.hxx" #include "bout/boundary_region.hxx" #include "bout/boutcomm.hxx" @@ -129,47 +131,15 @@ public: return parentSendY(g, handle); } - // Need to override this functions to trick mesh into communicating for - // FieldPerp type - void communicate(FieldPerp& f) override { - int nin = xstart; // Number of x points in inner guard cell - int nout = LocalNx - xend - 1; // Number of x points in outer guard cell - - if (registeredFieldPerps.count(&f) != 0) { - int id = registeredFieldPerps[&f]; - - if (xInMesh != nullptr && xInMesh->registeredFieldPerpIds.count(id) != 0) { - FieldPerp* xInField = xInMesh->registeredFieldPerpIds[id]; - for (int i = 0; i < nin * LocalNz; i++) { - IndPerp ind(i, 1, LocalNz); - f[ind] = (*xInField)[ind.xp(xend - xstart + 1)]; - } - } - - if (xOutMesh != nullptr && xOutMesh->registeredFieldPerpIds.count(id) != 0) { - FieldPerp* xOutField = xOutMesh->registeredFieldPerpIds[id]; - for (int i = 0; i < nout * LocalNz; i++) { - IndPerp ind((xend + 1) * LocalNz + i, 1, LocalNz); - f[ind] = (*xOutField)[ind.xm(xend - xstart + 1)]; - } - } - // No corner cells to communicate for FieldPerp - } - } - /// Use these methods to let the mesh know that this field has been /// created with it. It can then check in with its sibling meshes /// (representing other processors) to see if a corresponding field /// has been created for them which can be used to communicate guard /// cells with. - void registerField(FieldData& f, int id) { + void registerField(Field& f, int id) { registeredFields.emplace(&f, id); registeredFieldIds.emplace(id, &f); } - void registerField(FieldPerp& f, int id) { - registeredFieldPerps.emplace(&f, id); - registeredFieldPerpIds.emplace(id, &f); - } friend std::vector createFakeProcessors(int nx, int ny, int nz, int nxpe, int nype); @@ -257,10 +227,8 @@ public: private: FakeParallelMesh *yUpMesh, *yDownMesh, *xInMesh, *xOutMesh; bool communicatingX = false, communicatingY = false; - std::map registeredFields; - std::map registeredFieldIds; - std::map registeredFieldPerps; - std::map registeredFieldPerpIds; + std::map registeredFields; + std::map registeredFieldIds; std::unique_ptr mpiSmart; comm_handle parentSendX(FieldGroup& g, comm_handle handle, bool disable_corners) { @@ -270,10 +238,10 @@ private: return BoutMesh::sendY(g, handle); } - FieldGroup makeGroup(FakeParallelMesh* m, const std::vector ids) { + static FieldGroup makeGroup(FakeParallelMesh* m, const std::vector& ids) { FieldGroup g; - for (int i : ids) { - ASSERT1(m->registeredFieldIds.count(i) != 0); + for (const int i : ids) { + ASSERT1(m->registeredFieldIds.contains(i)); g.add(*m->registeredFieldIds[i]); } return g; diff --git a/tests/unit/field/test_field.cxx b/tests/unit/field/test_field.cxx index 40d675e5b3..2a3e600ab8 100644 --- a/tests/unit/field/test_field.cxx +++ b/tests/unit/field/test_field.cxx @@ -24,8 +24,8 @@ class FieldSubClass : public Field { FieldSubClass(Mesh* localmesh, CELL_LOC location_in, DirectionTypes directions_in) : Field(localmesh, location_in, directions_in) {} - bool is3D() const override { return false; } int size() const override { return 42; } + FieldType field_type() const override { return FieldType::field2d; } }; } // namespace diff --git a/tests/unit/field/test_vector2d.cxx b/tests/unit/field/test_vector2d.cxx index de662c7df8..674817abd2 100644 --- a/tests/unit/field/test_vector2d.cxx +++ b/tests/unit/field/test_vector2d.cxx @@ -97,12 +97,6 @@ TEST_F(Vector2DTest, ApplyBoundaryString) { EXPECT_DOUBLE_EQ(v.x(2, 2), 0.0); } -TEST_F(Vector2DTest, Is3D) { - Vector2D vector; - - EXPECT_FALSE(vector.is3D()); -} - TEST_F(Vector2DTest, BoutRealSize) { Vector2D vector; diff --git a/tests/unit/field/test_vector3d.cxx b/tests/unit/field/test_vector3d.cxx index 4db1c750a1..d8e9366574 100644 --- a/tests/unit/field/test_vector3d.cxx +++ b/tests/unit/field/test_vector3d.cxx @@ -100,12 +100,6 @@ TEST_F(Vector3DTest, ApplyBoundaryString) { EXPECT_DOUBLE_EQ(v.x(2, 2, 1), 0.0); } -TEST_F(Vector3DTest, Is3D) { - Vector3D vector; - - EXPECT_TRUE(vector.is3D()); -} - TEST_F(Vector3DTest, BoutRealSize) { Vector3D vector; diff --git a/tests/unit/mesh/test_mesh.cxx b/tests/unit/mesh/test_mesh.cxx index 4d0c9111dd..8d64f1d0a3 100644 --- a/tests/unit/mesh/test_mesh.cxx +++ b/tests/unit/mesh/test_mesh.cxx @@ -1,5 +1,6 @@ #include "gmock/gmock.h" #include "gtest/gtest.h" +#include #include "bout/boutexception.hxx" #include "bout/mesh.hxx" @@ -288,9 +289,9 @@ TEST_F(MeshTest, MsgLen) { Field2D f2D_1(0., &localmesh); Field2D f2D_2(0., &localmesh); - std::vector var_list{&f3D_1, &f2D_1, &f3D_2, &f2D_2}; + const std::vector var_list{&f3D_1, &f2D_1, &f3D_2, &f2D_2}; const int len = localmesh.msg_len(var_list, 0, nx, 0, ny); - EXPECT_EQ(len, 2 * (nx * ny * nz) + 2 * (nx * ny)); + EXPECT_EQ(len, (2 * (nx * ny * nz)) + (2 * (nx * ny))); }