-
-
/**
- 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 order derivatives in future
- versions of this library) are evaluated once and for all on the unit
- cell before doing the quadrature itself. Only the Jacobian matrix of
- the transformation from the unit cell to the real cell and the integration
- points in real space are calculated each time we move on to a new cell.
-
- The unit cell is defined to be the tensor product of the interval $[0,1]$
- in the present number of dimensions. In part of the literature, the convention
- is used that the unit cell be the tensor product of the interval $[-1,1]$,
- which is to distinguished properly.
-
- 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 dx_j} $$
- where the $\xi_i$ are the coordinates on the unit cell and the $x_i$ are
- the coordinates on the real cell.
- This 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.
-
- The #FEValues# object keeps track of those fields which really need to
- be computed, since the computation of the gradients of the ansatz functions
- on each real cell can be quite an expensive thing if it is not needed. The
- object knows about which fields are needed by the #UpdateFlags# object
- passed through the constructor. In debug mode, the accessor functions, which
- return values from the different fields, check whether the required field
- was initialized, thus avoiding use of unitialized data.
-
-
- {\bf Member functions}
-
- The functions of this class fall into different cathegories:
- \begin{itemize}
- \item #shape_value#, #shape_grad#, etc: return one of the values
- of this object at a time. In many cases you will want to get
- a whole bunch at a time for performance or convenience reasons,
- then use the #get_*# functions.
-
- \item #get_shape_values#, #get_shape_grads#, etc: these return
- a reference to a whole field. Usually these fields contain
- the values of all ansatz functions at all quadrature points.
-
- \item #get_function_values#, #get_function_gradients#: these
- two functions offer a simple way to avoid the detour of the
- ansatz functions, if you have a finite solution (resp. the
- vector of values associated with the different ansatz functions.)
- Then you may want to get information from the restriction of
- the finite element function to a certain cell, e.g. the values
- of the function at the quadrature points or the values of its
- gradient. These two functions provide the information needed:
- you pass it a vector holding the finite element solution and the
- functions return the values or gradients of the finite element
- function restricted to the cell which was given last time the
- #reinit# function was given.
-
- Though possible in principle, these functions do not call the
- #reinit# function, you have to do so yourself beforehand. On the
- other hand, a copy of the cell iterator is stored which was used
- last time the #reinit# function was called. This frees us from
- the need to pass the cell iterator again to these two functions,
- which guarantees that the cell used here is in sync with that used
- for the #reinit# function. You should, however, make sure that
- nothing substantial happens to the #DoFHandler# object or any
- other involved instance between the #reinit# and the #get_function_*#
- functions are called.
-
- \item #reinit#: initialize the #FEValues# object for a certain cell.
- See above for more information.
- \end{itemize}
- */
+ This class offers a multitude of arrays and other fields which are used by
+ the derived classes #FEValues# and #FEFaceValues#. In principle, it is the
+ back end of the front end for the unification of a certain finite element
+ and a quadrature formula which evaluates certain aspects of the finite
+ element at quadrature points.
+
+ 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 order derivatives in future
+ versions of this library) are evaluated once and for all on the unit
+ cell or face before doing the quadrature itself. Only the Jacobian matrix of
+ the transformation from the unit cell or face to the real cell or face and
+ the integration points in real space are calculated each time we move on
+ to a new face.
+
+ Actually, this class does none of the evaluations at startup itself; this is
+ all done by the derived classes. It only offers the basic functionality,
+ like providing those fields that are common to the derived classes and
+ access to these fields. Any computations are in the derived classes. See there
+ for more information.
+
+ @author Wolfgang Bangerth, 1998
+*/
template <int dim>
-class FEValues {
+class FEValuesBase {
public:
-
-
-
+
/**
* Number of quadrature points.
*/
const unsigned int n_quadrature_points;
/**
- * Total number of shape functions.
+ * Total number of shape functions
+ * per cell. If we use this base class
+ * to evaluate a finite element on
+ * faces of cells, this is still the
+ * number of degrees of freedom per
+ * cell, not per face.
*/
const unsigned int total_dofs;
-
- /**
- * Constructor. Fill all arrays with the
- * values of the shape functions of the
- * specified finite element using the
- * quadrature points of the given
- * quadrature rule.
- *
- * This function actually only fills
- * the fields related to the unit face,
- * the fields related to a real face (like
- * gradients, true quadrature points, etc.)
- * need to be initialized using the
- * #reinit# function.
- */
- FEValues (const FiniteElement<dim> &,
- const Quadrature<dim> &,
- const UpdateFlags);
-
- /**
- * Return the value of the #i#th shape
- * function at the #j# quadrature point.
- */
- double shape_value (const unsigned int i,
- const unsigned int j) const;
/**
- * Return a pointer to the matrix holding
- * all values of shape functions at all
- * integration points, on the present cell.
- * For the format of this matrix, see the
- * documentation for the matrix itself.
- */
- const dFMatrix & get_shape_values () const;
+ * Constructor. Set up the array sizes
+ * with #n_q_points# quadrature points,
+ * #n_ansatz_points# ansatz points (on
+ * the cell or face), #n_dof# ansatz
+ * functions per cell and with the
+ * given pattern to update the fields
+ * when the #reinit# function of the
+ * derived classes is called. The
+ * fields themselves are not set up,
+ * this must happen in the derived
+ * class's constructor, only the sizes
+ * are set correctly.
+ */
+ FEValuesBase (const unsigned int n_q_points,
+ const unsigned int n_ansatz_points,
+ const unsigned int n_dofs,
+ const UpdateFlags update_flags);
- /**
- * Return the values of the finite
- * element function characterized
- * by #fe_function# restricted to
- * #cell# at the quadrature points.
- *
- * The function assumes that the
- * #values# object already has the
- * right size.
- */
- void get_function_values (const dVector &fe_function,
- vector<double> &values) const;
-
- /**
+ /**
* Return the gradient of the #i#th shape
* function at the #j# quadrature point.
* If you want to get the derivative in
/**
* Return the position of the #i#th
* quadrature point in real space.
+ *
+ * If this object is used to evaluate
+ * finite elements on faces of cells,
+ * and for curved boundary cells, using
+ * biquadratic or higher mappings
+ * of the unit cell to the real cell,
+ * these points may not be on the
+ * plane submannifold on which the
+ * vertices of the face lie.
*/
const Point<dim> & quadrature_point (const unsigned int i) const;
* we have to take the continuous
* function's value at the ansatz function
* locations.
+ *
+ * For the evaluation of finite elements on
+ * faces of cells, #i# is the number
+ * of the ansatz function on the face, not
+ * on the cell.
*/
const Point<dim> & ansatz_point (const unsigned int i) const;
* Return the Jacobi determinant times
* the weight of the #i#th quadrature
* point.
+ *
+ * If faces of cells are concerned,
+ * the jacobi determinant is that of the
+ * transformation of the unit face to
+ * the present face, not of the unit
+ * cell to the real cell (unlike for
+ * the #jacobi_matrix# array of the
+ * derived classes which store the cell
+ * transformation's Jacobi matrix in
+ * all cases).
*/
double JxW (const unsigned int i) const;
* quadrature points.
*/
const vector<double> & get_JxW_values () const;
-
- /**
- * Reinitialize the gradients, Jacobi
- * determinants, etc for the given cell
- * and the given finite element.
- *
- * This function needs a boundary object
- * passed, since this class needs to know
- * how to handle faces which are located
- * on the boundary of the domain. In that
- * case, faces may be curved and the
- * calculation of quadrature points,
- * gradients and the like may need
- * additional effort, depending on the
- * mapping from the unit to the real cell
- * (linear mappings use straight boundary
- * segments, but higher order elements
- * may use other ways.)
- */
- void reinit (const typename DoFHandler<dim>::cell_iterator &,
- const FiniteElement<dim> &,
- const Boundary<dim> &);
/**
* Exception
int, int,
<< "Vector has wrong size " << arg1
<< ", expected size " << arg2);
-
- private:
/**
- * Store the values of the shape functions
- * at the quadrature points. Rows in this
- * matrix denote the values of a single
- * shape function at the different points,
- * columns are for a single point with the
- * different shape functions.
+ * Exception
*/
- dFMatrix shape_values;
-
+ DeclException0 (ExcInternalError);
/**
+ * Exception
+ */
+ DeclException0 (ExcNotImplemented);
+
+ protected:
+ /**
* Store the gradients of the shape
* functions at the quadrature points.
* Since unfortunately the full matrix
*/
vector<vector<Point<dim> > > shape_gradients;
- /**
- * Store the gradients of the shape
- * functions at the quadrature points on
- * the unit cell.
- * This field is set up upon construction
- * of the object and contains the gradients
- * on the reference element.
- */
- vector<vector<Point<dim> > > unit_shape_gradients;
-
/**
* Store an array of the weights of the
* quadrature points. This array is
* set up upon construction.
+ *
+ * If faces rather than cells are
+ * considered, the weights are stored
+ * only once still, since they are
+ * not transformed and are thus the same
+ * for all faces.
*/
vector<double> weights;
*/
vector<Point<dim> > quadrature_points;
- /**
- * Array of quadrature points in the unit
- * cell. This array is set up upon
- * construction and contains the quadrature
- * points on the reference element.
- */
- vector<Point<dim> > unit_quadrature_points;
-
- /**
+ /**
* Array of points denoting the off-point
* of the ansatz functions. In real space
* (no-one seems to need the off-point
* Store the jacobi matrices at the
* different quadrature points. This field
* is set each time #reinit# is called.
+ *
+ * If faces rather than cells are considered
+ * this is the Jacobi matrix of the
+ * transformation of the unit cell to the
+ * real cell, not of the unit face to the
+ * face. We need this full matrix for the
+ * transformation of the gradients to the
+ * real cell.
*/
vector<dFMatrix> jacobi_matrices;
+/**
+ Represent a finite element evaluated with a specific quadrature rule on
+ a cell.
+
+ The unit cell is defined to be the tensor product of the interval $[0,1]$
+ in the present number of dimensions. In part of the literature, the convention
+ is used that the unit cell be the tensor product of the interval $[-1,1]$,
+ which is to distinguished properly.
+
+ 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 dx_j} $$
+ where the $\xi_i$ are the coordinates on the unit cell and the $x_i$ are
+ the coordinates on the real cell.
+ This 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.
+
+ The #FEValues# object keeps track of those fields which really need to
+ be computed, since the computation of the gradients of the ansatz functions
+ on each real cell can be quite an expensive thing if it is not needed. The
+ object knows about which fields are needed by the #UpdateFlags# object
+ passed through the constructor. In debug mode, the accessor functions, which
+ return values from the different fields, check whether the required field
+ was initialized, thus avoiding use of unitialized data.
+
+
+ {\bf Member functions}
+
+ The functions of this class fall into different cathegories:
+ \begin{itemize}
+ \item #shape_value#, #shape_grad#, etc: return one of the values
+ of this object at a time. In many cases you will want to get
+ a whole bunch at a time for performance or convenience reasons,
+ then use the #get_*# functions.
+
+ \item #get_shape_values#, #get_shape_grads#, etc: these return
+ a reference to a whole field. Usually these fields contain
+ the values of all ansatz functions at all quadrature points.
+
+ \item #get_function_values#, #get_function_gradients#: these
+ two functions offer a simple way to avoid the detour of the
+ ansatz functions, if you have a finite solution (resp. the
+ vector of values associated with the different ansatz functions.)
+ Then you may want to get information from the restriction of
+ the finite element function to a certain cell, e.g. the values
+ of the function at the quadrature points or the values of its
+ gradient. These two functions provide the information needed:
+ you pass it a vector holding the finite element solution and the
+ functions return the values or gradients of the finite element
+ function restricted to the cell which was given last time the
+ #reinit# function was given.
+
+ Though possible in principle, these functions do not call the
+ #reinit# function, you have to do so yourself beforehand. On the
+ other hand, a copy of the cell iterator is stored which was used
+ last time the #reinit# function was called. This frees us from
+ the need to pass the cell iterator again to these two functions,
+ which guarantees that the cell used here is in sync with that used
+ for the #reinit# function. You should, however, make sure that
+ nothing substantial happens to the #DoFHandler# object or any
+ other involved instance between the #reinit# and the #get_function_*#
+ functions are called.
+
+ \item #reinit#: initialize the #FEValues# object for a certain cell.
+ See above for more information.
+ \end{itemize}
+
+ @author Wolfgang Bangerth, 1998
+ */
+template <int dim>
+class FEValues : public FEValuesBase<dim> {
+ public:
+
+
+
+ /**
+ * Constructor. Fill all arrays with the
+ * values of the shape functions of the
+ * specified finite element using the
+ * quadrature points of the given
+ * quadrature rule.
+ *
+ * This function actually only fills
+ * the fields related to the unit face,
+ * the fields related to a real face (like
+ * gradients, true quadrature points, etc.)
+ * need to be initialized using the
+ * #reinit# function.
+ */
+ FEValues (const FiniteElement<dim> &,
+ const Quadrature<dim> &,
+ const UpdateFlags);
+
+ /**
+ * Return the value of the #i#th shape
+ * function at the #j# quadrature point.
+ */
+ double shape_value (const unsigned int i,
+ const unsigned int j) const;
+
+ /**
+ * Return a pointer to the matrix holding
+ * all values of shape functions at all
+ * integration points, on the present cell.
+ * For the format of this matrix, see the
+ * documentation for the matrix itself.
+ */
+ const dFMatrix & get_shape_values () const;
+
+ /**
+ * Return the values of the finite
+ * element function characterized
+ * by #fe_function# restricted to
+ * #cell# at the quadrature points.
+ *
+ * The function assumes that the
+ * #values# object already has the
+ * right size.
+ */
+ void get_function_values (const dVector &fe_function,
+ vector<double> &values) const;
+
+
+ /**
+ * Reinitialize the gradients, Jacobi
+ * determinants, etc for the given cell
+ * and the given finite element.
+ *
+ * This function needs a boundary object
+ * passed, since this class needs to know
+ * how to handle faces which are located
+ * on the boundary of the domain. In that
+ * case, faces may be curved and the
+ * calculation of quadrature points,
+ * gradients and the like may need
+ * additional effort, depending on the
+ * mapping from the unit to the real cell
+ * (linear mappings use straight boundary
+ * segments, but higher order elements
+ * may use other ways.)
+ */
+ void reinit (const typename DoFHandler<dim>::cell_iterator &,
+ const FiniteElement<dim> &,
+ const Boundary<dim> &);
+
+ private:
+ /**
+ * Store the values of the shape functions
+ * at the quadrature points. Rows in this
+ * matrix denote the values of a single
+ * shape function at the different points,
+ * columns are for a single point with the
+ * different shape functions.
+ */
+ dFMatrix shape_values;
+
+ /**
+ * Store the gradients of the shape
+ * functions at the quadrature points on
+ * the unit cell.
+ * This field is set up upon construction
+ * of the object and contains the gradients
+ * on the reference element.
+ */
+ vector<vector<Point<dim> > > unit_shape_gradients;
+
+
+ /**
+ * Array of quadrature points in the unit
+ * cell. This array is set up upon
+ * construction and contains the quadrature
+ * points on the reference element.
+ */
+ vector<Point<dim> > unit_quadrature_points;
+
+};
+
+
+
/**
- 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 order derivatives in future
- versions of this library) are evaluated once and for all on the unit
- face before doing the quadrature itself. Only the Jacobian matrix of
- the transformation from the unit face to the real face and the integration
- points in real space are calculated each time we move on to a new face.
+ Represent a finite element evaluated with a specific quadrature rule on
+ the face of a cell.
The unit face is defined to be the tensor product of the interval $[0,1]$
in the present number of dimensions minus one. In part of the literature,
Finally, we will often need the outward normal to a cell at the quadrature
points. While this could in principle be easily done using the Jacobi
- matrices at the quadrature points and the normal vectors to the unit cell
- (also easily derived, since they have an appealingly easy form for the unit
- cell ;-), it is more efficiently done by the finite element class itself.
- For example for (bi-, tri-)linear mappings the normal vector is readily
- available without compicated matrix-vector-multiplications.
- */
-template <int dim>
-class FEFaceValues {
- public:
-
-
-
- /**
- * Number of quadrature points on
- * the face.
- */
- const unsigned int n_quadrature_points;
-
- /**
- * Total number of shape functions
- * on the cell adjacent to this face.
- * This number is not the same as the
- * number of shape functions of which
- * the center is located on the face.
- */
- const unsigned int total_dofs;
-
+ matrices at the quadrature points and the normal vectors to the unit cell
+ (also easily derived, since they have an appealingly easy form for the unit
+ cell ;-), it is more efficiently done by the finite element class itself.
+ For example for (bi-, tri-)linear mappings the normal vector is readily
+ available without compicated matrix-vector-multiplications.
+ */
+template <int dim>
+class FEFaceValues : public FEValuesBase<dim> {
+ public:
/**
* Constructor. Fill all arrays with the
* values of the shape functions of the
void get_function_values (const dVector &fe_function,
vector<double> &values) const;
- /**
- * Return the gradient of the #i#th shape
- * function at the #j# quadrature point.
- * If you want to get the derivative in
- * one of the coordinate directions, use
- * the appropriate function of the #Point#
- * class to extract one component. Since
- * only a reference to the gradient's value
- * is returned, there should be no major
- * performance drawback.
- * The function returns the gradient on the
- * real element, not the reference element.
- */
- const Point<dim> & shape_grad (const unsigned int i,
- const unsigned int j) const;
-
- /**
- * Return a pointer to the matrix holding
- * all gradients of shape functions at all
- * integration points, on the present cell.
- * For the format of this matrix, see the
- * documentation for the matrix itself.
- */
- const vector<vector<Point<dim> > > & get_shape_grads () const;
-
- /**
- * Return the gradients of the finite
- * element function characterized
- * by #fe_function# restricted to
- * #cell# at the quadrature points.
- *
- * The function assumes that the
- * #gradients# object already has the
- * right size.
- */
- void get_function_grads (const dVector &fe_function,
- vector<Point<dim> > &gradients) const;
-
- /**
- * Return the position of the #i#th
- * quadrature point in real space.
- *
- * For curved boundary cells, using
- * biquadratic or higher mappings
- * of the unit cell to the real cell,
- * these points may not be on the
- * plane submannifold on which the
- * vertices of the face lie.
- */
- const Point<dim> & quadrature_point (const unsigned int i) const;
-
- /**
- * Return a pointer to the vector of
- * quadrature points.
- */
- const vector<Point<dim> > & get_quadrature_points () const;
-
- /**
- * Return the point in real space where
- * the #i#th ansatz function is located
- * (location is in the sense of where it
- * assumes its nominal properties, e.g. at
- * the vertex of a cell, at the center of
- * a line, etc).
- *
- * This function is needed for the
- * interpolation problem: if we want to
- * transfer a continuous function to a
- * finite element function by interpolation
- * we have to take the continuous
- * function's value at the ansatz function
- * locations.
- */
- const Point<dim> & ansatz_point (const unsigned int i) const;
-
- /**
- * Return a pointer to the vector of points
- * denoting the location of the ansatz
- * functions.
- */
- const vector<Point<dim> > & get_ansatz_points () const;
-
- /**
- * Return the Jacobi determinant times
- * the weight of the #i#th quadrature
- * point. The Jacobi determinant is that
- * of the transformation of the unit
- * face to the real face, not of the
- * alike cells.
- */
- double JxW (const unsigned int i) const;
-
- /**
- * Return a pointer to the array holding
- * the JxW values at the different
- * quadrature points. The Jacobi
- * determinant is that
- * of the transformation of the unit
- * face to the real face, not of the
- * alike cells.
- */
- const vector<double> & get_JxW_values () const;
-
/**
* Return the outward normal vector to
* the cell at the #i#th quadrature
const FiniteElement<dim> &fe,
const Boundary<dim> &boundary);
- /**
- * Exception
- */
- DeclException2 (ExcInvalidIndex,
- int, int,
- << "The index " << arg1
- << " is out of range, it should be less than " << arg2);
- /**
- * Exception
- */
- DeclException0 (ExcAccessToUninitializedField);
- /**
- * Exception
- */
- DeclException0 (ExcCannotInitializeField);
- /**
- * Exception
- */
- DeclException0 (ExcInternalError);
- /**
- * Exception
- */
- DeclException0 (ExcNotImplemented);
- /**
- * Exception
- */
- DeclException2 (ExcWrongVectorSize,
- int, int,
- << "Vector has wrong size " << arg1
- << ", expected size " << arg2);
private:
/**
* Store the values of the shape functions
*/
dFMatrix shape_values[2*dim];
- /**
- * Store the gradients of the shape
- * functions at the quadrature points.
- * Since unfortunately the full matrix
- * classes of DEAL are not templated,
- * we have to store them in an
- * archetypic style.
- *
- * This field is reset each time
- * #reinit# is called and contains the
- * gradients on the real element, rather
- * than on the reference element. This
- * function does the transformation from
- * the unit cell to the real cell using
- * the #unit_shape_gradients# for the
- * selected face.
- */
- vector<vector<Point<dim> > > shape_gradients;
-
/**
* Store the gradients of the shape
* functions at the quadrature points on
*/
vector<vector<Point<dim> > > unit_shape_gradients[2*dim];
- /**
- * Store an array of the weights of the
- * quadrature points. This array is
- * set up upon construction.
- *
- * Since these weights are not transformed
- * they are the same for all faces.
- */
- vector<double> weights;
-
- /**
- * Store an array of weights times the
- * Jacobi determinant at the quadrature
- * points. This function is reset each time
- * #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;
-
- /**
- * Array of quadrature points. This array
- * is set up upon calling #reinit# and
- * contains the quadrature points on the
- * real element, rather than on the
- * reference element.
- */
- vector<Point<dim> > quadrature_points;
/**
* Array of quadrature points on the
*/
vector<Point<dim> > global_unit_quadrature_points[2*dim];
- /**
- * Array of points denoting the off-point
- * of the ansatz functions. In real space
- * (no-one seems to need the off-point
- * on the unit cell, so no function is
- * provided for this).
- */
- vector<Point<dim> > ansatz_points;
-
- /**
- * Store the jacobi matrices at the
- * different quadrature points. This field
- * is set each time #reinit# is called.
- * This is the Jacobi matrix of the
- * transformation of the unit cell to the
- * real cell, not of the unit face to the
- * face. We need this full matrix for the
- * transformation of the gradients to the
- * real cell.
- */
- vector<dFMatrix> jacobi_matrices;
-
/**
* List of values denoting the determinant
* of the transformation from the unit face
- * to the real face. Neede to actually
+ * to the real face. Needed to actually
* compute the JxW values.
*/
vector<double> face_jacobi_determinants;
*/
vector<Point<dim> > normal_vectors;
- /**
- * Store which fields are to be updated by
- * the reinit function.
- */
- UpdateFlags update_flags;
-
- /**
- * Store the cell selected last time
- * the #reinit# function was called
- * to make access
- * to the #get_function_*# functions
- * safer.
- */
- DoFHandler<dim>::cell_iterator present_cell;
-
/**
* Store the number of the face selected
* last time the #reinit# function was
-/*------------------------ Inline functions: FEValues ----------------------------*/
-
-
-
-template <int dim>
-inline
-const dFMatrix & FEValues<dim>::get_shape_values () const {
- return shape_values;
-};
+/*------------------------ Inline functions: FEValuesBase ------------------------*/
template <int dim>
inline
const vector<vector<Point<dim> > > &
-FEValues<dim>::get_shape_grads () const {
+FEValuesBase<dim>::get_shape_grads () const {
Assert (update_flags & update_gradients, ExcAccessToUninitializedField());
return shape_gradients;
};
template <int dim>
inline
const vector<Point<dim> > &
-FEValues<dim>::get_quadrature_points () const {
+FEValuesBase<dim>::get_quadrature_points () const {
Assert (update_flags & update_q_points, ExcAccessToUninitializedField());
return quadrature_points;
};
template <int dim>
inline
const vector<Point<dim> > &
-FEValues<dim>::get_ansatz_points () const {
+FEValuesBase<dim>::get_ansatz_points () const {
Assert (update_flags & update_ansatz_points, ExcAccessToUninitializedField());
return ansatz_points;
};
template <int dim>
inline
const vector<double> &
-FEValues<dim>::get_JxW_values () const {
+FEValuesBase<dim>::get_JxW_values () const {
Assert (update_flags & update_JxW_values, ExcAccessToUninitializedField());
return JxW_values;
};
-
-/*------------------------ Inline functions: FEFaceValues ------------------------*/
-
-
-template <int dim>
-inline
-const dFMatrix & FEFaceValues<dim>::get_shape_values () const {
- return shape_values[selected_face];
-};
-
-
-
-
-template <int dim>
-inline
-const vector<vector<Point<dim> > > &
-FEFaceValues<dim>::get_shape_grads () const {
- Assert (update_flags & update_gradients, ExcAccessToUninitializedField());
- return shape_gradients;
-};
-
+/*------------------------ Inline functions: FEValues ----------------------------*/
template <int dim>
inline
-const vector<Point<dim> > &
-FEFaceValues<dim>::get_quadrature_points () const {
- Assert (update_flags & update_q_points, ExcAccessToUninitializedField());
- return quadrature_points;
+const dFMatrix & FEValues<dim>::get_shape_values () const {
+ return shape_values;
};
-template <int dim>
-inline
-const vector<Point<dim> > &
-FEFaceValues<dim>::get_ansatz_points () const {
- Assert (update_flags & update_ansatz_points, ExcAccessToUninitializedField());
- return ansatz_points;
-};
-
+/*------------------------ Inline functions: FEFaceValues ------------------------*/
template <int dim>
inline
-const vector<double> &
-FEFaceValues<dim>::get_JxW_values () const {
- Assert (update_flags & update_JxW_values, ExcAccessToUninitializedField());
- return JxW_values;
+const dFMatrix & FEFaceValues<dim>::get_shape_values () const {
+ return shape_values[selected_face];
};
#include <grid/tria_iterator.h>
#include <grid/tria_accessor.h>
#include <grid/tria_boundary.h>
+#include <grid/dof_accessor.h>
-/*------------------------------- FEValues -------------------------------*/
-
-
-template <int dim>
-FEValues<dim>::FEValues (const FiniteElement<dim> &fe,
- const Quadrature<dim> &quadrature,
- const UpdateFlags update_flags) :
- n_quadrature_points(quadrature.n_quadrature_points),
- total_dofs(fe.total_dofs),
- shape_values(fe.total_dofs, quadrature.n_quadrature_points),
- shape_gradients(fe.total_dofs,
- vector<Point<dim> >(quadrature.n_quadrature_points)),
- unit_shape_gradients(fe.total_dofs,
- vector<Point<dim> >(quadrature.n_quadrature_points)),
- weights(quadrature.n_quadrature_points, 0),
- JxW_values(quadrature.n_quadrature_points, 0),
- quadrature_points(quadrature.n_quadrature_points, Point<dim>()),
- unit_quadrature_points(quadrature.get_quad_points()),
- ansatz_points (fe.total_dofs, Point<dim>()),
- jacobi_matrices (quadrature.n_quadrature_points,
- dFMatrix(dim,dim)),
- update_flags (update_flags)
-{
- for (unsigned int i=0; i<fe.total_dofs; ++i)
- for (unsigned int j=0; j<n_quadrature_points; ++j)
- {
- shape_values(i,j) = fe.shape_value(i, quadrature.quad_point(j));
- unit_shape_gradients[i][j]
- = fe.shape_grad(i, quadrature.quad_point(j));
- };
-
- for (unsigned int i=0; i<n_quadrature_points; ++i)
- {
- weights[i] = quadrature.weight(i);
- };
-};
-
-
-
-template <int dim>
-double FEValues<dim>::shape_value (const unsigned int i,
- const unsigned int j) const {
- Assert (i<shape_values.m(), ExcInvalidIndex (i, shape_values.m()));
- Assert (j<shape_values.n(), ExcInvalidIndex (j, shape_values.n()));
-
- return shape_values(i,j);
-};
-
+/*------------------------------- FEValuesBase ---------------------------*/
template <int dim>
-void FEValues<dim>::get_function_values (const dVector &fe_function,
- vector<double> &values) const {
- Assert (values.size() == n_quadrature_points,
- ExcWrongVectorSize(values.size(), n_quadrature_points));
-
- // get function values of dofs
- // on this cell
- vector<double> dof_values (total_dofs, 0);
- present_cell->get_dof_values (fe_function, dof_values);
-
- // initialize with zero
- fill_n (values.begin(), n_quadrature_points, 0);
+FEValuesBase<dim>::FEValuesBase (const unsigned int n_q_points,
+ const unsigned int n_ansatz_points,
+ const unsigned int n_dofs,
+ const UpdateFlags update_flags) :
+ n_quadrature_points (n_q_points),
+ total_dofs (n_dofs),
+ shape_gradients (n_dofs, vector<Point<dim> >(n_q_points)),
+ weights (n_q_points, 0),
+ JxW_values (n_q_points, 0),
+ quadrature_points (n_q_points, Point<dim>()),
+ ansatz_points (n_ansatz_points, Point<dim>()),
+ jacobi_matrices (n_q_points, dFMatrix(dim,dim)),
+ update_flags (update_flags) {};
- // add up contributions of ansatz
- // functions
- for (unsigned int point=0; point<n_quadrature_points; ++point)
- for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
- values[point] += (dof_values[shape_func] *
- shape_values(shape_func, point));
-};
template <int dim>
const Point<dim> &
-FEValues<dim>::shape_grad (const unsigned int i,
- const unsigned int j) const {
- Assert (i<shape_values.m(), ExcInvalidIndex (i, shape_values.m()));
- Assert (j<shape_values.n(), ExcInvalidIndex (j, shape_values.n()));
+FEValuesBase<dim>::shape_grad (const unsigned int i,
+ const unsigned int j) const {
+ Assert (i<shape_gradients.size(),
+ ExcInvalidIndex (i, shape_gradients.size()));
+ Assert (j<shape_gradients[i].size(),
+ ExcInvalidIndex (j, shape_gradients[i].size()));
Assert (update_flags & update_gradients, ExcAccessToUninitializedField());
return shape_gradients[i][j];
template <int dim>
-void FEValues<dim>::get_function_grads (const dVector &fe_function,
- vector<Point<dim> > &gradients) const {
+void FEValuesBase<dim>::get_function_grads (const dVector &fe_function,
+ vector<Point<dim> > &gradients) const {
Assert (gradients.size() == n_quadrature_points,
ExcWrongVectorSize(gradients.size(), n_quadrature_points));
template <int dim>
-const Point<dim> & FEValues<dim>::quadrature_point (const unsigned int i) const {
+const Point<dim> & FEValuesBase<dim>::quadrature_point (const unsigned int i) const {
Assert (i<n_quadrature_points, ExcInvalidIndex(i, n_quadrature_points));
Assert (update_flags & update_q_points, ExcAccessToUninitializedField());
template <int dim>
-const Point<dim> & FEValues<dim>::ansatz_point (const unsigned int i) const {
+const Point<dim> & FEValuesBase<dim>::ansatz_point (const unsigned int i) const {
Assert (i<ansatz_points.size(), ExcInvalidIndex(i, ansatz_points.size()));
Assert (update_flags & update_ansatz_points, ExcAccessToUninitializedField());
template <int dim>
-double FEValues<dim>::JxW (const unsigned int i) const {
+double FEValuesBase<dim>::JxW (const unsigned int i) const {
Assert (i<n_quadrature_points, ExcInvalidIndex(i, n_quadrature_points));
Assert (update_flags & update_JxW_values, ExcAccessToUninitializedField());
+
+/*------------------------------- FEValues -------------------------------*/
+
+template <int dim>
+FEValues<dim>::FEValues (const FiniteElement<dim> &fe,
+ const Quadrature<dim> &quadrature,
+ const UpdateFlags update_flags) :
+ FEValuesBase<dim> (quadrature.n_quadrature_points,
+ fe.total_dofs,
+ fe.total_dofs,
+ update_flags),
+ shape_values(fe.total_dofs, quadrature.n_quadrature_points),
+ unit_shape_gradients(fe.total_dofs,
+ vector<Point<dim> >(quadrature.n_quadrature_points)),
+ unit_quadrature_points(quadrature.get_quad_points())
+{
+ for (unsigned int i=0; i<fe.total_dofs; ++i)
+ for (unsigned int j=0; j<n_quadrature_points; ++j)
+ {
+ shape_values(i,j) = fe.shape_value(i, quadrature.quad_point(j));
+ unit_shape_gradients[i][j]
+ = fe.shape_grad(i, quadrature.quad_point(j));
+ };
+
+ for (unsigned int i=0; i<n_quadrature_points; ++i)
+ {
+ weights[i] = quadrature.weight(i);
+ };
+};
+
+
+
+template <int dim>
+double FEValues<dim>::shape_value (const unsigned int i,
+ const unsigned int j) const {
+ Assert (i<shape_values.m(), ExcInvalidIndex (i, shape_values.m()));
+ Assert (j<shape_values.n(), ExcInvalidIndex (j, shape_values.n()));
+
+ return shape_values(i,j);
+};
+
+
+
+template <int dim>
+void FEValues<dim>::get_function_values (const dVector &fe_function,
+ vector<double> &values) const {
+ Assert (values.size() == n_quadrature_points,
+ ExcWrongVectorSize(values.size(), n_quadrature_points));
+
+ // get function values of dofs
+ // on this cell
+ vector<double> dof_values (total_dofs, 0);
+ present_cell->get_dof_values (fe_function, dof_values);
+
+ // initialize with zero
+ fill_n (values.begin(), n_quadrature_points, 0);
+
+ // add up contributions of ansatz
+ // functions
+ for (unsigned int point=0; point<n_quadrature_points; ++point)
+ for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
+ values[point] += (dof_values[shape_func] *
+ shape_values(shape_func, point));
+};
+
+
+
+
+
template <int dim>
void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
const FiniteElement<dim> &fe,
FEFaceValues<dim>::FEFaceValues (const FiniteElement<dim> &fe,
const Quadrature<dim-1> &quadrature,
const UpdateFlags update_flags) :
- n_quadrature_points(quadrature.n_quadrature_points),
- total_dofs(fe.total_dofs),
- shape_gradients(fe.total_dofs,
- vector<Point<dim> >(quadrature.n_quadrature_points)),
- weights(quadrature.n_quadrature_points, 0),
- JxW_values(quadrature.n_quadrature_points, 0),
- quadrature_points(quadrature.n_quadrature_points, Point<dim>()),
+ FEValuesBase<dim> (quadrature.n_quadrature_points,
+ fe.dofs_per_face,
+ fe.total_dofs,
+ update_flags),
unit_quadrature_points(quadrature.get_quad_points()),
- ansatz_points (fe.dofs_per_face, Point<dim>()),
- jacobi_matrices (quadrature.n_quadrature_points,dFMatrix(dim,dim)),
face_jacobi_determinants (quadrature.n_quadrature_points,0),
normal_vectors (quadrature.n_quadrature_points,Point<dim>()),
- update_flags (update_flags),
selected_face(0)
{
for (unsigned int face=0; face<2*dim; ++face)
-template <int dim>
-const Point<dim> &
-FEFaceValues<dim>::shape_grad (const unsigned int i,
- const unsigned int j) const {
- Assert (i<shape_values[selected_face].m(),
- ExcInvalidIndex (i, shape_values[selected_face].m()));
- Assert (j<shape_values[selected_face].n(),
- ExcInvalidIndex (j, shape_values[selected_face].n()));
- Assert (update_flags & update_gradients,
- ExcAccessToUninitializedField());
-
- return shape_gradients[i][j];
-};
-
-
-
-template <int dim>
-void FEFaceValues<dim>::get_function_grads (const dVector &fe_function,
- vector<Point<dim> > &gradients) const {
- Assert (gradients.size() == n_quadrature_points,
- ExcWrongVectorSize(gradients.size(), n_quadrature_points));
-
- // get function values of dofs
- // on this cell
- vector<double> dof_values (total_dofs, 0);
- present_cell->get_dof_values (fe_function, dof_values);
-
- // initialize with zero
- fill_n (gradients.begin(), n_quadrature_points, Point<dim>());
-
- // add up contributions of ansatz
- // functions
- for (unsigned int point=0; point<n_quadrature_points; ++point)
- for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
- gradients[point] += (dof_values[shape_func] *
- shape_gradients[shape_func][point]);
-};
-
-
-
-template <int dim>
-const Point<dim> & FEFaceValues<dim>::quadrature_point (const unsigned int i) const {
- Assert (i<n_quadrature_points, ExcInvalidIndex(i, n_quadrature_points));
- Assert (update_flags & update_q_points,
- ExcAccessToUninitializedField());
-
- return quadrature_points[i];
-};
-
-
-
-template <int dim>
-const Point<dim> & FEFaceValues<dim>::ansatz_point (const unsigned int i) const {
- Assert (i<ansatz_points.size(), ExcInvalidIndex(i, ansatz_points.size()));
- Assert (update_flags & update_ansatz_points,
- ExcAccessToUninitializedField());
-
- return ansatz_points[i];
-};
-
-
-
template <int dim>
const Point<dim> & FEFaceValues<dim>::normal_vector (const unsigned int i) const {
Assert (i<normal_vectors.size(), ExcInvalidIndex(i, normal_vectors.size()));
-template <int dim>
-double FEFaceValues<dim>::JxW (const unsigned int i) const {
- Assert (i<n_quadrature_points, ExcInvalidIndex(i, n_quadrature_points));
- Assert (update_flags & update_JxW_values,
- ExcAccessToUninitializedField());
-
- return JxW_values[i];
-};
-
-
-#include <grid/dof_accessor.h>
-
template <int dim>
void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
const unsigned int face_no,
/*------------------------------- Explicit Instantiations -------------*/
+template class FEValuesBase<1>;
+template class FEValuesBase<2>;
+
template class FEValues<1>;
template class FEValues<2>;