From 6d286cb0ec7bfd658d53ad1e609460cf407dbb8d Mon Sep 17 00:00:00 2001 From: Pasquale Africa Date: Thu, 3 Sep 2020 12:53:09 +0000 Subject: [PATCH] Test for HDF5 now runs in parallel --- tests/simplex/data_out_write_hdf5_01.cc | 85 ++++----- ...ith_simplex_support=on.with_hdf5=on.output | 13 -- ...ith_simplex_support=on.with_hdf5=on.output | 7 + tests/simplex/data_out_write_hdf5_02.cc | 164 ++++++++++++++++++ ...ith_simplex_support=on.with_hdf5=on.output | 7 + 5 files changed, 216 insertions(+), 60 deletions(-) delete mode 100644 tests/simplex/data_out_write_hdf5_01.mpirun=4.with_simplex_support=on.with_hdf5=on.output create mode 100644 tests/simplex/data_out_write_hdf5_01.with_simplex_support=on.with_hdf5=on.output create mode 100644 tests/simplex/data_out_write_hdf5_02.cc create mode 100644 tests/simplex/data_out_write_hdf5_02.mpirun=4.with_simplex_support=on.with_hdf5=on.output diff --git a/tests/simplex/data_out_write_hdf5_01.cc b/tests/simplex/data_out_write_hdf5_01.cc index 48f03976fb..6decd2d490 100644 --- a/tests/simplex/data_out_write_hdf5_01.cc +++ b/tests/simplex/data_out_write_hdf5_01.cc @@ -73,47 +73,44 @@ test(const FiniteElement &fe, const unsigned int n_components) static unsigned int counter = 0; - for (unsigned int n_subdivisions = 1; n_subdivisions <= 2; ++n_subdivisions) - { - DataOut data_out; + DataOut data_out; - data_out.attach_dof_handler(dof_handler); - data_out.add_data_vector(solution, "solution"); + data_out.attach_dof_handler(dof_handler); + data_out.add_data_vector(solution, "solution"); - data_out.build_patches(mapping, n_subdivisions); + data_out.build_patches(mapping); - const std::string output_basename("test." + std::to_string(dim) + "." + - std::to_string(counter++)); + const std::string output_basename("test." + std::to_string(dim) + "." + + std::to_string(counter++)); - DataOutBase::DataOutFilter data_filter( - DataOutBase::DataOutFilterFlags(true, true)); - data_out.write_filtered_data(data_filter); - data_out.write_hdf5_parallel(data_filter, - output_basename + ".h5", - MPI_COMM_SELF); + DataOutBase::DataOutFilter data_filter( + DataOutBase::DataOutFilterFlags(true, true)); + data_out.write_filtered_data(data_filter); + data_out.write_hdf5_parallel(data_filter, + output_basename + ".h5", + MPI_COMM_SELF); - std::vector xdmf_entries({data_out.create_xdmf_entry( - data_filter, output_basename + ".h5", 0, MPI_COMM_SELF)}); + std::vector xdmf_entries({data_out.create_xdmf_entry( + data_filter, output_basename + ".h5", 0, MPI_COMM_SELF)}); - data_out.write_xdmf_file(xdmf_entries, - output_basename + ".xdmf", - MPI_COMM_SELF); + data_out.write_xdmf_file(xdmf_entries, + output_basename + ".xdmf", + MPI_COMM_SELF); - data_out.clear(); + data_out.clear(); - // Sadly hdf5 is binary and we can not use hd5dump because it might - // not be in the path. At least we can make sure that both the xdmf and - // the h5 file are created. - std::ifstream h5((output_basename + ".h5").c_str()); - AssertThrow(h5.good(), ExcIO()); + // Sadly hdf5 is binary and we can not use hd5dump because it might + // not be in the path. At least we can make sure that both the xdmf and + // the h5 file are created. + std::ifstream h5((output_basename + ".h5").c_str()); + AssertThrow(h5.good(), ExcIO()); - std::ifstream xdmf((output_basename + ".xdmf").c_str()); - AssertThrow(h5.good(), ExcIO()); + std::ifstream xdmf((output_basename + ".xdmf").c_str()); + AssertThrow(h5.good(), ExcIO()); - deallog << "Files " << output_basename + ".h5" - << " and " << output_basename + ".xdmf" - << " created succesfully!" << std::endl; - } + deallog << "Files " << output_basename + ".h5" + << " and " << output_basename + ".xdmf" + << " created succesfully!" << std::endl; } int @@ -121,26 +118,20 @@ main(int argc, char **argv) { initlog(); - Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); - { const unsigned int dim = 2; - test(Simplex::FE_P(2) /*=degree*/, 1); - test(FESystem(Simplex::FE_P(2 /*=degree*/), dim), dim); - test(FESystem(Simplex::FE_P(2 /*=degree*/), - dim, - Simplex::FE_P(1 /*=degree*/), - 1), - dim + 1); + test(Simplex::FE_P(2), 1); + test(FESystem(Simplex::FE_P(2), dim), dim); + test( + FESystem(Simplex::FE_P(2), dim, Simplex::FE_P(1), 1), + dim + 1); } { const unsigned int dim = 3; - test(Simplex::FE_P(2) /*=degree*/, 1); - test(FESystem(Simplex::FE_P(2 /*=degree*/), dim), dim); - test(FESystem(Simplex::FE_P(2 /*=degree*/), - dim, - Simplex::FE_P(1 /*=degree*/), - 1), - dim + 1); + test(Simplex::FE_P(2), 1); + test(FESystem(Simplex::FE_P(2), dim), dim); + test( + FESystem(Simplex::FE_P(2), dim, Simplex::FE_P(1), 1), + dim + 1); } } diff --git a/tests/simplex/data_out_write_hdf5_01.mpirun=4.with_simplex_support=on.with_hdf5=on.output b/tests/simplex/data_out_write_hdf5_01.mpirun=4.with_simplex_support=on.with_hdf5=on.output deleted file mode 100644 index 553a7f5878..0000000000 --- a/tests/simplex/data_out_write_hdf5_01.mpirun=4.with_simplex_support=on.with_hdf5=on.output +++ /dev/null @@ -1,13 +0,0 @@ - -DEAL::Files test.2.0.h5 and test.2.0.xdmf created succesfully! -DEAL::Files test.2.1.h5 and test.2.1.xdmf created succesfully! -DEAL::Files test.2.2.h5 and test.2.2.xdmf created succesfully! -DEAL::Files test.2.3.h5 and test.2.3.xdmf created succesfully! -DEAL::Files test.2.4.h5 and test.2.4.xdmf created succesfully! -DEAL::Files test.2.5.h5 and test.2.5.xdmf created succesfully! -DEAL::Files test.3.0.h5 and test.3.0.xdmf created succesfully! -DEAL::Files test.3.1.h5 and test.3.1.xdmf created succesfully! -DEAL::Files test.3.2.h5 and test.3.2.xdmf created succesfully! -DEAL::Files test.3.3.h5 and test.3.3.xdmf created succesfully! -DEAL::Files test.3.4.h5 and test.3.4.xdmf created succesfully! -DEAL::Files test.3.5.h5 and test.3.5.xdmf created succesfully! \ No newline at end of file diff --git a/tests/simplex/data_out_write_hdf5_01.with_simplex_support=on.with_hdf5=on.output b/tests/simplex/data_out_write_hdf5_01.with_simplex_support=on.with_hdf5=on.output new file mode 100644 index 0000000000..a967ed464e --- /dev/null +++ b/tests/simplex/data_out_write_hdf5_01.with_simplex_support=on.with_hdf5=on.output @@ -0,0 +1,7 @@ + +DEAL::Files test.2.0.h5 and test.2.0.xdmf created succesfully! +DEAL::Files test.2.1.h5 and test.2.1.xdmf created succesfully! +DEAL::Files test.2.2.h5 and test.2.2.xdmf created succesfully! +DEAL::Files test.3.0.h5 and test.3.0.xdmf created succesfully! +DEAL::Files test.3.1.h5 and test.3.1.xdmf created succesfully! +DEAL::Files test.3.2.h5 and test.3.2.xdmf created succesfully! \ No newline at end of file diff --git a/tests/simplex/data_out_write_hdf5_02.cc b/tests/simplex/data_out_write_hdf5_02.cc new file mode 100644 index 0000000000..9b17353106 --- /dev/null +++ b/tests/simplex/data_out_write_hdf5_02.cc @@ -0,0 +1,164 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2020 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// Test DataOut with HDF5 for simplex meshes in parallel. + +#include + +#include + +#include +#include + +#include +#include + +#include +#include +#include + +#include +#include + +#include +#include + +#include "../tests.h" + +using namespace dealii; + +template +class RightHandSideFunction : public Function +{ +public: + RightHandSideFunction(const unsigned int n_components) + : Function(n_components) + {} + + virtual double + value(const Point &p, const unsigned int component = 0) const + { + return p[component % dim] * p[component % dim]; + } +}; + +template +void +test(const FiniteElement &fe, const unsigned int n_components) +{ + MPI_Comm comm = MPI_COMM_WORLD; + + parallel::fullydistributed::Triangulation tria(comm); + + { + Triangulation tria_serial; + GridGenerator::subdivided_hyper_cube_with_simplices(tria_serial, + dim == 2 ? 4 : 2); + + GridTools::partition_triangulation(Utilities::MPI::n_mpi_processes(comm), + tria_serial); + + auto construction_data = TriangulationDescription::Utilities:: + create_description_from_triangulation(tria_serial, comm); + + tria.create_triangulation(construction_data); + } + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + + IndexSet owned_dofs = dof_handler.locally_owned_dofs(); + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs(dof_handler, locally_relevant_dofs); + + LinearAlgebra::distributed::Vector solution; + + solution.reinit(owned_dofs, locally_relevant_dofs, comm); + + MappingFE mapping(Simplex::FE_P(1)); + + VectorTools::interpolate(mapping, + dof_handler, + RightHandSideFunction(n_components), + solution); + solution.update_ghost_values(); + + static unsigned int counter = 0; + + DataOut data_out; + + data_out.attach_dof_handler(dof_handler); + data_out.add_data_vector(solution, "solution"); + + data_out.build_patches(mapping); + + const std::string output_basename("test." + std::to_string(dim) + "." + + std::to_string(counter++)); + + DataOutBase::DataOutFilter data_filter( + DataOutBase::DataOutFilterFlags(true, true)); + data_out.write_filtered_data(data_filter); + data_out.write_hdf5_parallel(data_filter, output_basename + ".h5", comm); + + std::vector xdmf_entries({data_out.create_xdmf_entry( + data_filter, output_basename + ".h5", 0, comm)}); + + data_out.write_xdmf_file(xdmf_entries, output_basename + ".xdmf", comm); + + data_out.clear(); + + // Sadly HDF5 is binary and we can not use hd5dump because it might not be + // in the path. At least we can make sure that both the .xdmf and the .h5 + // files are created. + if (Utilities::MPI::this_mpi_process(comm) == 0) + { + std::ifstream h5((output_basename + ".h5").c_str()); + AssertThrow(h5.good(), ExcIO()); + + std::ifstream xdmf((output_basename + ".xdmf").c_str()); + AssertThrow(h5.good(), ExcIO()); + + deallog << "Files " << output_basename + ".h5" + << " and " << output_basename + ".xdmf" + << " created succesfully!" << std::endl; + } +} + +int +main(int argc, char **argv) +{ + initlog(); + + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + { + const unsigned int dim = 2; + test(Simplex::FE_P(2), 1); + test(FESystem(Simplex::FE_P(2), dim), dim); + test( + FESystem(Simplex::FE_P(2), dim, Simplex::FE_P(1), 1), + dim + 1); + } + { + const unsigned int dim = 3; + test(Simplex::FE_P(2), 1); + test(FESystem(Simplex::FE_P(2), dim), dim); + test( + FESystem(Simplex::FE_P(2), dim, Simplex::FE_P(1), 1), + dim + 1); + } +} diff --git a/tests/simplex/data_out_write_hdf5_02.mpirun=4.with_simplex_support=on.with_hdf5=on.output b/tests/simplex/data_out_write_hdf5_02.mpirun=4.with_simplex_support=on.with_hdf5=on.output new file mode 100644 index 0000000000..a967ed464e --- /dev/null +++ b/tests/simplex/data_out_write_hdf5_02.mpirun=4.with_simplex_support=on.with_hdf5=on.output @@ -0,0 +1,7 @@ + +DEAL::Files test.2.0.h5 and test.2.0.xdmf created succesfully! +DEAL::Files test.2.1.h5 and test.2.1.xdmf created succesfully! +DEAL::Files test.2.2.h5 and test.2.2.xdmf created succesfully! +DEAL::Files test.3.0.h5 and test.3.0.xdmf created succesfully! +DEAL::Files test.3.1.h5 and test.3.1.xdmf created succesfully! +DEAL::Files test.3.2.h5 and test.3.2.xdmf created succesfully! \ No newline at end of file -- 2.39.5