From: Guido Kanschat Date: Wed, 25 Jan 2006 22:02:15 +0000 (+0000) Subject: block concept introduced into FiniteElement X-Git-Tag: v8.0.0~12511 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f17106d302dc7c15c9897295d265385251bc8cb0;p=dealii.git block concept introduced into FiniteElement git-svn-id: https://svn.dealii.org/trunk@12164 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe.h b/deal.II/deal.II/include/fe/fe.h index ef61097ec2..1431457b97 100644 --- a/deal.II/deal.II/include/fe/fe.h +++ b/deal.II/deal.II/include/fe/fe.h @@ -59,6 +59,28 @@ namespace hp * under consideration here is hanging nodes constraints and grid * transfer, respectively. * + *

Components and blocks

+ * + * 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. * *

Support points

* @@ -948,66 +970,6 @@ class FiniteElement : public Subscriptor, std::pair 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, 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, unsigned int> - face_system_to_base_index (const unsigned int index) const; - /** * Return in which of the vector * components of this finite @@ -1154,6 +1116,73 @@ class FiniteElement : public Subscriptor, 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, 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, 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 @@ -1180,6 +1209,15 @@ class FiniteElement : public Subscriptor, */ std::pair component_to_base (const unsigned int component) const; + + /** + * The vector block and the index + * inside the block for this + * shape function. + */ + std::pair + system_to_block_index (const unsigned int component) const; + //@} /** @@ -1615,6 +1653,170 @@ class FiniteElement : public Subscriptor, 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 &mapping, + const typename Triangulation::cell_iterator &cell, + const unsigned int offset, + typename Mapping::InternalDataBase &mapping_internal, + InternalDataBase &fe_internal, + FEValuesData &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 + compute_n_nonzero_components (const std::vector > &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 + * fill_*_values 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 *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 > 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 > unit_face_support_points; + + /** + * Support points used for + * interpolation functions of + * non-Lagrangian elements. + */ + std::vector > generalized_support_points; + + /** + * Face support points used for + * interpolation functions of + * non-Lagrangian elements. + */ + std::vector > generalized_face_support_points; + + private: /** * Store what * @p system_to_component_index @@ -1682,6 +1884,12 @@ class FiniteElement : public Subscriptor, */ std::vector,unsigned int> > face_system_to_base_table; + /** + * For each base element, store + * the first block in a block + * vector it will generate. + */ + std::vector first_block_of_base_table; /** * The base element establishing @@ -1785,51 +1993,6 @@ class FiniteElement : public Subscriptor, */ const std::vector 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 > 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 > unit_face_support_points; - - /** - * Support points used for - * interpolation functions of - * non-Lagrangian elements. - */ - std::vector > generalized_support_points; - - /** - * Face support points used for - * interpolation functions of - * non-Lagrangian elements. - */ - std::vector > generalized_face_support_points; - /** * For each shape function, give * a vector of bools (with size @@ -1861,128 +2024,6 @@ class FiniteElement : public Subscriptor, */ const std::vector 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 &mapping, - const typename Triangulation::cell_iterator &cell, - const unsigned int offset, - typename Mapping::InternalDataBase &mapping_internal, - InternalDataBase &fe_internal, - FEValuesData &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 - compute_n_nonzero_components (const std::vector > &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 - * fill_*_values 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 *clone() const = 0; - - private: /** * Second derivatives of shapes * functions are not computed @@ -2106,38 +2147,12 @@ class FiniteElement : public Subscriptor, typename Mapping::InternalDataBase &fe_internal, FEValuesData &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 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; friend class FEValues; friend class FEFaceValues; friend class FESubfaceValues; + template friend class FESystem; friend class hp::FECollection; }; @@ -2271,6 +2286,18 @@ FiniteElement::face_system_to_base_index (const unsigned int index) const +template +inline +unsigned int +FiniteElement::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 inline std::pair @@ -2283,6 +2310,23 @@ FiniteElement::component_to_base (const unsigned int index) const } +template +inline +std::pair +FiniteElement::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( + 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 inline bool diff --git a/deal.II/deal.II/include/fe/fe_base.h b/deal.II/deal.II/include/fe/fe_base.h index f6aef269fb..7413193744 100644 --- a/deal.II/deal.II/include/fe/fe_base.h +++ b/deal.II/deal.II/include/fe/fe_base.h @@ -2,7 +2,7 @@ // $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 @@ -30,9 +30,6 @@ template class FESystem; -/*!@addtogroup febase */ -/*@{*/ - /** * Dimension independent data for finite elements. See the derived * class FiniteElement class for information on its use. All @@ -45,6 +42,7 @@ template class FESystem; * 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 @@ -255,7 +253,19 @@ class FiniteElementData * 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 @@ -300,11 +310,14 @@ class FiniteElementData * @param conformity The finite * element space has continuity * of this Sobolev space. + * @param n_blocks Number of + * vector blocks. */ FiniteElementData (const std::vector &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. @@ -353,6 +366,11 @@ class FiniteElementData */ unsigned int n_components () const; + /** + * Number of blocks. + */ + unsigned int n_blocks () const; + /** * Maximal polynomial degree of a * shape function in a single @@ -384,15 +402,6 @@ class FiniteElementData }; - - -/** - * - * @author Wolfgang Bangerth, 1998, 2002; Ralf Hartmann, Guido Kanschat, 2001 - */ -/*@}*/ - - // --------- inline and template functions --------------- template @@ -459,6 +468,16 @@ FiniteElementData::n_components () const +template +inline +unsigned int +FiniteElementData::n_blocks () const +{ + return blocks; +} + + + template inline unsigned int diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index d2d21c565c..d880f8ee5b 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -2,7 +2,7 @@ // $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 @@ -203,7 +203,8 @@ DoFTools::compute_row_length_vector( for (cell = dofs.begin_active(); cell != end; ++cell) { const FiniteElement& fe = cell->get_fe(); - Assert (fe.is_primitive(), typename FiniteElement::ExcFENotPrimitive()); + Assert (fe.is_primitive(), + typename FiniteElement::ExcFENotPrimitive()); Assert (couplings.n_rows()==fe.n_components(), ExcDimensionMismatch(couplings.n_rows(), fe.n_components())); Assert (couplings.n_cols()==fe.n_components(), diff --git a/deal.II/deal.II/source/fe/fe.cc b/deal.II/deal.II/source/fe/fe.cc index 5572c0fa89..5a84e66d0a 100644 --- a/deal.II/deal.II/source/fe/fe.cc +++ b/deal.II/deal.II/source/fe/fe.cc @@ -2,7 +2,7 @@ // $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 @@ -112,15 +112,13 @@ FiniteElement::FiniteElement ( const std::vector > &nonzero_c) : FiniteElementData (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 @@ -183,18 +181,22 @@ FiniteElement::FiniteElement ( // 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 ; jdofs_per_cell ; ++j) + // the element is not primitive, + // leave these tables empty. + if (!cached_primitivity) { - system_to_component_table[j] = std::pair(0,j); - system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j); - } - for (unsigned int j=0 ; jdofs_per_face ; ++j) - { - face_system_to_component_table[j] = std::pair(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 ; jdofs_per_cell ; ++j) + { + system_to_component_table[j] = std::pair(0,j); + system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j); + } + for (unsigned int j=0 ; jdofs_per_face ; ++j) + { + face_system_to_component_table[j] = std::pair(0,j); + face_system_to_base_table[j] = std::make_pair(std::make_pair(0U,0U),j); + } } } diff --git a/deal.II/deal.II/source/fe/fe_data.cc b/deal.II/deal.II/source/fe/fe_data.cc index c57c985610..defac8ce05 100644 --- a/deal.II/deal.II/source/fe/fe_data.cc +++ b/deal.II/deal.II/source/fe/fe_data.cc @@ -2,7 +2,7 @@ // $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 @@ -31,6 +31,7 @@ FiniteElementData::FiniteElementData () dofs_per_face(0), dofs_per_cell (0), components(0), + blocks(0), degree(0), conforming_space(unknown) {} @@ -42,7 +43,8 @@ FiniteElementData:: FiniteElementData (const std::vector &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]), @@ -73,6 +75,8 @@ FiniteElementData (const std::vector &dofs_per_object, GeometryInfo::quads_per_cell * dofs_per_quad + GeometryInfo::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) { diff --git a/deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc b/deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc index b7c1ff6fb2..7c6b1fd63d 100644 --- a/deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc +++ b/deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc @@ -2,7 +2,7 @@ // $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 @@ -30,10 +30,14 @@ template FE_DGPNonparametric::FE_DGPNonparametric (const unsigned int degree) : - FiniteElement (FiniteElementData(get_dpo_vector(degree),1, FiniteElementData::L2), - std::vector(FiniteElementData(get_dpo_vector(degree),1).dofs_per_cell,true), - std::vector >(FiniteElementData(get_dpo_vector(degree),1).dofs_per_cell, - std::vector(1,true))), + FiniteElement ( + FiniteElementData(get_dpo_vector(degree), 1, degree, + FiniteElementData::L2), + std::vector( + FiniteElementData(get_dpo_vector(degree), 1, degree).dofs_per_cell,true), + std::vector >( + FiniteElementData(get_dpo_vector(degree),1, degree).dofs_per_cell, + std::vector(1,true))), degree(degree), polynomial_space (Polynomials::Legendre::generate_complete_basis(degree)) { diff --git a/deal.II/deal.II/source/fe/fe_dgq.cc b/deal.II/deal.II/source/fe/fe_dgq.cc index d88db932ce..c5694e7d7a 100644 --- a/deal.II/deal.II/source/fe/fe_dgq.cc +++ b/deal.II/deal.II/source/fe/fe_dgq.cc @@ -2,7 +2,7 @@ // $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 diff --git a/deal.II/deal.II/source/fe/fe_nedelec.cc b/deal.II/deal.II/source/fe/fe_nedelec.cc index 0cc239c62d..81915b73fa 100644 --- a/deal.II/deal.II/source/fe/fe_nedelec.cc +++ b/deal.II/deal.II/source/fe/fe_nedelec.cc @@ -34,10 +34,16 @@ template FE_Nedelec::FE_Nedelec (const unsigned int degree) : - FiniteElement (FiniteElementData(get_dpo_vector(degree), dim, degree+1, FiniteElementData::Hcurl), - std::vector (FiniteElementData(get_dpo_vector(degree),dim).dofs_per_cell,false), - std::vector >(FiniteElementData(get_dpo_vector(degree),dim).dofs_per_cell, - std::vector(dim,true))), + FiniteElement ( + FiniteElementData(get_dpo_vector(degree), dim, + degree+1, FiniteElementData::Hcurl, 1), + std::vector ( + FiniteElementData(get_dpo_vector(degree), dim, + degree+1).dofs_per_cell,false), + std::vector >( + FiniteElementData(get_dpo_vector(degree), dim, + degree+1).dofs_per_cell, + std::vector(dim,true))), degree(degree) { Assert (dim >= 2, ExcImpossibleInDim(dim)); @@ -54,14 +60,6 @@ FE_Nedelec::FE_Nedelec (const unsigned int degree) // 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 > tmp1, tmp2; - this->system_to_component_table.swap (tmp1); - this->face_system_to_component_table.swap (tmp2); } diff --git a/deal.II/deal.II/source/fe/fe_q.cc b/deal.II/deal.II/source/fe/fe_q.cc index 4febfb4e0d..80a56d2a68 100644 --- a/deal.II/deal.II/source/fe/fe_q.cc +++ b/deal.II/deal.II/source/fe/fe_q.cc @@ -2,7 +2,7 @@ // $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 @@ -435,7 +435,7 @@ template std::vector FE_Q::face_lexicographic_to_hierarchic_numbering (const unsigned int degree) { - const FiniteElementData face_data(FE_Q::get_dpo_vector(degree),1); + const FiniteElementData face_data(FE_Q::get_dpo_vector(degree),1,degree); std::vector face_renumber (face_data.dofs_per_cell); FETools::lexicographic_to_hierarchic_numbering (face_data, face_renumber); return face_renumber; diff --git a/deal.II/deal.II/source/fe/fe_q_hierarchical.cc b/deal.II/deal.II/source/fe/fe_q_hierarchical.cc index 4e998dcaff..931202cd22 100644 --- a/deal.II/deal.II/source/fe/fe_q_hierarchical.cc +++ b/deal.II/deal.II/source/fe/fe_q_hierarchical.cc @@ -2,7 +2,7 @@ // $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 @@ -755,7 +755,7 @@ std::vector FE_Q_Hierarchical:: face_fe_q_hierarchical_to_hierarchic_numbering (const unsigned int degree) { - FiniteElementData fe_data(FE_Q_Hierarchical::get_dpo_vector(degree),1); + FiniteElementData fe_data(FE_Q_Hierarchical::get_dpo_vector(degree),1,degree); return invert_numbering(FE_Q_Hierarchical:: hierarchic_to_fe_q_hierarchical_numbering (fe_data)); } diff --git a/deal.II/deal.II/source/fe/fe_raviart_thomas.cc b/deal.II/deal.II/source/fe/fe_raviart_thomas.cc index 9f90a22ec2..4c14f9a742 100644 --- a/deal.II/deal.II/source/fe/fe_raviart_thomas.cc +++ b/deal.II/deal.II/source/fe/fe_raviart_thomas.cc @@ -40,7 +40,7 @@ FE_RaviartThomas::FE_RaviartThomas (const unsigned int deg) FE_PolyTensor, dim> ( deg, FiniteElementData(get_dpo_vector(deg), - dim, deg+1, FiniteElementData::Hdiv), + dim, deg+1, FiniteElementData::Hdiv, 1), std::vector(PolynomialsRaviartThomas::compute_n_pols(deg), true), std::vector >(PolynomialsRaviartThomas::compute_n_pols(deg), std::vector(dim,true))), @@ -91,14 +91,6 @@ FE_RaviartThomas::FE_RaviartThomas (const unsigned int deg) 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 > tmp1, tmp2; - this->system_to_component_table.swap (tmp1); - this->face_system_to_component_table.swap (tmp2); } diff --git a/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc b/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc index a786989222..f01a405cbd 100644 --- a/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc +++ b/deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc @@ -2,7 +2,7 @@ // $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 @@ -35,12 +35,10 @@ FE_RaviartThomasNodal::FE_RaviartThomasNodal (const unsigned int deg) FE_PolyTensor, dim> ( deg, FiniteElementData(get_dpo_vector(deg), - dim, deg+1, FiniteElementData::Hdiv), + dim, deg+1, FiniteElementData::Hdiv, 1), get_ria_vector (deg), - std::vector >( - FiniteElementData(get_dpo_vector(deg), - dim,deg+1).dofs_per_cell, - std::vector(dim,true))) + std::vector >(PolynomialsRaviartThomas::compute_n_pols(deg), + std::vector(dim,true))) { Assert (dim >= 2, ExcImpossibleInDim(dim)); const unsigned int n_dofs = this->dofs_per_cell; @@ -290,27 +288,10 @@ template std::vector FE_RaviartThomasNodal::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(get_dpo_vector(deg),dim).dofs_per_cell == - dofs_per_cell, - ExcInternalError()); - Assert (FiniteElementData(get_dpo_vector(deg),dim).dofs_per_face == - dofs_per_face, - ExcInternalError()); - + const unsigned int dofs_per_cell = PolynomialsRaviartThomas::compute_n_pols(deg); + unsigned int dofs_per_face = deg+1; + for (unsigned int d=2;d::multiply_dof_numbers (const FiniteElementData &fe_data, 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 (dpo, fe_data.n_components() * N, fe_data.tensor_degree(), - fe_data.conforming_space); + return FiniteElementData (dpo, + fe_data.n_components() * N, + fe_data.tensor_degree(), + fe_data.conforming_space, + fe_data.n_blocks() * N); } @@ -1849,13 +1852,13 @@ FESystem::multiply_dof_numbers (const FiniteElementData &fe1, // makes the degree of the system // unknown. unsigned int degree = std::max(fe1.tensor_degree(), fe2.tensor_degree()); - return FiniteElementData (dpo, - fe1.n_components() * N1 + - fe2.n_components() * N2, - degree, - typename - FiniteElementData::Conformity(fe1.conforming_space - & fe2.conforming_space)); + return FiniteElementData ( + dpo, + fe1.n_components() * N1 + fe2.n_components() * N2, + degree, + typename FiniteElementData::Conformity(fe1.conforming_space + & fe2.conforming_space), + fe1.n_blocks() * N1 + fe2.n_blocks() * N2); } @@ -1885,14 +1888,14 @@ FESystem::multiply_dof_numbers (const FiniteElementData &fe1, // 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 (dpo, - fe1.n_components() * N1 + - fe2.n_components() * N2 + - fe3.n_components() * N3, degree, - typename - FiniteElementData::Conformity(fe1.conforming_space - & fe2.conforming_space - & fe3.conforming_space)); + return FiniteElementData ( + dpo, + fe1.n_components() * N1 + fe2.n_components() * N2 + fe3.n_components() * N3, + degree, + typename FiniteElementData::Conformity(fe1.conforming_space + & fe2.conforming_space + & fe3.conforming_space), + fe1.n_blocks() * N1 + fe2.n_blocks() * N2 + fe3.n_blocks() * N3); } diff --git a/deal.II/doc/doxygen/headers/glossary.h b/deal.II/doc/doxygen/headers/glossary.h index d972afc2c1..1da7ff5e2c 100644 --- a/deal.II/doc/doxygen/headers/glossary.h +++ b/deal.II/doc/doxygen/headers/glossary.h @@ -9,11 +9,41 @@ * *
* - *
@anchor GlossActive Active cells
+ *
@anchor GlossActive Active cells
*
Mesh cells not refined any further in the hierarchy.
* + *
@anchor GlossBlock block
+ *
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 dim components are the + * directional derivatives. Since the shape functions are linear + * combinations of those, they constitute only a single block. The + * primal function u would be in the second block, but in the + * dim+1st component. + *
+ * + *
@anchor GlossComponent component
+ * + *
For vector functions, component denotes the index in the + * vector. For instance, in the mixed Laplacian system, the first + * dim components are the derivatives in each coordinate + * direction and the last component is the primal function u. + * + * 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
* - *
@anchor GlossFaceOrientation Face orientation
+ *
@anchor GlossFaceOrientation Face orientation
*
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 @@ -51,7 +81,7 @@ * the QProjector class and its users. * * - *
@anchor GlossGeneralizedSupport Generalized support points
+ *
@anchor GlossGeneralizedSupport Generalized support points
*
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 @@ -70,7 +100,7 @@ *
* * - *
@anchor GlossInterpolation Interpolation with finite elements
+ *
@anchor GlossInterpolation Interpolation with finite elements
*
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 @@ -81,12 +111,12 @@ * vector. * * - *
@anchor GlossLagrange Lagrange elements
+ *
@anchor GlossLagrange Lagrange elements
*
Finite elements based on Lagrangian interpolation at * @ref GlossSupport "support points".
* * - *
@anchor GlossNodes Node values or node functionals
+ *
@anchor GlossNodes Node values or node functionals
* *
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 @@ -123,17 +153,20 @@ * Gauss points on edges(faces) and anisotropic Gauss points in the interior * * + *
@anchor GlossPrimitive Primitive finite elements
+ *
Finite element shape function sets with a unique relation from + * shape function number to vector @ref GlossComponent.
* - *
@anchor GlossReferenceCell Reference cell
+ *
@anchor GlossReferenceCell Reference cell
*
The hypercube [0,1]dim, on which all parametric finite * element shape functions are defined.
* * - *
@anchor GlossShape Shape functions
The restriction of + *
@anchor GlossShape Shape functions
The restriction of * the finite element basis functions to a single grid cell.
* * - *
@anchor GlossSupport Support points
Support points are + *
@anchor GlossSupport Support points
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 @@ -154,18 +187,18 @@ *
* * - *
@anchor GlossTargetComponent Target component
When + *
@anchor GlossTargetComponent Target component
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.
* * - *
@anchor GlossUnitCell Unit cell
+ *
@anchor GlossUnitCell Unit cell
*
See @ref GlossReferenceCell "Reference cell".
* * - *
@anchor GlossUnitSupport Unit support points
+ *
@anchor GlossUnitSupport Unit support points
*
@ref GlossSupport "Support points" on the reference cell, defined in * FiniteElementBase. For example, the usual Q1 element in 1d has support * points at x=0 and x=1 (and similarly, in higher diff --git a/deal.II/doc/doxygen/stylesheet.css b/deal.II/doc/doxygen/stylesheet.css index 0bfc78c445..8d6bd0ee15 100644 --- a/deal.II/doc/doxygen/stylesheet.css +++ b/deal.II/doc/doxygen/stylesheet.css @@ -306,3 +306,8 @@ INPUT.search { font-size: 75%; } TD.tiny { font-size: 75%; } + +dt.glossary { + font-size: large; + font-weight: bold; +}