]> https://gitweb.dealii.org/ - dealii.git/commitdiff
block concept introduced into FiniteElement
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 25 Jan 2006 22:02:15 +0000 (22:02 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 25 Jan 2006 22:02:15 +0000 (22:02 +0000)
git-svn-id: https://svn.dealii.org/trunk@12164 0785d39b-7218-0410-832d-ea1e28bc413d

15 files changed:
deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_base.h
deal.II/deal.II/source/dofs/dof_tools.cc
deal.II/deal.II/source/fe/fe.cc
deal.II/deal.II/source/fe/fe_data.cc
deal.II/deal.II/source/fe/fe_dgp_nonparametric.cc
deal.II/deal.II/source/fe/fe_dgq.cc
deal.II/deal.II/source/fe/fe_nedelec.cc
deal.II/deal.II/source/fe/fe_q.cc
deal.II/deal.II/source/fe/fe_q_hierarchical.cc
deal.II/deal.II/source/fe/fe_raviart_thomas.cc
deal.II/deal.II/source/fe/fe_raviart_thomas_nodal.cc
deal.II/deal.II/source/fe/fe_system.cc
deal.II/doc/doxygen/headers/glossary.h
deal.II/doc/doxygen/stylesheet.css

index ef61097ec24d7bba058abf6d4f428ce98b7e509e..1431457b97d79f07101d1ba6699c41ec7cd32587 100644 (file)
@@ -59,6 +59,28 @@ namespace hp
  * 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>
  *
@@ -948,66 +970,6 @@ class FiniteElement : public Subscriptor,
     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
@@ -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<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
@@ -1180,6 +1209,15 @@ class FiniteElement : public Subscriptor,
                                      */
     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;
+
                                     //@}
     
                                     /**
@@ -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<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
@@ -1682,6 +1884,12 @@ class FiniteElement : public Subscriptor,
                                      */
     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
@@ -1785,51 +1993,6 @@ class FiniteElement : public Subscriptor,
                                      */
     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
@@ -1861,128 +2024,6 @@ class FiniteElement : public Subscriptor,
                                      */
     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
@@ -2106,38 +2147,12 @@ class FiniteElement : public Subscriptor,
                            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>;
 };
 
@@ -2271,6 +2286,18 @@ FiniteElement<dim>::face_system_to_base_index (const unsigned int index) const
 
 
 
+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>
@@ -2283,6 +2310,23 @@ FiniteElement<dim>::component_to_base (const unsigned int index) const
 }
 
 
+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
index f6aef269fb04e5e627d5a138c5603181a32f2ab2..741319374469d912a9d2f36353c6b5f49058bc92 100644 (file)
@@ -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<int dim> 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<int dim> 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 <int dim>
@@ -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<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.
@@ -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 <int dim>
@@ -459,6 +468,16 @@ FiniteElementData<dim>::n_components () const
 
 
 
+template <int dim>
+inline
+unsigned int 
+FiniteElementData<dim>::n_blocks () const
+{
+  return blocks;
+}
+
+
+
 template <int dim>
 inline
 unsigned int 
index d2d21c565cc9a745df1617f804ac77531f275bcb..d880f8ee5b3f13640c13cb87064dbb69e8637514 100644 (file)
@@ -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<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(),
index 5572c0fa8939a24290d3cf96d672aa96d80cbe13..5a84e66d0a296f5fda88fd6a2fb32cdd63cc33df 100644 (file)
@@ -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<dim>::FiniteElement (
   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
@@ -183,18 +181,22 @@ FiniteElement<dim>::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 ; 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);      
+       }
     }
 }
 
index c57c9856105030b122917bd28e1ed844fe41749d..defac8ce053a54a976c9e3b93f6ec16592242f21 100644 (file)
@@ -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<dim>::FiniteElementData ()
                dofs_per_face(0),
                dofs_per_cell (0),
                components(0),
+               blocks(0),
                degree(0),
                conforming_space(unknown)
 {}
@@ -42,7 +43,8 @@ FiniteElementData<dim>::
 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]),
@@ -73,6 +75,8 @@ FiniteElementData (const std::vector<unsigned int> &dofs_per_object,
                               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)
 {
index b7c1ff6fb25690b5883dd9ede085904a229ac65e..7c6b1fd63d01110775e62e9ae4e27d66c1e8004b 100644 (file)
@@ -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
 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))
 {
index d88db932cef88fdc99d9b0fb54a3faa665d98970..c5694e7d7a061b02c7c7002695b1e8c45a0ffc1c 100644 (file)
@@ -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
index 0cc239c62de99f7ac3db6eb2c5d08872d22be2b1..81915b73fa2d2632048dd4e8964b5258c3888354 100644 (file)
 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));
@@ -54,14 +60,6 @@ FE_Nedelec<dim>::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<std::pair<unsigned,unsigned> > tmp1, tmp2;
-  this->system_to_component_table.swap (tmp1);
-  this->face_system_to_component_table.swap (tmp2);
 }
 
 
index 4febfb4e0de3b67ea33ee07e102e3cf914ed5af2..80a56d2a682922968b3ed0f6845c46b3f64e1f12 100644 (file)
@@ -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 <int dim>
 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;
index 4e998dcaffe99f3b3722f8999c58d2c7de5b47d4..931202cd22884b968fb6e209e4de07c639d3978a 100644 (file)
@@ -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<unsigned int>
 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));
 }
index 9f90a22ec2648c95f1d9564c633531b5033d812a..4c14f9a7428e97755e070624fc1eeebe15ee07ed 100644 (file)
@@ -40,7 +40,7 @@ FE_RaviartThomas<dim>::FE_RaviartThomas (const unsigned int deg)
                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))),
@@ -91,14 +91,6 @@ FE_RaviartThomas<dim>::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<std::pair<unsigned,unsigned> > tmp1, tmp2;
-  this->system_to_component_table.swap (tmp1);
-  this->face_system_to_component_table.swap (tmp2);
 }
 
 
index a78698922286d88718b6ade443065fae4eb414b3..f01a405cbd715b0f122664726e3cedbd4ad8cf80 100644 (file)
@@ -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<dim>::FE_RaviartThomasNodal (const unsigned int deg)
                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;
@@ -290,27 +288,10 @@ template <int dim>
 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.
index 9ed9b15973727c23ef409b47ee6a9a30ff1af0d8..f06ff678ea73ff221cab9a7a001cd2b8a9353849 100644 (file)
@@ -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
@@ -1730,8 +1730,11 @@ FESystem<dim>::multiply_dof_numbers (const FiniteElementData<dim> &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<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);
 }
 
 
@@ -1849,13 +1852,13 @@ FESystem<dim>::multiply_dof_numbers (const FiniteElementData<dim> &fe1,
                                   // 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);
 }
 
 
@@ -1885,14 +1888,14 @@ FESystem<dim>::multiply_dof_numbers (const FiniteElementData<dim> &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<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);
 }
 
 
index d972afc2c12db0921a3ad2dfc39b00841e6c1c90..1da7ff5e2c30e6b99d4b2faa953db479a9adf341 100644 (file)
@@ -9,11 +9,41 @@
  *
  * <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
@@ -51,7 +81,7 @@
  * 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
index 0bfc78c445fba90b7685e838e753a2c08169d2a9..8d6bd0ee154ed5d7a467987047edd7c03d5d01de 100644 (file)
@@ -306,3 +306,8 @@ INPUT.search { font-size: 75%;
 }
 TD.tiny      { font-size: 75%;
 }
+
+dt.glossary {
+              font-size: large;
+              font-weight: bold;
+}

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.