FiniteElementBase (const FiniteElementData<dim> &fe_data,
const std::vector<bool> &restriction_is_additive_flags);
+ /**
+ * Return the value of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element.
+ *
+ * An
+ * @p{ExcUnitShapeValuesDoNotExist}
+ * is thrown if the shape values
+ * of the @p{FiniteElement} under
+ * consideration depend on the
+ * shape of the cell in real
+ * space.
+ */
+ virtual double shape_value (const unsigned int i,
+ const Point<dim> &p) const;
+
+ /**
+ * Return the gradient of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element, and
+ * likewise the gradient is the
+ * gradient on the unit cell with
+ * respect to unit cell
+ * coordinates.
+ *
+ * An
+ * @p{ExcUnitShapeValuesDoNotExist}
+ * is thrown if the shape values
+ * of the @p{FiniteElement} under
+ * consideration depend on the
+ * shape of the cell in real
+ * space.
+ */
+ virtual Tensor<1,dim> shape_grad (const unsigned int i,
+ const Point<dim> &p) const;
+
+ /**
+ * Return the tensor of second
+ * derivatives of the @p{i}th
+ * shape function at point @p{p}
+ * on the unit cell. The
+ * derivatives are derivatives on
+ * the unit cell with respect to
+ * unit cell coordinates.
+ *
+ * An
+ * @p{ExcUnitShapeValuesDoNotExist}
+ * is thrown if the shape values
+ * of the @p{FiniteElement} under
+ * consideration depend on the
+ * shape of the cell in real
+ * space.
+ */
+ virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const;
+
/**
* Projection from a fine grid
* space onto a coarse grid
*/
unsigned int memory_consumption () const;
+ /**
+ * Exception
+ */
+ DeclException0 (ExcUnitShapeValuesDoNotExist);
+
/**
* Exception
*/
* Destructor.
*/
~FE_DGQ ();
+
+ /**
+ * Return the value of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element.
+ */
+ virtual double shape_value (const unsigned int i,
+ const Point<dim> &p) const;
+
+ /**
+ * Return the gradient of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element, and
+ * likewise the gradient is the
+ * gradient on the unit cell with
+ * respect to unit cell
+ * coordinates.
+ */
+ virtual Tensor<1,dim> shape_grad (const unsigned int i,
+ const Point<dim> &p) const;
+ /**
+ * Return the tensor of second
+ * derivatives of the @p{i}th
+ * shape function at point @p{p}
+ * on the unit cell. The
+ * derivatives are derivatives on
+ * the unit cell with respect to
+ * unit cell coordinates.
+ */
+ virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const;
+
/**
* Return the polynomial degree
* of this finite element,
* Destructor.
*/
~FE_Q ();
+
+ /**
+ * Return the value of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element.
+ */
+ virtual double shape_value (const unsigned int i,
+ const Point<dim> &p) const;
+
+ /**
+ * Return the gradient of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element, and
+ * likewise the gradient is the
+ * gradient on the unit cell with
+ * respect to unit cell
+ * coordinates.
+ */
+ virtual Tensor<1,dim> shape_grad (const unsigned int i,
+ const Point<dim> &p) const;
+
+ /**
+ * Return the tensor of second
+ * derivatives of the @p{i}th
+ * shape function at point @p{p}
+ * on the unit cell. The
+ * derivatives are derivatives on
+ * the unit cell with respect to
+ * unit cell coordinates.
+ */
+ virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const;
//TODO:[RH] make get_renumber private or move it some other place
/**
* Read-only access to the
* renumber vector.
+ *
+ * This function shouldn't be
+ * used, i.e. the users code
+ * shouldn't rely on the actual
+ * numbering of the degrees of
+ * dreedom on each cell. This,
+ * because this internal
+ * numbering might change in
+ * future.
*/
const std::vector<unsigned int> & get_renumber() const;
* shape function numbering.
*/
std::vector<unsigned int> renumber;
-
+
+ /**
+ * Inverse renumber
+ * vector. i.e. mapping from
+ * shape function numbering to
+ * lexicographic numbering.
+ */
+ std::vector<unsigned int> renumber_inverse;
+
/**
* Mapping from lexicographic to
* shape function numbering on first face.
*/
virtual ~FESystem ();
+ /**
+ * Return the value of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element.
+ *
+ * An
+ * @p{ExcUnitShapeValuesDoNotExist}
+ * is thrown if the shape values
+ * of the @p{FiniteElement}
+ * (corresponding to the @p{i}th
+ * shape function) depend on the
+ * shape of the cell in real
+ * space.
+ */
+ virtual double shape_value (const unsigned int i,
+ const Point<dim> &p) const;
+
+ /**
+ * Return the gradient of the
+ * @p{i}th shape function at the
+ * point @p{p}. @p{p} is a point
+ * on the reference element, and
+ * likewise the gradient is the
+ * gradient on the unit cell with
+ * respect to unit cell
+ * coordinates.
+ *
+ * An
+ * @p{ExcUnitShapeValuesDoNotExist}
+ * is thrown if the shape values
+ * of the @p{FiniteElement}
+ * (corresponding to the @p{i}th
+ * shape function) depend on the
+ * shape of the cell in real
+ * space.
+ */
+ virtual Tensor<1,dim> shape_grad (const unsigned int i,
+ const Point<dim> &p) const;
+
+ /**
+ * Return the tensor of second
+ * derivatives of the @p{i}th
+ * shape function at point @p{p}
+ * on the unit cell. The
+ * derivatives are derivatives on
+ * the unit cell with respect to
+ * unit cell coordinates.
+ *
+ * An
+ * @p{ExcUnitShapeValuesDoNotExist}
+ * is thrown if the shape values
+ * of the @p{FiniteElement}
+ * (corresponding to the @p{i}th
+ * shape function) depend on the
+ * shape of the cell in real
+ * space.
+ */
+ virtual Tensor<2,dim> shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const;
/**
* Number of different base
};
+template <int dim>
+double
+FiniteElementBase<dim>::shape_value (const unsigned int,
+ const Point<dim> &) const
+{
+ AssertThrow(false, ExcUnitShapeValuesDoNotExist());
+ return 0.;
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FiniteElementBase<dim>::shape_grad (const unsigned int,
+ const Point<dim> &) const
+{
+ AssertThrow(false, ExcUnitShapeValuesDoNotExist());
+ return Tensor<1,dim> ();
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FiniteElementBase<dim>::shape_grad_grad (const unsigned int,
+ const Point<dim> &) const
+{
+ AssertThrow(false, ExcUnitShapeValuesDoNotExist());
+ return Tensor<2,dim> ();
+}
+
template <int dim>
const FullMatrix<double> &
+template <int dim>
+double
+FE_DGQ<dim>::shape_value (const unsigned int i,
+ const Point<dim> &p) const
+{
+ return poly->compute_value(i, p);
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FE_DGQ<dim>::shape_grad (const unsigned int i,
+ const Point<dim> &p) const
+{
+ return poly->compute_grad(i, p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FE_DGQ<dim>::shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const
+{
+ return poly->compute_grad_grad(i, p);
+}
+
//----------------------------------------------------------------------
// Auxiliary functions
std::vector<bool> (1,false)),
degree(degree),
renumber(dofs_per_cell, 0),
+ renumber_inverse(dofs_per_cell, 0),
face_renumber(dofs_per_face, 0),
poly(0)
{
// face function is empty
build_renumbering (*this, degree, renumber);
build_face_renumbering (degree, face_renumber);
+
+ // build inverse of renumbering
+ // vector
+ for (unsigned int i=0; i<dofs_per_cell; ++i)
+ renumber_inverse[renumber[i]]=i;
// copy constraint matrices if they
// are defined. otherwise set them
+template <int dim>
+double
+FE_Q<dim>::shape_value (const unsigned int i,
+ const Point<dim> &p) const
+{
+ return poly->compute_value(renumber_inverse[i], p);
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FE_Q<dim>::shape_grad (const unsigned int i,
+ const Point<dim> &p) const
+{
+ return poly->compute_grad(renumber_inverse[i], p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FE_Q<dim>::shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const
+{
+ return poly->compute_grad_grad(renumber_inverse[i], p);
+}
+
+
//----------------------------------------------------------------------
// Auxiliary functions
//----------------------------------------------------------------------
+template <int dim>
+double
+FESystem<dim>::shape_value (const unsigned int i,
+ const Point<dim> &p) const
+{
+ Assert((i<dofs_per_cell), ExcIndexRange(i, 0, dofs_per_cell));
+
+ pair<unsigned,unsigned> comp = system_to_component_index(i);
+
+ return base_element(component_to_base_table[comp.first])
+ .shape_value(comp.second, p);
+}
+
+
+
+template <int dim>
+Tensor<1,dim>
+FESystem<dim>::shape_grad (const unsigned int i,
+ const Point<dim> &p) const
+{
+ Assert((i<dofs_per_cell), ExcIndexRange(i, 0, dofs_per_cell));
+
+ pair<unsigned,unsigned> comp = system_to_component_index(i);
+
+ return base_element(component_to_base_table[comp.first])
+ .shape_grad(comp.second, p);
+}
+
+
+
+template <int dim>
+Tensor<2,dim>
+FESystem<dim>::shape_grad_grad (const unsigned int i,
+ const Point<dim> &p) const
+{
+ Assert((i<dofs_per_cell), ExcIndexRange(i, 0, dofs_per_cell));
+
+
+ pair<unsigned,unsigned> comp = system_to_component_index(i);
+
+ return base_element(component_to_base_table[comp.first])
+ .shape_grad_grad(comp.second, p);
+}
+
+
+
template <int dim>
UpdateFlags
FESystem<dim>::update_once (const UpdateFlags flags) const