From: Wolfgang Bangerth Date: Thu, 12 Jan 2006 03:13:13 +0000 (+0000) Subject: Create a sub-module hpcollections, describe it, and link the collections classes... X-Git-Tag: v8.0.0~12677 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=db0923b77ac041e2c1f47e23e89b752b369d953a;p=dealii.git Create a sub-module hpcollections, describe it, and link the collections classes to it. git-svn-id: https://svn.dealii.org/trunk@11985 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_collection.h b/deal.II/deal.II/include/fe/fe_collection.h index 237ad4eabc..837e9bc85c 100644 --- a/deal.II/deal.II/include/fe/fe_collection.h +++ b/deal.II/deal.II/include/fe/fe_collection.h @@ -24,9 +24,10 @@ namespace hp /** * This class acts as a collection of finite element objects used in the * hp::DoFHandler. It is thus to a hp::DoFHandler what a - * FiniteElement is to a ::DoFHandler. This collection stores copies - * of the original elements added to it, and is therefore handling its memory - * manegement itself. + * FiniteElement is to a ::DoFHandler. + * + * It implements the concepts stated in the @ref hpcollection module described + * in the doxygen documentation. * * In addition to offering access to the elements of the collection, this * class provides access to the maximal number of degrees of freedom per @@ -34,7 +35,7 @@ namespace hp * the worst case when using the finite elements associated with the cells of * a triangulation. * - * @ingroup hp + * @ingroup hp hpcollection * * @author Wolfgang Bangerth, 2003 */ diff --git a/deal.II/deal.II/include/fe/mapping_collection.h b/deal.II/deal.II/include/fe/mapping_collection.h index 8446e44b58..4ef23db561 100644 --- a/deal.II/deal.II/include/fe/mapping_collection.h +++ b/deal.II/deal.II/include/fe/mapping_collection.h @@ -28,16 +28,19 @@ namespace hp * This class implements a collection of mapping objects in the same way as * the hp::FECollection implements a collection of finite element classes. * + * It implements the concepts stated in the @ref hpcollection module described + * in the doxygen documentation. + * * Although it is recommended to supply an appropriate mapping for each finite * element kind used in a hp-computation, the MappingCollection class * implements a conversion constructor from a single mapping. Therefore it is - * possible to offer only a single mapping to the hpFEValues class instead of - * a MappingCollection. This is for the convenience of the user, as many + * possible to offer only a single mapping to the hp::FEValues class instead of + * a hp::MappingCollection. This is for the convenience of the user, as many * simple geometries do not require different mappings along the boundary to * achieve optimal convergence rates. Hence providing a single mapping object * will usually suffice. * - * @ingroup hp + * @ingroup hp hpcollection * * @author Oliver Kayser-Herold, 2005 */ diff --git a/deal.II/deal.II/include/fe/q_collection.h b/deal.II/deal.II/include/fe/q_collection.h index 8fd68127ef..625fca557c 100644 --- a/deal.II/deal.II/include/fe/q_collection.h +++ b/deal.II/deal.II/include/fe/q_collection.h @@ -29,18 +29,10 @@ namespace hp * This class implements a collection of quadrature objects in the same way as * the hp::FECollection implements a collection of finite element classes. * - * Although it is strongly recommended to supply an appropriate quadrature - * for each finite element type used in a hp-computation, the QCollection - * class implements a conversion constructor from a single quadrature. - * Therefore it is possible to offer only a single quadrature to the - * hpFEValues class instead of a QCollection. The reason for this - * mechanism lies in the structure of the deal.II library. At many places - * throughout the library pseudo quadrature rules are constructed to obtain - * the values of the finite element function at certain points. An example - * is the DataOut class. Due to this conversion constructor, these lines - * of code continue to work without modifications. + * It implements the concepts stated in the @ref hpcollection module described + * in the doxygen documentation. * - * @ingroup hp + * @ingroup hp hpcollection * * @author Oliver Kayser-Herold, 2005 */ diff --git a/deal.II/doc/doxygen/headers/deal.II/dox.h b/deal.II/doc/doxygen/headers/deal.II/dox.h index aa13e5fc9b..36c7a786c1 100644 --- a/deal.II/doc/doxygen/headers/deal.II/dox.h +++ b/deal.II/doc/doxygen/headers/deal.II/dox.h @@ -56,4 +56,67 @@ * Classes and functions that have to do with hp finite elements. */ +/** + * @defgroup hpcollection hp Collections + * + * In the implementation of the hp finite element method, each cell might have + * a different finite element associated with it. To handle this, the + * hp::DoFHandler must have a whole set of finite element classes associated + * with it. This concept is represented by the hp::FECollection class: Objects + * of this type act as containers that hold a whole set of finite element + * objects. Instead of storing pointers to finite element objects on each + * cell, we then only store an index for each cell that identifies the finite + * element object within the collection that should be used by this cell. The + * DoFHandler object associated with the given cell can then assign degrees of + * freedom to each cell in accordance with the finite element used for it. + * + * A similar situation arises when integrating terms on a cell: one may want + * to use different quadrature formulas for different finite elements. For + * example, on cells where we use a Q1 element, a QGauss(2) object (i.e. a + * quadrature formula with two points in each space direction) may be + * sufficient, but on another cell where a Q3 element is used, this would lead + * to underintegration and we should use a QGauss(4) formula instead. Just as + * above, the exists a class hp::QCollection that acts as a collection of + * quadrature formulas + * + * Finally, one may want to use different orders for the boundary + * approximation for cells with different orders for the finite element. The + * hp::MappingCollection class allows to do this. + * + * All of these three classes, the hp::FECollection, hp::QCollection, and + * hp::MappingCollection classes, implement a similar interface. They have + * functions add_*() to add a finite element, quadrature formula, + * or mapping to the collection. They have a get_*(unsigned int) + * function that allows to retrieve a reference to a given element of the + * collection. And they have a n_*() function that returns the + * number of elements in the collection. Some of the classes, in particular + * that holding finite element objects, also implement other functions + * specific to their purpose. + * + * The similarity goes beyond the interface: When adding an element to the + * collection, all of the classes create a copy of the argument. This allows + * to pass a temporary object to the function adding the element. For example, + * the following works: + * @verbatim + * FECollection fe_collection; + * for (unsigned int degree=1; degree<5; ++degree) + * fe_collection.add_fe (FE_Q(degree)); + * @endverbatim + * + * This way, one can add elements of polynomial degree 1 through 4 to the + * collection. It is not necessary to retain the added object: the collection + * makes a coyp of it, it does not only store a pointer to the given finite + * element object. This same observation also holds for the other collection + * classes. + * + * It is customary that within an hp finite element program, one keeps + * collections of finite elements and quadrature formulas with the same number + * of elements, each element of the one collection matching the element in the + * other. This is not necessary, but it often makes coding a lot simpler. If a + * collection of mappings is used, the same holds for hp::MappingCollection + * objects as well. + * + * @ingroup hp + */ +