* under consideration here is hanging nodes constraints and grid
* transfer, respectively.
*
+ * <h3>Components and blocks</h3>
+ *
+ * For vector valued elements shape functions may have nonzero entries
+ * in one or several @ref GlossComponent "components" of the vector
+ * valued function. If the element is @ref GlossPrimitive "primitive",
+ * there is indeed a single component with a nonzero entry for each
+ * shape function. This component can be determined by
+ * system_to_component_index(), the number of components is
+ * FiniteElementData::n_components().
+ *
+ * Furthermore, you may want to split your linear system into @ref
+ * GlossBlock "blocks" for the use in BlockVector, BlockSparseMatrix,
+ * BlockMatrixArray and so on. If you use non-primitive elements, you
+ * cannot determine the block number by
+ * system_to_component_index(). Instead, you can use
+ * system_to_block_index(), which will automatically take care of the
+ * additional components occupied by vector valued elements. The
+ * number of generated blocks can be determined by
+ * FiniteElementData::n_blocks().
+ *
+ * If you decide to operate by base element and multiplicity, the
+ * function first_block_of_base() will be helpful.
*
* <h3>Support points</h3>
*
std::pair<unsigned int, unsigned int>
face_system_to_component_index (const unsigned int index) const;
- /**
- * Return for shape function
- * @p index the base element it
- * belongs to, the number of the
- * copy of this base element
- * (which is between zero and the
- * multiplicity of this element),
- * and the index of this shape
- * function within this base
- * element.
- *
- * If the element is not composed of
- * others, then base and instance
- * are always zero, and the index
- * is equal to the number of the
- * shape function. If the element
- * is composed of single
- * instances of other elements
- * (i.e. all with multiplicity
- * one) all of which are scalar,
- * then base values and dof
- * indices within this element
- * are equal to the
- * @p system_to_component_table. It
- * differs only in case the
- * element is composed of other
- * elements and at least one of
- * them is vector-valued itself.
- *
- * This function returns valid
- * values also in the case of
- * vector-valued
- * (i.e. non-primitive) shape
- * functions, in contrast to the
- * @p system_to_component_index
- * function.
- */
- std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
- system_to_base_index (const unsigned int index) const;
-
- /**
- * Same as
- * @p system_to_base_index, but
- * for degrees of freedom located
- * on a face. The range of allowed
- * indices is therefore
- * 0..dofs_per_face.
- *
- * You will rarely need this
- * function in application
- * programs, since almost all
- * application codes only need to
- * deal with cell indices, not
- * face indices. The function is
- * mainly there for use inside
- * the library.
- */
- std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
- face_system_to_base_index (const unsigned int index) const;
-
/**
* Return in which of the vector
* components of this finite
unsigned int
element_multiplicity (const unsigned int index) const = 0;
+ /**
+ * Return for shape function
+ * @p index the base element it
+ * belongs to, the number of the
+ * copy of this base element
+ * (which is between zero and the
+ * multiplicity of this element),
+ * and the index of this shape
+ * function within this base
+ * element.
+ *
+ * If the element is not composed of
+ * others, then base and instance
+ * are always zero, and the index
+ * is equal to the number of the
+ * shape function. If the element
+ * is composed of single
+ * instances of other elements
+ * (i.e. all with multiplicity
+ * one) all of which are scalar,
+ * then base values and dof
+ * indices within this element
+ * are equal to the
+ * @p system_to_component_table. It
+ * differs only in case the
+ * element is composed of other
+ * elements and at least one of
+ * them is vector-valued itself.
+ *
+ * This function returns valid
+ * values also in the case of
+ * vector-valued
+ * (i.e. non-primitive) shape
+ * functions, in contrast to the
+ * @p system_to_component_index
+ * function.
+ */
+ std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
+ system_to_base_index (const unsigned int index) const;
+
+ /**
+ * Same as
+ * @p system_to_base_index, but
+ * for degrees of freedom located
+ * on a face. The range of allowed
+ * indices is therefore
+ * 0..dofs_per_face.
+ *
+ * You will rarely need this
+ * function in application
+ * programs, since almost all
+ * application codes only need to
+ * deal with cell indices, not
+ * face indices. The function is
+ * mainly there for use inside
+ * the library.
+ */
+ std::pair<std::pair<unsigned int, unsigned int>, unsigned int>
+ face_system_to_base_index (const unsigned int index) const;
+
+ /**
+ * Given a base element number,
+ * return the first block of a
+ * BlockVector it would generate.
+ */
+ unsigned int first_block_of_base(unsigned int b) const;
+
/**
* Given a vector component,
* return an index which base
*/
std::pair<unsigned int,unsigned int>
component_to_base (const unsigned int component) const;
+
+ /**
+ * The vector block and the index
+ * inside the block for this
+ * shape function.
+ */
+ std::pair<unsigned int,unsigned int>
+ system_to_block_index (const unsigned int component) const;
+
//@}
/**
TableIndices<2>
interface_constraints_size () const;
+ /**
+ * Store whether all shape
+ * functions are primitive. Since
+ * finding this out is a very
+ * common operation, we cache the
+ * result, i.e. compute the value
+ * in the constructor for simpler
+ * access.
+ */
+ const bool cached_primitivity;
+
+ /**
+ * Compute second derivatives by
+ * finite differences of
+ * gradients.
+ */
+ void compute_2nd (const Mapping<dim> &mapping,
+ const typename Triangulation<dim>::cell_iterator &cell,
+ const unsigned int offset,
+ typename Mapping<dim>::InternalDataBase &mapping_internal,
+ InternalDataBase &fe_internal,
+ FEValuesData<dim> &data) const;
+
+ /**
+ * Given the pattern of nonzero
+ * components for each shape
+ * function, compute for each
+ * entry how many components are
+ * non-zero for each shape
+ * function. This function is
+ * used in the constructor of
+ * this class.
+ */
+ static
+ std::vector<unsigned int>
+ compute_n_nonzero_components (const std::vector<std::vector<bool> > &nonzero_components);
+
+ /**
+ * Exception
+ *
+ * @ingroup Exceptions
+ */
+ DeclException0 (ExcBoundaryFaceUsed);
+ /**
+ * Exception
+ *
+ * @ingroup Exceptions
+ */
+ DeclException0 (ExcJacobiDeterminantHasWrongSign);
+
+ /**
+ * Determine the values a finite
+ * element should compute on
+ * initialization of data for
+ * @p FEValues.
+ *
+ * Given a set of flags
+ * indicating what quantities are
+ * requested from a @p FEValues
+ * object, @p update_once and
+ * @p update_each compute which
+ * values must really be
+ * computed. Then, the
+ * <tt>fill_*_values</tt> functions
+ * are called with the result of
+ * these.
+ *
+ * Furthermore, values must be
+ * computed either on the unit
+ * cell or on the physical
+ * cell. For instance, the
+ * function values of @p FE_Q do
+ * only depend on the quadrature
+ * points on the unit
+ * cell. Therefore, this flags
+ * will be returned by
+ * @p update_once. The gradients
+ * require computation of the
+ * covariant transformation
+ * matrix. Therefore,
+ * @p update_covariant_transformation
+ * and @p update_gradients will
+ * be returned by
+ * @p update_each.
+ *
+ * For an example see the same
+ * function in the derived class
+ * @p FE_Q.
+ */
+ virtual UpdateFlags update_once (const UpdateFlags flags) const = 0;
+
+ /**
+ * Complementary function for
+ * @p update_once.
+ *
+ * While @p update_once returns
+ * the values to be computed on
+ * the unit cell for yielding the
+ * required data, this function
+ * determines the values that
+ * must be recomputed on each
+ * cell.
+ *
+ * Refer to @p update_once for
+ * more details.
+ */
+ virtual UpdateFlags update_each (const UpdateFlags flags) const = 0;
+
+ /**
+ * @p clone function instead of
+ * a copy constructor.
+ *
+ * This function is needed by the
+ * constructors of FESystem as well
+ * as by the hp::FECollection class.
+ */
+ virtual FiniteElement<dim> *clone() const = 0;
+
+ /**
+ * List of support points on the
+ * unit cell, in case the finite
+ * element has any. The
+ * constructor leaves this field
+ * empty, derived classes may
+ * write in some contents.
+ *
+ * Finite elements that allow
+ * some kind of interpolation
+ * operation usually have support
+ * points. On the other hand,
+ * elements that define their
+ * degrees of freedom by, for
+ * example, moments on faces, or
+ * as derivatives, don't have
+ * support points. In that case,
+ * this field remains empty.
+ */
+ std::vector<Point<dim> > unit_support_points;
+
+ /**
+ * Same for the faces. See the
+ * description of the
+ * @p get_unit_face_support_points
+ * function for a discussion of
+ * what contributes a face
+ * support point.
+ */
+ std::vector<Point<dim-1> > unit_face_support_points;
+
+ /**
+ * Support points used for
+ * interpolation functions of
+ * non-Lagrangian elements.
+ */
+ std::vector<Point<dim> > generalized_support_points;
+
+ /**
+ * Face support points used for
+ * interpolation functions of
+ * non-Lagrangian elements.
+ */
+ std::vector<Point<dim-1> > generalized_face_support_points;
+
+ private:
/**
* Store what
* @p system_to_component_index
*/
std::vector<std::pair<std::pair<unsigned int,unsigned int>,unsigned int> >
face_system_to_base_table;
+ /**
+ * For each base element, store
+ * the first block in a block
+ * vector it will generate.
+ */
+ std::vector<unsigned int> first_block_of_base_table;
/**
* The base element establishing
*/
const std::vector<bool> restriction_is_additive_flags;
- /**
- * List of support points on the
- * unit cell, in case the finite
- * element has any. The
- * constructor leaves this field
- * empty, derived classes may
- * write in some contents.
- *
- * Finite elements that allow
- * some kind of interpolation
- * operation usually have support
- * points. On the other hand,
- * elements that define their
- * degrees of freedom by, for
- * example, moments on faces, or
- * as derivatives, don't have
- * support points. In that case,
- * this field remains empty.
- */
- std::vector<Point<dim> > unit_support_points;
-
- /**
- * Same for the faces. See the
- * description of the
- * @p get_unit_face_support_points
- * function for a discussion of
- * what contributes a face
- * support point.
- */
- std::vector<Point<dim-1> > unit_face_support_points;
-
- /**
- * Support points used for
- * interpolation functions of
- * non-Lagrangian elements.
- */
- std::vector<Point<dim> > generalized_support_points;
-
- /**
- * Face support points used for
- * interpolation functions of
- * non-Lagrangian elements.
- */
- std::vector<Point<dim-1> > generalized_face_support_points;
-
/**
* For each shape function, give
* a vector of bools (with size
*/
const std::vector<unsigned int> n_nonzero_components_table;
- /**
- * Store whether all shape
- * functions are primitive. Since
- * finding this out is a very
- * common operation, we cache the
- * result, i.e. compute the value
- * in the constructor for simpler
- * access.
- */
- const bool cached_primitivity;
-
- /**
- * Compute second derivatives by
- * finite differences of
- * gradients.
- */
- void compute_2nd (const Mapping<dim> &mapping,
- const typename Triangulation<dim>::cell_iterator &cell,
- const unsigned int offset,
- typename Mapping<dim>::InternalDataBase &mapping_internal,
- InternalDataBase &fe_internal,
- FEValuesData<dim> &data) const;
-
- /**
- * Given the pattern of nonzero
- * components for each shape
- * function, compute for each
- * entry how many components are
- * non-zero for each shape
- * function. This function is
- * used in the constructor of
- * this class.
- */
- static
- std::vector<unsigned int>
- compute_n_nonzero_components (const std::vector<std::vector<bool> > &nonzero_components);
-
-
- /**
- * Exception
- *
- * @ingroup Exceptions
- */
- DeclException0 (ExcBoundaryFaceUsed);
- /**
- * Exception
- *
- * @ingroup Exceptions
- */
- DeclException0 (ExcJacobiDeterminantHasWrongSign);
-
- protected:
-
- /**
- * Determine the values a finite
- * element should compute on
- * initialization of data for
- * @p FEValues.
- *
- * Given a set of flags
- * indicating what quantities are
- * requested from a @p FEValues
- * object, @p update_once and
- * @p update_each compute which
- * values must really be
- * computed. Then, the
- * <tt>fill_*_values</tt> functions
- * are called with the result of
- * these.
- *
- * Furthermore, values must be
- * computed either on the unit
- * cell or on the physical
- * cell. For instance, the
- * function values of @p FE_Q do
- * only depend on the quadrature
- * points on the unit
- * cell. Therefore, this flags
- * will be returned by
- * @p update_once. The gradients
- * require computation of the
- * covariant transformation
- * matrix. Therefore,
- * @p update_covariant_transformation
- * and @p update_gradients will
- * be returned by
- * @p update_each.
- *
- * For an example see the same
- * function in the derived class
- * @p FE_Q.
- */
- virtual UpdateFlags update_once (const UpdateFlags flags) const = 0;
-
- /**
- * Complementary function for
- * @p update_once.
- *
- * While @p update_once returns
- * the values to be computed on
- * the unit cell for yielding the
- * required data, this function
- * determines the values that
- * must be recomputed on each
- * cell.
- *
- * Refer to @p update_once for
- * more details.
- */
- virtual UpdateFlags update_each (const UpdateFlags flags) const = 0;
-
- /**
- * @p clone function instead of
- * a copy constructor.
- *
- * This function is needed by the
- * constructors of FESystem as well
- * as by the hp::FECollection class.
- */
- virtual FiniteElement<dim> *clone() const = 0;
-
- private:
/**
* Second derivatives of shapes
* functions are not computed
typename Mapping<dim>::InternalDataBase &fe_internal,
FEValuesData<dim> &data) const = 0;
- /**
- * Allow the FESystem class to access the
- * restriction and prolongation matrices
- * directly. Hence, FESystem has the
- * possibility to see if these matrices
- * are initialized without accessing
- * these matrices through the
- * @p get_restriction_matrix and
- * @p get_prolongation_matrix
- * functions. This is important as these
- * functions include assertions that
- * throw if the matrices are not already
- * initialized.
- */
- template <int dim_> friend class FESystem;
-
- /**
- * Make the inner class a
- * friend. This is not strictly
- * necessary, but the Intel
- * compiler seems to want this.
- */
friend class InternalDataBase;
-
- /**
- * Declare some other classes as
- * friends of this class.
- */
friend class FEValuesBase<dim>;
friend class FEValues<dim>;
friend class FEFaceValues<dim>;
friend class FESubfaceValues<dim>;
+ template <int dim_> friend class FESystem;
friend class hp::FECollection<dim>;
};
+template <int dim>
+inline
+unsigned int
+FiniteElement<dim>::first_block_of_base (const unsigned int index) const
+{
+ Assert(index < first_block_of_base_table.size(),
+ ExcIndexRange(index, 0, first_block_of_base_table.size()));
+
+ return first_block_of_base_table[index];
+}
+
+
template <int dim>
inline
std::pair<unsigned int,unsigned int>
}
+template <int dim>
+inline
+std::pair<unsigned int,unsigned int>
+FiniteElement<dim>::system_to_block_index (const unsigned int index) const
+{
+ Assert (index < this->n_blocks(),
+ ExcIndexRange(index, 0, this->n_blocks()));
+ // The block is computed simply as
+ // first block of this base plus
+ // the index within the base blocks
+ return std::pair<unsigned int, unsigned int>(
+ first_block_of_base(system_to_base_table[index].first.first)
+ + system_to_base_table[index].first.second,
+ system_to_base_table[index].second);
+}
+
+
template <int dim>
inline
bool
// $Id$
// Version: $Name$
//
-// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
template<int dim> class FESystem;
-/*!@addtogroup febase */
-/*@{*/
-
/**
* Dimension independent data for finite elements. See the derived
* class FiniteElement class for information on its use. All
* version, FiniteElementData and FiniteElement should be private
* base classes of FiniteElement.
*
+ * @ingroup febase
* @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2001, 2003, 2005
*/
template <int dim>
* base element.
*/
const unsigned int components;
-
+ /**
+ * The number of vector blocks a
+ * BlockVector for this element
+ * should contain. For primitive
+ * elements, this is equal to
+ * #components, but for vector
+ * valued base elements it may
+ * differ. Actually, in systems
+ * this is the sum of the base
+ * element multiplicities.
+ */
+ const unsigned int blocks;
+
/**
* Maximal polynomial degree of a
* shape function in a single
* @param conformity The finite
* element space has continuity
* of this Sobolev space.
+ * @param n_blocks Number of
+ * vector blocks.
*/
FiniteElementData (const std::vector<unsigned int> &dofs_per_object,
const unsigned int n_components,
- const unsigned int degree = deal_II_numbers::invalid_unsigned_int,
- const Conformity conformity = unknown);
+ const unsigned int degree,
+ const Conformity conformity = unknown,
+ const unsigned int n_blocks = deal_II_numbers::invalid_unsigned_int);
/**
* Number of dofs per vertex.
*/
unsigned int n_components () const;
+ /**
+ * Number of blocks.
+ */
+ unsigned int n_blocks () const;
+
/**
* Maximal polynomial degree of a
* shape function in a single
};
-
-
-/**
- *
- * @author Wolfgang Bangerth, 1998, 2002; Ralf Hartmann, Guido Kanschat, 2001
- */
-/*@}*/
-
-
// --------- inline and template functions ---------------
template <int dim>
+template <int dim>
+inline
+unsigned int
+FiniteElementData<dim>::n_blocks () const
+{
+ return blocks;
+}
+
+
+
template <int dim>
inline
unsigned int
// $Id$
// Version: $Name$
//
-// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
for (cell = dofs.begin_active(); cell != end; ++cell)
{
const FiniteElement<DH::dimension>& fe = cell->get_fe();
- Assert (fe.is_primitive(), typename FiniteElement<DH::dimension>::ExcFENotPrimitive());
+ Assert (fe.is_primitive(),
+ typename FiniteElement<DH::dimension>::ExcFENotPrimitive());
Assert (couplings.n_rows()==fe.n_components(),
ExcDimensionMismatch(couplings.n_rows(), fe.n_components()));
Assert (couplings.n_cols()==fe.n_components(),
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
const std::vector<std::vector<bool> > &nonzero_c)
:
FiniteElementData<dim> (fe_data),
- system_to_component_table (this->dofs_per_cell),
- face_system_to_component_table(this->dofs_per_face),
+ cached_primitivity(false),
system_to_base_table(this->dofs_per_cell),
face_system_to_base_table(this->dofs_per_face),
component_to_base_table (this->components,
std::make_pair(0U, 0U)),
restriction_is_additive_flags(r_i_a_f),
- nonzero_components (nonzero_c),
- cached_primitivity(false)
+ nonzero_components (nonzero_c)
{
// Special handling of vectors of
// length one: in this case, we
// initialize some tables in the
// default way, i.e. if there is
// only one (vector-)component; if
- // there are several, then the
- // constructor of the derived class
- // needs to overwrite these arrays
- for (unsigned int j=0 ; j<this->dofs_per_cell ; ++j)
+ // the element is not primitive,
+ // leave these tables empty.
+ if (!cached_primitivity)
{
- system_to_component_table[j] = std::pair<unsigned,unsigned>(0,j);
- system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j);
- }
- for (unsigned int j=0 ; j<this->dofs_per_face ; ++j)
- {
- face_system_to_component_table[j] = std::pair<unsigned,unsigned>(0,j);
- face_system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j);
+ system_to_component_table.resize(this->dofs_per_cell);
+ face_system_to_component_table.resize(this->dofs_per_face);
+ for (unsigned int j=0 ; j<this->dofs_per_cell ; ++j)
+ {
+ system_to_component_table[j] = std::pair<unsigned,unsigned>(0,j);
+ system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j);
+ }
+ for (unsigned int j=0 ; j<this->dofs_per_face ; ++j)
+ {
+ face_system_to_component_table[j] = std::pair<unsigned,unsigned>(0,j);
+ face_system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j);
+ }
}
}
// $Id$
// Version: $Name$
//
-// Copyright (C) 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
dofs_per_face(0),
dofs_per_cell (0),
components(0),
+ blocks(0),
degree(0),
conforming_space(unknown)
{}
FiniteElementData (const std::vector<unsigned int> &dofs_per_object,
const unsigned int n_components,
const unsigned int degree,
- const Conformity conformity)
+ const Conformity conformity,
+ const unsigned int n_blocks)
:
dofs_per_vertex(dofs_per_object[0]),
dofs_per_line(dofs_per_object[1]),
GeometryInfo<dim>::quads_per_cell * dofs_per_quad +
GeometryInfo<dim>::hexes_per_cell * dofs_per_hex),
components(n_components),
+ blocks(n_blocks == deal_II_numbers::invalid_unsigned_int
+ ? n_components : n_blocks),
degree(degree),
conforming_space(conformity)
{
// $Id$
// Version: $Name$
//
-// Copyright (C) 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
template <int dim>
FE_DGPNonparametric<dim>::FE_DGPNonparametric (const unsigned int degree)
:
- FiniteElement<dim> (FiniteElementData<dim>(get_dpo_vector(degree),1, FiniteElementData<dim>::L2),
- std::vector<bool>(FiniteElementData<dim>(get_dpo_vector(degree),1).dofs_per_cell,true),
- std::vector<std::vector<bool> >(FiniteElementData<dim>(get_dpo_vector(degree),1).dofs_per_cell,
- std::vector<bool>(1,true))),
+ FiniteElement<dim> (
+ FiniteElementData<dim>(get_dpo_vector(degree), 1, degree,
+ FiniteElementData<dim>::L2),
+ std::vector<bool>(
+ FiniteElementData<dim>(get_dpo_vector(degree), 1, degree).dofs_per_cell,true),
+ std::vector<std::vector<bool> >(
+ FiniteElementData<dim>(get_dpo_vector(degree),1, degree).dofs_per_cell,
+ std::vector<bool>(1,true))),
degree(degree),
polynomial_space (Polynomials::Legendre::generate_complete_basis(degree))
{
// $Id$
// Version: $Name$
//
-// Copyright (C) 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
template <int dim>
FE_Nedelec<dim>::FE_Nedelec (const unsigned int degree)
:
- FiniteElement<dim> (FiniteElementData<dim>(get_dpo_vector(degree), dim, degree+1, FiniteElementData<dim>::Hcurl),
- std::vector<bool> (FiniteElementData<dim>(get_dpo_vector(degree),dim).dofs_per_cell,false),
- std::vector<std::vector<bool> >(FiniteElementData<dim>(get_dpo_vector(degree),dim).dofs_per_cell,
- std::vector<bool>(dim,true))),
+ FiniteElement<dim> (
+ FiniteElementData<dim>(get_dpo_vector(degree), dim,
+ degree+1, FiniteElementData<dim>::Hcurl, 1),
+ std::vector<bool> (
+ FiniteElementData<dim>(get_dpo_vector(degree), dim,
+ degree+1).dofs_per_cell,false),
+ std::vector<std::vector<bool> >(
+ FiniteElementData<dim>(get_dpo_vector(degree), dim,
+ degree+1).dofs_per_cell,
+ std::vector<bool>(dim,true))),
degree(degree)
{
Assert (dim >= 2, ExcImpossibleInDim(dim));
// on cell and face
initialize_unit_support_points ();
initialize_unit_face_support_points ();
-
- // then make
- // system_to_component_table
- // invalid, since this has no
- // meaning for the present element
- std::vector<std::pair<unsigned,unsigned> > tmp1, tmp2;
- this->system_to_component_table.swap (tmp1);
- this->face_system_to_component_table.swap (tmp2);
}
// $Id$
// Version: $Name$
//
-// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
std::vector<unsigned int>
FE_Q<dim>::face_lexicographic_to_hierarchic_numbering (const unsigned int degree)
{
- const FiniteElementData<dim-1> face_data(FE_Q<dim-1>::get_dpo_vector(degree),1);
+ const FiniteElementData<dim-1> face_data(FE_Q<dim-1>::get_dpo_vector(degree),1,degree);
std::vector<unsigned int> face_renumber (face_data.dofs_per_cell);
FETools::lexicographic_to_hierarchic_numbering (face_data, face_renumber);
return face_renumber;
// $Id$
// Version: $Name$
//
-// Copyright (C) 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
FE_Q_Hierarchical<dim>::
face_fe_q_hierarchical_to_hierarchic_numbering (const unsigned int degree)
{
- FiniteElementData<dim-1> fe_data(FE_Q_Hierarchical<dim-1>::get_dpo_vector(degree),1);
+ FiniteElementData<dim-1> fe_data(FE_Q_Hierarchical<dim-1>::get_dpo_vector(degree),1,degree);
return invert_numbering(FE_Q_Hierarchical<dim-1>::
hierarchic_to_fe_q_hierarchical_numbering (fe_data));
}
FE_PolyTensor<PolynomialsRaviartThomas<dim>, dim> (
deg,
FiniteElementData<dim>(get_dpo_vector(deg),
- dim, deg+1, FiniteElementData<dim>::Hdiv),
+ dim, deg+1, FiniteElementData<dim>::Hdiv, 1),
std::vector<bool>(PolynomialsRaviartThomas<dim>::compute_n_pols(deg), true),
std::vector<std::vector<bool> >(PolynomialsRaviartThomas<dim>::compute_n_pols(deg),
std::vector<bool>(dim,true))),
this->interface_constraints(target_row,j) = face_embeddings[d](i,j);
++target_row;
}
-//TODO:[WB] What is this?
- // then make
- // system_to_component_table
- // invalid, since this has no
- // meaning for the present element
- std::vector<std::pair<unsigned,unsigned> > tmp1, tmp2;
- this->system_to_component_table.swap (tmp1);
- this->face_system_to_component_table.swap (tmp2);
}
// $Id$
// Version: $Name$
//
-// Copyright (C) 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
FE_PolyTensor<PolynomialsRaviartThomas<dim>, dim> (
deg,
FiniteElementData<dim>(get_dpo_vector(deg),
- dim, deg+1, FiniteElementData<dim>::Hdiv),
+ dim, deg+1, FiniteElementData<dim>::Hdiv, 1),
get_ria_vector (deg),
- std::vector<std::vector<bool> >(
- FiniteElementData<dim>(get_dpo_vector(deg),
- dim,deg+1).dofs_per_cell,
- std::vector<bool>(dim,true)))
+ std::vector<std::vector<bool> >(PolynomialsRaviartThomas<dim>::compute_n_pols(deg),
+ std::vector<bool>(dim,true)))
{
Assert (dim >= 2, ExcImpossibleInDim(dim));
const unsigned int n_dofs = this->dofs_per_cell;
std::vector<bool>
FE_RaviartThomasNodal<dim>::get_ria_vector (const unsigned int deg)
{
- unsigned int dofs_per_cell, dofs_per_face;
- switch (dim)
- {
- case 2:
- dofs_per_face = deg+1;
- dofs_per_cell = 2*(deg+1)*(deg+2);
- break;
- case 3:
- dofs_per_face = (deg+1)*(deg+1);
- dofs_per_cell = 3*(deg+1)*(deg+1)*(deg+2);
- break;
- default:
- Assert (false, ExcNotImplemented());
- }
- Assert (FiniteElementData<dim>(get_dpo_vector(deg),dim).dofs_per_cell ==
- dofs_per_cell,
- ExcInternalError());
- Assert (FiniteElementData<dim>(get_dpo_vector(deg),dim).dofs_per_face ==
- dofs_per_face,
- ExcInternalError());
-
+ const unsigned int dofs_per_cell = PolynomialsRaviartThomas<dim>::compute_n_pols(deg);
+ unsigned int dofs_per_face = deg+1;
+ for (unsigned int d=2;d<dim;++d)
+ dofs_per_face *= deg+1;
// all face dofs need to be
// non-additive, since they have
// continuity requirements.
// $Id$
// Version: $Name$
//
-// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
if (dim>1) dpo.push_back(fe_data.dofs_per_quad * N);
if (dim>2) dpo.push_back(fe_data.dofs_per_hex * N);
- return FiniteElementData<dim> (dpo, fe_data.n_components() * N, fe_data.tensor_degree(),
- fe_data.conforming_space);
+ return FiniteElementData<dim> (dpo,
+ fe_data.n_components() * N,
+ fe_data.tensor_degree(),
+ fe_data.conforming_space,
+ fe_data.n_blocks() * N);
}
// makes the degree of the system
// unknown.
unsigned int degree = std::max(fe1.tensor_degree(), fe2.tensor_degree());
- return FiniteElementData<dim> (dpo,
- fe1.n_components() * N1 +
- fe2.n_components() * N2,
- degree,
- typename
- FiniteElementData<dim>::Conformity(fe1.conforming_space
- & fe2.conforming_space));
+ return FiniteElementData<dim> (
+ dpo,
+ fe1.n_components() * N1 + fe2.n_components() * N2,
+ degree,
+ typename FiniteElementData<dim>::Conformity(fe1.conforming_space
+ & fe2.conforming_space),
+ fe1.n_blocks() * N1 + fe2.n_blocks() * N2);
}
// degree is the maximal degree of the components
unsigned int degree = std::max(fe1.tensor_degree(), fe2.tensor_degree());
degree = std::max(degree, fe3.tensor_degree());
- return FiniteElementData<dim> (dpo,
- fe1.n_components() * N1 +
- fe2.n_components() * N2 +
- fe3.n_components() * N3, degree,
- typename
- FiniteElementData<dim>::Conformity(fe1.conforming_space
- & fe2.conforming_space
- & fe3.conforming_space));
+ return FiniteElementData<dim> (
+ dpo,
+ fe1.n_components() * N1 + fe2.n_components() * N2 + fe3.n_components() * N3,
+ degree,
+ typename FiniteElementData<dim>::Conformity(fe1.conforming_space
+ & fe2.conforming_space
+ & fe3.conforming_space),
+ fe1.n_blocks() * N1 + fe2.n_blocks() * N2 + fe3.n_blocks() * N3);
}
*
* <dl>
*
- * <dt>@anchor GlossActive <b>Active cells</b></dt>
+ * <dt class="glossary">@anchor GlossActive Active cells</dt>
* <dd>Mesh cells not refined any further in the hierarchy.</dd>
*
+ * <dt class="glossary">@anchor GlossBlock block</dt>
+ * <dd>Originally, blocks were introduced in BlockVector,
+ * BlockSparseMatrix and related classes. These are used to reflect the
+ * structure of a PDE system in linear algebra, in particular allowing
+ * for modular solvers. In DoFHandler, this block structure is
+ * prepared by DoFRenumbering::component_wise().
+ *
+ * Originally, this concept was intermixed with the idea of the vector
+ * @ref GlossComponent. Since the introduction of non-@ref GlossPrimitive
+ * "primitive" elements, they became different. Take for instance the
+ * solution of the mixed Laplacian system with
+ * FE_RaviartThomas. There, the first <tt>dim</dim> components are the
+ * directional derivatives. Since the shape functions are linear
+ * combinations of those, they constitute only a single block. The
+ * primal function <i>u</i> would be in the second block, but in the
+ * <tt>dim+1</tt>st component.
+ * </dd>
+ *
+ * <dt class="glossary">@anchor GlossComponent component</dt>
+ *
+ * <dd>For vector functions, component denotes the index in the
+ * vector. For instance, in the mixed Laplacian system, the first
+ * <tt>dim</tt> components are the derivatives in each coordinate
+ * direction and the last component is the primal function <i>u</i>.
+ *
+ * Originally, components were not distinguished from @ref GlossBlocks
+ * "blocks", but since the introduction of non-@ref GlossPrimitive
+ * "primitive" elements, they have to be distinguished. See
+ * FiniteElementData::n_components() and the documentation of
+ * FimiteElement</dd>
*
- * <dt>@anchor GlossFaceOrientation <b>Face orientation</b></dt>
+ * <dt class="glossary">@anchor GlossFaceOrientation Face orientation</dt>
* <dd>In a triangulation, the normal vector to a face
* can be deduced from the face orientation by
* applying the right hand side rule (x,y -> normal). We note, that
* the QProjector class and its users.
*
*
- * <dt>@anchor GlossGeneralizedSupport <b>Generalized support points</b></dt>
+ * <dt class="glossary">@anchor GlossGeneralizedSupport Generalized support points</dt>
* <dd>While @ref GlossSupport "support points" allow very simple interpolation
* into the finite element space, their concept is restricted to
* @ref GlossLagrange "Lagrange elements". For other elements, more general
* </dd>
*
*
- * <dt>@anchor GlossInterpolation <b>Interpolation with finite elements</b></dt>
+ * <dt class="glossary">@anchor GlossInterpolation Interpolation with finite elements</dt>
* <dd>The purpose of interpolation with finite elements is computing
* a vector of coefficients representing a finite element function,
* such that the @ref GlossNodes "node values" of the original
* vector.
*
*
- * <dt>@anchor GlossLagrange <b>Lagrange elements</b></dt>
+ * <dt class="glossary">@anchor GlossLagrange Lagrange elements</dt>
* <dd>Finite elements based on Lagrangian interpolation at
* @ref GlossSupport "support points".</dd>
*
*
- * <dt>@anchor GlossNodes <b>Node values or node functionals</b></dt>
+ * <dt class="glossary">@anchor GlossNodes Node values or node functionals</dt>
*
* <dd>It is customary to define a FiniteElement as a pair consisting
* of a local function space and a set of node values $N_i$ on the
* <td>Gauss points on edges(faces) and anisotropic Gauss points in the interior</td></tr>
* </table>
*
+ * <dt class="glossary">@anchor GlossPrimitive Primitive finite elements</dt>
+ * <dd>Finite element shape function sets with a unique relation from
+ * shape function number to vector @ref GlossComponent.</dd>
*
- * <dt>@anchor GlossReferenceCell <b>Reference cell</b></dt>
+ * <dt class="glossary">@anchor GlossReferenceCell Reference cell</dt>
* <dd>The hypercube [0,1]<sup>dim</sup>, on which all parametric finite
* element shape functions are defined.</dd>
*
*
- * <dt>@anchor GlossShape <b>Shape functions</b></dt> <dd>The restriction of
+ * <dt class="glossary">@anchor GlossShape Shape functions</dt> <dd>The restriction of
* the finite element basis functions to a single grid cell.</dd>
*
*
- * <dt>@anchor GlossSupport <b>Support points</b></dt> <dd>Support points are
+ * <dt class="glossary">@anchor GlossSupport Support points</dt> <dd>Support points are
* by definition those points $p_i$, such that for the shape functions
* $v_j$ holds $v_j(p_i) = \delta_{ij}$. Therefore, a finite element
* interpolation can be defined uniquely by the values in the support
* </dd>
*
*
- * <dt>@anchor GlossTargetComponent <b>Target component</b></dt> <dd>When
+ * <dt class="glossary">@anchor GlossTargetComponent Target component</dt> <dd>When
* vectors and matrices are grouped into blocks by component, it is
* often desirable to collect several of the original components into
* a single one. This could be for instance, grouping the velocities
* of a Stokes system into a single block.</dd>
*
*
- * <dt>@anchor GlossUnitCell <b>Unit cell</b></dt>
+ * <dt class="glossary">@anchor GlossUnitCell Unit cell</dt>
* <dd>See @ref GlossReferenceCell "Reference cell".</dd>
*
*
- * <dt>@anchor GlossUnitSupport <b>Unit support points</b></dt>
+ * <dt class="glossary">@anchor GlossUnitSupport Unit support points</dt>
* <dd>@ref GlossSupport "Support points" on the reference cell, defined in
* FiniteElementBase. For example, the usual Q1 element in 1d has support
* points at <tt>x=0</tt> and <tt>x=1</tt> (and similarly, in higher
}
TD.tiny { font-size: 75%;
}
+
+dt.glossary {
+ font-size: large;
+ font-weight: bold;
+}