// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
* Return the RefinementCase
* of this cell.
*/
- RefinementCase<structdim> refinement_case() const;
+ RefinementCase<structdim> refinement_case () const;
/**
* Index of the @p ith child.
* used.
*/
double measure () const;
+
+ /**
+ * Return true if the current object is a
+ * translation of the given argument.
+ *
+ * @note For the purpose of a
+ * triangulation, cells, faces, etc are
+ * only characterized by their
+ * vertices. The current function
+ * therefore only compares the locations
+ * of vertices. For many practical
+ * applications, however, it is not only
+ * the vertices that determine whether
+ * one cell is a translation of another,
+ * but also how the cell is mapped from
+ * the reference cell to its location in
+ * real space. For example, if we are
+ * using higher order mappings, then not
+ * only do the vertices have to be
+ * translations of each other, but also
+ * the points along edges. In these
+ * questions, therefore, it would be
+ * appropriate to ask the mapping, not
+ * the current function, whether two
+ * objects are translations of each
+ * other.
+ */
+ bool
+ is_translation_of (const TriaIterator<TriaAccessor<structdim,dim,spacedim> > &o) const;
+
/**
* @}
*/
// $Id$
// Version: $Name$
//
-// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008 by the deal.II authors
+// Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
}
+template <int structdim, int dim, int spacedim>
+bool
+TriaAccessor<structdim, dim, spacedim>::
+is_translation_of (const TriaIterator<TriaAccessor<structdim,dim,spacedim> > &o) const
+{
+ // go through the vertices and check... The
+ // cell is a translation of the previous
+ // one in case the distance between the
+ // individual vertices in the two cell is
+ // the same for all the vertices. So do the
+ // check by first getting the distance on
+ // the first vertex, and then checking
+ // whether all others have the same down to
+ // rounding errors (we have to be careful
+ // here because the calculation of the
+ // displacement between one cell and the
+ // next can already result in the loss of
+ // one or two digits), so we choose 1e-12
+ // times the distance between the zeroth
+ // vertices here.
+ bool is_translation = true;
+ const Point<spacedim> dist = o->vertex(0) - this->vertex(0);
+ const double tol_square = 1e-24 * dist.norm_square();
+ for (unsigned int i=1; i<GeometryInfo<structdim>::vertices_per_cell; ++i)
+ {
+ const Point<spacedim> dist_new = (o->vertex(i) - this->vertex(i)) - dist;
+ if (dist_new.norm_square() > tol_square)
+ {
+ is_translation = false;
+ break;
+ }
+ }
+ return is_translation;
+}
+
+
+
/*------------------------ Functions: CellAccessor<dim,spacedim> -----------------------*/
// case that there has not been any cell
// before
if (&*this->present_cell == 0)
- {
- cell_similarity = CellSimilarity::none;
- return;
- }
-
- // in MappingQ, data can have been
- // modified during the last call. Then,
- // we can't use that data on the new
- // cell.
- if (cell_similarity == CellSimilarity::invalid_next_cell)
- {
+ cell_similarity = CellSimilarity::none;
+ else
+ // in MappingQ, data can have been
+ // modified during the last call. Then,
+ // we can't use that data on the new
+ // cell.
+ if (cell_similarity == CellSimilarity::invalid_next_cell)
cell_similarity = CellSimilarity::none;
- return;
- }
-
- const typename Triangulation<dim,spacedim>::cell_iterator & present_cell =
- *this->present_cell;
-
- // test for translation
- {
- // otherwise, go through the vertices
- // and check... The cell is a
- // translation of the previous one in
- // case the distance between the
- // individual vertices in the two
- // cell is the same for all the
- // vertices. So do the check by first
- // getting the distance on the first
- // vertex, and then checking whether
- // all others have the same down to
- // rounding errors (we have to be
- // careful here because the
- // calculation of the displacement
- // between one cell and the next can
- // already result in the loss of one
- // or two digits), so we choose 1e-12
- // times the distance between the
- // zeroth vertices here.
- bool is_translation = true;
- const Point<spacedim> dist = cell->vertex(0) - present_cell->vertex(0);
- const double tolerance = 1e-24 * dist.norm_square();
- for (unsigned int i=1; i<GeometryInfo<dim>::vertices_per_cell; ++i)
- {
- Point<spacedim> dist_new = cell->vertex(i) - present_cell->vertex(i) - dist;
- if (dist_new.norm_square() > tolerance)
- {
- is_translation = false;
- break;
- }
- }
-
- cell_similarity = (is_translation == true
- ?
- CellSimilarity::translation
- :
- CellSimilarity::none);
- return;
- }
-
+ else
+ cell_similarity = (cell->is_translation_of
+ (static_cast<const typename Triangulation<dim,spacedim>::cell_iterator>(*this->present_cell))
+ ?
+ CellSimilarity::translation
+ :
+ CellSimilarity::none);
+
// TODO: here, one could implement other
// checks for similarity, e.g. for
// children of a parallelogram.
<h3>deal.II</h3>
<ol>
+ <li>
+ <p>
+ New: The new function TriaAccessor::is_translation_of computes
+ whether a cell, face, or edge is a translation of another.
+ <br>
+ (Martin Kronbichler, WB 2009/05/19)
+ </p>
+ </li>
+
<li>
<p>
New: The DoFTools::make_sparsity_pattern functions have acquired a
(WB 2009/04/29)
</p>
</li>
-
+
<li>
<p>
Fixed: The DoFRenumbering::component_wise function for MGDoFHandler objects