/**
* 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)
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
+ * <tt>boundary_dofs</tt> for
+ * non-essential boundary
+ * conditions.
+ *
+ * Since the patches are defined
+ * through refinement, th
+ *
+ * @arg <tt>block_list</tt>: the
+ * SparsityPattern into which the
+ * patches will be stored.
+ * @arg <tt>dof_handler</tt>: The
+ * multilevel dof handler
+ * providing the topology operated
+ * on.
+ * @arg
+ * <tt>interior_dofs_only</tt>:
+ * 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 <tt>boundary_dofs</tt>:
+ * include degrees of freedom,
+ * which would have excluded by
+ * <tt>interior_dofs_only</tt>,
+ * but are lying on the boundary
+ * of the domain, and thus need smoothing.
+ */
+ template <class DH>
+ 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 <tt>block_list</tt>: the
+ * SparsityPattern into which the
+ * patches will be stored.
+ * @arg <tt>dof_handler</tt>: The
+ * multilevel dof handler
+ * providing the topology operated
+ * on.
+ * @arg <tt>level</tt> The grid
+ * level used for building the list.
+ * @arg
+ * <tt>interior_dofs_only</tt>:
+ * if true, exclude degrees of freedom on
+ * the boundary of the domain.
+ */
+ template <class DH>
+ 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
}
}
+
+ template <class DH>
+ void make_single_patch(
+ SparsityPattern& block_list,
+ const DH& dof_handler,
+ const unsigned int level,
+ const bool interior_only)
+ {
+ const FiniteElement<DH::dimension>& 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<unsigned int> indices;
+ std::vector<bool> 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<GeometryInfo<DH::dimension>::faces_per_cell;++face)
+ if (cell->at_boundary(face) || cell->neighbor(face)->level() != cell->level())
+ for (unsigned int i=0;i<dpf;++i)
+ exclude[fe.face_to_cell_index(i,face)] = true;
+ for (unsigned int j=0;j<indices.size();++j)
+ if (!exclude[j])
+ block_list.add(0, indices[j]);
+ }
+ else
+ {
+ for (unsigned int j=0;j<indices.size();++j)
+ block_list.add(0, indices[j]);
+ }
+ }
+ }
+
+
+ template <class DH>
+ 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<unsigned int> indices;
+ std::vector<bool> exclude;
+
+ for (unsigned int block = 0;pcell != endc;++pcell)
+ {
+ if (pcell->has_children()) continue;
+ for (unsigned int child=0;child<pcell->n_children();++child)
+ {
+ cell = pcell->child(child);
+ // For hp, only this line
+ // here would have to be
+ // replaced.
+ const FiniteElement<DH::dimension>& 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<n_dofs;++i)
+ if (!exclude[i])
+ block_list.add(block, indices[i]);
+ }
+ ++block;
+ }
+ }
+
+
+
template <class DH>
void make_vertex_patches(
SparsityPattern& block_list,