]> https://gitweb.dealii.org/ - dealii.git/commitdiff
.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 6 Mar 1998 10:52:12 +0000 (10:52 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Fri, 6 Mar 1998 10:52:12 +0000 (10:52 +0000)
git-svn-id: https://svn.dealii.org/trunk@27 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/fe/fe.h
deal.II/deal.II/include/fe/fe_lib.lagrange.h

index a8d6e10854872c3e201adb010bebf539d1969df5..0e67b6713c2776ad5a7baf79b7348a7f9e45f732 100644 (file)
@@ -20,9 +20,25 @@ template <int dim> class Quadrature;
   Represent a finite element evaluated with a specific quadrature rule.
   This class is an optimization which avoids evaluating the shape functions
   at the quadrature points each time a quadrature takes place. Rather, the
-  values and gradients (and possibly higher oder derivatives in future
+  values and gradients (and possibly higher order derivatives in future
   versions of this library) are evaluated once and for all before doing the
   quadrature itself.
+
+  Objects of this class store a multitude of different values needed to
+  do the assemblage steps on real cells rather than on the unit cell. Among
+  these values are the values and gradients of the shape functions at the
+  quadrature points on the real and the unit cell, the location of the
+  quadrature points on the real and on the unit cell, the weights of the
+  quadrature points, the Jacobian matrices of the mapping from the unit to
+  the real cell at the quadrature points and so on.
+
+  The Jacobian matrix is defined to be
+  $$ J_{ij} = {d\xi_i \over d\x_j} $$
+  which is the form needed to compute the gradient on the real cell from
+  the gradient on the unit cell. If we want to transform the area element
+  $dx dy$ from the real to the unit cell, we have to take the determinant of
+  the inverse matrix, which is the reciprocal value of the determinant of the
+  matrix defined above.
   */
 template <int dim>
 class FEValues {
@@ -123,7 +139,11 @@ class FEValues {
                                      * Store an array of weights times the
                                      * Jacobi determinant at the quadrature
                                      * points. This function is reset each time
-                                     * #reinit# is called.
+                                     * #reinit# is called. The Jacobi determinant
+                                     * is actually the reciprocal value of the
+                                     * Jacobi matrices stored in this class,
+                                     * see the general documentation of this
+                                     * class for more information.
                                      */
     vector<double>       JxW_values;
 
@@ -258,15 +278,25 @@ class FiniteElementBase {
     const dFMatrix & constraints () const;
     
                                     /**
-                                     * Compute the jacobian matrix and the
+                                     * Compute the Jacobian matrix and the
                                      * quadrature points from the given cell
                                      * and the given quadrature points on the
-                                     * unit cell. The jacobian matrix is to
+                                     * unit cell. The Jacobian matrix is to
                                      * be computed at every quadrature point.
                                      * This function has to be in the finite
                                      * element class, since different finite
                                      * elements need different transformations
                                      * of the unit cell to a real cell.
+                                     *
+                                     * Refer to the documentation of the
+                                     * \Ref{FEValues} class for a definition
+                                     * of the Jacobi matrix.
+                                     *
+                                     * It is provided for the finite element
+                                     * class in one space dimension, but for
+                                     * higher dimensions, it depends on the
+                                     * present fe and needs reimplementation
+                                     * by the user.
                                      */
     virtual void fill_fe_values (const Triangulation<dim>::cell_iterator &cell,
                                 const vector<Point<dim> >               &unit_points,
@@ -431,6 +461,33 @@ class FiniteElement<1> : public FiniteElementBase<1> {
                                      * as for the base class.
                                      */
     bool operator == (const FiniteElement<1> &f) const;
+
+                                    /**
+                                     * Compute the Jacobian matrix and the
+                                     * quadrature points from the given cell
+                                     * and the given quadrature points on the
+                                     * unit cell. The Jacobian matrix is to
+                                     * be computed at every quadrature point.
+                                     *
+                                     * Refer to the documentation of the
+                                     * \Ref{FEValues} class for a definition
+                                     * of the Jacobi matrix.
+                                     *
+                                     * For one dimensional finite elements,
+                                     * these transformations are usually the
+                                     * same, linear ones, so we provide
+                                     * them in the FE<1> base class. You may,
+                                     * however override this implementation
+                                     * if you would like to use finite elements
+                                     * of higher than first order with
+                                     * non-equidistant integration points, e.g.
+                                     * with an exponential dependence from the
+                                     * distance to the origin.
+                                     */
+    virtual void fill_fe_values (const Triangulation<1>::cell_iterator &cell,
+                                const vector<Point<1> >               &unit_points,
+                                vector<dFMatrix>  &jacobians,
+                                vector<Point<1> > &points) const;
 };
 
 
@@ -526,6 +583,31 @@ class FiniteElement<2> : public FiniteElementBase<2> {
                                      * as for the base class.
                                      */
     bool operator == (const FiniteElement<2> &f) const;
+
+                                    /**
+                                     * Compute the Jacobian matrix and the
+                                     * quadrature points from the given cell
+                                     * and the given quadrature points on the
+                                     * unit cell. The Jacobian matrix is to
+                                     * be computed at every quadrature point.
+                                     *
+                                     * Refer to the documentation of the
+                                     * \Ref{FEValues} class for a definition
+                                     * of the Jacobi matrix.
+                                     *
+                                     * For two dimensional finite elements,
+                                     * these transformations are usually
+                                     * dependent on the actual finite element,
+                                     * which is expressed by the names
+                                     * sub- and isoparametric elements. This
+                                     * function is therefore not implemented
+                                     * by the FE<2> base class, but is made
+                                     * pure virtual.
+                                     */
+    virtual void fill_fe_values (const Triangulation<2>::cell_iterator &cell,
+                                const vector<Point<2> >               &unit_points,
+                                vector<dFMatrix>  &jacobians,
+                                vector<Point<2> > &points) const;
 };
 
 
index 28731171cae2d4bd92820a499ce632ba7cbc386f..e14f219e403f194fac65527726193356fc03daab 100644 (file)
 
 
 /**
-  Define a (bi-, tri-, etc)linear finite element in #dim# space dimensions.
+  Define a (bi-, tri-, etc)linear finite element in #dim# space dimensions,
+  along with (bi-, tri-)linear (therefore isoparametric) transforms from the
+  unit cell to the real cell.
+
+  The linear, isoparametric mapping from a point $\vec \xi$ on the unit cell
+  to a point $\vec x$ on the real cell is defined as
+  $$ \vec x(\vec \xi)  = \sum_j {\vec p_j} N_j(\xi) $$
+  where $\vec p_j$ is the vector to the $j$th corner point of the cell in
+  real space and $N_j(\vec \xi)$ is the value of the basis function associated
+  with the $j$th corner point, on the unit cell at point $\vec \xi$. The sum
+  over $j$ runs over all corner points.
   */
 template <int dim>
 class FELinear : public FiniteElement<dim> {
   public:
+                                    /**
+                                     * Constructor
+                                     */
     FELinear ();
+
+                                    /**
+                                     * Return the value of the #i#th shape
+                                     * function at point #p# on the unit cell.
+                                     */
     virtual double shape_value(const unsigned int i,
                               const Point<dim>& p) const;
+
+                                    /**
+                                     * Return the gradient of the #i#th shape
+                                     * function at point #p# on the unit cell.
+                                     */
     virtual Point<dim> shape_grad(const unsigned int i,
                                  const Point<dim>& p) const;
+
+                                    /**
+                                     * Compute the Jacobian matrix and the
+                                     * quadrature points from the given cell
+                                     * and the given quadrature points on the
+                                     * unit cell. The Jacobian matrix is to
+                                     * be computed at every quadrature point.
+                                     *
+                                     * Refer to the documentation of the
+                                     * \Ref{FEValues} class for a definition
+                                     * of the Jacobi matrix.
+                                     *
+                                     * For one dimensional elements, this
+                                     * function simply passes through to
+                                     * the one implemented in the base class.
+                                     * For two dimensional finite elements,
+                                     * these transformations are usually
+                                     * dependent on the actual finite element,
+                                     * which is expressed by the names
+                                     * sub- and isoparametric elements. This
+                                     * function is therefore not implemented
+                                     * by the FE<2> base class, but is made
+                                     * pure virtual.
+                                     */
+    virtual void fill_fe_values (const Triangulation<dim>::cell_iterator &cell,
+                                const vector<Point<dim> >               &unit_points,
+                                vector<dFMatrix>    &jacobians,
+                                vector<Point<dim> > &points) const;
 };
 
 
@@ -27,15 +78,58 @@ class FELinear : public FiniteElement<dim> {
 
 /**
   Define a (bi-, tri-, etc)quadratic finite element in #dim# space dimensions.
+  In one space dimension, a linear (subparametric) mapping from the unit cell
+  to the real cell is implemented.
   */
 template <int dim>
 class FEQuadratic : public FiniteElement<dim> {
   public:
+                                    /**
+                                     * Constructor
+                                     */
     FEQuadratic ();
+
+                                    /**
+                                     * Return the value of the #i#th shape
+                                     * function at point #p# on the unit cell.
+                                     */
     virtual double shape_value(const unsigned int i,
                               const Point<dim>& p) const;
+
+                                    /**
+                                     * Return the gradient of the #i#th shape
+                                     * function at point #p# on the unit cell.
+                                     */
     virtual Point<dim> shape_grad(const unsigned int i,
                                  const Point<dim>& p) const;
+
+                                    /**
+                                     * Compute the Jacobian matrix and the
+                                     * quadrature points from the given cell
+                                     * and the given quadrature points on the
+                                     * unit cell. The Jacobian matrix is to
+                                     * be computed at every quadrature point.
+                                     *
+                                     * Refer to the documentation of the
+                                     * \Ref{FEValues} class for a definition
+                                     * of the Jacobi matrix.
+                                     *
+                                     * For one dimensional elements, this
+                                     * function simply passes through to
+                                     * the one implemented in the base class.
+                                     * For two dimensional finite elements,
+                                     * these transformations are usually
+                                     * dependent on the actual finite element,
+                                     * which is expressed by the names
+                                     * sub- and isoparametric elements. This
+                                     * function is therefore not implemented
+                                     * by the FE<2> base class, but is made
+                                     * pure virtual.
+                                     */
+    virtual void fill_fe_values (const Triangulation<dim>::cell_iterator &cell,
+                                const vector<Point<dim> >               &unit_points,
+                                vector<dFMatrix>    &jacobians,
+                                vector<Point<dim> > &points) const;
 };
 
 
@@ -43,15 +137,58 @@ class FEQuadratic : public FiniteElement<dim> {
 
 /**
   Define a (bi-, tri-, etc)cubic finite element in #dim# space dimensions.
+  In one space dimension, a linear (subparametric) mapping from the unit cell
+  to the real cell is implemented.
   */
 template <int dim>
 class FECubic : public FiniteElement<dim> {
   public:
+                                    /**
+                                     * Constructor
+                                     */
     FECubic ();
+
+                                    /**
+                                     * Return the value of the #i#th shape
+                                     * function at point #p# on the unit cell.
+                                     */
     virtual double shape_value(const unsigned int i,
                               const Point<dim>& p) const;
+
+                                    /**
+                                     * Return the gradient of the #i#th shape
+                                     * function at point #p# on the unit cell.
+                                     */
     virtual Point<dim> shape_grad(const unsigned int i,
                                  const Point<dim>& p) const;
+
+                                    /**
+                                     * Compute the Jacobian matrix and the
+                                     * quadrature points from the given cell
+                                     * and the given quadrature points on the
+                                     * unit cell. The Jacobian matrix is to
+                                     * be computed at every quadrature point.
+                                     *
+                                     * Refer to the documentation of the
+                                     * \Ref{FEValues} class for a definition
+                                     * of the Jacobi matrix.
+                                     *
+                                     * For one dimensional elements, this
+                                     * function simply passes through to
+                                     * the one implemented in the base class.
+                                     * For two dimensional finite elements,
+                                     * these transformations are usually
+                                     * dependent on the actual finite element,
+                                     * which is expressed by the names
+                                     * sub- and isoparametric elements. This
+                                     * function is therefore not implemented
+                                     * by the FE<2> base class, but is made
+                                     * pure virtual.
+                                     */
+    virtual void fill_fe_values (const Triangulation<dim>::cell_iterator &cell,
+                                const vector<Point<dim> >               &unit_points,
+                                vector<dFMatrix>    &jacobians,
+                                vector<Point<dim> > &points) const;
 };
 
 

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.