From: bangerth fe.component_mask(block_mask)
where block_mask
+ * is a variable of type BlockMask. In other words, if you have a block
+ * mask but need to call a function that only accepts a component mask, this
+ * syntax can be used to obtain the necessary component mask.
+ *
+ * Creation of block masks:
+ * Block masks are typically created by asking the finite element
+ * to generate a block mask from certain selected vector components using
+ * code such as this where we create a mask that only denotes the
+ * velocity components of a Stokes element (see @ref vector_valued):
+ * @code
+ * FESystem[false, true]
. Similarly, using
+ * @code
+ * FEValuesExtractors::Vector velocities(0);
+ * BlockMask velocity_mask = stokes_fe.block_mask (velocities);
+ * @endcode
+ * would result in a mask [true, false]
in any dimension.
+ *
+ * Note, however, that if we had defined the finite element in the following
+ * way:
+ * @code
+ * FESystem[false, false, true]
because the element has
+ * dim+1
components and equally many blocks. See the
+ * discussion on what a block represents exactly in the
+ * @ref GlossBlock "block entry of this glossary".
+ * std::vector@
, i.e., an array with as many entries as the
- * finite element has components (e.g., in the Stokes case, there are
+ * ComponentMask object which one can think of as an array with
+ * as many entries as the finite element has components (e.g., in the Stokes case, there are
* dim+1
components) and where each entry is either true or false.
* In the example where we would like to interpolate boundary values only for
* the velocity components of the Stokes system, this component mask would then
@@ -296,12 +381,55 @@
* functions with these names but only some of them have a component mask
* argument.
*
- * Many of the functions that take a component mask accept a vector of length
- * zero to indicate all components, i.e., as if the vector had the
+ * Semantics of component masks:
+ * Many of the functions that take a component mask object that has been default
+ * constructed to indicate all components, i.e., as if the vector had the
* correct length and was filled with only true
values. The reason
- * is that the empty vector can be constructed in place using the code snippet
- * std::vector@
and can thus be used as a default
+ * is that default initialized objects can be constructed in place using the code snippet
+ * ComponentMask()
and can thus be used as a default
* argument in function signatures.
+ *
+ * In other words, ComponentMask objects can be in one of two states: They can have
+ * been initialized by a vector of booleans with a nonzero length; in that case,
+ * they represent a mask of a particular length where some elements may be true
+ * and others may be false. Or, the ComponentMask may have been default initialized
+ * (using the default constructor) in which case it represents an array of indefinite
+ * length (i.e., a length appropriate to the circumstances) in which every entry
+ * is true.
+ *
+ * Creation of component masks:
+ * Component masks are typically created by asking the finite element
+ * to generate a component mask from certain selected components using
+ * code such as this where we create a mask that only denotes the
+ * velocity components of a Stokes element (see @ref vector_valued):
+ * @code
+ * FESystem[false, false, true]
. Similarly, using
+ * @code
+ * FEValuesExtractors::Vector velocities(0);
+ * ComponentMask velocity_mask = stokes_fe.component_mask (velocities);
+ * @endcode
+ * would result in a mask [true, true, false]
in 2d. Of
+ * course, in 3d, the result would be [true, true, true, false]
.
+ *
+ * @note Just as one can think of composed elements as being made up of
+ * @ref GlossComponent "components" or @ref GlossBlock "blocks", there are
+ * component masks (represented by the ComponentMask class) and
+ * @ref GlossBlockMask "block masks" (represented by the BlockMask class).
+ * The FiniteElement class has functions that convert between the two kinds of
+ * objects.
+ *
+ * @note Not all component masks actually make sense. For example, if you have
+ * a FE_RaviartThomas object in 2d, then it doesn't make any sense to have a
+ * component mask of the form [true, false]
because you try to
+ * select individual vector components of a finite element where each shape
+ * function has both $x$ and $y$ velocities. In essence, while you can of
+ * course create such a component mask, there is nothing you can do with it.
*
*
*
diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h
index 6811424846..d935b2f623 100644
--- a/deal.II/doc/news/changes.h
+++ b/deal.II/doc/news/changes.h
@@ -24,12 +24,40 @@ inconvenience this causes.
operator[]
to query whether a particular block
+ * has been selected.
+ *
+ * The semantics of this class are the same as the related ComponentMask
+ * class, i.e., a default constructed mask represents all possible
+ * blocks. See there for more information about these semantics.
+ *
+ * Objects of this kind are used in many places where one wants to restrict
+ * operations to a certain subset of blocks, e.g. in DoFTools::extract_dofs.
+ * These objects can
+ * either be created by hand, or, simpler, by asking the finite element
+ * to generate a block mask from certain selected blocks using
+ * code such as this where we create a mask that only denotes the
+ * velocity block of a Stokes element (see @ref vector_valued):
+ * @code
+ * FESystem[false, true]
. (Compare this to the corresponding
+ * component mask discussed in the ComponentMask documentation.)
+ * Similarly, using
+ * @code
+ * FEValuesExtractors::Vector velocities(0);
+ * BlockMask velocity_mask = stokes_fe.block_mask (velocities);
+ * @endcode
+ * would result in a mask [true, false]
in both 2d and 3d.
+ *
+ * @ingroup fe
+ * @author Wolfgang Bangerth
+ * @date 2012
+ * @ingroup vector_valued
+ */
+class BlockMask
+{
+ public:
+ /**
+ * Initialize a block mask. The default is that a block
+ * mask represents a set of blocks that are all
+ * selected, i.e., calling this constructor results in
+ * a block mask that always returns true
+ * whenever asked whether a block is selected.
+ */
+ BlockMask ();
+
+ /**
+ * Initialize an object of this type with a set of selected
+ * blocks specified by the argument.
+ *
+ * @param block_mask A vector of true/false
+ * entries that determine which blocks of a finite element
+ * are selected. If the length of the given vector is zero,
+ * then this interpreted as the case where every block
+ * is selected.
+ */
+ BlockMask (const std::vectorn
blocks. This is true if either
+ * it was initilized with a vector with exactly n
+ * entries of type bool
(in this case, @p n must
+ * equal the result of size()) or if it was initialized
+ * with an empty vector (or using the default constructor) in
+ * which case it can represent a mask with an arbitrary number
+ * of blocks and will always say that a block is
+ * selected.
+ */
+ bool
+ represents_n_blocks (const unsigned int n) const;
+
+ /**
+ * Return the number of blocks that are selected by this mask.
+ *
+ * Since empty block masks represent a block mask that
+ * would return true
for every block, this
+ * function may not know the true size of the block
+ * mask and it therefore requires an argument that denotes the
+ * overall number of blocks.
+ *
+ * If the object has been initialized with a non-empty mask (i.e.,
+ * if the size() function returns something greater than zero,
+ * or equivalently if represents_the_all_selected_mask() returns
+ * false) then the argument can be omitted and the result of size()
+ * is taken.
+ */
+ unsigned int
+ n_selected_blocks (const unsigned int overall_number_of_blocks = numbers::invalid_unsigned_int) const;
+
+ /**
+ * Return the index of the first selected block. The argument
+ * is there for the same reason it exists with the
+ * n_selected_blocks() function.
+ *
+ * The function throws an exception if no block is selected at all.
+ */
+ unsigned int
+ first_selected_block (const unsigned int overall_number_of_blocks = numbers::invalid_unsigned_int) const;
+
+ /**
+ * Return true if this mask represents a default
+ * constructed mask that corresponds to one in which
+ * all blocks are selected. If true, then the size()
+ * function will return zero.
+ */
+ bool
+ represents_the_all_selected_mask () const;
+
+ /**
+ * Return a block mask that contains the union of the
+ * blocks selected by the current object and the one
+ * passed as an argument.
+ */
+ BlockMask operator | (const BlockMask &mask) const;
+
+ /**
+ * Return a block mask that has only those elements set that
+ * are set both in the current object as well as the one
+ * passed as an argument.
+ */
+ BlockMask operator & (const BlockMask &mask) const;
+
+ /**
+ * Return whether this object and the argument are identical.
+ */
+ bool operator== (const BlockMask &mask) const;
+
+ /**
+ * Determine an estimate for the
+ * memory consumption (in bytes)
+ * of this object.
+ */
+ std::size_t
+ memory_consumption () const;
+
+ private:
+ /**
+ * The actual block mask.
+ */
+ std::vector[all blocks selected]
to the
+ * stream. Otherwise, it prints the block mask in a form like
+ * [true,true,true,false]
.
+ *
+ * @param out The stream to write to.
+ * @param mask The mask to write.
+ * @return A reference to the first argument.
+ */
+std::ostream & operator << (std::ostream &out,
+ const BlockMask &mask);
+
+
+// -------------------- inline functions ---------------------
+
+inline
+BlockMask::BlockMask()
+{}
+
+
+inline
+BlockMask::BlockMask(const std::vectoroperator[]
to query whether a particular component
+ * has been selected.
+ *
+ * Objects of this kind are used in many places where one wants to restrict
+ * operations to a certain subset of components, e.g. in DoFTools::make_zero_boundary_values
+ * or VectorTools::interpolate_boundary_values. These objects can
+ * either be created by hand, or, simpler, by asking the finite element
+ * to generate a component mask from certain selected components using
+ * code such as this where we create a mask that only denotes the
+ * velocity components of a Stokes element (see @ref vector_valued):
+ * @code
+ * FESystem[false, false, true]
. Similarly, using
+ * @code
+ * FEValuesExtractors::Vector velocities(0);
+ * ComponentMask velocity_mask = stokes_fe.component_mask (velocities);
+ * @endcode
+ * would result in a mask [true, true, false]
in 2d. Of
+ * course, in 3d, the result would be [true, true, true, false]
.
+ *
+ * @ingroup fe
+ * @author Wolfgang Bangerth
+ * @date 2012
+ * @ingroup vector_valued
+ */
+class ComponentMask
+{
+ public:
+ /**
+ * Initialize a component mask. The default is that a component
+ * mask represents a set of components that are all
+ * selected, i.e., calling this constructor results in
+ * a component mask that always returns true
+ * whenever asked whether a component is selected.
+ */
+ ComponentMask ();
+
+ /**
+ * Initialize an object of this type with a set of selected
+ * components specified by the argument.
+ *
+ * @param component_mask A vector of true/false
+ * entries that determine which components of a finite element
+ * are selected. If the length of the given vector is zero,
+ * then this interpreted as the case where every component
+ * is selected.
+ */
+ ComponentMask (const std::vectorn
components. This is true if either
+ * it was initilized with a vector with exactly n
+ * entries of type bool
(in this case, @p n must
+ * equal the result of size()) or if it was initialized
+ * with an empty vector (or using the default constructor) in
+ * which case it can represent a mask with an arbitrary number
+ * of components and will always say that a component is
+ * selected.
+ */
+ bool
+ represents_n_components (const unsigned int n) const;
+
+ /**
+ * Return the number of components that are selected by this mask.
+ *
+ * Since empty component masks represent a component mask that
+ * would return true
for every component, this
+ * function may not know the true size of the component
+ * mask and it therefore requires an argument that denotes the
+ * overall number of components.
+ *
+ * If the object has been initialized with a non-empty mask (i.e.,
+ * if the size() function returns something greater than zero,
+ * or equivalently if represents_the_all_selected_mask() returns
+ * false) then the argument can be omitted and the result of size()
+ * is taken.
+ */
+ unsigned int
+ n_selected_components (const unsigned int overall_number_of_components = numbers::invalid_unsigned_int) const;
+
+ /**
+ * Return the index of the first selected component. The argument
+ * is there for the same reason it exists with the
+ * n_selected_components() function.
+ *
+ * The function throws an exception if no component is selected at all.
+ */
+ unsigned int
+ first_selected_component (const unsigned int overall_number_of_components = numbers::invalid_unsigned_int) const;
+
+ /**
+ * Return true if this mask represents a default
+ * constructed mask that corresponds to one in which
+ * all components are selected. If true, then the size()
+ * function will return zero.
+ */
+ bool
+ represents_the_all_selected_mask () const;
+
+ /**
+ * Return a component mask that contains the union of the
+ * components selected by the current object and the one
+ * passed as an argument.
+ */
+ ComponentMask operator | (const ComponentMask &mask) const;
+
+ /**
+ * Return a component mask that has only those elements set that
+ * are set both in the current object as well as the one
+ * passed as an argument.
+ */
+ ComponentMask operator & (const ComponentMask &mask) const;
+
+ /**
+ * Return whether this object and the argument are identical.
+ */
+ bool operator== (const ComponentMask &mask) const;
+
+ /**
+ * Determine an estimate for the
+ * memory consumption (in bytes)
+ * of this object.
+ */
+ std::size_t
+ memory_consumption () const;
+
+ private:
+ /**
+ * The actual component mask.
+ */
+ std::vector[all components selected]
to the
+ * stream. Otherwise, it prints the component mask in a form like
+ * [true,true,true,false]
.
+ *
+ * @param out The stream to write to.
+ * @param mask The mask to write.
+ * @return A reference to the first argument.
+ */
+std::ostream & operator << (std::ostream &out,
+ const ComponentMask &mask);
+
+
+// -------------------- inline functions ---------------------
+
+inline
+ComponentMask::ComponentMask()
+{}
+
+
+inline
+ComponentMask::ComponentMask(const std::vectordim
components are true that correspond to the given
+ * argument. See @ref GlossComponentMask "the glossary"
+ * for more information.
+ *
+ * @param vector An object that represents dim
+ * vector components of this finite element.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ ComponentMask
+ component_mask (const FEValuesExtractors::Vector &vector) const;
+
+ /**
+ * Return a component mask with as many elements as this
+ * object has vector components and of which exactly the
+ * dim*(dim+1)/2
components are true that
+ * correspond to the given argument. See @ref GlossComponentMask "the glossary"
+ * for more information.
+ *
+ * @param sym_tensor An object that represents dim*(dim+1)/2
+ * components of this finite element that are jointly to be
+ * interpreted as forming a symmetric tensor.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ ComponentMask
+ component_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
+
+ /**
+ * Given a block mask (see @ref GlossBlockMask "this glossary entry"),
+ * produce a component mask (see @ref GlossComponentMask "this glossary entry")
+ * that represents the components that correspond to the blocks selected in
+ * the input argument. This is essentially a conversion operator from
+ * BlockMask to ComponentMask.
+ *
+ * @param block_mask The mask that selects individual blocks of the finite
+ * element
+ * @return A mask that selects those components corresponding to the selected
+ * blocks of the input argument.
+ */
+ ComponentMask
+ component_mask (const BlockMask &block_mask) const;
+
+ /**
+ * Return a block mask with as many elements as this
+ * object has blocks and of which exactly the
+ * one component is true that corresponds to the given
+ * argument. See @ref GlossBlockMask "the glossary"
+ * for more information.
+ *
+ * @note This function will only succeed if the scalar referenced
+ * by the argument encompasses a complete block. In other words,
+ * if, for example, you pass an extractor for the single
+ * $x$ velocity and this object represents an FE_RaviartThomas
+ * object, then the single scalar object you selected is part
+ * of a larger block and consequently there is no block mask that
+ * would represent it. The function will then produce an exception.
+ *
+ * @param scalar An object that represents a single scalar
+ * vector component of this finite element.
+ * @return A component mask that is false in all components
+ * except for the one that corresponds to the argument.
+ */
+ BlockMask
+ block_mask (const FEValuesExtractors::Scalar &scalar) const;
+
+ /**
+ * Return a component mask with as many elements as this
+ * object has vector components and of which exactly the
+ * dim
components are true that correspond to the given
+ * argument. See @ref GlossBlockMask "the glossary"
+ * for more information.
+ *
+ * @note The same caveat applies as to the version of the function above:
+ * The extractor object passed as argument must be so that it corresponds
+ * to full blocks and does not split blocks of this element.
+ *
+ * @param vector An object that represents dim
+ * vector components of this finite element.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ BlockMask
+ block_mask (const FEValuesExtractors::Vector &vector) const;
+
+ /**
+ * Return a component mask with as many elements as this
+ * object has vector components and of which exactly the
+ * dim*(dim+1)/2
components are true that
+ * correspond to the given argument. See @ref GlossBlockMask "the glossary"
+ * for more information.
+ *
+ * @note The same caveat applies as to the version of the function above:
+ * The extractor object passed as argument must be so that it corresponds
+ * to full blocks and does not split blocks of this element.
+ *
+ * @param sym_tensor An object that represents dim*(dim+1)/2
+ * components of this finite element that are jointly to be
+ * interpreted as forming a symmetric tensor.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ BlockMask
+ block_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
+
+ /**
+ * Given a component mask (see @ref GlossComponentMask "this glossary entry"),
+ * produce a block mask (see @ref GlossBlockMask "this glossary entry")
+ * that represents the blocks that correspond to the components selected in
+ * the input argument. This is essentially a conversion operator from
+ * ComponentMask to BlockMask.
+ *
+ * @note This function will only succeed if the components referenced
+ * by the argument encompasses complete blocks. In other words,
+ * if, for example, you pass an component mask for the single
+ * $x$ velocity and this object represents an FE_RaviartThomas
+ * object, then the single component you selected is part
+ * of a larger block and consequently there is no block mask that
+ * would represent it. The function will then produce an exception.
+ *
+ * @param component_mask The mask that selects individual components of the finite
+ * element
+ * @return A mask that selects those blocks corresponding to the selected
+ * blocks of the input argument.
+ */
+ BlockMask
+ block_mask (const ComponentMask &component_mask) const;
+
//@}
/**
@@ -2213,7 +2369,7 @@ class FiniteElement : public Subscriptor,
*/
static
std::vectordim
elements, and second order symmetric tensors consisting of
- * (dim*dim + dim)/2
components
- *
- * See the description of the @ref vector_valued module for examples how to
- * use the features of this namespace.
- *
- * @ingroup feaccess vector_valued
- */
-namespace FEValuesExtractors
-{
- /**
- * Extractor for a single scalar component
- * of a vector-valued element. The result
- * of applying an object of this type to an
- * FEValues, FEFaceValues or
- * FESubfaceValues object is of type
- * FEValuesViews::Scalar. The concept of
- * extractors is defined in the
- * documentation of the namespace
- * FEValuesExtractors and in the @ref
- * vector_valued module.
- *
- * @ingroup feaccess vector_valued
- */
- struct Scalar
- {
- /**
- * The selected scalar component of the
- * vector.
- */
- unsigned int component;
-
- /**
- * Constructor. Take the selected
- * vector component as argument.
- */
- Scalar (const unsigned int component);
- };
-
-
- /**
- * Extractor for a vector of
- * spacedim
components of a
- * vector-valued element. The value of
- * spacedim
is defined by the
- * FEValues object the extractor is applied
- * to. The result of applying an object of
- * this type to an FEValues, FEFaceValues
- * or FESubfaceValues object is of type
- * FEValuesViews::Vector.
- *
- * The concept of
- * extractors is defined in the
- * documentation of the namespace
- * FEValuesExtractors and in the @ref
- * vector_valued module.
- *
- * Note that in the current context, a
- * vector is meant in the sense physics
- * uses it: it has spacedim
- * components that behave in specific ways
- * under coordinate system
- * transformations. Examples include
- * velocity or displacement fields. This is
- * opposed to how mathematics uses the word
- * "vector" (and how we use this word in
- * other contexts in the library, for
- * example in the Vector class), where it
- * really stands for a collection of
- * numbers. An example of this latter use
- * of the word could be the set of
- * concentrations of chemical species in a
- * flame; however, these are really just a
- * collection of scalar variables, since
- * they do not change if the coordinate
- * system is rotated, unlike the components
- * of a velocity vector, and consequently,
- * this class should not be used for this
- * context.
- *
- * @ingroup feaccess vector_valued
- */
- struct Vector
- {
- /**
- * The first component of the vector
- * view.
- */
- unsigned int first_vector_component;
-
- /**
- * Constructor. Take the first
- * component of the selected vector
- * inside the FEValues object as
- * argument.
- */
- Vector (const unsigned int first_vector_component);
- };
-
-
- /**
- * Extractor for a symmetric tensor of a
- * rank specified by the template
- * argument. For a second order symmetric
- * tensor, this represents a collection of
- * (dim*dim + dim)/2
- * components of a vector-valued
- * element. The value of dim
- * is defined by the FEValues object the
- * extractor is applied to. The result of
- * applying an object of this type to an
- * FEValues, FEFaceValues or
- * FESubfaceValues object is of type
- * FEValuesViews::SymmetricTensor.
- *
- * The concept of
- * extractors is defined in the
- * documentation of the namespace
- * FEValuesExtractors and in the @ref
- * vector_valued module.
- *
- * @ingroup feaccess vector_valued
- *
- * @author Andrew McBride, 2009
- */
- template i * n_components + c
.
+ * - The value stored at this
+ * position indicates the row
+ * in #shape_values and the
+ * other tables where the
+ * corresponding datum is stored
+ * for all the quadrature points.
+ *
+ * In the general, vector-valued
+ * context, the number of
+ * components is larger than one,
+ * but for a given shape
+ * function, not all vector
+ * components may be nonzero
+ * (e.g., if a shape function is
+ * primitive, then exactly one
+ * vector component is non-zero,
+ * while the others are all
+ * zero). For such zero
+ * components, #shape_values and
+ * friends do not have a
+ * row. Consequently, for vector
+ * components for which shape
+ * function i is zero, the entry
+ * in the current table is
+ * numbers::invalid_unsigned_int.
+ *
+ * On the other hand, the table
+ * is guaranteed to have at least
+ * one valid index for each shape
+ * function. In particular, for a
+ * primitive finite element, each
+ * shape function has exactly one
+ * nonzero component and so for
+ * each i, there is exactly one
+ * valid index within the range
+ * [i*n_components,
+ * (i+1)*n_components)
.
*/
std::vectordim
elements, and second order symmetric tensors consisting of
+ * (dim*dim + dim)/2
components
+ *
+ * See the description of the @ref vector_valued module for examples how to
+ * use the features of this namespace.
+ *
+ * @ingroup feaccess vector_valued
+ */
+namespace FEValuesExtractors
+{
+ /**
+ * Extractor for a single scalar component
+ * of a vector-valued element. The result
+ * of applying an object of this type to an
+ * FEValues, FEFaceValues or
+ * FESubfaceValues object is of type
+ * FEValuesViews::Scalar. The concept of
+ * extractors is defined in the
+ * documentation of the namespace
+ * FEValuesExtractors and in the @ref
+ * vector_valued module.
+ *
+ * @ingroup feaccess vector_valued
+ */
+ struct Scalar
+ {
+ /**
+ * The selected scalar component of the
+ * vector.
+ */
+ unsigned int component;
+
+ /**
+ * Constructor. Take the selected
+ * vector component as argument.
+ */
+ Scalar (const unsigned int component);
+ };
+
+
+ /**
+ * Extractor for a vector of
+ * spacedim
components of a
+ * vector-valued element. The value of
+ * spacedim
is defined by the
+ * FEValues object the extractor is applied
+ * to. The result of applying an object of
+ * this type to an FEValues, FEFaceValues
+ * or FESubfaceValues object is of type
+ * FEValuesViews::Vector.
+ *
+ * The concept of
+ * extractors is defined in the
+ * documentation of the namespace
+ * FEValuesExtractors and in the @ref
+ * vector_valued module.
+ *
+ * Note that in the current context, a
+ * vector is meant in the sense physics
+ * uses it: it has spacedim
+ * components that behave in specific ways
+ * under coordinate system
+ * transformations. Examples include
+ * velocity or displacement fields. This is
+ * opposed to how mathematics uses the word
+ * "vector" (and how we use this word in
+ * other contexts in the library, for
+ * example in the Vector class), where it
+ * really stands for a collection of
+ * numbers. An example of this latter use
+ * of the word could be the set of
+ * concentrations of chemical species in a
+ * flame; however, these are really just a
+ * collection of scalar variables, since
+ * they do not change if the coordinate
+ * system is rotated, unlike the components
+ * of a velocity vector, and consequently,
+ * this class should not be used for this
+ * context.
+ *
+ * @ingroup feaccess vector_valued
+ */
+ struct Vector
+ {
+ /**
+ * The first component of the vector
+ * view.
+ */
+ unsigned int first_vector_component;
+
+ /**
+ * Constructor. Take the first
+ * component of the selected vector
+ * inside the FEValues object as
+ * argument.
+ */
+ Vector (const unsigned int first_vector_component);
+ };
+
+
+ /**
+ * Extractor for a symmetric tensor of a
+ * rank specified by the template
+ * argument. For a second order symmetric
+ * tensor, this represents a collection of
+ * (dim*dim + dim)/2
+ * components of a vector-valued
+ * element. The value of dim
+ * is defined by the FEValues object the
+ * extractor is applied to. The result of
+ * applying an object of this type to an
+ * FEValues, FEFaceValues or
+ * FESubfaceValues object is of type
+ * FEValuesViews::SymmetricTensor.
+ *
+ * The concept of
+ * extractors is defined in the
+ * documentation of the namespace
+ * FEValuesExtractors and in the @ref
+ * vector_valued module.
+ *
+ * @ingroup feaccess vector_valued
+ *
+ * @author Andrew McBride, 2009
+ */
+ template dim
components are true that correspond to the given
+ * argument.
+ *
+ * @note This function is the equivalent of
+ * FiniteElement::component_mask() with the same arguments.
+ * It verifies that it gets the same result from every one
+ * of the elements that are stored in this FECollection. If
+ * this is not the case, it throws an exception.
+ *
+ * @param vector An object that represents dim
+ * vector components of this finite element.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ ComponentMask
+ component_mask (const FEValuesExtractors::Vector &vector) const;
+
+ /**
+ * Return a component mask with as many elements as this
+ * object has vector components and of which exactly the
+ * dim*(dim+1)/2
components are true that
+ * correspond to the given argument.
+ *
+ * @note This function is the equivalent of
+ * FiniteElement::component_mask() with the same arguments.
+ * It verifies that it gets the same result from every one
+ * of the elements that are stored in this FECollection. If
+ * this is not the case, it throws an exception.
+ *
+ * @param sym_tensor An object that represents dim*(dim+1)/2
+ * components of this finite element that are jointly to be
+ * interpreted as forming a symmetric tensor.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ ComponentMask
+ component_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
+
+ /**
+ * Given a block mask (see @ref GlossBlockMask "this glossary entry"),
+ * produce a component mask (see @ref GlossComponentMask "this glossary entry")
+ * that represents the components that correspond to the blocks selected in
+ * the input argument. This is essentially a conversion operator from
+ * BlockMask to ComponentMask.
+ *
+ * @note This function is the equivalent of
+ * FiniteElement::component_mask() with the same arguments.
+ * It verifies that it gets the same result from every one
+ * of the elements that are stored in this FECollection. If
+ * this is not the case, it throws an exception.
+ *
+ * @param block_mask The mask that selects individual blocks of the finite
+ * element
+ * @return A mask that selects those components corresponding to the selected
+ * blocks of the input argument.
+ */
+ ComponentMask
+ component_mask (const BlockMask &block_mask) const;
+
+ /**
+ * Return a block mask with as many elements as this
+ * object has blocks and of which exactly the
+ * one component is true that corresponds to the given
+ * argument. See @ref GlossBlockMask "the glossary"
+ * for more information.
+ *
+ * @note This function will only succeed if the scalar referenced
+ * by the argument encompasses a complete block. In other words,
+ * if, for example, you pass an extractor for the single
+ * $x$ velocity and this object represents an FE_RaviartThomas
+ * object, then the single scalar object you selected is part
+ * of a larger block and consequently there is no block mask that
+ * would represent it. The function will then produce an exception.
+ *
+ * @note This function is the equivalent of
+ * FiniteElement::component_mask() with the same arguments.
+ * It verifies that it gets the same result from every one
+ * of the elements that are stored in this FECollection. If
+ * this is not the case, it throws an exception.
+ *
+ * @param scalar An object that represents a single scalar
+ * vector component of this finite element.
+ * @return A component mask that is false in all components
+ * except for the one that corresponds to the argument.
+ */
+ BlockMask
+ block_mask (const FEValuesExtractors::Scalar &scalar) const;
+
+ /**
+ * Return a component mask with as many elements as this
+ * object has vector components and of which exactly the
+ * dim
components are true that correspond to the given
+ * argument. See @ref GlossBlockMask "the glossary"
+ * for more information.
+ *
+ * @note This function is the equivalent of
+ * FiniteElement::component_mask() with the same arguments.
+ * It verifies that it gets the same result from every one
+ * of the elements that are stored in this FECollection. If
+ * this is not the case, it throws an exception.
+ *
+ * @note The same caveat applies as to the version of the function above:
+ * The extractor object passed as argument must be so that it corresponds
+ * to full blocks and does not split blocks of this element.
+ *
+ * @param vector An object that represents dim
+ * vector components of this finite element.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ BlockMask
+ block_mask (const FEValuesExtractors::Vector &vector) const;
+
+ /**
+ * Return a component mask with as many elements as this
+ * object has vector components and of which exactly the
+ * dim*(dim+1)/2
components are true that
+ * correspond to the given argument. See @ref GlossBlockMask "the glossary"
+ * for more information.
+ *
+ * @note The same caveat applies as to the version of the function above:
+ * The extractor object passed as argument must be so that it corresponds
+ * to full blocks and does not split blocks of this element.
+ *
+ * @note This function is the equivalent of
+ * FiniteElement::component_mask() with the same arguments.
+ * It verifies that it gets the same result from every one
+ * of the elements that are stored in this FECollection. If
+ * this is not the case, it throws an exception.
+ *
+ * @param sym_tensor An object that represents dim*(dim+1)/2
+ * components of this finite element that are jointly to be
+ * interpreted as forming a symmetric tensor.
+ * @return A component mask that is false in all components
+ * except for the ones that corresponds to the argument.
+ */
+ BlockMask
+ block_mask (const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const;
+
+ /**
+ * Given a component mask (see @ref GlossComponentMask "this glossary entry"),
+ * produce a block mask (see @ref GlossBlockMask "this glossary entry")
+ * that represents the blocks that correspond to the components selected in
+ * the input argument. This is essentially a conversion operator from
+ * ComponentMask to BlockMask.
+ *
+ * @note This function will only succeed if the components referenced
+ * by the argument encompasses complete blocks. In other words,
+ * if, for example, you pass an component mask for the single
+ * $x$ velocity and this object represents an FE_RaviartThomas
+ * object, then the single component you selected is part
+ * of a larger block and consequently there is no block mask that
+ * would represent it. The function will then produce an exception.
+ *
+ * @note This function is the equivalent of
+ * FiniteElement::component_mask() with the same arguments.
+ * It verifies that it gets the same result from every one
+ * of the elements that are stored in this FECollection. If
+ * this is not the case, it throws an exception.
+ *
+ * @param component_mask The mask that selects individual components of the finite
+ * element
+ * @return A mask that selects those blocks corresponding to the selected
+ * blocks of the input argument.
+ */
+ BlockMask
+ block_mask (const ComponentMask &component_mask) const;
+
/**
* Exception
diff --git a/deal.II/include/deal.II/lac/constraint_matrix.templates.h b/deal.II/include/deal.II/lac/constraint_matrix.templates.h
index 970ee190aa..715527800e 100644
--- a/deal.II/include/deal.II/lac/constraint_matrix.templates.h
+++ b/deal.II/include/deal.II/lac/constraint_matrix.templates.h
@@ -1737,7 +1737,7 @@ namespace internals
// similar as before, now with shortcut for
- // deal.II sparse matrices. this lets use
+ // deal.II sparse matrices. this lets us
// avoid using extra arrays, and does all the
// operations just in place, i.e., in the
// respective matrix row
diff --git a/deal.II/include/deal.II/multigrid/mg_constrained_dofs.h b/deal.II/include/deal.II/multigrid/mg_constrained_dofs.h
index a5f929451c..e011893a8e 100644
--- a/deal.II/include/deal.II/multigrid/mg_constrained_dofs.h
+++ b/deal.II/include/deal.II/multigrid/mg_constrained_dofs.h
@@ -59,7 +59,7 @@ class MGConstrainedDoFs : public Subscriptor
template