From: Peter Munch Date: Sat, 14 May 2022 12:00:39 +0000 (+0200) Subject: Add test X-Git-Tag: v9.4.0-rc1~84^2~1 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=53a836f5b9e4941dcc7cec9c75c6263f7d814b30;p=dealii.git Add test --- diff --git a/tests/distributed_grids/solution_transfer_05.cc b/tests/distributed_grids/solution_transfer_05.cc new file mode 100644 index 0000000000..e021b80574 --- /dev/null +++ b/tests/distributed_grids/solution_transfer_05.cc @@ -0,0 +1,95 @@ +// --------------------------------------------------------------------- +// +// 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 distributed solution_transfer with averaging (without averaging the +// program fails) + +#include + +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +int +main(int argc, char **argv) +{ + constexpr unsigned int dim = 2; + constexpr unsigned int degree = 1; + using Number = double; + + Utilities::MPI::MPI_InitFinalize mpi_init(argc, argv, 1); + MPILogInitAll all; + + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + GridGenerator::subdivided_hyper_rectangle(tria, + {2, 1}, + {0., 0.0}, + {2.0, 1.0}); + if (tria.create_cell_iterator(CellId("0_0:"))->is_locally_owned()) + tria.create_cell_iterator(CellId("0_0:"))->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(FE_Q(degree)); + + LinearAlgebra::distributed::Vector v; + v.reinit(dof_handler.locally_owned_dofs(), + DoFTools::extract_locally_relevant_dofs(dof_handler), + MPI_COMM_WORLD); + + std::map> values; + values[CellId("1_0:")] = Vector{0.5, 0.5, 0.5, 0.5}; + values[CellId("0_1:0")] = Vector{0.5, 0.5, 0.5, 0.5}; + values[CellId("0_1:1")] = Vector{0.5, 0.5, 0.5, 0.1}; + values[CellId("0_1:2")] = Vector{0.5, 0.5, 0.5, 0.5}; + values[CellId("0_1:3")] = Vector{0.5, 0.1, 0.5, 0.5}; + + for (const auto cell : dof_handler.active_cell_iterators()) + if (cell->is_artificial() == false) + cell->set_dof_values(values[cell->id()], v); + + v.print(deallog.get_file_stream()); + + parallel::distributed:: + SolutionTransfer> + solution_trans(dof_handler, true /*enabling averaging*/); + + v.update_ghost_values(); + solution_trans.prepare_for_coarsening_and_refinement(v); + + if (tria.create_cell_iterator(CellId("1_0:"))->is_locally_owned()) + tria.create_cell_iterator(CellId("1_0:"))->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + + dof_handler.distribute_dofs(FE_Q(degree)); + v.reinit(dof_handler.locally_owned_dofs(), + DoFTools::extract_locally_relevant_dofs(dof_handler), + MPI_COMM_WORLD); + solution_trans.interpolate(v); + + v.print(deallog.get_file_stream()); +} diff --git a/tests/distributed_grids/solution_transfer_05.mpirun=1.with_p4est=true.output b/tests/distributed_grids/solution_transfer_05.mpirun=1.with_p4est=true.output new file mode 100644 index 0000000000..d332aacc79 --- /dev/null +++ b/tests/distributed_grids/solution_transfer_05.mpirun=1.with_p4est=true.output @@ -0,0 +1,9 @@ + +Process #0 +Local range: [0, 11), global size: 11 +Vector data: +5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 1.000e-01 5.000e-01 5.000e-01 +Process #0 +Local range: [0, 15), global size: 15 +Vector data: +5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 3.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 diff --git a/tests/distributed_grids/solution_transfer_05.mpirun=2.with_p4est=true.output b/tests/distributed_grids/solution_transfer_05.mpirun=2.with_p4est=true.output new file mode 100644 index 0000000000..728762d27f --- /dev/null +++ b/tests/distributed_grids/solution_transfer_05.mpirun=2.with_p4est=true.output @@ -0,0 +1,19 @@ + +Process #0 +Local range: [0, 9), global size: 11 +Vector data: +5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 1.000e-01 5.000e-01 5.000e-01 5.000e-01 +Process #0 +Local range: [0, 9), global size: 15 +Vector data: +5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 3.000e-01 5.000e-01 5.000e-01 5.000e-01 + +Process #1 +Local range: [9, 11), global size: 11 +Vector data: +5.000e-01 5.000e-01 +Process #1 +Local range: [9, 15), global size: 15 +Vector data: +5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 5.000e-01 +