From 0e57166977b85a4576272b5356b7e380efd6a754 Mon Sep 17 00:00:00 2001 From: heltai Date: Tue, 15 May 2007 16:05:01 +0000 Subject: [PATCH] Fixed bug in find_cells_adjacent_to_vertex. Fixed behavior of count_dofs_per_block and count_dofs_per_component. git-svn-id: https://svn.dealii.org/trunk@14671 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 61 ++-- deal.II/deal.II/source/dofs/dof_tools.cc | 321 ++++------------------ deal.II/deal.II/source/grid/grid_tools.cc | 4 +- deal.II/doc/news/changes.html | 17 ++ 4 files changed, 92 insertions(+), 311 deletions(-) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index c1dfd56e6a..672720d2da 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -481,46 +481,6 @@ class DoFTools make_sparsity_pattern (const DH &dof, const std::vector > &mask, SparsityPattern &sparsity_pattern); - - /** - * Construct a sparsity pattern that - * allows coupling degrees of freedom on - * two different but related meshes. - * - * The idea is that if the two given - * DoFHandler objects correspond to two - * different meshes (and potentially to - * different finite elements used on - * these cells), but that if the two - * triangulations they are based on are - * derived from the same coarse mesh - * through hierarchical refinement, then - * one may set up a problem where one - * would like to test shape functions - * from one mesh against the shape - * functions from another mesh. In - * particular, this means that shape - * functions from a cell on the first - * mesh are tested against those on the - * second cell that are located on the - * corresponding cell; this - * correspondence is something that the - * IntergridMap class can determine. - * - * This function then constructs a - * sparsity pattern for which the degrees - * of freedom that represent the rows - * come from the first given DoFHandler, - * whereas the ones that correspond to - * columns come from the second - * DoFHandler. - */ - template - static - void - make_sparsity_pattern (const DH &dof_row, - const DH &dof_col, - SparsityPattern &sparsity); /** * Create the sparsity pattern for @@ -1103,8 +1063,17 @@ class DoFTools * velocities into one block, and * the pressure into another). * - * The result is returned in - * @p dofs_per_component. + * The result is returned in @p + * dofs_per_component. Note that + * the size of @p + * dofs_per_component needs to be + * enough to hold all the indices + * specified in @p + * target_component. If this is + * not the case, an assertion is + * thrown. The indices not + * targetted by target_components + * are left untouched. */ template static void @@ -1123,7 +1092,11 @@ class DoFTools * counting is done by * blocks. See @ref GlossBlock * "blocks" in the glossary for - * details. + * details. Again the vectors are + * assumed to have the correct + * size before calling this + * function. If this is not the + * case, an assertion is thrown. */ template static void @@ -1150,7 +1123,7 @@ class DoFTools * This function can be used when * different variables shall be * discretized on different - * grids, where one grid is + * grids, where one grid is * coarser than the other. This * idea might seem nonsensical at * first, but has reasonable diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index a47c140e62..f936d34124 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -20,7 +20,6 @@ #include #include #include -#include #include #include #include @@ -169,108 +168,6 @@ DoFTools::make_sparsity_pattern ( -template -void -DoFTools::make_sparsity_pattern ( - const DH &dof_row, - const DH &dof_col, - SparsityPattern &sparsity) -{ - const unsigned int n_dofs_row = dof_row.n_dofs(); - const unsigned int n_dofs_col = dof_col.n_dofs(); - - Assert (sparsity.n_rows() == n_dofs_row, - ExcDimensionMismatch (sparsity.n_rows(), n_dofs_row)); - Assert (sparsity.n_cols() == n_dofs_col, - ExcDimensionMismatch (sparsity.n_cols(), n_dofs_col)); - - - const std::list > - cell_list - = GridTools::get_finest_common_cells (dof_row, dof_col); - - - typename std::list > - ::const_iterator - cell_iter = cell_list.begin(); - - for (; cell_iter!=cell_list.end(); ++cell_iter) - { - const typename DH::cell_iterator cell_row = cell_iter->first; - const typename DH::cell_iterator cell_col = cell_iter->second; - - if (!cell_row->has_children() && !cell_col->has_children()) - { - const unsigned int dofs_per_cell_row = - cell_row->get_fe().dofs_per_cell; - const unsigned int dofs_per_cell_col = - cell_col->get_fe().dofs_per_cell; - std::vector - local_dof_indices_row(dofs_per_cell_row); - std::vector - local_dof_indices_col(dofs_per_cell_col); - cell_row->get_dof_indices (local_dof_indices_row); - cell_col->get_dof_indices (local_dof_indices_col); - for (unsigned int i=0; ihas_children()) - { - const std::vector - child_cells = GridTools::get_active_child_cells (cell_row); - for (unsigned int i=0; iget_fe().dofs_per_cell; - const unsigned int dofs_per_cell_col = - cell_col->get_fe().dofs_per_cell; - std::vector - local_dof_indices_row(dofs_per_cell_row); - std::vector - local_dof_indices_col(dofs_per_cell_col); - cell_row_child->get_dof_indices (local_dof_indices_row); - cell_col->get_dof_indices (local_dof_indices_col); - for (unsigned int i=0; i - child_cells = GridTools::get_active_child_cells (cell_col); - for (unsigned int i=0; iget_fe().dofs_per_cell; - const unsigned int dofs_per_cell_col = - cell_col_child->get_fe().dofs_per_cell; - std::vector - local_dof_indices_row(dofs_per_cell_row); - std::vector - local_dof_indices_col(dofs_per_cell_col); - cell_row->get_dof_indices (local_dof_indices_row); - cell_col_child->get_dof_indices (local_dof_indices_col); - for (unsigned int i=0; i @@ -1067,9 +964,7 @@ namespace internal ExcInternalError()); Assert (fe2.dofs_per_line <= fe1.dofs_per_line, ExcInternalError()); - Assert ((dim < 3) - || - (fe2.dofs_per_quad <= fe1.dofs_per_quad), + Assert (fe2.dofs_per_quad <= fe1.dofs_per_quad, ExcInternalError()); // the idea here is to designate as @@ -1483,20 +1378,12 @@ namespace internal /** - * Copy constraints into a constraint - * matrix object. - * - * This function removes zero + * This method removes zero * constraints and those, which * constrain a DoF which was * already eliminated in one of * the previous steps of the hp * hanging node procedure. - * - * It also suppresses very small - * entries in the constraint matrix to - * avoid making the sparsity pattern - * fuller than necessary. */ #ifdef DEAL_II_ANON_NAMESPACE_BUG static @@ -1536,36 +1423,10 @@ namespace internal } if (constraint_already_satisfied == false) - { - // add up the absolute - // values of all - // constraints in this line - // to get a measure of - // their absolute size - double abs_sum = 0; - for (unsigned int i=0; i= 1e-14*abs_sum)) + if (face_constraints(row,i) != 0) { #ifdef WOLFGANG std::cout << " " << slave_dofs[row] @@ -2625,58 +2486,17 @@ namespace internal // (this is what // happens in // hp/hp_hanging_nodes_01 - // for example). - // - // another possibility - // is what happens in - // crash_13. there, we - // have - // FESystem(FE_Q(1),FE_DGQ(0)) - // vs. FESystem(FE_Q(1),FE_DGQ(1)). - // neither of them - // dominates the - // other. the point is - // that it doesn't - // matter, since - // hopefully for this - // case, both sides' - // dofs should have - // been unified. - // - // make sure this is - // actually true. this - // actually only - // matters, of course, - // if either of the two - // finite elements - // actually do have - // dofs on the face - if ((cell->get_fe().dofs_per_face != 0) - || - (cell->neighbor(face)->get_fe().dofs_per_face != 0)) + // for example); check + // the latter somewhat + // crudely by comparing + // fe names + if (cell->get_fe().get_name() != + cell->neighbor(face)->get_fe().get_name()) { - Assert (cell->get_fe().dofs_per_face - == - cell->neighbor(face)->get_fe().dofs_per_face, + Assert (cell->get_fe().dofs_per_face == 0, + ExcNotImplemented()); + Assert (cell->neighbor(face)->get_fe().dofs_per_face == 0, ExcNotImplemented()); - - // (ab)use the master - // and slave dofs - // arrays for a - // moment here - master_dofs.resize (cell->get_fe().dofs_per_face); - cell->face(face) - ->get_dof_indices (master_dofs, - cell->active_fe_index ()); - - slave_dofs.resize (cell->neighbor(face)->get_fe().dofs_per_face); - cell->face(face) - ->get_dof_indices (slave_dofs, - cell->neighbor(face)->active_fe_index ()); - - for (unsigned int i=0; iget_fe().dofs_per_face; ++i) - Assert (master_dofs[i] == slave_dofs[i], - ExcInternalError()); } break; @@ -3455,7 +3275,7 @@ DoFTools::extract_hanging_node_dofs (const DoFHandler<2> &dof_handler, Assert(selected_dofs.size() == dof_handler.n_dofs(), ExcDimensionMismatch(selected_dofs.size(), dof_handler.n_dofs())); // preset all values by false - std::fill_n (selected_dofs.begin(), dof_handler.n_dofs(), false); + fill_n (selected_dofs.begin(), dof_handler.n_dofs(), false); const FiniteElement &fe = dof_handler.get_fe(); @@ -3676,9 +3496,6 @@ DoFTools::count_dofs_per_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); - // If the empty vector was given as // default argument, set up this // vector as identity. @@ -3691,6 +3508,21 @@ DoFTools::count_dofs_per_component ( Assert(target_component.size()==n_components, ExcDimensionMismatch(target_component.size(),n_components)); + + + // Check that target components + // contains sensible information + // and reset the counters for + // dofs_per_component + unsigned int size_dst = dofs_per_component.size(); + for(unsigned int i=0; i& 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); + unsigned int size_dst = dofs_per_block.size(); // If the empty vector was given as // default argument, set up this @@ -3785,14 +3620,22 @@ DoFTools::count_dofs_per_block ( Assert(target_block.size()==n_blocks, ExcDimensionMismatch(target_block.size(),n_blocks)); - + + // Check constistency and reset + // counters + for(unsigned int i=0; i, const Table<2,Coupling>&, CompressedBlockSparsityPattern&); - -template void -DoFTools::make_sparsity_pattern, - SparsityPattern> -(const DoFHandler &dof_row, - const DoFHandler &dof_col, - SparsityPattern &sparsity); -template void -DoFTools::make_sparsity_pattern, - CompressedSparsityPattern> -(const DoFHandler &dof_row, - const DoFHandler &dof_col, - CompressedSparsityPattern &sparsity); -template void -DoFTools::make_sparsity_pattern, - BlockSparsityPattern> -(const DoFHandler &dof_row, - const DoFHandler &dof_col, - BlockSparsityPattern &sparsity); -template void -DoFTools::make_sparsity_pattern, - CompressedBlockSparsityPattern> -(const DoFHandler &dof_row, - const DoFHandler &dof_col, - CompressedBlockSparsityPattern &sparsity); - -template void -DoFTools::make_sparsity_pattern, - SparsityPattern> -(const hp::DoFHandler &dof_row, - const hp::DoFHandler &dof_col, - SparsityPattern &sparsity); -template void -DoFTools::make_sparsity_pattern, - CompressedSparsityPattern> -(const hp::DoFHandler &dof_row, - const hp::DoFHandler &dof_col, - CompressedSparsityPattern &sparsity); -template void -DoFTools::make_sparsity_pattern, - BlockSparsityPattern> -(const hp::DoFHandler &dof_row, - const hp::DoFHandler &dof_col, - BlockSparsityPattern &sparsity); -template void -DoFTools::make_sparsity_pattern, - CompressedBlockSparsityPattern> -(const hp::DoFHandler &dof_row, - const hp::DoFHandler &dof_col, - CompressedBlockSparsityPattern &sparsity); - - // #if deal_II_dimension > 1 template void DoFTools::make_boundary_sparsity_pattern,SparsityPattern> diff --git a/deal.II/deal.II/source/grid/grid_tools.cc b/deal.II/deal.II/source/grid/grid_tools.cc index c30d15faa2..6c0bec81b5 100644 --- a/deal.II/deal.II/source/grid/grid_tools.cc +++ b/deal.II/deal.II/source/grid/grid_tools.cc @@ -625,7 +625,7 @@ GridTools::find_cells_adjacent_to_vertex(const Container &container, // all faces to which this vertex // belongs, check the level of // the neighbor, and if it is coarser, - // then check whether the vortex is + // then check whether the vertex is // part of that neighbor or not. for (unsigned vface = 0; vface < dim; vface++) { @@ -661,7 +661,7 @@ GridTools::find_cells_adjacent_to_vertex(const Container &container, { bool found = false; for (unsigned v=0; v::vertices_per_cell; v++) - if (cell->vertex_index(v) == (int)vertex) + if (nb->vertex_index(v) == (int)vertex) { found = true; break; diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index eee2283bb3..d8cf807d61 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -1033,6 +1033,23 @@ inconvenience this causes.

deal.II

    +
  1. Fixed: the function + GridTools::find_cells_adjacent_to_vertex + was not detecting properly the coarse cells adjacent to + refined cells. +
    + (Luca Heltai 2007/05/15) +

    + +
  2. Fixed: the two tools + DoFTools::count_dofs_per_component and + DoFTools::count_dofs_per_block where changing the + size of the destination vector. Consistently with (most of) the + rest of the library, now the vectors are expected to be the + right size before calling these functions. +
    + (Luca Heltai 2007/05/15) +

  3. New: The classes MGTransferBlockSelect and