From 25729e423dc1d6fc6cd0d7dfc0d251d68330b5be Mon Sep 17 00:00:00 2001 From: bangerth Date: Mon, 14 Nov 2011 05:14:46 +0000 Subject: [PATCH] Implement DoFTools::count_dofs_per_component for hp::DoFHandler as well. git-svn-id: https://svn.dealii.org/trunk@24747 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 5 + deal.II/include/deal.II/dofs/dof_tools.h | 18 +- deal.II/source/dofs/dof_tools.cc | 168 ++++++++++--------- deal.II/source/dofs/dof_tools.inst.in | 27 ++- tests/hp/count_dofs_per_block_01.cc | 154 +++++++++++++++++ tests/hp/count_dofs_per_block_01/cmp/generic | 19 +++ tests/hp/count_dofs_per_block_02.cc | 154 +++++++++++++++++ tests/hp/count_dofs_per_block_02/cmp/generic | 19 +++ 8 files changed, 472 insertions(+), 92 deletions(-) create mode 100644 tests/hp/count_dofs_per_block_01.cc create mode 100644 tests/hp/count_dofs_per_block_01/cmp/generic create mode 100644 tests/hp/count_dofs_per_block_02.cc create mode 100644 tests/hp/count_dofs_per_block_02/cmp/generic diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 5023a1fc9b..cc9bcb4456 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -34,6 +34,11 @@ inconvenience this causes.

General

    +
  1. New: The DoFTools::count_dofs_per_block function now also works +for objects of type hp::DoFHandler. +
    +(Jason Sheldon, 2011/11/13) +
  2. Fixed: The Intel compiler complains that it can't copy Trilinos vector reference objects, preventing the compiling of step-32. This is now fixed.
    diff --git a/deal.II/include/deal.II/dofs/dof_tools.h b/deal.II/include/deal.II/dofs/dof_tools.h index bc9870f6c1..e5202f3d36 100644 --- a/deal.II/include/deal.II/dofs/dof_tools.h +++ b/deal.II/include/deal.II/dofs/dof_tools.h @@ -1491,7 +1491,7 @@ namespace DoFTools const unsigned int level, const std::vector& selected_dofs = std::vector(), unsigned int offset = 0); - + /** * Create a sparsity pattern which * in each row lists the degrees @@ -1558,7 +1558,7 @@ namespace DoFTools const bool boundary_patches = false, const bool level_boundary_patches = false, const bool single_cell_patches = false); - + /** * Create a sparsity pattern which * in each row lists the degrees of @@ -1587,7 +1587,7 @@ namespace DoFTools * conditions. * * Since the patches are defined - * through refinement, th + * through refinement, th * * @arg block_list: the * SparsityPattern into which the @@ -1654,7 +1654,7 @@ namespace DoFTools const DH& dof_handler, const unsigned int level, const bool interior_dofs_only = false); - + //@} /** * Extract a vector that represents the @@ -1843,12 +1843,12 @@ namespace DoFTools * in the target_blocks argument * if given. */ - template + template void - count_dofs_per_block (const DoFHandler& dof_handler, - std::vector& dofs_per_block, - std::vector target_blocks - = std::vector()); + count_dofs_per_block (const DH &dof, + std::vector &dofs_per_block, + std::vector target_block + = std::vector()); /** * @deprecated See the previous diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc index 560e752a66..cf6181d11c 100644 --- a/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/source/dofs/dof_tools.cc @@ -4671,78 +4671,84 @@ namespace DoFTools - template + template void - count_dofs_per_block (const DoFHandler& dof_handler, - std::vector &dofs_per_block, - std::vector target_block) + count_dofs_per_block (const DH &dof_handler, + std::vector &dofs_per_block, + std::vector target_block) { - const FiniteElement& fe = dof_handler.get_fe(); - - std::fill (dofs_per_block.begin(), dofs_per_block.end(), 0U); + const dealii::hp::FECollection + fe_collection (dof_handler.get_fe()); + Assert (fe_collection.size() < 256, ExcNotImplemented()); - // If the empty vector was given as - // default argument, set up this - // vector as identity. - if (target_block.size()==0) + for (unsigned int this_fe=0; this_fe& fe = fe_collection[this_fe]; + std::fill (dofs_per_block.begin(), dofs_per_block.end(), 0U); + // If the empty vector was given as + // default argument, set up this + // vector as identity. + if (target_block.size()==0) + { + target_block.resize(fe.n_blocks()); + for (unsigned int i=0; i dofs_by_block (dof_handler.n_locally_owned_dofs()); - internal::extract_dofs_by_component (dof_handler, std::vector(), - true, dofs_by_block); + AssertDimension (dofs_per_block.size(), n_target_blocks); - // next count what we got - for (unsigned int block=0; block dofs_by_block (dof_handler.n_locally_owned_dofs()); + internal::extract_dofs_by_component (dof_handler, std::vector(), + true, dofs_by_block); + + // next count what we got + for (unsigned int block=0; block * tria - = (dynamic_cast*> - (&dof_handler.get_tria()))) - { - std::vector local_dof_count = dofs_per_block; - MPI_Allreduce ( &local_dof_count[0], &dofs_per_block[0], n_target_blocks, - MPI_UNSIGNED, MPI_SUM, tria->get_communicator()); - } + // if we are working on a parallel + // mesh, we now need to collect + // this information from all + // processors + if (const parallel::distributed::Triangulation * tria + = (dynamic_cast*> + (&dof_handler.get_tria()))) + { + std::vector local_dof_count = dofs_per_block; + MPI_Allreduce ( &local_dof_count[0], &dofs_per_block[0], n_target_blocks, + MPI_UNSIGNED, MPI_SUM, tria->get_communicator()); + } #endif #endif + } } @@ -5871,7 +5877,7 @@ namespace DoFTools typename DH::cell_iterator cell; typename DH::cell_iterator endc = dof_handler.end(level); std::vector indices; - + unsigned int i=0; for (cell=dof_handler.begin(level); cell != endc; ++i, ++cell) { @@ -5893,8 +5899,8 @@ namespace DoFTools } } } - - + + template void make_single_patch( SparsityPattern& block_list, @@ -5906,15 +5912,15 @@ namespace DoFTools block_list.reinit(1, dof_handler.n_dofs(level), dof_handler.n_dofs(level)); typename DH::cell_iterator cell; typename DH::cell_iterator endc = dof_handler.end(level); - + std::vector indices; std::vector exclude; - + for (cell=dof_handler.begin(level); cell != endc;++cell) { indices.resize(cell->get_fe().dofs_per_cell); cell->get_mg_dof_indices(indices); - + if (interior_only) { // Exclude degrees of @@ -5923,7 +5929,7 @@ namespace DoFTools exclude.resize(fe.dofs_per_cell); std::fill(exclude.begin(), exclude.end(), false); const unsigned int dpf = fe.dofs_per_face; - + for (unsigned int face=0;face::faces_per_cell;++face) if (cell->at_boundary(face) || cell->neighbor(face)->level() != cell->level()) for (unsigned int i=0;i void make_child_patches( SparsityPattern& block_list, @@ -5955,10 +5961,10 @@ namespace DoFTools typename DH::cell_iterator cell; typename DH::cell_iterator pcell = dof_handler.begin(level-1); typename DH::cell_iterator endc = dof_handler.end(level-1); - + std::vector indices; std::vector exclude; - + for (unsigned int block = 0;pcell != endc;++pcell) { if (!pcell->has_children()) continue; @@ -5982,7 +5988,7 @@ namespace DoFTools // which are on faces // of the parent const unsigned int dpf = fe.dofs_per_face; - + for (unsigned int d=0;d::vertex_to_face[child][d]; @@ -6000,7 +6006,7 @@ namespace DoFTools for (unsigned int i=0;i void make_vertex_patches( SparsityPattern& block_list, @@ -6037,7 +6043,7 @@ namespace DoFTools // Estimate for the number of // dofs at this point std::vector vertex_dof_count(dof_handler.get_tria().n_vertices(), 0); - + // Identify all vertices active // on this level and remember // some data about them @@ -6061,7 +6067,7 @@ namespace DoFTools // From now on, only vertices // with positive dof count are // "in". - + // Remove vertices at boundaries // or in corners for (unsigned int vg=0;vg indices; std::vector exclude; - + for (cell=dof_handler.begin(level); cell != endc; ++cell) { const FiniteElement& fe = cell->get_fe(); indices.resize(fe.dofs_per_cell); cell->get_mg_dof_indices(indices); - + for (unsigned int v=0;v::vertices_per_cell;++v) { const unsigned int vg = cell->vertex_index(v); const unsigned int block = vertex_mapping[vg]; if (block == numbers::invalid_unsigned_int) continue; - + if (interior_only) { // Exclude degrees of @@ -6122,7 +6128,7 @@ namespace DoFTools exclude.resize(fe.dofs_per_cell); std::fill(exclude.begin(), exclude.end(), false); const unsigned int dpf = fe.dofs_per_face; - + for (unsigned int d=0;d::vertex_to_face[v][d]; diff --git a/deal.II/source/dofs/dof_tools.inst.in b/deal.II/source/dofs/dof_tools.inst.in index 866832b649..996ba002a9 100644 --- a/deal.II/source/dofs/dof_tools.inst.in +++ b/deal.II/source/dofs/dof_tools.inst.in @@ -543,10 +543,21 @@ DoFTools::extract_level_dofs template void -DoFTools::count_dofs_per_block ( +DoFTools::count_dofs_per_block > ( const DoFHandler&, std::vector&, std::vector); +template +void +DoFTools::count_dofs_per_block > ( + const hp::DoFHandler&, + std::vector&, std::vector); + +template +void +DoFTools::count_dofs_per_block > ( + const MGDoFHandler&, + std::vector&, std::vector); template void @@ -619,10 +630,22 @@ DoFTools::map_dofs_to_support_points template void -DoFTools::count_dofs_per_block ( +DoFTools::count_dofs_per_block > ( const DoFHandler&, std::vector&, std::vector); +template +void +DoFTools::count_dofs_per_block > ( + const hp::DoFHandler&, + std::vector&, std::vector); + +template +void +DoFTools::count_dofs_per_block > ( + const MGDoFHandler&, + std::vector&, std::vector); + #endif template diff --git a/tests/hp/count_dofs_per_block_01.cc b/tests/hp/count_dofs_per_block_01.cc new file mode 100644 index 0000000000..a921f8e04f --- /dev/null +++ b/tests/hp/count_dofs_per_block_01.cc @@ -0,0 +1,154 @@ +//---------------------------- count_dofs_per_block_01.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2009, 2011 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- count_dofs_per_block_01.cc --------------------------- + +// like bits/count_dofs_per_block_01 but for hp::DoFHandler + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + + +std::string output_file_name = "count_dofs_per_block_01/output"; + + + +void print (const std::vector &v) +{ + deallog << v.size(); + for (unsigned int i=0; i +void +check () +{ + // create tria and dofhandler + // objects. set different boundary + // and sub-domain ids + Triangulation tria; + GridGenerator::hyper_cube(tria, 0., 1.); + tria.refine_global (1); + for (int i=0; i<2; ++i) + { + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement (); + } + + // produce two copies of the Taylor-Hood + // element + hp::FECollection fe; + fe.push_back (FESystem (FE_Q(1), dim, FE_DGQ(0), 1)); + fe.push_back (FESystem (FE_Q(1), dim, FE_DGQ(0), 1)); + + hp::DoFHandler dof_handler (tria); + + typename hp::DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + cell->set_active_fe_index (std::rand() % fe.size()); + + dof_handler.distribute_dofs (fe); + + // no grouping + { + std::vector dpc(dim+1); + DoFTools::count_dofs_per_component (dof_handler, dpc); + print (dpc); + } + + { + std::vector dpc(dim+1); + DoFTools::count_dofs_per_block (dof_handler, dpc); + print (dpc); + } + + + // grouping into less groups than + // components + { + std::vector group(dim+1, 0); + group[dim] = 1; + std::vector dpc(2); + DoFTools::count_dofs_per_component (dof_handler, dpc, false, group); + Assert (dpc.size() == 2, ExcInternalError()); + print (dpc); + } + + { + std::vector group(dim+1, 0); + group[dim] = 1; + std::vector dpc(2); + DoFTools::count_dofs_per_block (dof_handler, dpc, group); + Assert (dpc.size() == 2, ExcInternalError()); + print (dpc); + } + + // grouping into more groups than + // components + { + std::vector group(dim+1, 2*dim); + group[dim] = 0; + std::vector dpc(2*dim+1); + DoFTools::count_dofs_per_component (dof_handler, dpc, false, group); + Assert (dpc.size() == 2*dim+1, ExcInternalError()); + print (dpc); + } + + { + std::vector group(dim+1, 2*dim); + group[dim] = 0; + std::vector dpc(2*dim+1); + DoFTools::count_dofs_per_block (dof_handler, dpc, group); + Assert (dpc.size() == 2*dim+1, ExcInternalError()); + print (dpc); + } +} + + + +int main() +{ + std::ofstream logfile(output_file_name.c_str()); + logfile << std::setprecision (2); + deallog << std::setprecision (2); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + check<1> (); + check<2> (); + check<3> (); +} diff --git a/tests/hp/count_dofs_per_block_01/cmp/generic b/tests/hp/count_dofs_per_block_01/cmp/generic new file mode 100644 index 0000000000..1728b32c62 --- /dev/null +++ b/tests/hp/count_dofs_per_block_01/cmp/generic @@ -0,0 +1,19 @@ + +DEAL::2 5 4 +DEAL::2 5 4 +DEAL::2 5 4 +DEAL::2 5 4 +DEAL::3 4 0 5 +DEAL::3 4 0 5 +DEAL::3 18 18 10 +DEAL::3 18 18 10 +DEAL::2 36 10 +DEAL::2 36 10 +DEAL::5 10 0 0 0 36 +DEAL::5 10 0 0 0 36 +DEAL::4 60 60 60 22 +DEAL::4 60 60 60 22 +DEAL::2 180 22 +DEAL::2 180 22 +DEAL::7 22 0 0 0 0 0 180 +DEAL::7 22 0 0 0 0 0 180 diff --git a/tests/hp/count_dofs_per_block_02.cc b/tests/hp/count_dofs_per_block_02.cc new file mode 100644 index 0000000000..9a2bb0d3e9 --- /dev/null +++ b/tests/hp/count_dofs_per_block_02.cc @@ -0,0 +1,154 @@ +//---------------------------- count_dofs_per_block_02.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2009, 2011 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------- count_dofs_per_block_02.cc --------------------------- + +// like _01 but with different elements + + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + + +std::string output_file_name = "count_dofs_per_block_02/output"; + + + +void print (const std::vector &v) +{ + deallog << v.size(); + for (unsigned int i=0; i +void +check () +{ + // create tria and dofhandler + // objects. set different boundary + // and sub-domain ids + Triangulation tria; + GridGenerator::hyper_cube(tria, 0., 1.); + tria.refine_global (1); + for (int i=0; i<2; ++i) + { + tria.begin_active()->set_refine_flag(); + tria.execute_coarsening_and_refinement (); + } + + // produce two copies of the Taylor-Hood + // element + hp::FECollection fe; + fe.push_back (FESystem (FE_Q(1), dim, FE_DGQ(0), 1)); + fe.push_back (FESystem (FE_Q(2), dim, FE_DGQ(1), 1)); + + hp::DoFHandler dof_handler (tria); + + typename hp::DoFHandler::active_cell_iterator + cell = dof_handler.begin_active(), + endc = dof_handler.end(); + for (; cell!=endc; ++cell) + cell->set_active_fe_index (std::rand() % fe.size()); + + dof_handler.distribute_dofs (fe); + + // no grouping + { + std::vector dpc(dim+1); + DoFTools::count_dofs_per_component (dof_handler, dpc); + print (dpc); + } + + { + std::vector dpc(dim+1); + DoFTools::count_dofs_per_block (dof_handler, dpc); + print (dpc); + } + + + // grouping into less groups than + // components + { + std::vector group(dim+1, 0); + group[dim] = 1; + std::vector dpc(2); + DoFTools::count_dofs_per_component (dof_handler, dpc, false, group); + Assert (dpc.size() == 2, ExcInternalError()); + print (dpc); + } + + { + std::vector group(dim+1, 0); + group[dim] = 1; + std::vector dpc(2); + DoFTools::count_dofs_per_block (dof_handler, dpc, group); + Assert (dpc.size() == 2, ExcInternalError()); + print (dpc); + } + + // grouping into more groups than + // components + { + std::vector group(dim+1, 2*dim); + group[dim] = 0; + std::vector dpc(2*dim+1); + DoFTools::count_dofs_per_component (dof_handler, dpc, false, group); + Assert (dpc.size() == 2*dim+1, ExcInternalError()); + print (dpc); + } + + { + std::vector group(dim+1, 2*dim); + group[dim] = 0; + std::vector dpc(2*dim+1); + DoFTools::count_dofs_per_block (dof_handler, dpc, group); + Assert (dpc.size() == 2*dim+1, ExcInternalError()); + print (dpc); + } +} + + + +int main() +{ + std::ofstream logfile(output_file_name.c_str()); + logfile << std::setprecision (2); + deallog << std::setprecision (2); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + check<1> (); + check<2> (); + check<3> (); +} diff --git a/tests/hp/count_dofs_per_block_02/cmp/generic b/tests/hp/count_dofs_per_block_02/cmp/generic new file mode 100644 index 0000000000..d5d6551f28 --- /dev/null +++ b/tests/hp/count_dofs_per_block_02/cmp/generic @@ -0,0 +1,19 @@ + +DEAL::2 8 7 +DEAL::2 8 7 +DEAL::2 8 7 +DEAL::2 8 7 +DEAL::3 7 0 8 +DEAL::3 7 0 8 +DEAL::3 45 45 28 +DEAL::3 45 45 28 +DEAL::2 90 28 +DEAL::2 90 28 +DEAL::5 28 0 0 0 90 +DEAL::5 28 0 0 0 90 +DEAL::4 205 205 205 92 +DEAL::4 205 205 205 92 +DEAL::2 615 92 +DEAL::2 615 92 +DEAL::7 92 0 0 0 0 0 615 +DEAL::7 92 0 0 0 0 0 615 -- 2.39.5