From 9013c02f05dc010fcf444690e0b91dfaba54fb88 Mon Sep 17 00:00:00 2001 From: wolf Date: Wed, 19 Sep 2001 14:02:42 +0000 Subject: [PATCH] Slightly extend the infrastructure around base elements and components. Significantly extend the docs which were not very clear at some places. git-svn-id: https://svn.dealii.org/trunk@5027 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/fe/fe.h | 40 +++++++++++++------- deal.II/deal.II/include/fe/fe_base.h | 13 ++++++- deal.II/deal.II/include/fe/fe_dgq.h | 17 +++++++++ deal.II/deal.II/include/fe/fe_q.h | 17 +++++++++ deal.II/deal.II/include/fe/fe_system.h | 52 +++++++------------------- deal.II/deal.II/source/fe/fe.cc | 13 +------ deal.II/deal.II/source/fe/fe_dgq.cc | 20 ++++++++++ deal.II/deal.II/source/fe/fe_q.cc | 19 ++++++++++ deal.II/deal.II/source/fe/fe_system.cc | 29 +++++++++++--- 9 files changed, 150 insertions(+), 70 deletions(-) diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index 34ae55012c..ecfbc594db 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -69,24 +69,38 @@ class FiniteElement : public FiniteElementBase virtual ~FiniteElement (); /** - * Number of base elements in a mixed - * discretization. This function returns - * 1 for primitive elements. + * Number of base elements in a + * mixed discretization. + * + * Note that even for vector + * valued finite elements, the + * number of components needs not + * coincide with the number of + * base elements, since they may + * be reused. For example, if you + * create a @ref{FESystem} with + * three identical finite element + * classes by using the + * constructor that takes one + * finite element and a + * multiplicity, then the number + * of base elements is still one, + * although the number of + * components of the finite + * element is equal to the + * multiplicity. */ - virtual unsigned int n_base_elements () const; + virtual unsigned int n_base_elements () const = 0; /** * Access to base element - * objects. By default, - * #base_element(0)# is #this#. - * This function is overloaded by - * system elements to allow - * access to the different - * components of mixed - * discretizations. + * objects. If the element is + * scalar, then + * @p{base_element(0)} is + * @p{this}. */ - virtual const FiniteElement & base_element (const unsigned int index) const; - + virtual const FiniteElement & base_element (const unsigned int index) const = 0; + /** * Determine an estimate for the * memory consumption (in bytes) diff --git a/deal.II/deal.II/include/fe/fe_base.h b/deal.II/deal.II/include/fe/fe_base.h index f3935c6aa8..aa533bd5b3 100644 --- a/deal.II/deal.II/include/fe/fe_base.h +++ b/deal.II/deal.II/include/fe/fe_base.h @@ -625,7 +625,7 @@ class FiniteElementBase : public Subscriptor, * shape functions of the base * element. */ - unsigned int component_to_base(unsigned int index) const; + unsigned int component_to_base (unsigned int index) const; /** * Access the @p{restriction_is_additive_flag} @@ -881,6 +881,17 @@ class FiniteElementBase : public Subscriptor, * the result allows access to * shape functions of the base * element. + * + * This variable is set to the + * correct size by the + * constructor of this class, but + * needs to be initialized by + * derived classes, unless its + * size is one and the only entry + * is a zero, which is the case + * for scalar elements. In that + * case, the initialization by + * the base class is sufficient. */ std::vector component_to_base_table; diff --git a/deal.II/deal.II/include/fe/fe_dgq.h b/deal.II/deal.II/include/fe/fe_dgq.h index 7f22635c3f..8a89d33ae9 100644 --- a/deal.II/deal.II/include/fe/fe_dgq.h +++ b/deal.II/deal.II/include/fe/fe_dgq.h @@ -86,6 +86,23 @@ class FE_DGQ : public FiniteElement */ unsigned int get_degree () const; + /** + * Number of base elements in a + * mixed discretization. Since + * this is a scalar element, + * return one. + */ + virtual unsigned int n_base_elements () const; + + /** + * Access to base element + * objects. Since this element is + * scalar, @p{base_element(0)} is + * @p{this}, and all other + * indices throw an error. + */ + virtual const FiniteElement & base_element (const unsigned int index) const; + /** * Determine an estimate for the * memory consumption (in bytes) diff --git a/deal.II/deal.II/include/fe/fe_q.h b/deal.II/deal.II/include/fe/fe_q.h index 40bfd075ff..02e892b418 100644 --- a/deal.II/deal.II/include/fe/fe_q.h +++ b/deal.II/deal.II/include/fe/fe_q.h @@ -291,6 +291,23 @@ class FE_Q : public FiniteElement */ unsigned int get_degree () const; + /** + * Number of base elements in a + * mixed discretization. Since + * this is a scalar element, + * return one. + */ + virtual unsigned int n_base_elements () const; + + /** + * Access to base element + * objects. Since this element is + * scalar, @p{base_element(0)} is + * @p{this}, and all other + * indices throw an error. + */ + virtual const FiniteElement & base_element (const unsigned int index) const; + /** * Determine an estimate for the * memory consumption (in bytes) diff --git a/deal.II/deal.II/include/fe/fe_system.h b/deal.II/deal.II/include/fe/fe_system.h index 2b2b0b25b8..b92787b709 100644 --- a/deal.II/deal.II/include/fe/fe_system.h +++ b/deal.II/deal.II/include/fe/fe_system.h @@ -197,28 +197,25 @@ class FESystem : public FiniteElement virtual unsigned int n_base_elements () const; /** - * How often is a composing element used. - * + * How often is a composing + * element used. */ unsigned int element_multiplicity (const unsigned int index) const; /** - * Access to a composing element. - * - * If you assemble your system - * matrix, you usually will not - * want to have an FEValues object - * with a lot of equal entries. Ok, - * so initialize your FEValues with - * the @p{base_element} you get by - * this function. In a mixed - * discretization, you can choose - * the different base element types - * by index. - * + * Access to a composing + * element. The index needs to be + * smaller than the number of + * base elements. Note that the + * number of base elements may in + * turn be smaller than the + * number of components of the + * system element, if the + * multiplicities are greater + * than one. */ - virtual const FiniteElement & base_element(const unsigned int index) const; - + virtual const FiniteElement & base_element (const unsigned int index) const; + /** * Determine an estimate for the * memory consumption (in bytes) @@ -626,14 +623,6 @@ template <> void FESystem<1>::initialize_unit_face_support_points (); /* ------------------------- inline functions ------------------------- */ -template -inline unsigned int -FESystem::n_base_elements() const -{ - return base_elements.size(); -}; - - template inline unsigned int @@ -646,19 +635,6 @@ FESystem::element_multiplicity (const unsigned int index) const -template -inline const FiniteElement & -FESystem::base_element (const unsigned int index) const -{ - Assert (index < base_elements.size(), - ExcIndexRange(index, 0, base_elements.size())); - return *base_elements[index].first; -}; - - - - - /*---------------------------- fe_system.h ---------------------------*/ #endif /*---------------------------- fe_system.h ---------------------------*/ diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index d77a428651..f3d8892814 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -113,7 +113,7 @@ FiniteElementBase::FiniteElementBase (const FiniteElementData &fe_data face_system_to_component_table(dofs_per_face), component_to_system_table(components, std::vector(dofs_per_cell)), face_component_to_system_table(components, std::vector(dofs_per_face)), - component_to_base_table(components), + component_to_base_table (components, 0), restriction_is_additive_flags(restriction_is_additive_flags) { Assert(restriction_is_additive_flags.size()==fe_data.components, @@ -444,17 +444,6 @@ FiniteElement::get_subface_data (const UpdateFlags flags, -template -unsigned int -FiniteElement::n_base_elements() const -{ - // default implementation - return 1; -} - - - - template unsigned int FiniteElement::memory_consumption () const diff --git a/deal.II/deal.II/source/fe/fe_dgq.cc b/deal.II/deal.II/source/fe/fe_dgq.cc index 6574b23e89..77b6340ab7 100644 --- a/deal.II/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/deal.II/source/fe/fe_dgq.cc @@ -563,6 +563,26 @@ FE_DGQ::fill_fe_subface_values (const Mapping &mappi } + +template +unsigned int +FE_DGQ::n_base_elements () const +{ + return 1; +}; + + + +template +const FiniteElement & +FE_DGQ::base_element (const unsigned int index) const +{ + Assert (index==0, ExcIndexRange(index, 0, 1)); + return *this; +}; + + + template unsigned int FE_DGQ::memory_consumption () const diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index ccb2cb0504..8de210c466 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -1191,6 +1191,25 @@ FE_Q::fill_fe_subface_values (const Mapping &mapping +template +unsigned int +FE_Q::n_base_elements () const +{ + return 1; +}; + + + +template +const FiniteElement & +FE_Q::base_element (const unsigned int index) const +{ + Assert (index==0, ExcIndexRange(index, 0, 1)); + return *this; +}; + + + template unsigned int FE_Q::memory_consumption () const diff --git a/deal.II/deal.II/source/fe/fe_system.cc b/deal.II/deal.II/source/fe/fe_system.cc index 1e6c8eed27..191c766102 100644 --- a/deal.II/deal.II/source/fe/fe_system.cc +++ b/deal.II/deal.II/source/fe/fe_system.cc @@ -664,6 +664,8 @@ FESystem::build_cell_table() for (unsigned int base=0 ; base < n_base_elements() ; ++base) for (unsigned int m = 0; m < element_multiplicity(base); ++m) component_to_base_table[total_index++] = base; + Assert (total_index == component_to_base_table.size(), + ExcInternalError()); // Initialize index table // Multi-component base elements have @@ -793,17 +795,12 @@ template void FESystem::build_face_table() { - unsigned total_index = 0; - for (unsigned int base=0 ; base < n_base_elements() ; ++base) - for (unsigned int m = 0; m < element_multiplicity(base); ++m) - component_to_base_table[total_index++] = base; - // Initialize index table // Multi-component base elements have // to be thought of. // 1. Vertices - total_index = 0; + unsigned int total_index = 0; for (unsigned int vertex_number= 0 ; vertex_number < GeometryInfo::vertices_per_face ; ++vertex_number) { @@ -1423,6 +1420,26 @@ FESystem::compute_restriction_is_additive_flags (const FiniteElement & +template +const FiniteElement & +FESystem::base_element (const unsigned int index) const +{ + Assert (index < base_elements.size(), + ExcIndexRange(index, 0, base_elements.size())); + return *base_elements[index].first; +}; + + + +template +unsigned int +FESystem::n_base_elements () const +{ + return base_elements.size(); +}; + + + template unsigned int FESystem::memory_consumption () const -- 2.39.5