virtual
Mapping<dim,spacedim> * clone () const = 0;
+ /**
+ * Returns whether the mapping preserves
+ * vertex locations. Returns @p true for
+ * MappingQ, MappingQ1, MappingCartesian,
+ * but @p false for MappingQEulerian,
+ * MappingQ1Eulerian.
+ */
+ virtual
+ bool preserves_vertex_locations () const = 0;
+
/**
* Exception
*/
* along the coordinate directions. It is specifically developed for
* cartesian meshes. Apply this mapping to a general mesh to get
* strange results.
- *
+ *
* For more information about the <tt>spacedim</tt> template parameter
* check the documentation of FiniteElement or the one of
* Triangulation.
transform (const VectorSlice<const std::vector<Tensor<2,dim> > > input,
VectorSlice<std::vector<Tensor<2,spacedim> > > output,
const typename Mapping<dim,spacedim>::InternalDataBase &internal,
- const MappingType type) const;
-
+ const MappingType type) const;
+
virtual Point<spacedim>
transform_unit_to_real_cell (
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
transform_real_to_unit_cell (
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const Point<spacedim> &p) const;
-
+
/**
* Return a pointer to a copy of the
*/
virtual
Mapping<dim, spacedim> * clone () const;
-
+
+ /**
+ * Always returns @p true because
+ * MappingCartesian preserves vertex
+ * locations.
+ */
+ bool preserves_vertex_locations () const;
+
protected:
- /**
+ /**
* Storage for internal data of
* the scaling.
*/
* The volume element
*/
double volume_element;
-
+
/**
* Vector of all quadrature
* points. Especially, all
*/
std::vector<Point<dim> > quadrature_points;
};
-
+
/**
* Do the computation for the
* <tt>fill_*</tt> functions.
std::vector<Point<dim> >& normal_vectors) const;
private:
- virtual UpdateFlags update_once (const UpdateFlags) const;
+ virtual UpdateFlags update_once (const UpdateFlags) const;
virtual UpdateFlags update_each (const UpdateFlags) const;
-
+
/**
* Value to indicate that a given
* face or subface number is
* invalid.
*/
- static const unsigned int invalid_face_number = numbers::invalid_unsigned_int;
+ static const unsigned int invalid_face_number = numbers::invalid_unsigned_int;
};
/*@}*/
/* -------------- declaration of explicit specializations ------------- */
+#ifndef DOXYGEN
+
+template <int dim, int spacedim>
+inline
+bool
+MappingCartesian<dim,spacedim>::preserves_vertex_locations () const
+{
+ return true;
+}
+
+#endif // DOXYGEN
+
DEAL_II_NAMESPACE_CLOSE
#endif
* For more details about Qp-mappings, see the `mapping' report at
* <tt>deal.II/doc/reports/mapping_q/index.html</tt> in the `Reports'
* section of `Documentation'.
- *
+ *
* For more information about the <tt>spacedim</tt> template parameter
* check the documentation of FiniteElement or the one of
* Triangulation.
* Destructor.
*/
virtual ~MappingQ ();
-
+
/**
* Transforms the point @p p on
* the unit cell to the point
transform_unit_to_real_cell (
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const Point<dim> &p) const;
-
+
/**
* Transforms the point @p p on
* the real cell to the point
VectorSlice<std::vector<Tensor<2,spacedim> > > output,
const typename Mapping<dim,spacedim>::InternalDataBase &internal,
const MappingType type) const;
-
+
/**
* Return the degree of the
* mapping, i.e. the value which
*/
virtual
Mapping<dim,spacedim> * clone () const;
-
- /**
+
+ /**
* Storage for internal data of
* Q_degree transformation.
*/
* Constructor.
*/
InternalData (const unsigned int n_shape_functions);
-
+
/**
* Return an estimate (in
* used.
*/
bool use_mapping_q1_on_current_cell;
-
+
/**
* On interior cells
* @p MappingQ1 is used.
virtual void
add_quad_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell,
std::vector<Point<dim> > &a) const;
-
+
private:
-
+
virtual
typename Mapping<dim,spacedim>::InternalDataBase *
get_data (const UpdateFlags,
typename Mapping<dim,spacedim>::InternalDataBase *
get_subface_data (const UpdateFlags flags,
const Quadrature<dim-1>& quadrature) const;
-
+
/**
* Compute shape values and/or
* derivatives.
*/
void
set_laplace_on_quad_vector(Table<2,double> &loqvs) const;
-
+
/**
* This function is needed by the
* constructor of <tt>MappingQ<3></tt>.
* of the `mapping' report.
*/
void set_laplace_on_hex_vector(Table<2,double> &lohvs) const;
-
+
/**
* Computes the
* <tt>laplace_on_quad(hex)_vector</tt>.
*/
void apply_laplace_vector(const Table<2,double> &lvs,
std::vector<Point<dim> > &a) const;
-
+
/**
* Computes the support points of
* the mapping.
void compute_support_points_laplace(
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
std::vector<Point<dim> > &a) const;
-
+
/**
* Needed by the
* @p laplace_on_quad function
* `mapping' report.
*/
Table<2,double> laplace_on_quad_vector;
-
+
/**
* Needed by the
* @p laplace_on_hex function
DeclException1 (ExcLaplaceVectorNotSet,
int,
<< "laplace_vector not set for degree=" << arg1 << ".");
-
+
/**
* Degree @p p of the
* polynomials used as shape
* functions for the Qp mapping
* of cells at the boundary.
- */
+ */
const unsigned int degree;
/**
* functions on the boundary.
*/
const unsigned int n_outer;
-
+
/**
* Pointer to the
* @p dim-dimensional tensor
* boundary.
*/
const TensorProductPolynomials<dim> *tensor_pols;
-
+
/**
* Number of the Qp tensor
* product shape functions.
* Shape function for this mapping are the same as for the finite
* element FE_Q of order 1. Therefore, coupling these two yields
* an isoparametric element.
- *
+ *
* For more information about the <tt>spacedim</tt> template parameter
* check the documentation of FiniteElement or the one of
* Triangulation.
* Default constructor.
*/
MappingQ1 ();
-
+
virtual Point<spacedim>
transform_unit_to_real_cell (
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const Point<dim> &p) const;
-
+
/**
* Transforms the point @p p on
* the real cell to the point
transform_real_to_unit_cell (
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const Point<spacedim> &p) const;
-
+
virtual void
transform (const VectorSlice<const std::vector<Tensor<1,dim> > > input,
VectorSlice<std::vector<Tensor<1,spacedim> > > output,
VectorSlice<std::vector<Tensor<2,spacedim> > > output,
const typename Mapping<dim,spacedim>::InternalDataBase &internal,
const MappingType type) const;
-
+
/**
* Return a pointer to a copy of the
* present object. The caller of this
virtual
Mapping<dim,spacedim> * clone () const;
- /**
+ /**
* Storage for internal data of
* d-linear transformation.
*/
*/
double shape (const unsigned int qpoint,
const unsigned int shape_nr) const;
-
+
/**
* Shape function at quadrature
* point. See above.
*/
double &shape (const unsigned int qpoint,
const unsigned int shape_nr);
-
+
/**
* Gradient of shape function
* in quadrature point. See
* object.
*/
virtual unsigned int memory_consumption () const;
-
+
/**
* Values of shape
* functions. Access by
* Computed once.
*/
std::vector<double> shape_values;
-
+
/**
* Values of shape function
* derivatives. Access by
* Computed once.
*/
std::vector<Tensor<1,dim> > shape_derivatives;
-
+
/**
* Values of shape function
* second derivatives. Access
* Computed once.
*/
std::vector<Tensor<2,dim> > shape_second_derivatives;
-
+
/**
* Tensors of covariant
* transformation at each of
* Computed on each cell.
*/
std::vector<Tensor<2,spacedim> > covariant;
-
+
/**
* Tensors of contravariant
* transformation at each of
* Computed on each cell.
*/
std::vector<Tensor<2,spacedim> > contravariant;
-
+
/**
* Unit tangential vectors. Used
* for the computation of
* Filled once.
*/
std::vector<std::vector<Tensor<1,spacedim> > > unit_tangentials;
-
+
/**
* Auxiliary vectors for internal use.
*/
* the @p cell_of_current_support_points.
*/
std::vector<Point<spacedim> > mapping_support_points;
-
+
/**
* Stores the cell of which the
* @p mapping_support_points are
* stored.
*/
typename Triangulation<dim,spacedim>::cell_iterator cell_of_current_support_points;
-
+
/**
* Default value of this flag
* is @p true. If <tt>*this</tt>
typedef
typename QProjector<dim>::DataSetDescriptor
DataSetDescriptor;
-
+
/**
* Implementation of the interface in
* Mapping.
std::vector<double> &JxW_values,
typename std::vector<Tensor<1,dim> > &boundary_form,
typename std::vector<Point<spacedim> > &normal_vectors) const ;
-
+
/**
* Compute shape values and/or
* derivatives.
const Quadrature<dim> &quadrature,
const unsigned int n_orig_q_points,
InternalData &data) const;
-
+
/**
* Do the computation for the
* <tt>fill_*</tt> functions.
const CellSimilarity::Similarity cell_similarity,
InternalData &data,
std::vector<Point<spacedim> > &quadrature_points) const;
-
+
/**
* Do the computation for the
* <tt>fill_*</tt> functions.
std::vector<double> &JxW_values,
std::vector<Tensor<1,dim> > &boundary_form,
std::vector<Point<spacedim> > &normal_vectors) const;
-
+
/**
* Compute shape values and/or
* derivatives.
* support points.
*/
Point<spacedim> transform_unit_to_real_cell_internal (const InternalData &mdata) const;
-
+
/**
* Transforms the point @p p on
* the real cell to the point
InternalData &mdata,
Point<dim> &p_unit) const;
+ /**
+ * Always returns @p true because
+ * MappingQ1 preserves vertex locations.
+ */
+ virtual
+ bool preserves_vertex_locations () const;
+
private:
/**
* Implementation of the interface in
* </ul>
*/
virtual UpdateFlags update_once (const UpdateFlags flags) const;
-
+
/**
* Implementation of the interface in
* Mapping.
}
+
+template <int dim, int spacedim>
+inline
+bool
+MappingQ1<dim,spacedim>::preserves_vertex_locations () const
+{
+ return true;
+}
+
#endif // DOXYGEN
/* -------------- declaration of explicit specializations ------------- */
#include <base/config.h>
#include <base/smartpointer.h>
-#include <fe/mapping_q1.h>
+#include <fe/mapping_q1.h>
DEAL_II_NAMESPACE_OPEN
* functions. Each cell is thus shifted in space by values given to
* the mapping through a finite element field.
*
- * <h3>Usage</h3>
+ * <h3>Usage</h3>
*
* The constructor of this class takes two arguments: a reference to
* the vector that defines the mapping from the reference
* configuration to the current configuration and a reference to the
* DoFHandler. The vector should then represent a (flattened out
- * version of a) vector valued field defined at nodes defined by the
+ * version of a) vector valued field defined at nodes defined by the
* the DoFHandler, where the number of components of the vector
- * field equals the number of space dimensions. Thus, the
- * DoFHandler shall operate on a finite element that has as many
- * components as space dimensions. As an additional requirement, we
+ * field equals the number of space dimensions. Thus, the
+ * DoFHandler shall operate on a finite element that has as many
+ * components as space dimensions. As an additional requirement, we
* impose that it have as many degree of freedom per vertex as there
- * are space dimensions; since this object only evaluates the finite
+ * are space dimensions; since this object only evaluates the finite
* element field at the vertices, the values
* of all other degrees of freedom (not associated to vertices) are
* ignored. These requirements are met if the finite element which the
* Not specifying this template argument in applications using the PETSc
* vector classes leads to the construction of a copy of the vector
* which is not acccessible afterwards!
- *
+ *
* For more information about the <tt>spacedim</tt> template parameter
* check the documentation of FiniteElement or the one of
* Triangulation.
*/
MappingQ1Eulerian (const VECTOR &euler_transform_vectors,
const DoFHandler<dim,spacedim> &shiftmap_dof_handler);
-
+
/**
* Return a pointer to a copy of the
* present object. The caller of this
virtual
Mapping<dim,spacedim> * clone () const;
-
+ /**
+ * Always returns @p false because
+ * MappingQ1Eulerian does not preserve
+ * vertex locations.
+ */
+ bool preserves_vertex_locations () const;
+
+
/**
* Exception
*/
DeclException0 (ExcWrongNoOfComponents);
-
+
/**
* Exception.
*/
* shifts.
*/
const VECTOR &euler_transform_vectors;
-
+
/**
* Pointer to the DoFHandler to
* which the mapping vector is
* associated.
*/
const SmartPointer<const DoFHandler<dim,spacedim>,MappingQ1Eulerian<dim,VECTOR,spacedim> > shiftmap_dof_handler;
-
- private:
+
+ private:
/**
* Computes the support points of
* the mapping. For
virtual void compute_mapping_support_points(
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
std::vector<Point<spacedim> > &a) const;
-
+
};
/*@}*/
+/*----------------------------------------------------------------------*/
+
+#ifndef DOXYGEN
+
+template <int dim, class VECTOR, int spacedim>
+inline
+bool
+MappingQ1Eulerian<dim,VECTOR,spacedim>::preserves_vertex_locations () const
+{
+ return false;
+}
+
+#endif // DOXYGEN
+
DEAL_II_NAMESPACE_CLOSE
#endif
// Actually, we compute the
// inverse of the reordering
// vector, called reverse here.
- // Later, irs inverse is computed
+ // Later, its inverse is computed
// into new_indices, which is the
// return argument.
DataOutBase::Patch<DH::dimension, DH::space_dimension> &patch,
const CurvedCellRegion curved_cell_region)
{
- // use ucd_to_deal map as patch
- // vertices are in the old,
- // unnatural ordering
+ // use ucd_to_deal map as patch vertices
+ // are in the old, unnatural ordering. if
+ // the mapping does not preserve locations
+ // (e.g. MappingQEulerian), we need to
+ // compute the offset of the vertex for the
+ // graphical output. Otherwise, we can just
+ // use the vertex info.
for (unsigned int vertex=0; vertex<GeometryInfo<DH::dimension>::vertices_per_cell; ++vertex)
- patch.vertices[vertex] = data.mapping_collection[0].transform_unit_to_real_cell
- (cell_and_index->first,
- GeometryInfo<DH::dimension>::unit_cell_vertex (vertex));
+ if (data.mapping_collection[0].preserves_vertex_locations())
+ patch.vertices[vertex] = cell_and_index->first->vertex(vertex);
+ else
+ patch.vertices[vertex] = data.mapping_collection[0].transform_unit_to_real_cell
+ (cell_and_index->first,
+ GeometryInfo<DH::dimension>::unit_cell_vertex (vertex));
if (data.n_datasets > 0)
{