From a4e40562ffc4782c27301ae6bf22d8358737fb82 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Thu, 21 Dec 2017 18:27:59 +0100 Subject: [PATCH] Simplify the logic used and extract all locally relevant DoFs --- doc/news/changes/minor/20171209DanielArndt | 4 +- include/deal.II/dofs/dof_tools.h | 10 +- source/dofs/dof_tools.cc | 73 ++--- tests/dofs/dof_tools_04.cc | 15 +- tests/dofs/dof_tools_04.output | 12 +- tests/dofs/dof_tools_04a.cc | 79 ++--- ..._tools_04a.with_p4est=true.mpirun=1.output | 284 +----------------- ..._tools_04a.with_p4est=true.mpirun=3.output | 284 +----------------- 8 files changed, 85 insertions(+), 676 deletions(-) diff --git a/doc/news/changes/minor/20171209DanielArndt b/doc/news/changes/minor/20171209DanielArndt index 17850e1935..7fb7fa7815 100644 --- a/doc/news/changes/minor/20171209DanielArndt +++ b/doc/news/changes/minor/20171209DanielArndt @@ -1,5 +1,5 @@ Improved: DoFTools::extract_hanging_node_dofs works for parallel::shared::Triangualtion and parallel::distributed::Triangulation -and reports DoFs belonging to locally owned cells. +and reports locally relevant DoFs.
-(Daniel Arndt, 2017/12/09) +(Daniel Arndt, 2017/12/21) diff --git a/include/deal.II/dofs/dof_tools.h b/include/deal.II/dofs/dof_tools.h index 520ed7c066..066f8ff2e8 100644 --- a/include/deal.II/dofs/dof_tools.h +++ b/include/deal.II/dofs/dof_tools.h @@ -1249,9 +1249,13 @@ namespace DoFTools * The size of @p selected_dofs shall equal dof_handler.n_dofs(). * Previous contents of this array or overwritten. * - * In case of a parallel::shared::Triangualtion or a - * parallel::distributed::Triangulation only dofs belonging to locally owned - * cells are reported. + * In case of a parallel::shared::Triangulation or a + * parallel::distributed::Triangulation only locally relevant dofs are + * considered. Note that the vector returned through the second argument still + * has size dof_handler.n_dofs(). Consequently, it can be very large + * for large parallel computations -- in fact, it may be too large to store on + * each processor. In that case, you may want to choose the variant of this + * function that returns an IndexSet object. */ template void diff --git a/source/dofs/dof_tools.cc b/source/dofs/dof_tools.cc index d498af038a..b175999346 100644 --- a/source/dofs/dof_tools.cc +++ b/source/dofs/dof_tools.cc @@ -890,7 +890,7 @@ namespace DoFTools // this function is similar to the make_sparsity_pattern function, // see there for more information for (auto cell : dof_handler.active_cell_iterators()) - if (cell->is_locally_owned()) + if (!cell->is_artificial()) { for (unsigned int face=0; face::faces_per_cell; ++face) if (cell->face(face)->has_children()) @@ -902,8 +902,12 @@ namespace DoFTools selected_dofs.add_index(line->child(0)->vertex_dof_index(1,dof)); for (unsigned int child=0; child<2; ++child) - for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) - selected_dofs.add_index(line->child(child)->dof_index(dof)); + { + if (cell->neighbor_child_on_subface (face, child)->is_artificial()) + continue; + for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) + selected_dofs.add_index(line->child(child)->dof_index(dof)); + } } } @@ -918,52 +922,35 @@ namespace DoFTools const unsigned int dim = 3; IndexSet selected_dofs (dof_handler.n_dofs()); + IndexSet unconstrained_dofs (dof_handler.n_dofs()); const FiniteElement &fe = dof_handler.get_fe(); - // this function is similar to the make_sparsity_pattern function, - // see there for more information for (auto cell : dof_handler.active_cell_iterators()) - if (cell->is_locally_owned()) + if (!cell->is_artificial()) for (unsigned int f=0; f::faces_per_cell; ++f) - if (cell->face(f)->has_children()) - { - const typename dealii::DoFHandler::face_iterator - face = cell->face(f); - - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - selected_dofs.add_index(face->child(0)->vertex_dof_index(2,dof)); - - // dof numbers on the centers of the lines bounding this - // face - for (unsigned int line=0; line<4; ++line) - for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) - selected_dofs.add_index(face->line(line)->child(0)->vertex_dof_index(1,dof)); - - // next the dofs on the lines interior to the face; the - // order of these lines is laid down in the FiniteElement - // class documentation - for (unsigned int dof=0; dofchild(0)->line(1)->dof_index(dof)); - for (unsigned int dof=0; dofchild(1)->line(2)->dof_index(dof)); - for (unsigned int dof=0; dofchild(2)->line(3)->dof_index(dof)); - for (unsigned int dof=0; dofchild(3)->line(0)->dof_index(dof)); - - // dofs on the bordering lines - for (unsigned int line=0; line<4; ++line) - for (unsigned int child=0; child<2; ++child) - for (unsigned int dof=0; dof!=fe.dofs_per_line; ++dof) - selected_dofs.add_index(face->line(line)->child(child)->dof_index(dof)); + { + const typename dealii::DoFHandler::face_iterator + face = cell->face(f); + if (cell->face(f)->has_children()) + { + for (unsigned int child=0; child<4; ++child) + if (!cell->neighbor_child_on_subface (f, child)->is_artificial()) + { + // simply take all DoFs that live on this subface + std::vector ldi(fe.dofs_per_face); + face->child(child)->get_dof_indices(ldi); + selected_dofs.add_indices(ldi.begin(), ldi.end()); + } - // finally, for the dofs interior to the four child faces - for (unsigned int child=0; child<4; ++child) - for (unsigned int dof=0; dof!=fe.dofs_per_quad; ++dof) - selected_dofs.add_index(face->child(child)->dof_index(dof)); - } - selected_dofs.compress(); + // and subtract (in the end) all the indices which a shared + // between this face and its subfaces + for (unsigned int vertex=0; vertex<4; ++vertex) + for (unsigned int dof=0; dof!=fe.dofs_per_vertex; ++dof) + unconstrained_dofs.add_index(face->vertex_dof_index(vertex,dof)); + } + } + selected_dofs.subtract_set(unconstrained_dofs); return selected_dofs; } } diff --git a/tests/dofs/dof_tools_04.cc b/tests/dofs/dof_tools_04.cc index 78400ff1df..ee89ead282 100644 --- a/tests/dofs/dof_tools_04.cc +++ b/tests/dofs/dof_tools_04.cc @@ -27,8 +27,21 @@ template void check_this (const DoFHandler &dof_handler) { - std::vector hanging_node_dofs (dof_handler.n_dofs()); + const types::global_dof_index n_dofs = dof_handler.n_dofs(); + + std::vector hanging_node_dofs (n_dofs); DoFTools::extract_hanging_node_dofs (dof_handler, hanging_node_dofs); + + ConstraintMatrix constraints; + DoFTools::make_hanging_node_constraints(dof_handler, constraints); + constraints.close(); + + for (types::global_dof_index dof=0; dof(3)1 in 2d: DEAL::0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 DEAL::Checking FE_Q<3>(1)3 in 3d: -DEAL::000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000111111000111111111000111111000000000000000000000000111000111111000111111 +DEAL::000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000111111111111111111111111111000000000000000000000000111111111111111111111 DEAL::Checking FE_DGQ<3>(2)2 in 3d: DEAL::000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 DEAL::Checking FE_DGP<3>(3)1 in 3d: @@ -118,7 +118,7 @@ DEAL::00000000000000000000000000000000000000000000000000000000000000000000000000 DEAL::Checking FE_DGP<2>(3)1FE_DGP<2>(3)1FE_Q<2>(2)3 in 2d: DEAL::00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000111000000111000000000000000000000000001110000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001110001110000000000000000000000000011100000000000000000000000 DEAL::Checking FE_DGQ<3>(1)3FE_DGP<3>(3)1FE_Q<3>(1)3 in 3d: -DEAL::00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001111110000000000000000000000000000000000000000000000011100000000000000000000000000000000000000000000111111000000000000000000000000000000000000000000000001110000000000000000000000000000000000000000000011100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000111000000000000000000000000000000000000000000000001110000000000000000000000000000000000000000000011100000000000000000000000000000000000000000000000111000000000000000000000000000000000000000000001110000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 +DEAL::00000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000001111111110000000000000000000000000000000000000000000011100000000000000000000000000000000000000000000111111111000000000000000000000000000000000000000000001110000000000000000000000000000000000000000000011100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000111111000000000000000000000000000000000000000000001110000000000000000000000000000000000000000000011111100000000000000000000000000000000000000000000111000000000000000000000000000000000000000000001110000000000000000000000000000000000000000000000000000000000000000000000000000000000000000 DEAL::Checking (FESystem<2>(FE_Q<2>(1),3))3FE_DGQ<2>(0)1FE_Q<2>(1)3 in 2d: DEAL::0000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000111111111111000000000000000000000000000000000000000000000000000011111111111100 DEAL::Checking FE_DGQ<2>(3)1FESystem<2>(FE_DGQ<2>(0),3)1FESystem<2>(FE_Q<2>(2),1, FE_DGQ<2>(0),1)2 in 2d: diff --git a/tests/dofs/dof_tools_04a.cc b/tests/dofs/dof_tools_04a.cc index 481daf79a7..aa0421b7c1 100644 --- a/tests/dofs/dof_tools_04a.cc +++ b/tests/dofs/dof_tools_04a.cc @@ -14,11 +14,11 @@ // --------------------------------------------------------------------- -// The same as dof_tools_04 but for a parallel::distributed::parallel::distributed::Triangulation -// instead of a parallel::distributed::Triangulation: +// Similar to dof_tools_04 but for a parallel::distributed::Triangulation +// instead of a (serial) Triangulation: // check // DoFTools::extract_hanging_node_constraints -// uses a slightly different refinement and less different FiniteElements +// using a slightly different refinement and less different FiniteElements #include "../tests.h" #include @@ -44,72 +44,27 @@ template void check_this (const DoFHandler &dof_handler) { - types::global_dof_index n_dofs = dof_handler.n_dofs(); + const types::global_dof_index n_dofs = dof_handler.n_dofs(); + + IndexSet locally_relevant_dofs; + DoFTools::extract_locally_relevant_dofs (dof_handler, locally_relevant_dofs); std::vector is_hanging_node_constrained (n_dofs); DoFTools::extract_hanging_node_dofs (dof_handler, is_hanging_node_constrained); - MappingQGeneric mapping(1); - std::map< types::global_dof_index, Point< dim > > support_point_of_dof; - DoFTools::map_dofs_to_support_points (mapping, dof_handler, support_point_of_dof); - - std::vector> constrained_points; - - for (const auto &pair : support_point_of_dof) - if (is_hanging_node_constrained[pair.first]) - constrained_points.push_back(pair.second); - - const int root = 0; - const int mylen = constrained_points.size()*dim; - const unsigned int n_processes = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); - const unsigned int my_id = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); - std::vector recvcounts (n_processes); - - MPI_Gather(&mylen, 1, MPI_INT, recvcounts.data(), 1, MPI_INT, 0, MPI_COMM_WORLD); - - int totlen = 0; - std::vector displs; - std::vector> all_points; + ConstraintMatrix constraints (locally_relevant_dofs); + DoFTools::make_hanging_node_constraints(dof_handler, constraints); + constraints.close(); - if (my_id==0) - { - displs.resize(n_processes); - - displs[0] = 0; - totlen += recvcounts[0]; - - for (unsigned int i=1; i &a, const Point &b) - { - for (unsigned int i=0; i b(i)) - return false; - } - return false; - }); - all_points.erase(std::unique(all_points.begin(), all_points.end()), all_points.end()); - - deallog << all_points.size() << std::endl; - for (const auto &point: all_points) - deallog << point << std::endl; - deallog << std::endl; + deallog << "OK" << std::endl; } @@ -137,6 +92,7 @@ check (const FiniteElement &fe, cell->set_refine_flag(); tria.execute_coarsening_and_refinement (); } + DoFHandler dof_handler (tria); dof_handler.distribute_dofs (fe); @@ -157,6 +113,7 @@ int main (int argc, char **argv) Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); mpi_initlog(); CHECK_ALL(Q,1); + CHECK_ALL(Q,2); return 0; diff --git a/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=1.output b/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=1.output index 6feb20d8b2..c3d49df7e0 100644 --- a/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=1.output +++ b/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=1.output @@ -1,285 +1,9 @@ DEAL::Checking Q1 in 2d: -DEAL::4 -DEAL::0.125000 0.500000 -DEAL::0.375000 0.500000 -DEAL::0.500000 0.125000 -DEAL::0.500000 0.375000 -DEAL:: +DEAL::OK DEAL::Checking Q1 in 3d: -DEAL::40 -DEAL::0.00000 0.500000 0.125000 -DEAL::0.00000 0.500000 0.375000 -DEAL::0.00000 0.500000 0.625000 -DEAL::0.00000 0.500000 0.875000 -DEAL::0.125000 0.500000 0.00000 -DEAL::0.125000 0.500000 0.250000 -DEAL::0.125000 0.500000 0.500000 -DEAL::0.125000 0.500000 0.750000 -DEAL::0.125000 0.500000 1.00000 -DEAL::0.250000 0.500000 0.125000 -DEAL::0.250000 0.500000 0.375000 -DEAL::0.250000 0.500000 0.625000 -DEAL::0.250000 0.500000 0.875000 -DEAL::0.375000 0.500000 0.00000 -DEAL::0.375000 0.500000 0.250000 -DEAL::0.375000 0.500000 0.500000 -DEAL::0.375000 0.500000 0.750000 -DEAL::0.375000 0.500000 1.00000 -DEAL::0.500000 0.00000 0.125000 -DEAL::0.500000 0.00000 0.375000 -DEAL::0.500000 0.00000 0.625000 -DEAL::0.500000 0.00000 0.875000 -DEAL::0.500000 0.125000 0.00000 -DEAL::0.500000 0.125000 0.250000 -DEAL::0.500000 0.125000 0.750000 -DEAL::0.500000 0.125000 1.00000 -DEAL::0.500000 0.125000 0.500000 -DEAL::0.500000 0.250000 0.125000 -DEAL::0.500000 0.250000 0.375000 -DEAL::0.500000 0.250000 0.625000 -DEAL::0.500000 0.250000 0.875000 -DEAL::0.500000 0.375000 0.00000 -DEAL::0.500000 0.375000 0.250000 -DEAL::0.500000 0.375000 0.750000 -DEAL::0.500000 0.375000 1.00000 -DEAL::0.500000 0.375000 0.500000 -DEAL::0.500000 0.500000 0.125000 -DEAL::0.500000 0.500000 0.375000 -DEAL::0.500000 0.500000 0.625000 -DEAL::0.500000 0.500000 0.875000 -DEAL:: +DEAL::OK DEAL::Checking Q2 in 2d: -DEAL::12 -DEAL::0.0625000 0.500000 -DEAL::0.125000 0.500000 -DEAL::0.187500 0.500000 -DEAL::0.312500 0.500000 -DEAL::0.375000 0.500000 -DEAL::0.437500 0.500000 -DEAL::0.500000 0.0625000 -DEAL::0.500000 0.125000 -DEAL::0.500000 0.187500 -DEAL::0.500000 0.312500 -DEAL::0.500000 0.375000 -DEAL::0.500000 0.437500 -DEAL:: +DEAL::OK DEAL::Checking Q2 in 3d: -DEAL::216 -DEAL::0.00000 0.500000 0.0625000 -DEAL::0.00000 0.500000 0.125000 -DEAL::0.00000 0.500000 0.187500 -DEAL::0.00000 0.500000 0.312500 -DEAL::0.00000 0.500000 0.375000 -DEAL::0.00000 0.500000 0.437500 -DEAL::0.00000 0.500000 0.562500 -DEAL::0.00000 0.500000 0.625000 -DEAL::0.00000 0.500000 0.687500 -DEAL::0.00000 0.500000 0.812500 -DEAL::0.00000 0.500000 0.875000 -DEAL::0.00000 0.500000 0.937500 -DEAL::0.0625000 0.500000 0.00000 -DEAL::0.0625000 0.500000 0.0625000 -DEAL::0.0625000 0.500000 0.125000 -DEAL::0.0625000 0.500000 0.187500 -DEAL::0.0625000 0.500000 0.250000 -DEAL::0.0625000 0.500000 0.312500 -DEAL::0.0625000 0.500000 0.375000 -DEAL::0.0625000 0.500000 0.437500 -DEAL::0.0625000 0.500000 0.500000 -DEAL::0.0625000 0.500000 0.562500 -DEAL::0.0625000 0.500000 0.625000 -DEAL::0.0625000 0.500000 0.687500 -DEAL::0.0625000 0.500000 0.750000 -DEAL::0.0625000 0.500000 0.812500 -DEAL::0.0625000 0.500000 0.875000 -DEAL::0.0625000 0.500000 0.937500 -DEAL::0.0625000 0.500000 1.00000 -DEAL::0.125000 0.500000 0.00000 -DEAL::0.125000 0.500000 0.250000 -DEAL::0.125000 0.500000 0.500000 -DEAL::0.125000 0.500000 0.750000 -DEAL::0.125000 0.500000 1.00000 -DEAL::0.187500 0.500000 0.00000 -DEAL::0.187500 0.500000 0.0625000 -DEAL::0.187500 0.500000 0.125000 -DEAL::0.187500 0.500000 0.187500 -DEAL::0.187500 0.500000 0.250000 -DEAL::0.187500 0.500000 0.312500 -DEAL::0.187500 0.500000 0.375000 -DEAL::0.187500 0.500000 0.437500 -DEAL::0.187500 0.500000 0.500000 -DEAL::0.187500 0.500000 0.562500 -DEAL::0.187500 0.500000 0.625000 -DEAL::0.187500 0.500000 0.687500 -DEAL::0.187500 0.500000 0.750000 -DEAL::0.187500 0.500000 0.812500 -DEAL::0.187500 0.500000 0.875000 -DEAL::0.187500 0.500000 0.937500 -DEAL::0.187500 0.500000 1.00000 -DEAL::0.250000 0.500000 0.0625000 -DEAL::0.250000 0.500000 0.125000 -DEAL::0.250000 0.500000 0.187500 -DEAL::0.250000 0.500000 0.312500 -DEAL::0.250000 0.500000 0.375000 -DEAL::0.250000 0.500000 0.437500 -DEAL::0.250000 0.500000 0.562500 -DEAL::0.250000 0.500000 0.625000 -DEAL::0.250000 0.500000 0.687500 -DEAL::0.250000 0.500000 0.812500 -DEAL::0.250000 0.500000 0.875000 -DEAL::0.250000 0.500000 0.937500 -DEAL::0.312500 0.500000 0.00000 -DEAL::0.312500 0.500000 0.0625000 -DEAL::0.312500 0.500000 0.125000 -DEAL::0.312500 0.500000 0.187500 -DEAL::0.312500 0.500000 0.250000 -DEAL::0.312500 0.500000 0.312500 -DEAL::0.312500 0.500000 0.375000 -DEAL::0.312500 0.500000 0.437500 -DEAL::0.312500 0.500000 0.500000 -DEAL::0.312500 0.500000 0.562500 -DEAL::0.312500 0.500000 0.625000 -DEAL::0.312500 0.500000 0.687500 -DEAL::0.312500 0.500000 0.750000 -DEAL::0.312500 0.500000 0.812500 -DEAL::0.312500 0.500000 0.875000 -DEAL::0.312500 0.500000 0.937500 -DEAL::0.312500 0.500000 1.00000 -DEAL::0.375000 0.500000 0.00000 -DEAL::0.375000 0.500000 0.250000 -DEAL::0.375000 0.500000 0.500000 -DEAL::0.375000 0.500000 0.750000 -DEAL::0.375000 0.500000 1.00000 -DEAL::0.437500 0.500000 0.00000 -DEAL::0.437500 0.500000 0.0625000 -DEAL::0.437500 0.500000 0.125000 -DEAL::0.437500 0.500000 0.187500 -DEAL::0.437500 0.500000 0.250000 -DEAL::0.437500 0.500000 0.312500 -DEAL::0.437500 0.500000 0.375000 -DEAL::0.437500 0.500000 0.437500 -DEAL::0.437500 0.500000 0.500000 -DEAL::0.437500 0.500000 0.562500 -DEAL::0.437500 0.500000 0.625000 -DEAL::0.437500 0.500000 0.687500 -DEAL::0.437500 0.500000 0.750000 -DEAL::0.437500 0.500000 0.812500 -DEAL::0.437500 0.500000 0.875000 -DEAL::0.437500 0.500000 0.937500 -DEAL::0.437500 0.500000 1.00000 -DEAL::0.500000 0.00000 0.0625000 -DEAL::0.500000 0.00000 0.125000 -DEAL::0.500000 0.00000 0.187500 -DEAL::0.500000 0.00000 0.312500 -DEAL::0.500000 0.00000 0.375000 -DEAL::0.500000 0.00000 0.437500 -DEAL::0.500000 0.00000 0.562500 -DEAL::0.500000 0.00000 0.625000 -DEAL::0.500000 0.00000 0.687500 -DEAL::0.500000 0.00000 0.812500 -DEAL::0.500000 0.00000 0.875000 -DEAL::0.500000 0.00000 0.937500 -DEAL::0.500000 0.0625000 0.00000 -DEAL::0.500000 0.0625000 0.0625000 -DEAL::0.500000 0.0625000 0.187500 -DEAL::0.500000 0.0625000 0.250000 -DEAL::0.500000 0.0625000 0.312500 -DEAL::0.500000 0.0625000 0.437500 -DEAL::0.500000 0.0625000 0.562500 -DEAL::0.500000 0.0625000 0.687500 -DEAL::0.500000 0.0625000 0.750000 -DEAL::0.500000 0.0625000 0.812500 -DEAL::0.500000 0.0625000 0.937500 -DEAL::0.500000 0.0625000 1.00000 -DEAL::0.500000 0.0625000 0.500000 -DEAL::0.500000 0.125000 0.00000 -DEAL::0.500000 0.125000 0.0625000 -DEAL::0.500000 0.125000 0.187500 -DEAL::0.500000 0.125000 0.250000 -DEAL::0.500000 0.125000 0.312500 -DEAL::0.500000 0.125000 0.437500 -DEAL::0.500000 0.125000 0.562500 -DEAL::0.500000 0.125000 0.687500 -DEAL::0.500000 0.125000 0.750000 -DEAL::0.500000 0.125000 0.812500 -DEAL::0.500000 0.125000 0.937500 -DEAL::0.500000 0.125000 1.00000 -DEAL::0.500000 0.125000 0.500000 -DEAL::0.500000 0.187500 0.00000 -DEAL::0.500000 0.187500 0.0625000 -DEAL::0.500000 0.187500 0.187500 -DEAL::0.500000 0.187500 0.250000 -DEAL::0.500000 0.187500 0.312500 -DEAL::0.500000 0.187500 0.437500 -DEAL::0.500000 0.187500 0.687500 -DEAL::0.500000 0.187500 0.750000 -DEAL::0.500000 0.187500 0.812500 -DEAL::0.500000 0.187500 0.937500 -DEAL::0.500000 0.187500 1.00000 -DEAL::0.500000 0.187500 0.562500 -DEAL::0.500000 0.187500 0.500000 -DEAL::0.500000 0.250000 0.0625000 -DEAL::0.500000 0.250000 0.125000 -DEAL::0.500000 0.250000 0.187500 -DEAL::0.500000 0.250000 0.312500 -DEAL::0.500000 0.250000 0.375000 -DEAL::0.500000 0.250000 0.437500 -DEAL::0.500000 0.250000 0.562500 -DEAL::0.500000 0.250000 0.625000 -DEAL::0.500000 0.250000 0.687500 -DEAL::0.500000 0.250000 0.812500 -DEAL::0.500000 0.250000 0.875000 -DEAL::0.500000 0.250000 0.937500 -DEAL::0.500000 0.312500 0.00000 -DEAL::0.500000 0.312500 0.0625000 -DEAL::0.500000 0.312500 0.187500 -DEAL::0.500000 0.312500 0.250000 -DEAL::0.500000 0.312500 0.312500 -DEAL::0.500000 0.312500 0.437500 -DEAL::0.500000 0.312500 0.687500 -DEAL::0.500000 0.312500 0.750000 -DEAL::0.500000 0.312500 0.812500 -DEAL::0.500000 0.312500 0.937500 -DEAL::0.500000 0.312500 1.00000 -DEAL::0.500000 0.312500 0.562500 -DEAL::0.500000 0.312500 0.500000 -DEAL::0.500000 0.375000 0.00000 -DEAL::0.500000 0.375000 0.0625000 -DEAL::0.500000 0.375000 0.187500 -DEAL::0.500000 0.375000 0.250000 -DEAL::0.500000 0.375000 0.312500 -DEAL::0.500000 0.375000 0.687500 -DEAL::0.500000 0.375000 0.750000 -DEAL::0.500000 0.375000 0.812500 -DEAL::0.500000 0.375000 0.937500 -DEAL::0.500000 0.375000 1.00000 -DEAL::0.500000 0.375000 0.437500 -DEAL::0.500000 0.375000 0.562500 -DEAL::0.500000 0.375000 0.500000 -DEAL::0.500000 0.437500 0.00000 -DEAL::0.500000 0.437500 0.0625000 -DEAL::0.500000 0.437500 0.187500 -DEAL::0.500000 0.437500 0.250000 -DEAL::0.500000 0.437500 0.312500 -DEAL::0.500000 0.437500 0.437500 -DEAL::0.500000 0.437500 0.687500 -DEAL::0.500000 0.437500 0.750000 -DEAL::0.500000 0.437500 0.812500 -DEAL::0.500000 0.437500 0.937500 -DEAL::0.500000 0.437500 1.00000 -DEAL::0.500000 0.437500 0.562500 -DEAL::0.500000 0.437500 0.500000 -DEAL::0.500000 0.500000 0.0625000 -DEAL::0.500000 0.500000 0.125000 -DEAL::0.500000 0.500000 0.187500 -DEAL::0.500000 0.500000 0.312500 -DEAL::0.500000 0.500000 0.375000 -DEAL::0.500000 0.500000 0.437500 -DEAL::0.500000 0.500000 0.562500 -DEAL::0.500000 0.500000 0.625000 -DEAL::0.500000 0.500000 0.687500 -DEAL::0.500000 0.500000 0.812500 -DEAL::0.500000 0.500000 0.875000 -DEAL::0.500000 0.500000 0.937500 -DEAL:: +DEAL::OK diff --git a/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=3.output b/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=3.output index 6feb20d8b2..c3d49df7e0 100644 --- a/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=3.output +++ b/tests/dofs/dof_tools_04a.with_p4est=true.mpirun=3.output @@ -1,285 +1,9 @@ DEAL::Checking Q1 in 2d: -DEAL::4 -DEAL::0.125000 0.500000 -DEAL::0.375000 0.500000 -DEAL::0.500000 0.125000 -DEAL::0.500000 0.375000 -DEAL:: +DEAL::OK DEAL::Checking Q1 in 3d: -DEAL::40 -DEAL::0.00000 0.500000 0.125000 -DEAL::0.00000 0.500000 0.375000 -DEAL::0.00000 0.500000 0.625000 -DEAL::0.00000 0.500000 0.875000 -DEAL::0.125000 0.500000 0.00000 -DEAL::0.125000 0.500000 0.250000 -DEAL::0.125000 0.500000 0.500000 -DEAL::0.125000 0.500000 0.750000 -DEAL::0.125000 0.500000 1.00000 -DEAL::0.250000 0.500000 0.125000 -DEAL::0.250000 0.500000 0.375000 -DEAL::0.250000 0.500000 0.625000 -DEAL::0.250000 0.500000 0.875000 -DEAL::0.375000 0.500000 0.00000 -DEAL::0.375000 0.500000 0.250000 -DEAL::0.375000 0.500000 0.500000 -DEAL::0.375000 0.500000 0.750000 -DEAL::0.375000 0.500000 1.00000 -DEAL::0.500000 0.00000 0.125000 -DEAL::0.500000 0.00000 0.375000 -DEAL::0.500000 0.00000 0.625000 -DEAL::0.500000 0.00000 0.875000 -DEAL::0.500000 0.125000 0.00000 -DEAL::0.500000 0.125000 0.250000 -DEAL::0.500000 0.125000 0.750000 -DEAL::0.500000 0.125000 1.00000 -DEAL::0.500000 0.125000 0.500000 -DEAL::0.500000 0.250000 0.125000 -DEAL::0.500000 0.250000 0.375000 -DEAL::0.500000 0.250000 0.625000 -DEAL::0.500000 0.250000 0.875000 -DEAL::0.500000 0.375000 0.00000 -DEAL::0.500000 0.375000 0.250000 -DEAL::0.500000 0.375000 0.750000 -DEAL::0.500000 0.375000 1.00000 -DEAL::0.500000 0.375000 0.500000 -DEAL::0.500000 0.500000 0.125000 -DEAL::0.500000 0.500000 0.375000 -DEAL::0.500000 0.500000 0.625000 -DEAL::0.500000 0.500000 0.875000 -DEAL:: +DEAL::OK DEAL::Checking Q2 in 2d: -DEAL::12 -DEAL::0.0625000 0.500000 -DEAL::0.125000 0.500000 -DEAL::0.187500 0.500000 -DEAL::0.312500 0.500000 -DEAL::0.375000 0.500000 -DEAL::0.437500 0.500000 -DEAL::0.500000 0.0625000 -DEAL::0.500000 0.125000 -DEAL::0.500000 0.187500 -DEAL::0.500000 0.312500 -DEAL::0.500000 0.375000 -DEAL::0.500000 0.437500 -DEAL:: +DEAL::OK DEAL::Checking Q2 in 3d: -DEAL::216 -DEAL::0.00000 0.500000 0.0625000 -DEAL::0.00000 0.500000 0.125000 -DEAL::0.00000 0.500000 0.187500 -DEAL::0.00000 0.500000 0.312500 -DEAL::0.00000 0.500000 0.375000 -DEAL::0.00000 0.500000 0.437500 -DEAL::0.00000 0.500000 0.562500 -DEAL::0.00000 0.500000 0.625000 -DEAL::0.00000 0.500000 0.687500 -DEAL::0.00000 0.500000 0.812500 -DEAL::0.00000 0.500000 0.875000 -DEAL::0.00000 0.500000 0.937500 -DEAL::0.0625000 0.500000 0.00000 -DEAL::0.0625000 0.500000 0.0625000 -DEAL::0.0625000 0.500000 0.125000 -DEAL::0.0625000 0.500000 0.187500 -DEAL::0.0625000 0.500000 0.250000 -DEAL::0.0625000 0.500000 0.312500 -DEAL::0.0625000 0.500000 0.375000 -DEAL::0.0625000 0.500000 0.437500 -DEAL::0.0625000 0.500000 0.500000 -DEAL::0.0625000 0.500000 0.562500 -DEAL::0.0625000 0.500000 0.625000 -DEAL::0.0625000 0.500000 0.687500 -DEAL::0.0625000 0.500000 0.750000 -DEAL::0.0625000 0.500000 0.812500 -DEAL::0.0625000 0.500000 0.875000 -DEAL::0.0625000 0.500000 0.937500 -DEAL::0.0625000 0.500000 1.00000 -DEAL::0.125000 0.500000 0.00000 -DEAL::0.125000 0.500000 0.250000 -DEAL::0.125000 0.500000 0.500000 -DEAL::0.125000 0.500000 0.750000 -DEAL::0.125000 0.500000 1.00000 -DEAL::0.187500 0.500000 0.00000 -DEAL::0.187500 0.500000 0.0625000 -DEAL::0.187500 0.500000 0.125000 -DEAL::0.187500 0.500000 0.187500 -DEAL::0.187500 0.500000 0.250000 -DEAL::0.187500 0.500000 0.312500 -DEAL::0.187500 0.500000 0.375000 -DEAL::0.187500 0.500000 0.437500 -DEAL::0.187500 0.500000 0.500000 -DEAL::0.187500 0.500000 0.562500 -DEAL::0.187500 0.500000 0.625000 -DEAL::0.187500 0.500000 0.687500 -DEAL::0.187500 0.500000 0.750000 -DEAL::0.187500 0.500000 0.812500 -DEAL::0.187500 0.500000 0.875000 -DEAL::0.187500 0.500000 0.937500 -DEAL::0.187500 0.500000 1.00000 -DEAL::0.250000 0.500000 0.0625000 -DEAL::0.250000 0.500000 0.125000 -DEAL::0.250000 0.500000 0.187500 -DEAL::0.250000 0.500000 0.312500 -DEAL::0.250000 0.500000 0.375000 -DEAL::0.250000 0.500000 0.437500 -DEAL::0.250000 0.500000 0.562500 -DEAL::0.250000 0.500000 0.625000 -DEAL::0.250000 0.500000 0.687500 -DEAL::0.250000 0.500000 0.812500 -DEAL::0.250000 0.500000 0.875000 -DEAL::0.250000 0.500000 0.937500 -DEAL::0.312500 0.500000 0.00000 -DEAL::0.312500 0.500000 0.0625000 -DEAL::0.312500 0.500000 0.125000 -DEAL::0.312500 0.500000 0.187500 -DEAL::0.312500 0.500000 0.250000 -DEAL::0.312500 0.500000 0.312500 -DEAL::0.312500 0.500000 0.375000 -DEAL::0.312500 0.500000 0.437500 -DEAL::0.312500 0.500000 0.500000 -DEAL::0.312500 0.500000 0.562500 -DEAL::0.312500 0.500000 0.625000 -DEAL::0.312500 0.500000 0.687500 -DEAL::0.312500 0.500000 0.750000 -DEAL::0.312500 0.500000 0.812500 -DEAL::0.312500 0.500000 0.875000 -DEAL::0.312500 0.500000 0.937500 -DEAL::0.312500 0.500000 1.00000 -DEAL::0.375000 0.500000 0.00000 -DEAL::0.375000 0.500000 0.250000 -DEAL::0.375000 0.500000 0.500000 -DEAL::0.375000 0.500000 0.750000 -DEAL::0.375000 0.500000 1.00000 -DEAL::0.437500 0.500000 0.00000 -DEAL::0.437500 0.500000 0.0625000 -DEAL::0.437500 0.500000 0.125000 -DEAL::0.437500 0.500000 0.187500 -DEAL::0.437500 0.500000 0.250000 -DEAL::0.437500 0.500000 0.312500 -DEAL::0.437500 0.500000 0.375000 -DEAL::0.437500 0.500000 0.437500 -DEAL::0.437500 0.500000 0.500000 -DEAL::0.437500 0.500000 0.562500 -DEAL::0.437500 0.500000 0.625000 -DEAL::0.437500 0.500000 0.687500 -DEAL::0.437500 0.500000 0.750000 -DEAL::0.437500 0.500000 0.812500 -DEAL::0.437500 0.500000 0.875000 -DEAL::0.437500 0.500000 0.937500 -DEAL::0.437500 0.500000 1.00000 -DEAL::0.500000 0.00000 0.0625000 -DEAL::0.500000 0.00000 0.125000 -DEAL::0.500000 0.00000 0.187500 -DEAL::0.500000 0.00000 0.312500 -DEAL::0.500000 0.00000 0.375000 -DEAL::0.500000 0.00000 0.437500 -DEAL::0.500000 0.00000 0.562500 -DEAL::0.500000 0.00000 0.625000 -DEAL::0.500000 0.00000 0.687500 -DEAL::0.500000 0.00000 0.812500 -DEAL::0.500000 0.00000 0.875000 -DEAL::0.500000 0.00000 0.937500 -DEAL::0.500000 0.0625000 0.00000 -DEAL::0.500000 0.0625000 0.0625000 -DEAL::0.500000 0.0625000 0.187500 -DEAL::0.500000 0.0625000 0.250000 -DEAL::0.500000 0.0625000 0.312500 -DEAL::0.500000 0.0625000 0.437500 -DEAL::0.500000 0.0625000 0.562500 -DEAL::0.500000 0.0625000 0.687500 -DEAL::0.500000 0.0625000 0.750000 -DEAL::0.500000 0.0625000 0.812500 -DEAL::0.500000 0.0625000 0.937500 -DEAL::0.500000 0.0625000 1.00000 -DEAL::0.500000 0.0625000 0.500000 -DEAL::0.500000 0.125000 0.00000 -DEAL::0.500000 0.125000 0.0625000 -DEAL::0.500000 0.125000 0.187500 -DEAL::0.500000 0.125000 0.250000 -DEAL::0.500000 0.125000 0.312500 -DEAL::0.500000 0.125000 0.437500 -DEAL::0.500000 0.125000 0.562500 -DEAL::0.500000 0.125000 0.687500 -DEAL::0.500000 0.125000 0.750000 -DEAL::0.500000 0.125000 0.812500 -DEAL::0.500000 0.125000 0.937500 -DEAL::0.500000 0.125000 1.00000 -DEAL::0.500000 0.125000 0.500000 -DEAL::0.500000 0.187500 0.00000 -DEAL::0.500000 0.187500 0.0625000 -DEAL::0.500000 0.187500 0.187500 -DEAL::0.500000 0.187500 0.250000 -DEAL::0.500000 0.187500 0.312500 -DEAL::0.500000 0.187500 0.437500 -DEAL::0.500000 0.187500 0.687500 -DEAL::0.500000 0.187500 0.750000 -DEAL::0.500000 0.187500 0.812500 -DEAL::0.500000 0.187500 0.937500 -DEAL::0.500000 0.187500 1.00000 -DEAL::0.500000 0.187500 0.562500 -DEAL::0.500000 0.187500 0.500000 -DEAL::0.500000 0.250000 0.0625000 -DEAL::0.500000 0.250000 0.125000 -DEAL::0.500000 0.250000 0.187500 -DEAL::0.500000 0.250000 0.312500 -DEAL::0.500000 0.250000 0.375000 -DEAL::0.500000 0.250000 0.437500 -DEAL::0.500000 0.250000 0.562500 -DEAL::0.500000 0.250000 0.625000 -DEAL::0.500000 0.250000 0.687500 -DEAL::0.500000 0.250000 0.812500 -DEAL::0.500000 0.250000 0.875000 -DEAL::0.500000 0.250000 0.937500 -DEAL::0.500000 0.312500 0.00000 -DEAL::0.500000 0.312500 0.0625000 -DEAL::0.500000 0.312500 0.187500 -DEAL::0.500000 0.312500 0.250000 -DEAL::0.500000 0.312500 0.312500 -DEAL::0.500000 0.312500 0.437500 -DEAL::0.500000 0.312500 0.687500 -DEAL::0.500000 0.312500 0.750000 -DEAL::0.500000 0.312500 0.812500 -DEAL::0.500000 0.312500 0.937500 -DEAL::0.500000 0.312500 1.00000 -DEAL::0.500000 0.312500 0.562500 -DEAL::0.500000 0.312500 0.500000 -DEAL::0.500000 0.375000 0.00000 -DEAL::0.500000 0.375000 0.0625000 -DEAL::0.500000 0.375000 0.187500 -DEAL::0.500000 0.375000 0.250000 -DEAL::0.500000 0.375000 0.312500 -DEAL::0.500000 0.375000 0.687500 -DEAL::0.500000 0.375000 0.750000 -DEAL::0.500000 0.375000 0.812500 -DEAL::0.500000 0.375000 0.937500 -DEAL::0.500000 0.375000 1.00000 -DEAL::0.500000 0.375000 0.437500 -DEAL::0.500000 0.375000 0.562500 -DEAL::0.500000 0.375000 0.500000 -DEAL::0.500000 0.437500 0.00000 -DEAL::0.500000 0.437500 0.0625000 -DEAL::0.500000 0.437500 0.187500 -DEAL::0.500000 0.437500 0.250000 -DEAL::0.500000 0.437500 0.312500 -DEAL::0.500000 0.437500 0.437500 -DEAL::0.500000 0.437500 0.687500 -DEAL::0.500000 0.437500 0.750000 -DEAL::0.500000 0.437500 0.812500 -DEAL::0.500000 0.437500 0.937500 -DEAL::0.500000 0.437500 1.00000 -DEAL::0.500000 0.437500 0.562500 -DEAL::0.500000 0.437500 0.500000 -DEAL::0.500000 0.500000 0.0625000 -DEAL::0.500000 0.500000 0.125000 -DEAL::0.500000 0.500000 0.187500 -DEAL::0.500000 0.500000 0.312500 -DEAL::0.500000 0.500000 0.375000 -DEAL::0.500000 0.500000 0.437500 -DEAL::0.500000 0.500000 0.562500 -DEAL::0.500000 0.500000 0.625000 -DEAL::0.500000 0.500000 0.687500 -DEAL::0.500000 0.500000 0.812500 -DEAL::0.500000 0.500000 0.875000 -DEAL::0.500000 0.500000 0.937500 -DEAL:: +DEAL::OK -- 2.39.5