From: Wolfgang Bangerth Date: Fri, 6 Feb 2009 04:02:49 +0000 (+0000) Subject: Fix erroneous resizing of output argument in count_dofs_per_component/block. X-Git-Tag: v8.0.0~8061 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=6ef60b276bbee04643a4dafc8fe1403a6481922a;p=dealii.git Fix erroneous resizing of output argument in count_dofs_per_component/block. git-svn-id: https://svn.dealii.org/trunk@18339 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index bc720ee484..e4402d1949 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -3966,23 +3966,34 @@ DoFTools::count_dofs_per_component ( std::vector target_component) { const FiniteElement& fe = dof_handler.get_fe(); - const unsigned int n_components = fe.n_components(); - dofs_per_component.resize (n_components); - std::fill (dofs_per_component.begin(), dofs_per_component.end(), 0U); + std::fill (dofs_per_component.begin(), dofs_per_component.end(), 0U); + // If the empty vector was given as // default argument, set up this // vector as identity. if (target_component.size()==0) { - target_component.resize(n_components); - for (unsigned int i=0;i target_block) { const FiniteElement& fe = dof_handler.get_fe(); - const unsigned int n_blocks = fe.n_blocks(); - dofs_per_block.resize (n_blocks); + std::fill (dofs_per_block.begin(), dofs_per_block.end(), 0U); // If the empty vector was given as @@ -4069,13 +4079,25 @@ DoFTools::count_dofs_per_block ( // vector as identity. if (target_block.size()==0) { - target_block.resize(n_blocks); - for (unsigned int i=0;ideal.II
    +
  1. +

    + Changed: The DoFTools::count_dofs_per_component and + DoFTools::count_dofs_per_block erroneously resized the output argument + to the number of components or blocks in the finite element, respectively, + even if the target component/block list given as an additional argument + needed a different number of output elements. This is now fixed. +
    + (WB 2009/02/05) +

    +
  2. New: The GridGenerator::half_hyper_shell function now also exists in 3d. diff --git a/tests/bits/Makefile b/tests/bits/Makefile index 06e71ed3a2..88b61d211e 100644 --- a/tests/bits/Makefile +++ b/tests/bits/Makefile @@ -1,6 +1,6 @@ ############################################################ # Makefile,v 1.15 2002/06/13 12:51:13 hartmann Exp -# Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors +# Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors ############################################################ ############################################################ @@ -99,7 +99,8 @@ tests_x = grid_generator_?? \ joa_* \ neighboring_cells_at_two_faces \ solution_transfer \ - make_boundary_constraints_* + make_boundary_constraints_* \ + count_dofs_per_block* # tests for the hp branch: # fe_collection_* diff --git a/tests/bits/count_dofs_per_block_01.cc b/tests/bits/count_dofs_per_block_01.cc new file mode 100644 index 0000000000..c59bf926d9 --- /dev/null +++ b/tests/bits/count_dofs_per_block_01.cc @@ -0,0 +1,145 @@ +//---------------------------- count_dofs_per_block_01.cc --------------------------- +// $Id$ +// Version: $Name$ +// +// Copyright (C) 2009 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 --------------------------- + +// the DoFTools::count_dofs_per_{component,block} functions resized the output +// array to fe.n_components or fe.n_blocks even if a grouping argument was +// given. this would appear wrong and can lead to all sorts of interesting +// behavior if not caught by an assertion early enough + + +#include "../tests.h" +#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 (); + } + + // Taylor-Hood element + FESystem fe (FE_Q(1), dim, FE_DGQ(0), 1); + DoFHandler dof_handler (tria); + 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/bits/count_dofs_per_block_01/cmp/generic b/tests/bits/count_dofs_per_block_01/cmp/generic new file mode 100644 index 0000000000..1728b32c62 --- /dev/null +++ b/tests/bits/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