* @p{n.ds}.
*/
const Tensor<1,dim> & boundary_form (const unsigned int i) const;
+
+ /**
+ * Get the natural orientation of
+ * the face.
+ *
+ * If this is true, then the
+ * natural normal vector of the
+ * face points outwards (faces
+ * 0,1 in two dimensions, faces
+ * 0,2,4 in three). On the other
+ * faces, the natural normal
+ * vector points inwards.
+ *
+ * With this function, you have
+ * access to a consistent and
+ * unique orientation of all
+ * faces.
+ */
+ bool orientation () const;
/**
* Return the list of outward normal
* formula here.
*/
const Quadrature<dim-1> quadrature;
+
+ /**
+ * Value for orientation function.
+ */
+ bool my_orientation;
+ /**
+ * Data table for orientation
+ * values. Since the faces are
+ * numbered in a creative
+ * fashion, we stores these
+ * values instead of computing
+ * them.
+ */
+ static const bool orientation_table[];
};
* number @p{face_no} of @p{cell}
* and the given finite element.
*/
- void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
- const unsigned int face_no);
+ void reinit (const typename DoFHandler<dim>::cell_iterator& cell,
+ const unsigned int face_no);
private:
* number @p{face_no} of @p{cell}
* and the given finite element.
*/
- void reinit (const typename DoFHandler<dim>::cell_iterator &cell,
- const unsigned int face_no,
- const unsigned int subface_no);
+ void reinit (const typename DoFHandler<dim>::cell_iterator& cell,
+ const unsigned int face_no,
+ const unsigned int subface_no);
/**
* Exception
}
+template <int dim>
+bool
+FEFaceValuesBase<dim>::orientation () const
+{
+ return my_orientation;
+}
+
+
template <int dim>
const Quadrature<dim-1> &
FEFaceValuesBase<dim>::get_quadrature () const
/*------------------------------- FEFaceValuesBase --------------------------*/
+#if deal_II_dimension == 1
+template <>
+const bool
+FEFaceValuesBase<deal_II_dimension>::orientation_table[] =
+{
+ false, true
+};
+#endif
+
+
+#if deal_II_dimension == 2
+template <>
+const bool
+FEFaceValuesBase<deal_II_dimension>::orientation_table[] =
+{
+ true, true, false, false
+};
+#endif
+
+
+#if deal_II_dimension == 3
+template <>
+const bool
+FEFaceValuesBase<deal_II_dimension>::orientation_table[] =
+{
+ true, false, true, false, true, false
+};
+#endif
+
+
template <int dim>
FEFaceValuesBase<dim>::FEFaceValuesBase (const unsigned int n_q_points,
const unsigned int dofs_per_cell,
typename FEValuesBase<dim>::ExcFEDontMatch());
this->present_cell = cell;
-
+ this->my_orientation = this->orientation_table[face_no];
+
this->get_mapping().fill_fe_face_values(cell, face_no,
this->quadrature,
*this->mapping_data,
ExcIndexRange (subface_no, 0, GeometryInfo<dim>::subfaces_per_face));
this->present_cell = cell;
+ this->my_orientation = this->orientation_table[face_no];
this->get_mapping().fill_fe_subface_values(cell, face_no, subface_no,
this->quadrature,
<h3>deal.II</h3>
<ol>
+ <li> <p>
+ Improved: <code class="class">FEFaceValuesBase</code> has a new
+ function <code class="member">orientation</code> accessing a unique
+ and consistent orientation for each face.
+ <br>
+ (GK 2003/06/23)
+
<li> <p>
Changed: Embedding and restriction matrices for intergrid transfer are
now computed in the constructor of most elements, rather than taken from
elements are available since the old tables were only precomputed for
certain polynomial degrees and are now available for all.
<br>
- (WB 2002/06/09)
+ (WB 2003/06/09)
</p>
<li> <p>
this function if it implements the desired interpolation, and falls back
to the old algorithm if this is not the case.
<br>
- (WB 2002/06/09)
+ (WB 2003/06/09)
</p>
<li> <p>
class="member">get_name</class>, which can be used to identify a finite
element by its name.
<br>
- (WB 2002/06/09)
+ (WB 2003/06/09)
</p>
<li> <p>
New: Raviart-Thomas elements are now implemented in the <code
class="class">FE_RaviartThomas</class> class.
<br>
- (WB 2002/06/09)
+ (WB 2003/06/09)
</p>
</ol>