From: Timo Heister Date: Mon, 11 Jul 2022 18:30:30 +0000 (-0400) Subject: parallel intermediate out: compression levels X-Git-Tag: v9.5.0-rc1~1072^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=003613576b4527a8de6cfaeb37ef17ed32a50abb;p=dealii.git parallel intermediate out: compression levels - support different compression levels - add a test --- diff --git a/source/base/data_out_base.cc b/source/base/data_out_base.cc index 67166f65b3..9caf6a6033 100644 --- a/source/base/data_out_base.cc +++ b/source/base/data_out_base.cc @@ -113,6 +113,30 @@ namespace } } + /** + * Convert between the enum specified inside VtkFlags and the preprocessor + * constant defined by boost::iostreams::zlib. + */ + int + get_boost_zlib_compression_level( + const DataOutBase::VtkFlags::ZlibCompressionLevel level) + { + switch (level) + { + case (DataOutBase::VtkFlags::no_compression): + return boost::iostreams::zlib::no_compression; + case (DataOutBase::VtkFlags::best_speed): + return boost::iostreams::zlib::best_speed; + case (DataOutBase::VtkFlags::best_compression): + return boost::iostreams::zlib::best_compression; + case (DataOutBase::VtkFlags::default_compression): + return boost::iostreams::zlib::default_compression; + default: + Assert(false, ExcNotImplemented()); + return boost::iostreams::zlib::no_compression; + } + } + /** * Do a zlib compression followed by a base64 encoding of the given data. The * result is then written to the given stream. @@ -7619,7 +7643,8 @@ namespace DataOutBase if (compression != VtkFlags::no_compression) # ifdef DEAL_II_WITH_ZLIB - f.push(boost::iostreams::zlib_compressor()); + f.push(boost::iostreams::zlib_compressor( + get_boost_zlib_compression_level(compression))); # else AssertThrow( false, diff --git a/tests/data_out/data_out_intermediate_parallel_02.cc b/tests/data_out/data_out_intermediate_parallel_02.cc new file mode 100644 index 0000000000..1b67a6ca1c --- /dev/null +++ b/tests/data_out/data_out_intermediate_parallel_02.cc @@ -0,0 +1,93 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2022 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::write_deal_II_intermediate_in_parallel() and +// DataOutReader::read_whole_parallel_file() with compression + +#include + +#include +#include +#include + +#include + +#include +#include +#include +#include + +#include + +#include + +#include "../tests.h" + + +template +void +check() +{ + Triangulation tria; + GridGenerator::hyper_cube(tria, 0., 1.); + tria.refine_global(1); + + Vector cell_data(tria.n_active_cells()); + for (unsigned int i = 0; i < tria.n_active_cells(); ++i) + cell_data(i) = i * 1.0; + + FE_Q fe(1); + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + + Vector x(dof_handler.n_dofs()); + + DataOut data_out; + data_out.attach_dof_handler(dof_handler); + data_out.add_data_vector(cell_data, + "cell_data", + DataOut::type_cell_data); + data_out.add_data_vector(x, "solution"); + + data_out.build_patches(); + + data_out.write_deal_II_intermediate_in_parallel( + "test.pd2", MPI_COMM_WORLD, DataOutBase::VtkFlags::best_compression); + + const unsigned int my_rank = + dealii::Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + if (my_rank == 0) + { + // Read the data back in and dump it into the deallog: + std::ifstream in("test.pd2"); + Assert(in, dealii::ExcIO()); + DataOutReader reader; + reader.read_whole_parallel_file(in); + reader.write_deal_II_intermediate(deallog.get_file_stream()); + } + + deallog << "OK" << std::endl; +} + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll all; + + check<2>(); + + return 0; +} diff --git a/tests/data_out/data_out_intermediate_parallel_02.with_zlib=on.mpirun=2.output b/tests/data_out/data_out_intermediate_parallel_02.with_zlib=on.mpirun=2.output new file mode 100644 index 0000000000..5bd7169e4e --- /dev/null +++ b/tests/data_out/data_out_intermediate_parallel_02.with_zlib=on.mpirun=2.output @@ -0,0 +1,95 @@ + +2 2 +[deal.II intermediate format graphics data] +[written by deal.II x.y.z] +[Version: 4] +2 +solution +cell_data +8 +[deal.II intermediate Patch<2,2>] +3 +0.00000 0.00000 0.500000 0.00000 0.00000 0.500000 0.500000 0.500000 +4294967295 1 4294967295 2 +0 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 + + +[deal.II intermediate Patch<2,2>] +3 +0.500000 0.00000 1.00000 0.00000 0.500000 0.500000 1.00000 0.500000 +0 4294967295 4294967295 3 +1 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 + + +[deal.II intermediate Patch<2,2>] +3 +0.00000 0.500000 0.500000 0.500000 0.00000 1.00000 0.500000 1.00000 +4294967295 3 0 4294967295 +2 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 2.00000 2.00000 2.00000 2.00000 + + +[deal.II intermediate Patch<2,2>] +3 +0.500000 0.500000 1.00000 0.500000 0.500000 1.00000 1.00000 1.00000 +2 4294967295 1 4294967295 +3 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 3.00000 3.00000 3.00000 3.00000 + + +[deal.II intermediate Patch<2,2>] +3 +0.00000 0.00000 0.500000 0.00000 0.00000 0.500000 0.500000 0.500000 +4294967295 5 4294967295 6 +4 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000 + + +[deal.II intermediate Patch<2,2>] +3 +0.500000 0.00000 1.00000 0.00000 0.500000 0.500000 1.00000 0.500000 +4 4294967295 4294967295 7 +5 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 1.00000 1.00000 1.00000 1.00000 + + +[deal.II intermediate Patch<2,2>] +3 +0.00000 0.500000 0.500000 0.500000 0.00000 1.00000 0.500000 1.00000 +4294967295 7 4 4294967295 +6 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 2.00000 2.00000 2.00000 2.00000 + + +[deal.II intermediate Patch<2,2>] +3 +0.500000 0.500000 1.00000 0.500000 0.500000 1.00000 1.00000 1.00000 +6 4294967295 5 4294967295 +7 1 +0 +2 4 +0.00000 0.00000 0.00000 0.00000 3.00000 3.00000 3.00000 3.00000 + + +0 + +DEAL:0::OK + +DEAL:1::OK +