From a59943d93e11b94ba63b6ecad550fc728c4084b6 Mon Sep 17 00:00:00 2001 From: kanschat Date: Sun, 30 Oct 2011 19:06:25 +0000 Subject: [PATCH] new patches git-svn-id: https://svn.dealii.org/trunk@24703 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/include/deal.II/dofs/dof_tools.h | 105 ++++++++++++++++++++++- deal.II/source/dofs/dof_tools.cc | 94 ++++++++++++++++++++ deal.II/source/dofs/dof_tools.inst.in | 8 ++ 3 files changed, 203 insertions(+), 4 deletions(-) diff --git a/deal.II/include/deal.II/dofs/dof_tools.h b/deal.II/include/deal.II/dofs/dof_tools.h index e387e694c4..540158a2d7 100644 --- a/deal.II/include/deal.II/dofs/dof_tools.h +++ b/deal.II/include/deal.II/dofs/dof_tools.h @@ -1492,10 +1492,14 @@ namespace DoFTools /** * Create a sparsity pattern which - * in each row lists the degrees of - * freedom associated to the - * corresponding cells - * corresponding to a vertex. + * in each row lists the degrees + * of freedom associated to the + * cells corresponding to a + * vertex. The sparsity pattern + * may be empty when entering this + * function and will be + * reinitialized to the correct + * size. * * The function has some boolean * arguments (lsited below) @@ -1553,6 +1557,99 @@ namespace DoFTools const bool level_boundary_patches = false, const bool single_cell_patches = false); + /** + * Create a sparsity pattern which + * in each row lists the degrees of + * freedom associated to the + * cells which are the children of + * the same cell. The + * sparsity pattern may be empty + * when entering this function and + * will be reinitialized to the + * correct size. + * + * The function has some boolean + * arguments (lsited below) + * controlling details of the + * generated patches. The default + * settings are those for + * Arnold-Falk-Winther type + * smoothers for divergence and + * curl conforming finite elements + * with essential boundary + * conditions. Other applications + * are possible, in particular + * changing + * boundary_dofs for + * non-essential boundary + * conditions. + * + * Since the patches are defined + * through refinement, th + * + * @arg block_list: the + * SparsityPattern into which the + * patches will be stored. + * @arg dof_handler: The + * multilevel dof handler + * providing the topology operated + * on. + * @arg + * interior_dofs_only: + * for each patch of cells around + * a vertex, collect only the + * interior degrees of freedom of + * the patch and disregard those + * on the boundary of the + * patch. This is for instance the + * setting for smoothers of + * Arnold-Falk-Winther type. + * @arg boundary_dofs: + * include degrees of freedom, + * which would have excluded by + * interior_dofs_only, + * but are lying on the boundary + * of the domain, and thus need smoothing. + */ + template + void make_child_patches(SparsityPattern& block_list, + const DH& dof_handler, + const unsigned int level, + const bool interior_dofs_only, + const bool boundary_dofs = false); + + /** + * Create a block list with only a + * single patch, which in turn + * contains all degrees of freedom + * on the given level. + * + * This function is mostly a + * closure on level 0 for functions + * like make_child_patches() and + * make_vertex_patches(), which may + * produce an empty patch list. + * + * @arg block_list: the + * SparsityPattern into which the + * patches will be stored. + * @arg dof_handler: The + * multilevel dof handler + * providing the topology operated + * on. + * @arg level The grid + * level used for building the list. + * @arg + * interior_dofs_only: + * if true, exclude degrees of freedom on + * the boundary of the domain. + */ + template + void make_single_patch(SparsityPattern& block_list, + const DH& dof_handler, + const unsigned int level, + const bool interior_dofs_only = false); + //@} /** * Extract a vector that represents the diff --git a/deal.II/source/dofs/dof_tools.cc b/deal.II/source/dofs/dof_tools.cc index 2de1b6646e..2a89f8b8c4 100644 --- a/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/source/dofs/dof_tools.cc @@ -5877,6 +5877,100 @@ namespace DoFTools } } + + template + void make_single_patch( + SparsityPattern& block_list, + const DH& dof_handler, + const unsigned int level, + const bool interior_only) + { + const FiniteElement& fe = dof_handler.get_fe(); + 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 + // freedom on faces + // opposite to the vertex + 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, + const DH& dof_handler, + const unsigned int level, + const bool interior_dofs_only, + const bool boundary_dofs) + { + Assert(level > 0 && level < dof_handler.get_tria().n_levels(), + ExcIndexRange(level, 1, dof_handler.get_tria().n_levels())); + + 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; + for (unsigned int child=0;childn_children();++child) + { + cell = pcell->child(child); + // For hp, only this line + // here would have to be + // replaced. + const FiniteElement& fe = dof_handler.get_fe(); + const unsigned int n_dofs = fe.dofs_per_cell; + indices.resize(n_dofs); + exclude.resize(n_dofs); + std::fill(exclude.begin(), exclude.end(), false); + cell->get_mg_dof_indices(indices); + + Assert(!interior_dofs_only, ExcNotImplemented()); + Assert(!boundary_dofs, ExcNotImplemented()); + + for (unsigned int i=0;i void make_vertex_patches( SparsityPattern& block_list, diff --git a/deal.II/source/dofs/dof_tools.inst.in b/deal.II/source/dofs/dof_tools.inst.in index d1d29ef355..c70c279e45 100644 --- a/deal.II/source/dofs/dof_tools.inst.in +++ b/deal.II/source/dofs/dof_tools.inst.in @@ -244,6 +244,14 @@ for (deal_II_dimension : DIMENSIONS) void DoFTools::make_vertex_patches (SparsityPattern&, const MGDoFHandler&, unsigned int, bool, bool, bool, bool); + template + void DoFTools::make_single_patch (SparsityPattern&, const MGDoFHandler&, + unsigned int, bool); + + template + void DoFTools::make_child_patches(SparsityPattern&, const MGDoFHandler&, + unsigned int, bool, bool); + #if deal_II_dimension < 3 template void -- 2.39.5