/*---------------------------- fe.h ---------------------------*/
#include <base/exceptions.h>
+#include <base/subscriptor.h>
#include <grid/point.h>
#include <grid/dof.h>
#include <grid/geometry_info.h>
* @author Wolfgang Bangerth, 1998
*/
template <int dim>
-struct FiniteElementBase : public FiniteElementData<dim> {
+struct FiniteElementBase :
+ public Subscriptor,
+ public FiniteElementData<dim>
+
+{
public:
/**
* Construct an object of this type.
/*---------------------------- fe_values.h ---------------------------*/
+#include <base/exceptions.h>
+#include <base/subscriptor.h>
#include <lac/dfmatrix.h>
#include <grid/point.h>
-#include <base/exceptions.h>
#include <grid/tria.h>
#include <fe/fe_update_flags.h>
* may use other ways.)
*/
void reinit (const typename DoFHandler<dim>::cell_iterator &,
- const FiniteElement<dim> &,
const Boundary<dim> &);
private:
+ /**
+ * Store the finite element for later use.
+ */
+ SmartPointer<const FiniteElement<dim> > fe;
/**
* Store the gradients of the shape
* functions at the quadrature points on
*/
void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
const unsigned int face_no,
- const FiniteElement<dim> &fe,
const Boundary<dim> &boundary);
+private:
+ /**
+ * Store the finite element for later use.
+ */
+ SmartPointer<const FiniteElement<dim> > fe;
};
void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
const unsigned int face_no,
const unsigned int subface_no,
- const FiniteElement<dim> &fe,
const Boundary<dim> &boundary);
/**
* Exception
*/
DeclException0 (ExcReinitCalledWithBoundaryFace);
+private:
+ /**
+ * Store the finite element for later use.
+ */
+ SmartPointer<const FiniteElement<dim> > fe;
};
fe.n_transform_functions,
1,
update_flags),
+ fe(&fe),
unit_shape_gradients(fe.total_dofs,
vector<Point<dim> >(quadrature.n_quadrature_points)),
unit_shape_gradients_transform(fe.n_transform_functions,
template <int dim>
void FEValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
- const FiniteElement<dim> &fe,
const Boundary<dim> &boundary) {
present_cell = cell;
// fill jacobi matrices and real
(update_flags & update_q_points) ||
(update_flags & update_gradients) ||
(update_flags & update_ansatz_points))
- fe.fill_fe_values (cell,
+ fe->fill_fe_values (cell,
unit_quadrature_points,
jacobi_matrices,
update_flags & (update_jacobians |
// compute gradients on real element if
// requested
if (update_flags & update_gradients)
- for (unsigned int i=0; i<fe.total_dofs; ++i)
+ for (unsigned int i=0; i<fe->total_dofs; ++i)
for (unsigned int j=0; j<n_quadrature_points; ++j)
for (unsigned int s=0; s<dim; ++s)
{
fe.total_dofs,
fe.n_transform_functions,
GeometryInfo<dim>::faces_per_cell,
- update_flags)
+ update_flags),
+ fe(&fe)
{
unit_face_quadrature_points = quadrature.get_quad_points();
weights = quadrature.get_weights ();
template <int dim>
void FEFaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
const unsigned int face_no,
- const FiniteElement<dim> &fe,
const Boundary<dim> &boundary) {
present_cell = cell;
selected_dataset = face_no;
(update_flags & update_gradients) ||
(update_flags & update_ansatz_points) ||
(update_flags & update_JxW_values))
- fe.fill_fe_face_values (cell,
+ fe->fill_fe_face_values (cell,
face_no,
unit_face_quadrature_points,
unit_quadrature_points[face_no],
// compute gradients on real element if
// requested
if (update_flags & update_gradients)
- for (unsigned int i=0; i<fe.total_dofs; ++i)
+ for (unsigned int i=0; i<fe->total_dofs; ++i)
{
fill_n (shape_gradients[i].begin(),
n_quadrature_points,
fe.total_dofs,
fe.n_transform_functions,
GeometryInfo<dim>::faces_per_cell * GeometryInfo<dim>::subfaces_per_face,
- update_flags)
+ update_flags),
+ fe(&fe)
{
Assert ((update_flags & update_ansatz_points) == false,
ExcInvalidUpdateFlag());
void FESubfaceValues<dim>::reinit (const typename DoFHandler<dim>::cell_iterator &cell,
const unsigned int face_no,
const unsigned int subface_no,
- const FiniteElement<dim> &fe,
const Boundary<dim> &boundary) {
Assert (cell->face(face_no)->at_boundary() == false,
ExcReinitCalledWithBoundaryFace());
(update_flags & update_q_points) ||
(update_flags & update_gradients) ||
(update_flags & update_JxW_values))
- fe.fill_fe_subface_values (cell,
+ fe->fill_fe_subface_values (cell,
face_no,
subface_no,
unit_face_quadrature_points,
// compute gradients on real element if
// requested
if (update_flags & update_gradients)
- for (unsigned int i=0; i<fe.total_dofs; ++i)
+ for (unsigned int i=0; i<fe->total_dofs; ++i)
{
fill_n (shape_gradients[i].begin(),
n_quadrature_points,
present_level,
present_index,
dof_handler),
- fe,
boundary);
const unsigned int n_dofs = dof_handler->get_selected_fe().total_dofs;
void KellyErrorEstimator<dim>::
integrate_over_regular_face (const active_cell_iterator &cell,
const unsigned int face_no,
- const FiniteElement<dim> &fe,
+ const FiniteElement<dim> &,
const Boundary<dim> &boundary,
const FunctionMap &neumann_bc,
const unsigned int n_q_points,
// initialize data of the restriction
// of this cell to the present face
- fe_face_values_cell.reinit (cell, face_no, fe, boundary);
+ fe_face_values_cell.reinit (cell, face_no, boundary);
// set up a vector of the gradients
// of the finite element function
// get restriction of finite element
// function of #neighbor# to the
// common face.
- fe_face_values_neighbor.reinit (neighbor, neighbor_neighbor,
- fe, boundary);
+ fe_face_values_neighbor.reinit (neighbor, neighbor_neighbor, boundary);
// get a list of the gradients of
// the finite element solution
void KellyErrorEstimator<dim>::
integrate_over_irregular_face (const active_cell_iterator &cell,
const unsigned int face_no,
- const FiniteElement<dim> &fe,
+ const FiniteElement<dim> &,
const Boundary<dim> &boundary,
const unsigned int n_q_points,
FEFaceValues<dim> &fe_face_values,
// present cell to the subface and
// store the gradient of the solution
// in psi
- fe_subface_values.reinit (cell, face_no, subface_no,
- fe, boundary);
+ fe_subface_values.reinit (cell, face_no, subface_no, boundary);
fe_subface_values.get_function_grads (solution, psi);
// restrict the finite element on the
// neighbor cell to the common #subface#.
// store the gradient in #neighbor_psi#
vector<Point<dim> > neighbor_psi (n_q_points);
- fe_face_values.reinit (neighbor_child, neighbor_neighbor,
- fe, boundary);
+ fe_face_values.reinit (neighbor_child, neighbor_neighbor, boundary);
fe_face_values.get_function_grads (solution, neighbor_psi);
// compute the jump in the gradients
cell_matrix.clear ();
cell_vector.clear ();
- fe_values.reinit (cell, face, fe, boundary);
+ fe_values.reinit (cell, face, boundary);
const dFMatrix &values = fe_values.get_shape_values ();
const vector<double> &weights = fe_values.get_JxW_values ();
{
double diff=0;
// initialize for this cell
- fe_values.reinit (cell, fe, boundary);
+ fe_values.reinit (cell, boundary);
switch (norm)
{