Skip to content
Snippets Groups Projects
Select Git revision
  • 71a0343dcd773af5a221e28e8dc75bc7bcedfd74
  • develop default protected
  • feature/gmsh-reader
  • origin/stage/bouguettaia
  • feature/kinetic-schemes
  • feature/reconstruction
  • feature/local-dt-fsi
  • feature/composite-scheme-sources
  • feature/composite-scheme-other-fluxes
  • feature/serraille
  • feature/variational-hydro
  • feature/composite-scheme
  • hyperplastic
  • feature/polynomials
  • feature/gks
  • feature/implicit-solver-o2
  • feature/coupling_module
  • feature/implicit-solver
  • feature/merge-local-dt-fsi
  • master protected
  • feature/escobar-smoother
  • v0.5.0 protected
  • v0.4.1 protected
  • v0.4.0 protected
  • v0.3.0 protected
  • v0.2.0 protected
  • v0.1.0 protected
  • Kidder
  • v0.0.4 protected
  • v0.0.3 protected
  • v0.0.2 protected
  • v0 protected
  • v0.0.1 protected
33 results

ConnectivityComputer.hpp

Blame
  • test_Messenger.cpp 14.51 KiB
    #include <catch2/catch.hpp>
    
    #include <utils/Array.hpp>
    #include <utils/Messenger.hpp>
    
    #include <utils/pugs_config.hpp>
    
    #include <fstream>
    
    #ifdef PUGS_HAS_MPI
    #include <mpi.h>
    #define IF_MPI(INSTRUCTION) INSTRUCTION
    #else
    #define IF_MPI(INSTRUCTION)
    #endif   // PUGS_HAS_MPI
    
    namespace mpi_check
    {
    struct integer
    {
      int m_int;
    
      operator int&()
      {
        return m_int;
      }
    
      operator const int&() const
      {
        return m_int;
      }
    
      integer&
      operator=(int i)
      {
        m_int = i;
        return *this;
      }
    
      template <typename T>
      bool
      operator==(const T& t) const
      {
        return m_int == static_cast<int>(t);
      }
    };
    
    struct tri_int
    {
      int m_int_0;
      int m_int_1;
      int m_int_2;
    
      bool
      operator==(tri_int t) const
      {
        return ((m_int_0 == t.m_int_0) and (m_int_1 == t.m_int_1) and (m_int_2 == t.m_int_2));
      }
    };
    
    template <typename T>
    void
    test_allToAll()
    {
      Array<T> data_array(parallel::size());
      for (size_t i = 0; i < data_array.size(); ++i) {
        data_array[i] = parallel::rank();
      }
      auto exchanged_array = parallel::allToAll(data_array);
    
      for (size_t i = 0; i < data_array.size(); ++i) {
        const size_t value = exchanged_array[i];
        REQUIRE(value == i);
      }
    }
    
    template <>
    void
    test_allToAll<bool>()
    {
      Array<bool> data_array(parallel::size());
      for (size_t i = 0; i < data_array.size(); ++i) {
        data_array[i] = ((parallel::rank() % 2) == 0);
      }
      auto exchanged_array = parallel::allToAll(data_array);
    
      for (size_t i = 0; i < data_array.size(); ++i) {
        REQUIRE(exchanged_array[i] == ((i % 2) == 0));
      }
    }
    
    template <>
    void
    test_allToAll<tri_int>()
    {
      Array<tri_int> data_array(parallel::size());
      for (size_t i = 0; i < data_array.size(); ++i) {
        const int val = 1 + parallel::rank();
        data_array[i] = tri_int{val, 2 * val, val + 3};
      }
      auto exchanged_array = parallel::allToAll(data_array);
    
      for (size_t i = 0; i < data_array.size(); ++i) {
        const int val = 1 + i;
        REQUIRE(exchanged_array[i] == tri_int{val, 2 * val, val + 3});
      }
    }
    
    }   // namespace mpi_check
    
    // clazy:excludeall=non-pod-global-static
    
    TEST_CASE("Messenger", "[mpi]")
    {
      SECTION("communication info")
      {
        int rank = 0;
        IF_MPI(MPI_Comm_rank(MPI_COMM_WORLD, &rank));
        REQUIRE(rank == static_cast<int>(parallel::rank()));
    
        int size = 1;
        IF_MPI(MPI_Comm_size(MPI_COMM_WORLD, &size));
        REQUIRE(size == static_cast<int>(parallel::size()));
      }
    
      SECTION("reduction")
      {
        const int min_value = parallel::allReduceMin(parallel::rank() + 3);
        REQUIRE(min_value == 3);
    
        const int max_value = parallel::allReduceMax(parallel::rank() + 3);
        REQUIRE(max_value == static_cast<int>((parallel::size() - 1) + 3));
    
        const bool and_value = parallel::allReduceAnd(true);
        REQUIRE(and_value == true);
    
        const bool and_value_2 = parallel::allReduceAnd(parallel::rank() > 0);
        REQUIRE(and_value_2 == false);
    
        const bool or_value = parallel::allReduceOr(false);
        REQUIRE(or_value == false);
    
        const bool or_value_2 = parallel::allReduceOr(parallel::rank() > 0);
        REQUIRE(or_value_2 == (parallel::size() > 1));
      }
    
      SECTION("all to all")
      {
        // chars
        mpi_check::test_allToAll<char>();
        mpi_check::test_allToAll<wchar_t>();
    
        // integers
        mpi_check::test_allToAll<int8_t>();
        mpi_check::test_allToAll<int16_t>();
        mpi_check::test_allToAll<int32_t>();
        mpi_check::test_allToAll<int64_t>();
        mpi_check::test_allToAll<uint8_t>();
        mpi_check::test_allToAll<uint16_t>();
        mpi_check::test_allToAll<uint32_t>();
        mpi_check::test_allToAll<uint64_t>();
        mpi_check::test_allToAll<signed long long int>();
        mpi_check::test_allToAll<unsigned long long int>();
    
        // floats
        mpi_check::test_allToAll<float>();
        mpi_check::test_allToAll<double>();
        mpi_check::test_allToAll<long double>();
    
        // bools
        mpi_check::test_allToAll<bool>();
    
        // trivial simple type
        mpi_check::test_allToAll<mpi_check::integer>();
    
        // compound trivial type
        mpi_check::test_allToAll<mpi_check::tri_int>();
    
    #ifndef NDEBUG
        SECTION("checking invalid all to all")
        {
          if (parallel::size() > 1) {
            Array<int> invalid_all_to_all(parallel::size() + 1);
            REQUIRE_THROWS_AS(parallel::allToAll(invalid_all_to_all), AssertError);
    
            Array<int> different_size_all_to_all(parallel::size() * (parallel::rank() + 1));
            REQUIRE_THROWS_AS(parallel::allToAll(different_size_all_to_all), AssertError);
          }
        }
    #endif   // NDEBUG
      }
    
      SECTION("broadcast value")
      {
        {
          // simple type
          size_t value{(3 + parallel::rank()) * 2};
          parallel::broadcast(value, 0);
          REQUIRE(value == 6);
        }
    
        {
          // trivial simple type
          mpi_check::integer value{static_cast<int>((3 + parallel::rank()) * 2)};
          parallel::broadcast(value, 0);
          REQUIRE((value == 6));
        }
    
        {
          // compound trivial type
          mpi_check::tri_int value{static_cast<int>((3 + parallel::rank()) * 2), static_cast<int>(2 + parallel::rank()),
                                   static_cast<int>(4 - parallel::rank())};
          parallel::broadcast(value, 0);
          REQUIRE((value == mpi_check::tri_int{6, 2, 4}));
        }
      }
    
      SECTION("broadcast array")
      {
        {
          // simple type
          Array<size_t> array(3);
          array[0] = (3 + parallel::rank()) * 2;
          array[1] = 2 + parallel::rank();
          array[2] = 4 - parallel::rank();
          parallel::broadcast(array, 0);
          REQUIRE(((array[0] == 6) and (array[1] == 2) and (array[2] == 4)));
        }
    
        {
          // trivial simple type
          Array<mpi_check::integer> array(3);
          array[0] = static_cast<int>((3 + parallel::rank()) * 2);
          array[1] = static_cast<int>(2 + parallel::rank());
          array[2] = static_cast<int>(4 - parallel::rank());
          parallel::broadcast(array, 0);
          REQUIRE(((array[0] == 6) and (array[1] == 2) and (array[2] == 4)));
        }
    
        {
          // compound trivial type
          Array<mpi_check::tri_int> array(3);
          array[0] = mpi_check::tri_int{static_cast<int>((3 + parallel::rank()) * 2),
                                        static_cast<int>(2 + parallel::rank()), static_cast<int>(4 - parallel::rank())};
          array[1] = mpi_check::tri_int{static_cast<int>((2 + parallel::rank()) * 4),
                                        static_cast<int>(3 + parallel::rank()), static_cast<int>(1 - parallel::rank())};
          array[2] = mpi_check::tri_int{static_cast<int>((5 + parallel::rank())), static_cast<int>(-3 + parallel::rank()),
                                        static_cast<int>(parallel::rank())};
          parallel::broadcast(array, 0);
          REQUIRE(((array[0] == mpi_check::tri_int{6, 2, 4}) and (array[1] == mpi_check::tri_int{8, 3, 1}) and
                   (array[2] == mpi_check::tri_int{5, -3, 0})));
        }
      }
    
      SECTION("all gather value")
      {
        {
          // simple type
          size_t value{(3 + parallel::rank()) * 2};
          Array<size_t> gather_array = parallel::allGather(value);
          REQUIRE(gather_array.size() == parallel::size());
    
          for (size_t i = 0; i < gather_array.size(); ++i) {
            REQUIRE((gather_array[i] == (3 + i) * 2));
          }
        }
    
        {
          // trivial simple type
          mpi_check::integer value{static_cast<int>((3 + parallel::rank()) * 2)};
          Array<mpi_check::integer> gather_array = parallel::allGather(value);
          REQUIRE(gather_array.size() == parallel::size());
    
          for (size_t i = 0; i < gather_array.size(); ++i) {
            const int expected_value = (3 + i) * 2;
            REQUIRE(gather_array[i] == expected_value);
          }
        }
    
        {
          // compound trivial type
          mpi_check::tri_int value{static_cast<int>((3 + parallel::rank()) * 2), static_cast<int>(2 + parallel::rank()),
                                   static_cast<int>(4 - parallel::rank())};
          Array<mpi_check::tri_int> gather_array = parallel::allGather(value);
    
          REQUIRE(gather_array.size() == parallel::size());
          for (size_t i = 0; i < gather_array.size(); ++i) {
            mpi_check::tri_int expected_value{static_cast<int>((3 + i) * 2), static_cast<int>(2 + i),
                                              static_cast<int>(4 - i)};
            REQUIRE((gather_array[i] == expected_value));
          }
        }
      }
    
      SECTION("all gather array")
      {
        {
          // simple type
          Array<int> array(3);
          for (size_t i = 0; i < array.size(); ++i) {
            array[i] = (3 + parallel::rank()) * 2 + i;
          }
          Array<int> gather_array = parallel::allGather(array);
          REQUIRE(gather_array.size() == array.size() * parallel::size());
    
          for (size_t i = 0; i < gather_array.size(); ++i) {
            const int expected_value = (3 + i / array.size()) * 2 + (i % array.size());
            REQUIRE((gather_array[i] == expected_value));
          }
        }
    
        {
          // trivial simple type
          Array<mpi_check::integer> array(3);
          for (size_t i = 0; i < array.size(); ++i) {
            array[i] = (3 + parallel::rank()) * 2 + i;
          }
          Array<mpi_check::integer> gather_array = parallel::allGather(array);
          REQUIRE(gather_array.size() == array.size() * parallel::size());
    
          for (size_t i = 0; i < gather_array.size(); ++i) {
            const int expected_value = (3 + i / array.size()) * 2 + (i % array.size());
            REQUIRE((gather_array[i] == expected_value));
          }
        }
    
        {
          // compound trivial type
          Array<mpi_check::tri_int> array(3);
          for (size_t i = 0; i < array.size(); ++i) {
            array[i] =
              mpi_check::tri_int{static_cast<int>((3 + parallel::rank()) * 2), static_cast<int>(2 + parallel::rank() + i),
                                 static_cast<int>(4 - parallel::rank() - i)};
          }
          Array<mpi_check::tri_int> gather_array = parallel::allGather(array);
    
          REQUIRE(gather_array.size() == array.size() * parallel::size());
          for (size_t i = 0; i < gather_array.size(); ++i) {
            mpi_check::tri_int expected_value{static_cast<int>((3 + i / array.size()) * 2),
                                              static_cast<int>(2 + i / array.size() + (i % array.size())),
                                              static_cast<int>(4 - i / array.size() - (i % array.size()))};
            REQUIRE((gather_array[i] == expected_value));
          }
        }
      }
    
      SECTION("all array exchanges")
      {
        {   // simple type
          std::vector<Array<const size_t>> send_array_list(parallel::size());
          for (size_t i = 0; i < send_array_list.size(); ++i) {
            Array<size_t> send_array(i + 1);
            for (size_t j = 0; j < send_array.size(); ++j) {
              send_array[j] = (parallel::rank() + 1) * j;
            }
            send_array_list[i] = send_array;
          }
    
          std::vector<Array<size_t>> recv_array_list(parallel::size());
          for (size_t i = 0; i < recv_array_list.size(); ++i) {
            recv_array_list[i] = Array<size_t>(parallel::rank() + 1);
          }
          parallel::exchange(send_array_list, recv_array_list);
    
          for (size_t i = 0; i < parallel::size(); ++i) {
            const Array<const size_t> recv_array = recv_array_list[i];
            for (size_t j = 0; j < recv_array.size(); ++j) {
              REQUIRE(recv_array[j] == (i + 1) * j);
            }
          }
        }
    
        {   //  trivial simple type
          std::vector<Array<mpi_check::integer>> send_array_list(parallel::size());
          for (size_t i = 0; i < send_array_list.size(); ++i) {
            Array<mpi_check::integer> send_array(i + 1);
            for (size_t j = 0; j < send_array.size(); ++j) {
              send_array[j] = static_cast<int>((parallel::rank() + 1) * j);
            }
            send_array_list[i] = send_array;
          }
    
          std::vector<Array<mpi_check::integer>> recv_array_list(parallel::size());
          for (size_t i = 0; i < recv_array_list.size(); ++i) {
            recv_array_list[i] = Array<mpi_check::integer>(parallel::rank() + 1);
          }
          parallel::exchange(send_array_list, recv_array_list);
    
          for (size_t i = 0; i < parallel::size(); ++i) {
            const Array<const mpi_check::integer> recv_array = recv_array_list[i];
            for (size_t j = 0; j < recv_array.size(); ++j) {
              REQUIRE(recv_array[j] == (i + 1) * j);
            }
          }
        }
    
        {
          // compound trivial type
          std::vector<Array<mpi_check::tri_int>> send_array_list(parallel::size());
          for (size_t i = 0; i < send_array_list.size(); ++i) {
            Array<mpi_check::tri_int> send_array(i + 1);
            for (size_t j = 0; j < send_array.size(); ++j) {
              send_array[j] = mpi_check::tri_int{static_cast<int>((parallel::rank() + 1) * j),
                                                 static_cast<int>(parallel::rank()), static_cast<int>(j)};
            }
            send_array_list[i] = send_array;
          }
    
          std::vector<Array<mpi_check::tri_int>> recv_array_list(parallel::size());
          for (size_t i = 0; i < recv_array_list.size(); ++i) {
            recv_array_list[i] = Array<mpi_check::tri_int>(parallel::rank() + 1);
          }
          parallel::exchange(send_array_list, recv_array_list);
    
          for (size_t i = 0; i < parallel::size(); ++i) {
            const Array<const mpi_check::tri_int> recv_array = recv_array_list[i];
            for (size_t j = 0; j < recv_array.size(); ++j) {
              mpi_check::tri_int expected_value{static_cast<int>((i + 1) * j), static_cast<int>(i), static_cast<int>(j)};
              REQUIRE((recv_array[j] == expected_value));
            }
          }
        }
      }
    
    #ifndef NDEBUG
      SECTION("checking all array exchange invalid sizes")
      {
        std::vector<Array<const int>> send_array_list(parallel::size());
        for (size_t i = 0; i < send_array_list.size(); ++i) {
          Array<int> send_array(i + 1);
          send_array.fill(parallel::rank());
          send_array_list[i] = send_array;
        }
    
        std::vector<Array<int>> recv_array_list(parallel::size());
        REQUIRE_THROWS_AS(parallel::exchange(send_array_list, recv_array_list), AssertError);
      }
    #endif   // NDEBUG
    
      SECTION("checking barrier")
      {
        for (size_t i = 0; i < parallel::size(); ++i) {
          if (i == parallel::rank()) {
            std::ofstream file;
            if (i == 0) {
              file.open("barrier_test", std::ios_base::out);
            } else {
              file.open("barrier_test", std::ios_base::app);
            }
            file << i << "\n" << std::flush;
          }
          parallel::barrier();
        }
    
        {   // reading produced file
          std::ifstream file("barrier_test");
          std::vector<size_t> number_list;
          while (file) {
            size_t value;
            file >> value;
            if (file) {
              number_list.push_back(value);
            }
          }
          REQUIRE(number_list.size() == parallel::size());
          for (size_t i = 0; i < number_list.size(); ++i) {
            REQUIRE(number_list[i] == i);
          }
        }
        parallel::barrier();
    
        std::remove("barrier_test");
      }
    
      SECTION("errors")
      {
        int argc    = 0;
        char** argv = nullptr;
        REQUIRE_THROWS_WITH((parallel::Messenger::create(argc, argv)), "unexpected error: Messenger already created");
      }
    }