From 5b7f897e0bd017e2a16d7d8ad3bf08cbb9eb700d Mon Sep 17 00:00:00 2001 From: Timo Heister Date: Mon, 4 Aug 2014 15:23:41 -0400 Subject: [PATCH] add test --- tests/mpi/attach_data_01.cc | 173 +++++++++++++++++++++++ tests/mpi/attach_data_01.mpirun=1.output | 48 +++++++ tests/mpi/attach_data_01.mpirun=3.output | 62 ++++++++ 3 files changed, 283 insertions(+) create mode 100644 tests/mpi/attach_data_01.cc create mode 100644 tests/mpi/attach_data_01.mpirun=1.output create mode 100644 tests/mpi/attach_data_01.mpirun=3.output diff --git a/tests/mpi/attach_data_01.cc b/tests/mpi/attach_data_01.cc new file mode 100644 index 0000000000..995973c3a7 --- /dev/null +++ b/tests/mpi/attach_data_01.cc @@ -0,0 +1,173 @@ +// --------------------------------------------------------------------- +// $Id$ +// +// Copyright (C) 2009 - 2013 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 at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + + +// check register_data_attach and notify_ready_to_unpack + +#include "../tests.h" +#include "coarse_grid_common.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +#include + +std::string id_to_string(const CellId &id) +{ + std::ostringstream ss; + ss << id; + return ss.str(); +} + +template +void pack_function (const typename parallel::distributed::Triangulation::cell_iterator &cell, + const typename parallel::distributed::Triangulation::CellStatus status, + void *data) +{ + static int some_number = 0; + deallog << "packing cell " << cell->id() << " with data=" << some_number << " status="; + if (status==parallel::distributed::Triangulation::CELL_PERSIST) + deallog << "PERSIST"; + else if (status==parallel::distributed::Triangulation::CELL_REFINE) + deallog << "REFINE"; + else if (status==parallel::distributed::Triangulation::CELL_COARSEN) + deallog << "COARSEN"; + deallog << std::endl; + + if (status==parallel::distributed::Triangulation::CELL_COARSEN) + { + Assert(cell->has_children(), ExcInternalError()); + } + else + { + Assert(!cell->has_children(), ExcInternalError()); + } + + int * intdata = reinterpret_cast(data); + *intdata = some_number; + + ++some_number; +} + +template +void unpack_function (const typename parallel::distributed::Triangulation::cell_iterator &cell, + const typename parallel::distributed::Triangulation::CellStatus status, + const void *data) +{ + const int * intdata = reinterpret_cast(data); + + deallog << "unpacking cell " << cell->id() << " with data=" << (*intdata) << " status="; + if (status==parallel::distributed::Triangulation::CELL_PERSIST) + deallog << "PERSIST"; + else if (status==parallel::distributed::Triangulation::CELL_REFINE) + deallog << "REFINE"; + else if (status==parallel::distributed::Triangulation::CELL_COARSEN) + deallog << "COARSEN"; + deallog << std::endl; + + if (status==parallel::distributed::Triangulation::CELL_REFINE) + { + Assert(cell->has_children(), ExcInternalError()); + } + else + { + Assert(!cell->has_children(), ExcInternalError()); + } +} + + +template +void test() +{ + unsigned int myid = Utilities::MPI::this_mpi_process (MPI_COMM_WORLD); + unsigned int numprocs = Utilities::MPI::n_mpi_processes (MPI_COMM_WORLD); + + if (true) + { + parallel::distributed::Triangulation tr(MPI_COMM_WORLD); + + GridGenerator::subdivided_hyper_cube(tr, 2); + tr.refine_global(1); + deallog << "cells before: " << tr.n_global_active_cells() << std::endl; + + typename Triangulation::active_cell_iterator cell; + + for (cell = tr.begin_active(); + cell != tr.end(); + ++cell) + { + if (id_to_string(cell->id())=="0_1:0") + { + cell->set_refine_flag(); + } + else if (id_to_string(cell->parent()->id())=="3_0:") + cell->set_coarsen_flag(); + + if (cell->is_locally_owned()) + { + deallog << "myid=" << myid << " cellid=" << cell->id(); + if (cell->coarsen_flag_set()) + deallog << " coarsening" << std::endl; + else if (cell->refine_flag_set()) + deallog << " refining" << std::endl; + else + deallog << std::endl; + } + } + + unsigned int offset = tr.register_data_attach(sizeof(int), pack_function); + + deallog << "offset=" << offset << std::endl; + tr.execute_coarsening_and_refinement(); + + deallog << "cells after: " << tr.n_global_active_cells() << std::endl; + + /* + for (cell = tr.begin_active(); + cell != tr.end(); + ++cell) + { + if (cell->is_locally_owned()) + deallog << "myid=" << myid << " cellid=" << cell->id() << std::endl; + }*/ + + tr.notify_ready_to_unpack(offset, unpack_function); + + const unsigned int checksum = tr.get_checksum (); + deallog << "Checksum: " + << checksum + << std::endl; + } + + deallog << "OK" << std::endl; +} + + +int main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + + test<2>(); +} diff --git a/tests/mpi/attach_data_01.mpirun=1.output b/tests/mpi/attach_data_01.mpirun=1.output new file mode 100644 index 0000000000..c13bf0b12b --- /dev/null +++ b/tests/mpi/attach_data_01.mpirun=1.output @@ -0,0 +1,48 @@ + +DEAL:0::cells before: 16 +DEAL:0::myid=0 cellid=0_1:0 refining +DEAL:0::myid=0 cellid=0_1:1 +DEAL:0::myid=0 cellid=0_1:2 +DEAL:0::myid=0 cellid=0_1:3 +DEAL:0::myid=0 cellid=1_1:0 +DEAL:0::myid=0 cellid=1_1:1 +DEAL:0::myid=0 cellid=1_1:2 +DEAL:0::myid=0 cellid=1_1:3 +DEAL:0::myid=0 cellid=2_1:0 +DEAL:0::myid=0 cellid=2_1:1 +DEAL:0::myid=0 cellid=2_1:2 +DEAL:0::myid=0 cellid=2_1:3 +DEAL:0::myid=0 cellid=3_1:0 coarsening +DEAL:0::myid=0 cellid=3_1:1 coarsening +DEAL:0::myid=0 cellid=3_1:2 coarsening +DEAL:0::myid=0 cellid=3_1:3 coarsening +DEAL:0::offset=4 +DEAL:0::packing cell 0_1:0 with data=0 status=REFINE +DEAL:0::packing cell 0_1:1 with data=1 status=PERSIST +DEAL:0::packing cell 0_1:2 with data=2 status=PERSIST +DEAL:0::packing cell 0_1:3 with data=3 status=PERSIST +DEAL:0::packing cell 1_1:0 with data=4 status=PERSIST +DEAL:0::packing cell 1_1:1 with data=5 status=PERSIST +DEAL:0::packing cell 1_1:2 with data=6 status=PERSIST +DEAL:0::packing cell 1_1:3 with data=7 status=PERSIST +DEAL:0::packing cell 2_1:0 with data=8 status=PERSIST +DEAL:0::packing cell 2_1:1 with data=9 status=PERSIST +DEAL:0::packing cell 2_1:2 with data=10 status=PERSIST +DEAL:0::packing cell 2_1:3 with data=11 status=PERSIST +DEAL:0::packing cell 3_0: with data=12 status=COARSEN +DEAL:0::cells after: 16 +DEAL:0::unpacking cell 0_1:0 with data=0 status=REFINE +DEAL:0::unpacking cell 0_1:1 with data=1 status=PERSIST +DEAL:0::unpacking cell 0_1:2 with data=2 status=PERSIST +DEAL:0::unpacking cell 0_1:3 with data=3 status=PERSIST +DEAL:0::unpacking cell 1_1:0 with data=4 status=PERSIST +DEAL:0::unpacking cell 1_1:1 with data=5 status=PERSIST +DEAL:0::unpacking cell 1_1:2 with data=6 status=PERSIST +DEAL:0::unpacking cell 1_1:3 with data=7 status=PERSIST +DEAL:0::unpacking cell 2_1:0 with data=8 status=PERSIST +DEAL:0::unpacking cell 2_1:1 with data=9 status=PERSIST +DEAL:0::unpacking cell 2_1:2 with data=10 status=PERSIST +DEAL:0::unpacking cell 2_1:3 with data=11 status=PERSIST +DEAL:0::unpacking cell 3_0: with data=12 status=COARSEN +DEAL:0::Checksum: 2822439380 +DEAL:0::OK diff --git a/tests/mpi/attach_data_01.mpirun=3.output b/tests/mpi/attach_data_01.mpirun=3.output new file mode 100644 index 0000000000..ca1c6de254 --- /dev/null +++ b/tests/mpi/attach_data_01.mpirun=3.output @@ -0,0 +1,62 @@ + +DEAL:0::cells before: 16 +DEAL:0::myid=0 cellid=0_1:0 refining +DEAL:0::myid=0 cellid=0_1:1 +DEAL:0::myid=0 cellid=0_1:2 +DEAL:0::myid=0 cellid=0_1:3 +DEAL:0::offset=4 +DEAL:0::packing cell 0_1:0 with data=0 status=REFINE +DEAL:0::packing cell 0_1:1 with data=1 status=PERSIST +DEAL:0::packing cell 0_1:2 with data=2 status=PERSIST +DEAL:0::packing cell 0_1:3 with data=3 status=PERSIST +DEAL:0::cells after: 16 +DEAL:0::unpacking cell 0_1:0 with data=0 status=REFINE +DEAL:0::unpacking cell 0_1:1 with data=1 status=PERSIST +DEAL:0::Checksum: 2822439380 +DEAL:0::OK + +DEAL:1::cells before: 16 +DEAL:1::myid=1 cellid=1_1:0 +DEAL:1::myid=1 cellid=1_1:1 +DEAL:1::myid=1 cellid=1_1:2 +DEAL:1::myid=1 cellid=1_1:3 +DEAL:1::myid=1 cellid=2_1:0 +DEAL:1::myid=1 cellid=2_1:1 +DEAL:1::myid=1 cellid=2_1:2 +DEAL:1::myid=1 cellid=2_1:3 +DEAL:1::offset=4 +DEAL:1::packing cell 1_1:0 with data=0 status=PERSIST +DEAL:1::packing cell 1_1:1 with data=1 status=PERSIST +DEAL:1::packing cell 1_1:2 with data=2 status=PERSIST +DEAL:1::packing cell 1_1:3 with data=3 status=PERSIST +DEAL:1::packing cell 2_1:0 with data=4 status=PERSIST +DEAL:1::packing cell 2_1:1 with data=5 status=PERSIST +DEAL:1::packing cell 2_1:2 with data=6 status=PERSIST +DEAL:1::packing cell 2_1:3 with data=7 status=PERSIST +DEAL:1::cells after: 16 +DEAL:1::unpacking cell 0_1:2 with data=2 status=PERSIST +DEAL:1::unpacking cell 0_1:3 with data=3 status=PERSIST +DEAL:1::unpacking cell 1_1:0 with data=0 status=PERSIST +DEAL:1::unpacking cell 1_1:1 with data=1 status=PERSIST +DEAL:1::unpacking cell 1_1:2 with data=2 status=PERSIST +DEAL:1::unpacking cell 1_1:3 with data=3 status=PERSIST +DEAL:1::Checksum: 0 +DEAL:1::OK + + +DEAL:2::cells before: 16 +DEAL:2::myid=2 cellid=3_1:0 coarsening +DEAL:2::myid=2 cellid=3_1:1 coarsening +DEAL:2::myid=2 cellid=3_1:2 coarsening +DEAL:2::myid=2 cellid=3_1:3 coarsening +DEAL:2::offset=4 +DEAL:2::packing cell 3_0: with data=0 status=COARSEN +DEAL:2::cells after: 16 +DEAL:2::unpacking cell 2_1:0 with data=4 status=PERSIST +DEAL:2::unpacking cell 2_1:1 with data=5 status=PERSIST +DEAL:2::unpacking cell 2_1:2 with data=6 status=PERSIST +DEAL:2::unpacking cell 2_1:3 with data=7 status=PERSIST +DEAL:2::unpacking cell 3_0: with data=0 status=COARSEN +DEAL:2::Checksum: 0 +DEAL:2::OK + -- 2.39.5