From 2e8411cf1b864a661effb3db9e355df3c014b124 Mon Sep 17 00:00:00 2001 From: young Date: Thu, 4 Dec 2008 11:21:23 +0000 Subject: [PATCH] Added function DoFTools::make_zero_boundary_constraints git-svn-id: https://svn.dealii.org/trunk@17841 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/dofs/dof_tools.h | 11 ++ deal.II/deal.II/source/dofs/dof_tools.cc | 219 +++++++++++++++++++++++ deal.II/doc/news/changes.h | 10 ++ 3 files changed, 240 insertions(+) diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index 3e9a12263b..40e901b079 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -1752,6 +1752,17 @@ class DoFTools convert_couplings_to_blocks (const hp::DoFHandler& dof_handler, const Table<2, Coupling>& table_by_component, std::vector >& tables_by_block); + + /** + * Make a constraint matrix with + * zero boundary values. + */ + template class DH> + static void + make_zero_boundary_constraints (const DH &dof, + ConstraintMatrix &zero_boundary_constraints, + const std::vector &component_mask_=std::vector()); + /** * Map a coupling table from the * user friendly organization by diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index f1189af9a5..014275230f 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -26,6 +26,7 @@ #include #include #include +#include #include #include #include @@ -34,6 +35,7 @@ #include #include #include +#include #include #include @@ -5132,6 +5134,217 @@ DoFTools::convert_couplings_to_blocks ( } } +#if deal_II_dimension == 1 + +template class DH> +void +DoFTools::make_zero_boundary_constraints (const DH &dof, + ConstraintMatrix &zero_boundary_constraints, + const std::vector &component_mask_) +{ + Assert ((component_mask_.size() == 0) || + (component_mask_.size() == dof.get_fe().n_components()), + ExcMessage ("The number of components in the mask has to be either " + "zero or equal to the number of components in the finite " + "element.")); + + // check for boundary cells in both + // directions to secure indicators + // for the entire boundary, i.e. 0 + // is left boundary and 1 is right + // boundary + for (unsigned int direction = 0; direction < 2; ++direction) + { + + // first find the outermost active + // cell by traversing the coarse + // grid to its end and then looking + // for the children + typename DH::cell_iterator outermost_cell = dof.begin(0); + while (outermost_cell->neighbor(direction).state() == IteratorState::valid) + outermost_cell = outermost_cell->neighbor(direction); + + while (outermost_cell->has_children()) + outermost_cell = outermost_cell->child(direction); + + // then get the FE corresponding to + // this cell + const FiniteElement &fe = outermost_cell->get_fe(); + + // set the component mask to either + // the original value or a vector + // of trues + const std::vector component_mask ((component_mask_.size() == 0) ? + std::vector (fe.n_components(), true) : + component_mask_); + Assert (std::count(component_mask.begin(), component_mask.end(), true) > 0, + VectorTools::ExcNoComponentSelected()); + + // cast zero boundary constraints + // onto a matrix if component_mask + // == true + for (unsigned int i=0; ivertex_dof_index (direction, i)); + + } +} + +#else + +template class DH> +void +DoFTools::make_zero_boundary_constraints (const DH &dof, + ConstraintMatrix &zero_boundary_constraints, + const std::vector &component_mask_) +{ + Assert ((component_mask_.size() == 0) || + (component_mask_.size() == dof.get_fe().n_components()), + ExcMessage ("The number of components in the mask has to be either " + "zero or equal to the number of components in the finite " + "element.")); + + const unsigned int n_components = DoFTools::n_components(dof); + const bool fe_is_system = (n_components != 1); + + // set the component mask to either + // the original value or a vector + // of trues + const std::vector component_mask ((component_mask_.size() == 0) ? + std::vector (n_components, true) : + component_mask_); + Assert (std::count(component_mask.begin(), component_mask.end(), true) > 0, + VectorTools::ExcNoComponentSelected()); + + // a field to store the indices + std::vector face_dofs; + face_dofs.reserve (DoFTools::max_dofs_per_face(dof)); + + typename DH::active_cell_iterator cell = dof.begin_active(), + endc = dof.end(); + for (; cell!=endc; ++cell) + for (unsigned int face_no = 0; face_no < GeometryInfo::faces_per_cell; + ++face_no) + { + const FiniteElement &fe = cell->get_fe(); + + // we can presently deal only with + // primitive elements for boundary + // values. make sure that all shape + // functions that are non-zero for + // the components we are interested + // in, are in fact primitive + for (unsigned int i=0; iget_fe().dofs_per_cell; ++i) + { + const std::vector &nonzero_component_array + = cell->get_fe().get_nonzero_components (i); + for (unsigned int c=0; cget_fe().is_primitive (i), + ExcMessage ("This function can only deal with requested boundary " + "values that correspond to primitive (scalar) base " + "elements")); + } + + typename DH::face_iterator face = cell->face(face_no); + if (face->boundary_indicator () == 0) + // face is of the right component + { + // get indices and physical + // location on this face + face_dofs.resize (fe.dofs_per_face); + face->get_dof_indices (face_dofs, cell->active_fe_index()); + + if (fe_is_system) + { + // enter those dofs into the list + // that match the component + // signature. + for (unsigned int i=0; i&, const Table<2, Coupling>&, std::vector >&); +template +void +DoFTools::make_zero_boundary_constraints +(const DoFHandler &, + ConstraintMatrix &, + const std::vector &); DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 4cf5ff5eea..d3266aada9 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -503,6 +503,16 @@ inconvenience this causes.

deal.II

    +
  1. +

    + New: The function DoFTools::make_zero_boundary_constraints() applies + zero boundary constraints to a (symmetric preserved) matrix without + accessing the sparsity pattern. Use cases are for when the sparsity pattern + of the matrix is unknown / inaccessible. +
    + (Toby D. Young 2008/12/04) +

    +
  2. Updated: The function ConstraintMatrix::distribute_local_to_global() for -- 2.39.5