* (@p false). Which is the
* standard direction is
* documented with the
- * Triangulation class. In
+ * GeometryInfo class. In
* 1d and 2d, this is always
* @p true, but in 3d it may be
* different, see the respective
* discussion in the
* documentation of the
- * Triangulation class.
+ * GeometryInfo class.
*
* This function is really only
* for internal use in the
* know what this is all about.
*/
bool face_orientation (const unsigned int face) const;
+
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 180 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_flip (const unsigned int face) const;
+
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 90 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_rotation (const unsigned int face) const;
+
+ /**
+ * Return whether the line with index @p
+ * line is oriented in standard
+ * direction. @p true indicates, that the
+ * line is oriented from vertex 0 to vertex
+ * 1, whereas it is the other way around
+ * otherwise. In 1d and 2d, this is always
+ * @p true, but in 3d it may be different,
+ * see the respective discussion in the
+ * documentation of the
+ * GeometryInfo classe.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool line_orientation (const unsigned int line) const;
+ /**
+ * Set whether the quad with
+ * index @p face has its normal
+ * pointing in the standard
+ * direction (@p true) or
+ * whether it is the opposite
+ * (@p false). Which is the
+ * standard direction is
+ * documented with the
+ * GeometryInfo class.
+ *
+ * This function is only for
+ * internal use in the
+ * library. Setting this flag to
+ * any other value than the one
+ * that the triangulation has
+ * already set is bound to bring
+ * you desaster.
+ *
+ * It is only possible to set the
+ * face_orientation of cells in 3d
+ * (i.e. <code>celldim==3 &&
+ * dim==3</code>).
+ */
+ void set_face_orientation (const unsigned int face,
+ const bool orientation) const;
+
+ /**
+ * Set the flag indicating, what
+ * <code>face_flip()</code> will
+ * return.
+ *
+ * It is only possible to set the
+ * face_flip of cells in 3d
+ * (i.e. <code>celldim==3 &&
+ * dim==3</code>).
+ */
+ void set_face_flip (const unsigned int face,
+ const bool flip) const;
+
+ /**
+ * Set the flag indicating, what
+ * <code>face_rotation()</code> will
+ * return.
+ *
+ * It is only possible to set the
+ * face_rotation of cells in 3d
+ * (i.e. <code>celldim==3 &&
+ * dim==3</code>).
+ */
+ void set_face_rotation (const unsigned int face,
+ const bool rotation) const;
+
+ /**
+ * Set the flag indicating, what
+ * <code>line_orientation()</code> will
+ * return.
+ *
+ * It is only possible to set the
+ * line_orientation of faces in 3d
+ * (i.e. <code>celldim==2 &&
+ * dim==3</code>).
+ */
+ void set_line_orientation (const unsigned int line,
+ const bool orientation) const;
+
private:
/**
* Copy operator. This is normally
* (@p false). Which is the
* standard direction is
* documented with the
- * Triangulation class. In
+ * GeometryInfo class. In
* 1d and 2d, this is always
* @p true, but in 3d it may be
* different, see the respective
* discussion in the
* documentation of the
- * Triangulation class.
+ * GeometryInfo class.
*
* This function is really only
* for internal use in the
*/
bool face_orientation (const unsigned int face) const;
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 180 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_flip (const unsigned int face) const;
+
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 90 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_rotation (const unsigned int face) const;
+
+ /**
+ * Return whether the line with index @p
+ * line is oriented in standard
+ * direction. @p true indicates, that the
+ * line is oriented from vertex 0 to vertex
+ * 1, whereas it is the other way around
+ * otherwise. In 1d and 2d, this is always
+ * @p true, but in 3d it may be different,
+ * see the respective discussion in the
+ * documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool line_orientation (const unsigned int line) const;
+
+ /**
+ * Set the flag indicating, what
+ * <code>line_orientation()</code> will
+ * return.
+ *
+ * It is only possible to set the
+ * line_orientation of faces in 3d
+ * (i.e. <code>celldim==2 &&
+ * dim==3</code>).
+ */
+ void set_line_orientation (const unsigned int line,
+ const bool orientation) const;
+
private:
/**
* Copy operator. This is normally
* (@p false). Which is the
* standard direction is
* documented with the
- * Triangulation class. In
+ * GeometryInfo class. In
* 1d and 2d, this is always
* @p true, but in 3d it may be
* different, see the respective
* discussion in the
* documentation of the
- * Triangulation class.
+ * GeometryInfo class.
*
* This function is really only
* for internal use in the
*/
bool face_orientation (const unsigned int face) const;
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 180 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_flip (const unsigned int face) const;
+
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 90 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_rotation (const unsigned int face) const;
+
+ /**
+ * Return whether the line with index @p
+ * line is oriented in standard
+ * direction. @p true indicates, that the
+ * line is oriented from vertex 0 to vertex
+ * 1, whereas it is the other way around
+ * otherwise. In 1d and 2d, this is always
+ * @p true, but in 3d it may be different,
+ * see the respective discussion in the
+ * documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool line_orientation (const unsigned int line) const;
+
+ /**
+ * Set the flag indicating, what
+ * <code>line_orientation()</code> will
+ * return.
+ *
+ * It is only possible to set the
+ * line_orientation of faces in 3d
+ * (i.e. <code>celldim==2 &&
+ * dim==3</code>).
+ */
+ void set_line_orientation (const unsigned int line,
+ const bool orientation) const;
+
private:
/**
* Copy operator. This is
* (@p false). Which is the
* standard direction is
* documented with the
- * Triangulation class. In
+ * GeometryInfo class. In
* 1d and 2d, this is always
* @p true, but in 3d it may be
* different, see the respective
* discussion in the
* documentation of the
- * Triangulation class.
+ * GeometryInfo class.
*
* This function is really only
* for internal use in the
* know what this is all about.
*/
bool face_orientation (const unsigned int face) const;
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 180 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_flip (const unsigned int face) const;
+ /**
+ * Return whether the face with index @p
+ * face is rotated by 90 degrees (@p true)
+ * or or not (@p false). In 1d and 2d, this
+ * is always @p false, but in 3d it may be
+ * different, see the respective discussion
+ * in the documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool face_rotation (const unsigned int face) const;
+
+ /**
+ * Return whether the line with index @p
+ * line is oriented in standard
+ * direction. @p true indicates, that the
+ * line is oriented from vertex 0 to vertex
+ * 1, whereas it is the other way around
+ * otherwise. In 1d and 2d, this is always
+ * @p true, but in 3d it may be different,
+ * see the respective discussion in the
+ * documentation of the
+ * GeometryInfo class.
+ *
+ * This function is really only
+ * for internal use in the
+ * library unless you absolutely
+ * know what this is all about.
+ */
+ bool line_orientation (const unsigned int line) const;
+
/**
* Set whether the quad with
* index @p face has its normal
* (@p false). Which is the
* standard direction is
* documented with the
- * Triangulation class.
+ * GeometryInfo class.
*
* This function is only for
* internal use in the
*/
void set_face_orientation (const unsigned int face,
const bool orientation) const;
+
+ /**
+ * Set the flag indicating, what
+ * <code>face_flip()</code> will
+ * return.
+ *
+ * It is only possible to set the
+ * face_orientation of cells in 3d
+ * (i.e. <code>celldim==3 &&
+ * dim==3</code>).
+ */
+ void set_face_flip (const unsigned int face,
+ const bool flip) const;
+
+ /**
+ * Set the flag indicating, what
+ * <code>face_rotation()</code> will
+ * return.
+ *
+ * It is only possible to set the
+ * face_orientation of cells in 3d
+ * (i.e. <code>celldim==3 &&
+ * dim==3</code>).
+ */
+ void set_face_rotation (const unsigned int face,
+ const bool rotation) const;
private:
/**
* for the index of the
* child.
*
- * However, the function
- * is more complicated in 3d,
- * since there faces may have
- * more than one orientation, and
- * we have to use
- * @p face_orientation for both
- * this and the neighbor cell to
- * figure out which cell we want
- * to have.
+ * However, the function is more
+ * complicated in 3d, since there faces may
+ * have more than one orientation, and we
+ * have to use @p face_orientation, @p
+ * face_flip and @p face_rotation for both
+ * this and the neighbor cell to figure out
+ * which cell we want to have.
*
* This can lead to surprising
* results: if we are sitting on
* orderings are only identical
* if the face orientation is @p
* true, and reversed otherwise.
+ *
+ * Similarly, effects of
+ * <tt>face_flip()==true</tt> and
+ * <tt>face_rotation()==true()</tt>, both
+ * of which indicate a non-standard face
+ * have to be considered.
*
* Fortunately, this is only very rarely
* of concern, since one is usually only
}
+template <int dim>
+inline
+bool
+TriaObjectAccessor<1,dim>::face_flip (const unsigned int) const
+{
+ return false;
+}
+
+
+template <int dim>
+inline
+bool
+TriaObjectAccessor<1,dim>::face_rotation (const unsigned int) const
+{
+ return false;
+}
+
+
+template <int dim>
+inline
+bool
+TriaObjectAccessor<1,dim>::line_orientation (const unsigned int) const
+{
+ return true;
+}
+
+
template <int dim>
inline
}
+template <int dim>
+inline
+bool
+TriaObjectAccessor<2,dim>::face_flip (const unsigned int) const
+{
+ return false;
+}
+
+
+template <int dim>
+inline
+bool
+TriaObjectAccessor<2,dim>::face_rotation (const unsigned int) const
+{
+ return false;
+}
+
+
+template <>
+inline
+bool
+TriaObjectAccessor<2,3>::line_orientation (const unsigned int line) const
+{
+ // we cannot usae the quads() function here,
+ // since it returns a reference to
+ // TriaObjects<Quad>, but we need a
+ // (reference to) TriaObjectsQuad3D
+ const int dim=3;
+ Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
+ Assert (line<GeometryInfo<2>::lines_per_cell,
+ ExcIndexRange (line, 0, GeometryInfo<2>::lines_per_cell));
+ Assert (this->present_index * GeometryInfo<2>::lines_per_cell + line
+ < this->tria->faces->quads.line_orientations.size(),
+ ExcInternalError());
+ return this->tria->faces->quads.
+ line_orientations[this->present_index * GeometryInfo<2>::lines_per_cell + line];
+}
+
+
+template <int dim>
+inline
+bool
+TriaObjectAccessor<2,dim>::line_orientation (const unsigned int) const
+{
+ return true;
+}
+
template <int dim>
inline
Assert (i<12, ExcIndexRange(i,0,12));
// get the line index by asking the
- // quads. make sure we handle
- // reverted faces correctly
+ // quads. first assume standard orientation
//
- // so set up a table that for each
+ // set up a table that for each
// line describes a) from which
// quad to take it, b) which line
// therein it is if the face is
- // oriented correctly, and c) if in
- // the opposite direction
- static const unsigned int lookup_table[12][3] =
- { { 4, 0, 2 }, // take first four lines from bottom face
- { 4, 1, 3 },
- { 4, 2, 0 },
- { 4, 3, 1 },
-
- { 5, 0, 2 }, // second four lines from top face
- { 5, 1, 3 },
- { 5, 2, 0 },
- { 5, 3, 1 },
-
- { 0, 0, 2 }, // the rest randomly
- { 1, 0, 2 },
- { 0, 1, 3 },
- { 1, 1, 3 }};
+ // oriented correctly
+ static const unsigned int lookup_table[12][2] =
+ { { 4, 0 }, // take first four lines from bottom face
+ { 4, 1 },
+ { 4, 2 },
+ { 4, 3 },
+
+ { 5, 0 }, // second four lines from top face
+ { 5, 1 },
+ { 5, 2 },
+ { 5, 3 },
+
+ { 0, 0 }, // the rest randomly
+ { 1, 0 },
+ { 0, 1 },
+ { 1, 1 }};
+
+ // respect non-standard faces by calling the
+ // reordering function from GeometryInfo
+
+ const unsigned int quad_index=lookup_table[i][0];
+ const unsigned int std_line_index=lookup_table[i][1];
+
+ const unsigned int line_index=GeometryInfo<dim>::standard_to_real_face_line(
+ std_line_index,
+ face_orientation(quad_index),
+ face_flip(quad_index),
+ face_rotation(quad_index));
- return (this->quad(lookup_table[i][0])
- ->line_index(face_orientation(lookup_table[i][0]) ?
- lookup_table[i][1] :
- lookup_table[i][2]));
+ return (this->quad(quad_index)
+ ->line_index(line_index));
}
Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
Assert (face<GeometryInfo<3>::faces_per_cell,
ExcIndexRange (face, 0, GeometryInfo<3>::faces_per_cell));
- Assert (this->present_index * GeometryInfo<3>::faces_per_cell + face
- < this->tria->levels[this->present_level]
- ->hexes.face_orientations.size(),
- ExcInternalError());
+ Assert (this->present_index * GeometryInfo<3>::faces_per_cell + face
+ < this->tria->levels[this->present_level]
+ ->hexes.face_orientations.size(),
+ ExcInternalError());
+
+ return (this->tria->levels[this->present_level]
+ ->hexes.face_orientations[this->present_index *
+ GeometryInfo<3>::faces_per_cell
+ + face]);
+}
+
+
+template <>
+inline
+bool
+TriaObjectAccessor<3, 3>::
+face_flip (const unsigned int face) const
+{
+ const int dim=3;
+ Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
+ Assert (face<GeometryInfo<3>::faces_per_cell,
+ ExcIndexRange (face, 0, GeometryInfo<3>::faces_per_cell));
+ Assert (this->present_index * GeometryInfo<3>::faces_per_cell + face
+ < this->tria->levels[this->present_level]
+ ->hexes.face_flips.size(),
+ ExcInternalError());
+
+ return (this->tria->levels[this->present_level]
+ ->hexes.face_flips[this->present_index *
+ GeometryInfo<3>::faces_per_cell
+ + face]);
+}
+
+
+template <>
+inline
+bool
+TriaObjectAccessor<3, 3>::
+face_rotation (const unsigned int face) const
+{
+ const int dim=3;
+ Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
+ Assert (face<GeometryInfo<3>::faces_per_cell,
+ ExcIndexRange (face, 0, GeometryInfo<3>::faces_per_cell));
+ Assert (this->present_index * GeometryInfo<3>::faces_per_cell + face
+ < this->tria->levels[this->present_level]
+ ->hexes.face_rotations.size(),
+ ExcInternalError());
- return (this->tria->levels[this->present_level]
- ->hexes.face_orientations[this->present_index *
- GeometryInfo<3>::faces_per_cell
- + face]);
+ return (this->tria->levels[this->present_level]
+ ->hexes.face_rotations[this->present_index *
+ GeometryInfo<3>::faces_per_cell
+ + face]);
+}
+
+
+template <>
+inline
+bool
+TriaObjectAccessor<3, 3>::
+line_orientation (const unsigned int line) const
+{
+ const int dim=3;
+ Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
+ Assert (line<GeometryInfo<3>::lines_per_cell,
+ ExcIndexRange (line, 0, GeometryInfo<3>::lines_per_cell));
+
+ // get the line index by asking the
+ // quads. first assume standard orientation
+ //
+ // set up a table that for each
+ // line describes a) from which
+ // quad to take it, b) which line
+ // therein it is if the face is
+ // oriented correctly
+ static const unsigned int lookup_table[12][2] =
+ { { 4, 0 }, // take first four lines from bottom face
+ { 4, 1 },
+ { 4, 2 },
+ { 4, 3 },
+
+ { 5, 0 }, // second four lines from top face
+ { 5, 1 },
+ { 5, 2 },
+ { 5, 3 },
+
+ { 0, 0 }, // the rest randomly
+ { 1, 0 },
+ { 0, 1 },
+ { 1, 1 }};
+
+ const unsigned int quad_index=lookup_table[line][0];
+ const unsigned int std_line_index=lookup_table[line][1];
+
+ const unsigned int line_index=GeometryInfo<dim>::standard_to_real_face_line(
+ std_line_index,
+ face_orientation(quad_index),
+ face_flip(quad_index),
+ face_rotation(quad_index));
+
+ // now we got to the correct line and ask
+ // the quad for its line_orientation. however, if
+ // the face is rotated, it might be possible,
+ // that a standard orientation of the line
+ // with respect to the face corrsponds to a
+ // non-standard orientation for the line with
+ // respect to the cell.
+ //
+ // set up a table indicating if the two
+ // standard orientations coincide
+ //
+ // first index: two pairs of lines 0(lines
+ // 0/1) and 1(lines 2/3)
+ //
+ // second index: face_orientation; 0:
+ // opposite normal, 1: standard
+ //
+ // third index: face_flip; 0: standard, 1:
+ // face rotated by 180 degrees
+ //
+ // forth index: face_rotation: 0: standard,
+ // 1: face rotated by 90 degrees
+
+ static const bool bool_table[2][2][2][2] =
+ { { { { true, false }, // lines 0/1, face_orientation=false, face_flip=false, face_rotation=false and true
+ { false, true }}, // lines 0/1, face_orientation=false, face_flip=true, face_rotation=false and true
+ { { true, true }, // lines 0/1, face_orientation=true, face_flip=false, face_rotation=false and true
+ { false, false }}},// linea 0/1, face_orientation=true, face_flip=true, face_rotation=false and true
+
+ { { { true, true }, // lines 2/3 ...
+ { false, false }},
+ { { true, false },
+ { false, true }}}};
+
+
+ return (this->quad(quad_index)
+ ->line_orientation(line_index)
+ == bool_table[std_line_index/2]
+ [face_orientation(quad_index)]
+ [face_flip(quad_index)]
+ [face_rotation(quad_index)]);
}
* The TriaObject containing the data of quads.
*/
- TriaObjects<Quad> quads;
+ TriaObjectsQuad3D quads;
/**
* The TriaObject containing the data of lines.
DEAL_II_NAMESPACE_OPEN
+// TODO: See if we can unify the class hierarchy a bit using partial
+// specialization of classes here, e.g. declare a general class
+// TriaObjects<dim,G>.
+// To consider for this: in that case we would have to duplicate quite a few
+// things, e.g. TriaObjectsHex is derived from TriaObjects<Hexahedron> and
+// declares mainly additional data. This would have to be changed in case of a
+// partial specialization.
+
namespace internal
{
namespace Triangulation
* has.
*/
std::vector<bool> face_orientations;
+ /**
+ * flip = rotation by 180 degrees
+ */
+ std::vector<bool> face_flips;
+ /**
+ * rotation by 90 degrees
+ */
+ std::vector<bool> face_rotations;
/**
* Assert that enough space is
*/
unsigned int memory_consumption () const;
};
+
+
+/**
+ * For quadrilaterals in 3D the data of TriaObjects needs to be extended, as we
+ * can obtain faces (quads) with lines in non-standard-orientation, therefore we
+ * declare a class TriaObjectsQuad3D, which additionaly contains a bool-vector
+ * of the line-orientations.
+ * @ingroup grid
+ */
+
+ class TriaObjectsQuad3D: public TriaObjects<Quad>
+ {
+ public:
+
+ /**
+ * In effect, this field has
+ * <tt>4*n_quads</tt> elements,
+ * being the number of quads
+ * times the four lines each
+ * has.
+ */
+ std::vector<bool> line_orientations;
+
+ /**
+ * Assert that enough space is
+ * allocated to accomodate
+ * <tt>new_objs</tt> new objects.
+ * This function does not only call
+ * <tt>vector::reserve()</tt>, but
+ * does really append the needed
+ * elements.
+ */
+ void reserve_space (const unsigned int new_objs);
+
+ /**
+ * Clear all the data contained in this object.
+ */
+ void clear();
+
+ /**
+ * Check the memory consistency of the
+ * different containers. Should only be
+ * called with the prepro flag @p DEBUG
+ * set. The function should be called from
+ * the functions of the higher
+ * TriaLevel classes.
+ */
+ void monitor_memory (const unsigned int true_dimension) const;
+
+ /**
+ * Determine an estimate for the
+ * memory consumption (in bytes)
+ * of this object.
+ */
+ unsigned int memory_consumption () const;
+ };
+
}
}
{
Assert (corner<4, ExcIndexRange(corner,0,4));
- return line(corner%2)->vertex_index(corner/2);
+ // table used to switch the vertices, if the
+ // line orientation is wrong,
+ //
+ // first index: line orientation 0: false or
+ // 1: true=standard
+ //
+ // second index: vertex index to be switched
+ // (or not)
+
+ static const unsigned int switch_table[2][2]={{1,0},{0,1}};
+
+ return line(corner%2)->vertex_index(switch_table[line_orientation(corner%2)][corner/2]);
}
}
+template <int dim>
+void
+TriaObjectAccessor<2, dim>::
+set_line_orientation (const unsigned int ,
+ const bool ) const
+{
+ Assert(false,
+ ExcMessage("The line orientation of Quads can only be set in 3D. In 2D the lines are alway in standard orientation."));}
+
+
+
+#if deal_II_dimension == 3
+
+template <>
+void
+TriaObjectAccessor<2, 3>::
+set_line_orientation (const unsigned int line,
+ const bool orientation) const
+{
+ const int dim=3;
+ Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
+ Assert (line<GeometryInfo<dim>::lines_per_face,
+ ExcIndexRange (line, 0, GeometryInfo<3>::lines_per_face));
+ Assert (this->present_index * GeometryInfo<dim>::lines_per_face + line
+ < this->tria->faces->quads.line_orientations.size(),
+ ExcInternalError());
+
+ this->tria->faces
+ ->quads.line_orientations[this->present_index *
+ GeometryInfo<dim>::lines_per_face
+ +
+ line] = orientation;
+}
+
+#endif
/*------------------------ Functions: TriaObjectAccessor ---------------------------*/
{
Assert (corner<8, ExcIndexRange(corner,0,8));
- // get the corner indices by asking
- // either the bottom or the top
- // face for its vertices. make sure
- // we take into account that the
- // face might have non-standard
- // orientation; if this is not the
- // case, transpose vertices 1 and 2
- // of this face
- static const unsigned int vertex_translation[4] =
- { 0, 2, 1, 3 };
- if (corner<4)
- {
- if (face_orientation(4) == true)
- return quad(4)->vertex_index(corner);
- else
- return quad(4)->vertex_index(vertex_translation[corner]);
- }
- else
- {
- if (face_orientation(5) == true)
- return quad(5)->vertex_index(corner-4);
- else
- return quad(5)->vertex_index(vertex_translation[corner-4]);
- }
+ // get the corner indices by asking either
+ // the bottom or the top face for its
+ // vertices. handle non-standard faces by
+ // calling the vertex reordering function
+ // from GeometryInfo
+
+ // bottom face (4) for first four vertices,
+ // top face (5) for the rest
+ const unsigned int face_index=4+corner/4;
+
+ return quad(face_index)->vertex_index(GeometryInfo<dim>
+ ::standard_to_real_face_vertex(corner%4,
+ face_orientation(face_index),
+ face_flip(face_index),
+ face_rotation(face_index)));
}
face] = orientation;
}
+
+template <>
+void
+TriaObjectAccessor<3, 3>::
+set_face_flip (const unsigned int face,
+ const bool flip) const
+{
+ const int dim=3;
+ Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
+ Assert (face<GeometryInfo<3>::faces_per_cell,
+ ExcIndexRange (face, 0, GeometryInfo<3>::faces_per_cell));
+ Assert (this->present_index * GeometryInfo<3>::faces_per_cell + face
+ < this->tria->levels[this->present_level]
+ ->hexes.face_flips.size(),
+ ExcInternalError());
+
+ this->tria->levels[this->present_level]
+ ->hexes.face_flips[this->present_index *
+ GeometryInfo<3>::faces_per_cell
+ +
+ face] = flip;
+}
+
+
+template <>
+void
+TriaObjectAccessor<3, 3>::
+set_face_rotation (const unsigned int face,
+ const bool rotation) const
+{
+ const int dim=3;
+ Assert (used(), TriaAccessor<dim>::ExcCellNotUsed());
+ Assert (face<GeometryInfo<3>::faces_per_cell,
+ ExcIndexRange (face, 0, GeometryInfo<3>::faces_per_cell));
+ Assert (this->present_index * GeometryInfo<3>::faces_per_cell + face
+ < this->tria->levels[this->present_level]
+ ->hexes.face_rotations.size(),
+ ExcInternalError());
+
+ this->tria->levels[this->present_level]
+ ->hexes.face_rotations[this->present_index *
+ GeometryInfo<3>::faces_per_cell
+ +
+ face] = rotation;
+}
+
#endif
// Remark: The following explicit instantiations needed to be moved to
// this place here to work around a problem with gcc3.3 on Apple MacOSX.
ExcMessage ("The neighbor must be on the same level as this cell!"));
Assert (this->neighbor(face)->has_children() == true,
ExcMessage ("The neighbor must have children!"));
-
- static const unsigned int subface_translation[4]
- = { 0, 2, 1, 3 };
- // see whether face and
- // the neighbor's
- // counterface share the
- // same indexing of
- // children. if not so,
- // translate child
- // indices
+ // this function returns they neighbor's
+ // child on a given face and
+ // subface. However, irrespective of the
+ // actual face orientation, flip and rotation
+ // of the cells face it returns that child,
+ // which is neighboring the subface which
+ // would be subface number @p subface, if the
+ // face was in standard orientation. thus,
+ // our first task is to extract the
+ // information, which real_subface_no
+ // corresponds to that for the actual face,
+ // which might have non-standard orientation,
+ // flip and rotation.
+
+ // to this end: the subfaces are numbered
+ // according to the vertices of the
+ // face. thus we can use the
+ // standard_to_real_face_vertex function in
+ // GeometryInfo<dim> to translate the given
+ // subface to the real_subface_no
+ unsigned int real_subface_no=GeometryInfo<3>::
+ standard_to_real_face_vertex(subface,
+ this->face_orientation(face),
+ this->face_flip(face),
+ this->face_rotation(face));
+ // get the neighbor's number for the given
+ // face and the neighbor
const unsigned int neighbor_neighbor
= this->neighbor_of_neighbor (face);
- const bool face_orientations_match
- = (this->neighbor(face)->face_orientation(neighbor_neighbor) ==
- this->face_orientation(face));
+ const TriaIterator<3,CellAccessor<3> > neighbor=this->neighbor(face);
+
+ // now use the info provided by GeometryInfo
+ // to extract the neighbors child number
const unsigned int neighbor_child_index
- = GeometryInfo<3>::child_cell_on_face(neighbor_neighbor, subface,
- face_orientations_match);
+ = GeometryInfo<3>::child_cell_on_face(neighbor_neighbor, real_subface_no,
+ neighbor->face_orientation(neighbor_neighbor),
+ neighbor->face_flip(neighbor_neighbor),
+ neighbor->face_rotation(neighbor_neighbor));
const TriaIterator<3,CellAccessor<3> > neighbor_child=
- this->neighbor(face)->child(neighbor_child_index);
+ neighbor->child(neighbor_child_index);
- // if the face on the side of the present
- // cell is in the correct order, then make
- // sure that the neighbor child cell we
+ // make sure that the neighbor child cell we
// have found shares the desired subface.
- //
- // otherwise: if the face is turned
- // inside out when viewed from the
- // present cell, then the subface we are
- // interested in is not
- // this->face(face)->child(subface), but
- // instead
- // this->face(face)->child(subface'),
- // where subface' is the translated
- // subface number. this is so because we
- // only store the face only once, so its
- // children are defined in terms of its
- // own (circular) orientation, not in
- // terms of face_orientation as viewed
- // from one of the adjacent cells. in
- // that case, cell->face(f)->subface(sf)
- // may yield unexpected results; in that
- // case, a caller may need to adjust
- // according to the face_orientation
- // flag, though one in general only wants
- // to loop over all subfaces, and not
- // pick a particular one
- Assert(((this->face_orientation(face) == true)
- &&
- (this->face(face)->child(subface) ==
- neighbor_child->face(neighbor_neighbor)))
- ||
- ((this->face_orientation(face) == false)
- &&
- (this->face(face)->child(subface_translation[subface]) ==
- neighbor_child->face(neighbor_neighbor))),
+ Assert((this->face(face)->child(real_subface_no) ==
+ neighbor_child->face(neighbor_neighbor)),
ExcInternalError());
return neighbor_child;
new_size * GeometryInfo<3>::faces_per_cell
- face_orientations.size(),
true);
+ face_flips.reserve (new_size * GeometryInfo<3>::faces_per_cell);
+ face_flips.insert (face_flips.end(),
+ new_size * GeometryInfo<3>::faces_per_cell
+ - face_flips.size(),
+ false);
+ face_rotations.reserve (new_size * GeometryInfo<3>::faces_per_cell);
+ face_rotations.insert (face_rotations.end(),
+ new_size * GeometryInfo<3>::faces_per_cell
+ - face_rotations.size(),
+ false);
};
}
+ void
+ TriaObjectsQuad3D::reserve_space (const unsigned int new_quads)
+ {
+ const unsigned int new_size = new_quads +
+ std::count_if (used.begin(),
+ used.end(),
+ std::bind2nd (std::equal_to<bool>(), true));
+
+ // see above...
+ if (new_size>cells.size())
+ {
+ TriaObjects<Quad>::reserve_space(new_quads);
+ line_orientations.reserve (new_size * GeometryInfo<2>::lines_per_cell);
+ line_orientations.insert (line_orientations.end(),
+ new_size * GeometryInfo<2>::lines_per_cell
+ - line_orientations.size(),
+ true);
+ };
+ }
+
+
template<>
void
TriaObjects<Line>::monitor_memory (const unsigned int) const
== face_orientations.size(),
ExcMemoryInexact (cells.size() * GeometryInfo<3>::faces_per_cell,
face_orientations.size()));
+ Assert (cells.size() * GeometryInfo<3>::faces_per_cell
+ == face_flips.size(),
+ ExcMemoryInexact (cells.size() * GeometryInfo<3>::faces_per_cell,
+ face_flips.size()));
+ Assert (cells.size() * GeometryInfo<3>::faces_per_cell
+ == face_rotations.size(),
+ ExcMemoryInexact (cells.size() * GeometryInfo<3>::faces_per_cell,
+ face_rotations.size()));
+ }
+
+
+ void
+ TriaObjectsQuad3D::monitor_memory (const unsigned int) const
+ {
+ // check that we have not allocated
+ // too much memory. note that bool
+ // vectors allocate their memory in
+ // chunks of whole integers, so
+ // they may over-allocate by up to
+ // as many elements as an integer
+ // has bits
+ Assert (cells.size() * GeometryInfo<2>::lines_per_cell
+ == line_orientations.size(),
+ ExcMemoryInexact (cells.size() * GeometryInfo<2>::lines_per_cell,
+ line_orientations.size()));
+ TriaObjects<Quad>::monitor_memory (3);
+
}
{
TriaObjects<Hexahedron>::clear();
face_orientations.clear();
+ face_flips.clear();
+ face_rotations.clear();
+ }
+
+
+ void
+ TriaObjectsQuad3D::clear()
+ {
+ TriaObjects<Quad>::clear();
+ line_orientations.clear();
}
TriaObjectsHex::memory_consumption () const
{
return (MemoryConsumption::memory_consumption (face_orientations) +
- this->TriaObjects<Hexahedron>::memory_consumption() );
+ MemoryConsumption::memory_consumption (face_flips) +
+ MemoryConsumption::memory_consumption (face_rotations) +
+ TriaObjects<Hexahedron>::memory_consumption() );
+ }
+
+
+ unsigned int
+ TriaObjectsQuad3D::memory_consumption () const
+ {
+ return (MemoryConsumption::memory_consumption (line_orientations) +
+ this->TriaObjects<Quad>::memory_consumption() );
}
// explicit instantiations