From: Guido Kanschat Date: Mon, 24 Oct 2005 04:33:03 +0000 (+0000) Subject: attempt on documentation X-Git-Tag: v8.0.0~12951 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f8a8f0efcfedf18176f1ace27291fc7c2fe48386;p=dealii.git attempt on documentation git-svn-id: https://svn.dealii.org/trunk@11652 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/fe/fe_update_flags.h b/deal.II/deal.II/include/fe/fe_update_flags.h index 8426632196..ed2d3efd85 100644 --- a/deal.II/deal.II/include/fe/fe_update_flags.h +++ b/deal.II/deal.II/include/fe/fe_update_flags.h @@ -21,98 +21,184 @@ /*@{*/ /** - * Provide a set of flags which tells the FEValues::reinit - * function, which fields are to be updated for each cell. E.g. if you - * do not need the gradients since you want to assemble the mass - * matrix, you can switch that off. By default, all flags are off, - * i.e. no reinitialization will be done. + * The enum type given to the constructors of FEValues, FEFaceValues + * and FESubfaceValues, telling those objects which data will be + * needed on each mesh cell. * - * A variable of this type has to be passed to the constructor of the - * FEValues object. You can select more than one flag by concatenation - * using the | (bitwise or) operator. + * Selecting these flags in a restrictive way is crucial for the + * efficiency of FEValues::reinit(), FEFaceValues::reinit() and + * FESubfaceValues::reinit(). Therefore, only the flags actually + * needed should be selected. It is the responsibility of the involved + * Mapping and FiniteElement to add additional flags according to + * their own requirements. For instance, most finite elements will + * add #update_covariant_transformation if #update_gradients is + * selected. + + * By default, all flags are off, i.e. no reinitialization will be + * done. + * + * You can select more than one flag by concatenation + * using the bitwise or operator|(UpdateFlags,UpdateFlags). + * + *

Generating the actual flags

+ * + * When given a set of UpdateFlags @pflags, the FEValues object must + * determine, which values will have to be computed once only for the + * reference cell and which values will have to be updated for each + * cell. Here, it is important to note that in many cases, the + * FiniteElement will require additional updates from the Mapping. To + * this end, several auxiliary functions have been implemented: * + * FiniteElement::update_once(flags) and + * FiniteElement::update_each(flags) determine the values required by + * the FiniteElement once or on each cell. The same functions exist in Mapping. * - *

Description of Flags

+ * Since the FiniteElement does not know if a value reuired from + * Mapping should be computed once or for each cell, + * FEValuesBase::compute_update_flags() is used to compute the union + * of all values to be computed ever. It does this by first adding to + * the flags set by the user all flags (once and each) added by the + * FiniteElement. This new set of flags is then given to the Mapping + * and all flags required there are added, again once and each. * - * The following flags are declared: - * + * This union of all flags is given to Mapping::fill_fe_values() and + * FiniteElement::fill_fe_values, where it is split again into the + * information generated only once and the information that must be + * updated on each cell. * - * @author Wolfgang Bangerth, Guido Kanschat, 1998, 1999, 2000, 2001, Ralf Hartmann 2004 + * The flags finally stored in FEValues then are the union of all the + * flags required by the user, by FiniteElement and by Mapping, for + * computation once or on each cell. Subsequent calls to the functions + * @p update_once and @p update_each should just select among these + * flags, but should not add new flags. */ enum UpdateFlags { //! No update update_default = 0, //! Shape function values + /** + * Compute the values of the + * shape functions at the + * quadrature points on the + * real space cell. For the + * usual Lagrange elements, + * these values are equal to + * the values of the shape + * functions at the quadrature + * points on the unit cell, but + * they are different for more + * complicated elements, such + * as FE_RaviartThomas + * elements. + */ update_values = 0x0001, //! Shape function gradients + /** + * Compute the gradients of the + * shape functions in + * coordinates of the real + * cell. + */ update_gradients = 0x0002, //! Second derivatives of shape functions + /** + * Compute the second + * derivatives of the shape + * functions in coordinates of + * the real cell. + */ update_second_derivatives = 0x0004, - /*! Normal vector times surface + //! Outter normal vector, not normalized + /** + * Vector product of tangential + * vectors, yielding a normal + * vector with a length + * corresponding to the surface * element; may be more - * efficient than computing both. + * efficient than computing + * both. */ update_boundary_forms = 0x0008, //! Transformed quadrature points + /** + * Compute the quadrature + * points transformed into real + * cell coordinates. + */ update_q_points = 0x0010, //! Transformed quadrature weights + /** + * Compute the quadrature + * weights on the real cell, + * i.e. the weights of the + * quadrature rule multiplied + * with the determinant of the + * Jacoian of the + * transformation from + * reference to realcell. + */ update_JxW_values = 0x0020, - //! Transformed normal vectors + //! Normal vectors + /** + * Compute the unit outer + * normal vector on the face of + * a cell. + */ update_normal_vectors = 0x0040, //! Volume element + /** + * Compute the Jacobian of the + * transformation from the + * reference cell to the real + * cell. + */ update_jacobians = 0x0080, //! Gradient of volume element + /** + * Compute the dervatives of + * the Jacobian of the + * transformation. + */ update_jacobian_grads = 0x0100, //! Covariant transformation + /** + * Compute all values the + * Mapping needs to perform a + * contravariant transformation of + * vectors. For special + * mappings like + * MappingCartesian this may be + * simpler than + * #update_jacobians. + */ update_covariant_transformation = 0x0200, //! Contravariant transformation + /** + * Compute all values the + * Mapping needs to perform a + * contravariant transformation of + * vectors. For special + * mappings like + * MappingCartesian this may be + * simpler than + * #update_jacobians. + */ update_contravariant_transformation = 0x0400, //! Shape function values of transformation + /** + * Compute the shape function + * values of the transformation + * defined by the Mapping. + */ update_transformation_values = 0x0800, //! Shape function gradients of transformation + /** + * Compute the shape function + * gradients of the + * transformation defined by + * the Mapping. + */ update_transformation_gradients = 0x1000 };