* then a call to the
* @p{get_unit_support_points}
* yields a non-empty array.
+ *
+ * The result may be false if an
+ * element is not defined by
+ * interpolating shape functions,
+ * for example by P-elements on
+ * quadrilaterals. It will
+ * usually only be true if the
+ * element constructs its shape
+ * functions by the requirement
+ * that they be one at a certain
+ * point and zero at all the
+ * points associated with the
+ * other shape functions.
+ *
+ * In composed elements (i.e. for
+ * the @ref{FESystem} class, the
+ * result will be true if all all
+ * the base elements have defined
+ * support points.
*/
bool has_support_points () const;
+
+ /**
+ * Return the position of the
+ * support point of the
+ * @p{index}th shape function. If
+ * it does not exist, raise an
+ * exception.
+ *
+ * The default implementation
+ * simply returns the respective
+ * element from the array you get
+ * from
+ * @ref{get_unit_support_points},
+ * but derived elements may
+ * overload this function. In
+ * particular, note that the
+ * @ref{FESystem} class overloads
+ * it so that it can return the
+ * support points of individual
+ * base elements, of not all the
+ * base elements define support
+ * points. In this way, you can
+ * still ask for certain support
+ * points, even if
+ * @p{get_unit_support_points}
+ * only returns an empty array.
+ */
+ virtual
+ Point<dim>
+ unit_support_point (const unsigned int index) const;
/**
* Return the support points of
* is true, then a call to the
* @p{get_unit_support_points}
* yields a non-empty array.
+ *
+ * For more information, see the
+ * documentation for the
+ * @ref{has_support_points}
+ * function.
*/
bool has_face_support_points () const;
+ /**
+ * The function corresponding to
+ * the @p{unit_support_point}
+ * function, but for faces. See
+ * there for more information.
+ */
+ virtual
+ Point<dim-1>
+ unit_face_support_point (const unsigned int index) const;
+
/**
* Return in which of the vector
* components of this finite
base_element (const unsigned int index) const;
/**
- * Check for non-zero values on a face.
+ * Check for non-zero values on a
+ * face.
*
* This function returns
* @p{true}, if the shape
virtual bool has_support_on_face (const unsigned int shape_index,
const unsigned int face_index) const;
+ /**
+ * Implementation of the
+ * respective function in the
+ * base class.
+ */
+ virtual
+ Point<dim>
+ unit_support_point (const unsigned int index) const;
+
+ /**
+ * Implementation of the
+ * respective function in the
+ * base class.
+ */
+ virtual
+ Point<dim-1>
+ unit_face_support_point (const unsigned int index) const;
+
/**
* Determine an estimate for the
* memory consumption (in bytes)
+template <int dim>
+Point<dim>
+FiniteElementBase<dim>::unit_support_point (const unsigned index) const
+{
+ Assert (index < this->dofs_per_cell,
+ ExcIndexRange (index, 0, this->dofs_per_cell));
+ Assert (unit_support_points.size() == this->dofs_per_cell,
+ ExcFEHasNoSupportPoints ());
+ return unit_support_points[index];
+};
+
+
+
template <int dim>
const std::vector<Point<dim-1> > &
FiniteElementBase<dim>::get_unit_face_support_points () const
+template <int dim>
+Point<dim-1>
+FiniteElementBase<dim>::unit_face_support_point (const unsigned index) const
+{
+ Assert (index < this->dofs_per_face,
+ ExcIndexRange (index, 0, this->dofs_per_face));
+ Assert (unit_face_support_points.size() == this->dofs_per_face,
+ ExcFEHasNoSupportPoints ());
+ return unit_face_support_points[index];
+};
+
+
+
template <int dim>
unsigned int
FiniteElementBase<dim>::memory_consumption () const
+template <int dim>
+Point<dim>
+FESystem<dim>::unit_support_point (const unsigned index) const
+{
+ Assert (index < this->dofs_per_cell,
+ ExcIndexRange (index, 0, this->dofs_per_cell));
+ Assert ((unit_support_points.size() == this->dofs_per_cell) ||
+ (unit_support_points.size() == 0),
+ typename FiniteElementBase<dim>::ExcFEHasNoSupportPoints ());
+
+ // let's see whether we have the
+ // information pre-computed
+ if (this->unit_support_points.size() != 0)
+ return this->unit_support_points[index];
+ else
+ // no. ask the base element
+ // whether it would like to
+ // provide this information
+ return (base_element(system_to_base_index(index).first.first)
+ .unit_support_point(system_to_base_index(index).second));
+};
+
+
+
+template <int dim>
+Point<dim-1>
+FESystem<dim>::unit_face_support_point (const unsigned index) const
+{
+ Assert (index < this->dofs_per_face,
+ ExcIndexRange (index, 0, this->dofs_per_face));
+ Assert ((unit_face_support_points.size() == this->dofs_per_face) ||
+ (unit_face_support_points.size() == 0),
+ typename FiniteElementBase<dim>::ExcFEHasNoSupportPoints ());
+
+ // let's see whether we have the
+ // information pre-computed
+ if (this->unit_face_support_points.size() != 0)
+ return this->unit_face_support_points[index];
+ else
+ // no. ask the base element
+ // whether it would like to
+ // provide this information
+ return (base_element(face_system_to_base_index(index).first.first)
+ .unit_face_support_point(face_system_to_base_index(index).second));
+};
+
+
+
template <int dim>
unsigned int
template <int dim>
void
-VectorTools::interpolate_boundary_values (const Mapping<dim> &mapping,
- const DoFHandler<dim> &dof,
- const typename FunctionMap<dim>::type &function_map,
- std::map<unsigned int,double> &boundary_values,
- const std::vector<bool> &component_mask_)
+VectorTools::
+interpolate_boundary_values (const Mapping<dim> &mapping,
+ const DoFHandler<dim> &dof,
+ const typename FunctionMap<dim>::type &function_map,
+ std::map<unsigned int,double> &boundary_values,
+ const std::vector<bool> &component_mask_)
{
// if for whatever reason we were
// passed an empty map, return
// support points. this wil be used
// to obtain the quadrature points
// on the real cell's face
- const std::vector<Point<dim-1> >
- & unit_support_points = fe.get_unit_face_support_points();
+ std::vector<Point<dim-1> >
+ unit_support_points = fe.get_unit_face_support_points();
// check whether there are support
- // points on the face, if not, then
- // this FE does not allow to be
- // used in this function
- Assert (unit_support_points.size() != 0, ExcNonInterpolatingFE());
+ // points on the face. if not, then
+ // we should try a more clever
+ // way. the idea is that a finite
+ // element may not offer support
+ // points for all its shape
+ // functions, but maybe only
+ // some. if it offers support
+ // points for the components we are
+ // interested in in this function,
+ // then that's fine. if not, the
+ // function we call in the finite
+ // element will raise an
+ // exception. the support points
+ // for the other shape functions
+ // are left uninitialized (well,
+ // initialized by the default
+ // constructor), since we don't
+ // need them anyway.
+ if (unit_support_points.size() == 0)
+ {
+ unit_support_points.resize (fe.dofs_per_face);
+ for (unsigned int i=0; i<fe.dofs_per_face; ++i)
+ if (component_mask[fe.face_system_to_component_index(i).first]
+ == true)
+ unit_support_points[i] = fe.unit_face_support_point(i);
+ };
Quadrature<dim-1> aux_quad (unit_support_points);
FEFaceValues<dim> fe_values (mapping, fe, aux_quad, update_q_points);
<h3>deal.II</h3>
<ol>
+ <li> <p>
+ New: For finite element classes, the functions
+ <code class="member">unit_support_point</code> and
+ <code class="member">unit_face_support_point</code> return the position
+ of an individual support point. This is necessary when you want to get
+ information about the support points of certain components in a composed
+ finite element, where not all components provide support points, and the
+ composed element thus does not fill the array the
+ <code class="member">get_unit_support_points</code> function returns.
+ <br>
+ (WB 2002/11/23)
+ </p>
+
<li> <p>
Fixed: Vectors could not be given as input and output vectors to the
<code class="class">SolutionTransfer</code> class at the same time, but