* class chooses them such that the discretized boundary is globally
* continuously differentiable.
*
- * To use this class, make sure that the Boundary::@p get_normals_at_vertices
+ * To use this class, make sure that the Boundary::get_normals_at_vertices()
* function is implemented for the user's boundary object.
*
- * For more information about the <tt>spacedim</tt> template parameter check
- * the documentation of FiniteElement or the one of Triangulation.
- *
- * @author Wolfgang Bangerth, 2001
+ * @author Wolfgang Bangerth, 2001, 2015
*/
template<int dim, int spacedim=dim>
class MappingC1 : public MappingQ<dim,spacedim>
Mapping<dim,spacedim> *clone () const;
protected:
- /**
- * For <tt>dim=2,3</tt>. Append the support points of all shape functions
- * located on bounding lines to the vector @p a. Points located on the line
- * but on vertices are not included.
- *
- * Needed by the <tt>compute_support_points_simple(laplace)</tt> functions.
- * For <tt>dim=1</tt> this function is empty.
- *
- * This function chooses the respective points not such that they are
- * interpolating the boundary (as does the base class), but rather such that
- * the resulting cubic mapping is a continuous one.
- */
- virtual void
- add_line_support_points (const typename Triangulation<dim>::cell_iterator &cell,
- std::vector<Point<dim> > &a) const;
/**
- * For <tt>dim=3</tt>. Append the support points of all shape functions
- * located on bounding faces (quads in 3d) to the vector @p a. Points
- * located on the line but on vertices are not included.
- *
- * Needed by the @p compute_support_points_laplace function. For
- * <tt>dim=1</tt> and 2 this function is empty.
- *
- * This function chooses the respective points not such that they are
- * interpolating the boundary (as does the base class), but rather such that
- * the resulting cubic mapping is a continuous one.
+ * A class derived from MappingQGeneric that provides the generic
+ * mapping with support points on boundary objects so that the
+ * corresponding Q3 mapping ends up being C1.
*/
- virtual void
- add_quad_support_points(const typename Triangulation<dim>::cell_iterator &cell,
- std::vector<Point<dim> > &a) const;
+ class MappingC1Generic : public MappingQGeneric<dim,spacedim>
+ {
+ public:
+
+ /**
+ * Constructor.
+ */
+ MappingC1Generic ();
+
+ /**
+ * For <tt>dim=2,3</tt>. Append the support points of all shape functions
+ * located on bounding lines to the vector @p a. Points located on the line
+ * but on vertices are not included.
+ *
+ * Needed by the <tt>compute_support_points_simple(laplace)</tt> functions.
+ * For <tt>dim=1</tt> this function is empty.
+ *
+ * This function chooses the respective points not such that they are
+ * interpolating the boundary (as does the base class), but rather such that
+ * the resulting cubic mapping is a continuous one.
+ */
+ virtual void
+ add_line_support_points (const typename Triangulation<dim>::cell_iterator &cell,
+ std::vector<Point<dim> > &a) const;
+
+ /**
+ * For <tt>dim=3</tt>. Append the support points of all shape functions
+ * located on bounding faces (quads in 3d) to the vector @p a. Points
+ * located on the line but on vertices are not included.
+ *
+ * Needed by the @p compute_support_points_laplace function. For
+ * <tt>dim=1</tt> and 2 this function is empty.
+ *
+ * This function chooses the respective points not such that they are
+ * interpolating the boundary (as does the base class), but rather such that
+ * the resulting cubic mapping is a continuous one.
+ */
+ virtual void
+ add_quad_support_points(const typename Triangulation<dim>::cell_iterator &cell,
+ std::vector<Point<dim> > &a) const;
+ };
};
/*@}*/
#ifndef DOXYGEN
-template <> void MappingC1<1>::add_line_support_points (
+template <> void MappingC1<1>::MappingC1Generic::add_line_support_points (
const Triangulation<1>::cell_iterator &,
std::vector<Point<1> > &) const;
-template <> void MappingC1<2>::add_line_support_points (
+template <> void MappingC1<2>::MappingC1Generic::add_line_support_points (
const Triangulation<2>::cell_iterator &cell,
std::vector<Point<2> > &a) const;
-template <> void MappingC1<1>::add_quad_support_points (
+template <> void MappingC1<1>::MappingC1Generic::add_quad_support_points (
const Triangulation<1>::cell_iterator &,
std::vector<Point<1> > &) const;
-template <> void MappingC1<2>::add_quad_support_points (
+template <> void MappingC1<2>::MappingC1Generic::add_quad_support_points (
const Triangulation<2>::cell_iterator &,
std::vector<Point<2> > &) const;
* @author Ralf Hartmann, 2000, 2001, 2005; Guido Kanschat 2000, 2001, Wolfgang Bangerth, 2015
*/
template <int dim, int spacedim=dim>
-class MappingQ : public MappingQGeneric<dim,spacedim>
+class MappingQ : public Mapping<dim,spacedim>
{
public:
/**
*/
MappingQ (const MappingQ<dim,spacedim> &mapping);
+ /**
+ * Return the degree of the mapping, i.e. the value which was passed to the
+ * constructor.
+ */
+ unsigned int get_degree () const;
+
+ /**
+ * Always returns @p true because the default implementation of
+ * functions in this class preserves vertex locations.
+ */
+ virtual
+ bool preserves_vertex_locations () const;
+
/**
* Transforms the point @p p on the unit cell to the point @p p_real on the
* real cell @p cell and returns @p p_real.
* also store a pointer to an InternalData object that pertains
* to a $Q_1$ mapping.
*/
- class InternalData : public MappingQGeneric<dim,spacedim>::InternalData
+ class InternalData : public Mapping<dim,spacedim>::InternalDataBase
{
public:
/**
* Constructor.
*/
- InternalData (const unsigned int polynomial_degree);
+ InternalData ();
/**
* $Q_1$ mapping that is, by default, used on all interior cells.
*/
std_cxx11::unique_ptr<typename MappingQGeneric<dim,spacedim>::InternalData> mapping_q1_data;
+
+ /**
+ * A pointer to a structure to store the information for the full
+ * $Q_p$ mapping that is, by default, used on all boundary cells.
+ */
+ std_cxx11::unique_ptr<typename MappingQGeneric<dim,spacedim>::InternalData> mapping_qp_data;
};
protected:
+ // documentation can be found in Mapping::requires_update_flags()
+ virtual
+ UpdateFlags
+ requires_update_flags (const UpdateFlags update_flags) const;
+
// documentation can be found in Mapping::get_data()
virtual
InternalData *
protected:
+ /**
+ * The polynomial degree of the cells to be used on all cells at the
+ * boundary of the domain, or everywhere if so specified.
+ */
+ const unsigned int polynomial_degree;
+
/**
* If this flag is set @p true then @p MappingQ is used on all cells, not
* only on boundary cells.
* displacements. This also means that we really need to store
* our own Q1 mapping here, rather than simply resorting to
* StaticMappingQ1::mapping.
+ *
+ * @note If the polynomial degree used for the current object is one,
+ * then the qp_mapping and q1_mapping variables point to the same
+ * underlying object.
*/
- std_cxx11::unique_ptr<const MappingQGeneric<dim,spacedim> > q1_mapping;
+ std_cxx11::shared_ptr<const MappingQGeneric<dim,spacedim> > q1_mapping;
/**
- * Declare other MappingQ classes friends.
+ * Pointer to a Q_p mapping. This mapping is used on boundary cells unless
+ * use_mapping_q_on_all_cells was set in the call to the
+ * constructor (in which case it is used for all cells).
+ *
+ * @note MappingQEulerian and MappingC1 reset this pointer to an object of
+ * their own implementation to ensure that the Q_p mapping also knows
+ * about the proper shifts and transformations of the Eulerian
+ * displacements (Eulerian case) and proper choice of support
+ * points (C1 case).
+ *
+ * @note If the polynomial degree used for the current object is one,
+ * then the qp_mapping and q1_mapping variables point to the same
+ * underlying object.
*/
- template <int,int> friend class MappingQ;
+ std_cxx11::shared_ptr<const MappingQGeneric<dim,spacedim> > qp_mapping;
};
/*@}*/
mutable Threads::Mutex fe_values_mutex;
/**
- * Compute the positions of the support points in the current configuration.
- * See the documentation of MappingQGeneric::compute_mapping_support_points()
- * for more information.
+ * A class derived from MappingQGeneric that provides the generic
+ * mapping with support points on boundary objects so that the
+ * corresponding Q3 mapping ends up being C1.
*/
- virtual
- std::vector<Point<spacedim> >
- compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const;
+ class MappingQEulerianGeneric : public MappingQGeneric<dim,spacedim>
+ {
+ public:
+
+ /**
+ * Constructor.
+ */
+ MappingQEulerianGeneric (const unsigned int degree,
+ const MappingQEulerian<dim,VECTOR,spacedim> &mapping_q_eulerian);
+
+ /**
+ * Return the mapped vertices of the cell. For the current class, this function does
+ * not use the support points from the geometry of the current cell but
+ * instead evaluates an externally given displacement field in addition to
+ * the geometry of the cell.
+ */
+ virtual
+ std_cxx11::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell>
+ get_vertices (const typename Triangulation<dim,spacedim>::cell_iterator &cell) const;
+
+ /**
+ * Compute the positions of the support points in the current configuration.
+ * See the documentation of MappingQGeneric::compute_mapping_support_points()
+ * for more information.
+ */
+ virtual
+ std::vector<Point<spacedim> >
+ compute_mapping_support_points(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const;
+
+ private:
+ /**
+ * Reference to the surrounding object off of which we live.
+ */
+ const MappingQEulerian<dim,VECTOR,spacedim> &mapping_q_eulerian;
+ };
+
};
/*@}*/
+template <int dim, int spacedim>
+MappingC1<dim,spacedim>::MappingC1Generic::MappingC1Generic ()
+ :
+ MappingQGeneric<dim,spacedim> (3)
+{}
+
+
+
template <int dim, int spacedim>
MappingC1<dim,spacedim>::MappingC1 ()
:
MappingQ<dim,spacedim> (3)
{
Assert (dim > 1, ExcImpossibleInDim(dim));
+
+ // replace the mapping_qp objects of the base class by something
+ // that knows about generating data points based on the geometry
+ //
+ // we only need to replace the Qp mapping because that's the one that's
+ // used on boundary cells where it matters
+ this->qp_mapping.reset (new MappingC1<dim,spacedim>::MappingC1Generic());
}
template <>
void
-MappingC1<1>::add_line_support_points (const Triangulation<1>::cell_iterator &,
- std::vector<Point<1> > &) const
+MappingC1<1>::MappingC1Generic::add_line_support_points (const Triangulation<1>::cell_iterator &,
+
std::vector<Point<1> > &) const
{
const unsigned int dim = 1;
(void)dim;
template <>
void
-MappingC1<2>::add_line_support_points (const Triangulation<2>::cell_iterator &cell,
- std::vector<Point<2> > &a) const
+MappingC1<2>::MappingC1Generic::add_line_support_points (const Triangulation<2>::cell_iterator &cell,
+
std::vector<Point<2> > &a) const
{
const unsigned int dim = 2;
std::vector<Point<dim> > line_points (2);
static const StraightBoundary<dim> straight_boundary;
straight_boundary.get_intermediate_points_on_line (line, line_points);
a.insert (a.end(), line_points.begin(), line_points.end());
- };
- };
+ }
+ }
}
template<int dim, int spacedim>
void
-MappingC1<dim,spacedim>::add_line_support_points (const typename Triangulation<dim>::cell_iterator &,
-
std::vector<Point<dim> > &) const
+MappingC1<dim,spacedim>::MappingC1Generic::add_line_support_points (const typename Triangulation<dim>::cell_iterator &,
+ std::vector<Point<dim> > &) const
{
Assert (false, ExcNotImplemented());
}
template <>
void
-MappingC1<1>::add_quad_support_points (const Triangulation<1>::cell_iterator &,
- std::vector<Point<1> > &) const
+MappingC1<1>::MappingC1Generic::add_quad_support_points (const Triangulation<1>::cell_iterator &,
+
std::vector<Point<1> > &) const
{
const unsigned int dim = 1;
(void)dim;
template <>
void
-MappingC1<2>::add_quad_support_points (const Triangulation<2>::cell_iterator &,
- std::vector<Point<2> > &) const
+MappingC1<2>::MappingC1Generic::add_quad_support_points (const Triangulation<2>::cell_iterator &,
+
std::vector<Point<2> > &) const
{
const unsigned int dim = 2;
(void)dim;
template<int dim, int spacedim>
void
-MappingC1<dim,spacedim>::add_quad_support_points (const typename Triangulation<dim>::cell_iterator &,
-
std::vector<Point<dim> > &) const
+MappingC1<dim,spacedim>::MappingC1Generic::add_quad_support_points (const typename Triangulation<dim>::cell_iterator &,
+ std::vector<Point<dim> > &) const
{
Assert (false, ExcNotImplemented());
}
template<int dim, int spacedim>
-MappingQ<dim,spacedim>::InternalData::InternalData (const unsigned int polynomial_degree)
+MappingQ<dim,spacedim>::InternalData::InternalData ()
:
- MappingQGeneric<dim,spacedim>::InternalData(polynomial_degree),
use_mapping_q1_on_current_cell(false)
{}
std::size_t
MappingQ<dim,spacedim>::InternalData::memory_consumption () const
{
- return (MappingQGeneric<dim,spacedim>::InternalData::memory_consumption () +
+ return (Mapping<dim,spacedim>::InternalDataBase::memory_consumption () +
MemoryConsumption::memory_consumption (use_mapping_q1_on_current_cell) +
- MemoryConsumption::memory_consumption (mapping_q1_data));
+ MemoryConsumption::memory_consumption (mapping_q1_data) +
+ MemoryConsumption::memory_consumption (mapping_qp_data));
}
MappingQ<dim,spacedim>::MappingQ (const unsigned int degree,
const bool use_mapping_q_on_all_cells)
:
- MappingQGeneric<dim,spacedim>(degree),
+ polynomial_degree (degree),
// see whether we want to use *this* mapping objects on *all* cells,
// or defer to an explicit Q1 mapping on interior cells. if
(dim != spacedim)),
// create a Q1 mapping for use on interior cells (if necessary)
// or to create a good initial guess in transform_real_to_unit_cell()
- q1_mapping (new MappingQGeneric<dim,spacedim>(1))
+ q1_mapping (new MappingQGeneric<dim,spacedim>(1)),
+
+ // create a Q_p mapping; if p=1, simply share the Q_1 mapping already
+ // created via the shared_ptr objects
+ qp_mapping (this->polynomial_degree>1
+ ?
+ std_cxx11::shared_ptr<const MappingQGeneric<dim,spacedim> >(new MappingQGeneric<dim,spacedim>(degree))
+ :
+ q1_mapping)
{}
template<int dim, int spacedim>
MappingQ<dim,spacedim>::MappingQ (const MappingQ<dim,spacedim> &mapping)
:
- MappingQGeneric<dim,spacedim>(mapping.get_degree()),
+ polynomial_degree (mapping.polynomial_degree),
use_mapping_q_on_all_cells (mapping.use_mapping_q_on_all_cells),
// clone the Q1 mapping for use on interior cells (if necessary)
// or to create a good initial guess in transform_real_to_unit_cell()
- q1_mapping (dynamic_cast<MappingQGeneric<dim,spacedim>*>(mapping.q1_mapping->clone()))
+ q1_mapping (dynamic_cast<MappingQGeneric<dim,spacedim>*>(mapping.q1_mapping->clone())),
+ // create a Q_p mapping; if p=1, simply share the Q_1 mapping already
+ // created via the shared_ptr objects
+ qp_mapping (this->polynomial_degree>1
+ ?
+ std_cxx11::shared_ptr<const MappingQGeneric<dim,spacedim> >(dynamic_cast<MappingQGeneric<dim,spacedim>*>(mapping.qp_mapping->clone()))
+ :
+ q1_mapping)
{}
+template<int dim, int spacedim>
+unsigned int
+MappingQ<dim,spacedim>::get_degree() const
+{
+ return polynomial_degree;
+}
+
+
+
+template <int dim, int spacedim>
+inline
+bool
+MappingQ<dim,spacedim>::preserves_vertex_locations () const
+{
+ return true;
+}
+
+
+
+template<int dim, int spacedim>
+UpdateFlags
+MappingQ<dim,spacedim>::requires_update_flags (const UpdateFlags in) const
+{
+ return (q1_mapping->requires_update_flags(in)
+ |
+ qp_mapping->requires_update_flags(in));
+}
+
+
+
template<int dim, int spacedim>
typename MappingQ<dim,spacedim>::InternalData *
MappingQ<dim,spacedim>::get_data (const UpdateFlags update_flags,
const Quadrature<dim> &quadrature) const
{
- InternalData *data = new InternalData(this->polynomial_degree);
-
- // fill the data of both the Q_p and the Q_1 objects in parallel
- Threads::Task<>
- task = Threads::new_task (std_cxx11::function<void()>
- (std_cxx11::bind(&MappingQGeneric<dim,spacedim>::InternalData::initialize,
- data,
- update_flags,
- std_cxx11::cref(quadrature),
- quadrature.size())));
+ InternalData *data = new InternalData;
+
+ // build the Q1 and Qp internal data objects in parallel
+ Threads::Task<typename MappingQGeneric<dim,spacedim>::InternalData *>
+ do_get_data = Threads::new_task (&MappingQGeneric<dim,spacedim>::get_data,
+ *qp_mapping,
+ update_flags,
+ quadrature);
+
if (!use_mapping_q_on_all_cells)
data->mapping_q1_data.reset (q1_mapping->get_data (update_flags, quadrature));
- task.join ();
+ // wait for the task above to finish and use returned value
+ data->mapping_qp_data.reset (do_get_data.return_value());
return data;
}
MappingQ<dim,spacedim>::get_face_data (const UpdateFlags update_flags,
const Quadrature<dim-1>& quadrature) const
{
- InternalData *data = new InternalData(this->polynomial_degree);
- const Quadrature<dim> q (QProjector<dim>::project_to_all_faces(quadrature));
-
- // fill the data of both the Q_p and the Q_1 objects in parallel
- Threads::Task<>
- task = Threads::new_task (std_cxx11::function<void()>
- (std_cxx11::bind(&MappingQGeneric<dim,spacedim>::InternalData::initialize_face,
- data,
- update_flags,
- std_cxx11::cref(q),
- quadrature.size())));
+ InternalData *data = new InternalData;
+
+ // build the Q1 and Qp internal data objects in parallel
+ Threads::Task<typename MappingQGeneric<dim,spacedim>::InternalData *>
+ do_get_data = Threads::new_task (&MappingQGeneric<dim,spacedim>::get_face_data,
+ *qp_mapping,
+ update_flags,
+ quadrature);
+
if (!use_mapping_q_on_all_cells)
data->mapping_q1_data.reset (q1_mapping->get_face_data (update_flags, quadrature));
- task.join ();
+ // wait for the task above to finish and use returned value
+ data->mapping_qp_data.reset (do_get_data.return_value());
return data;
}
MappingQ<dim,spacedim>::get_subface_data (const UpdateFlags update_flags,
const Quadrature<dim-1>& quadrature) const
{
- InternalData *data = new InternalData(this->polynomial_degree);
- const Quadrature<dim> q (QProjector<dim>::project_to_all_subfaces(quadrature));
-
- // fill the data of both the Q_p and the Q_1 objects in parallel
- Threads::Task<>
- task = Threads::new_task (std_cxx11::function<void()>
- (std_cxx11::bind(&MappingQGeneric<dim,spacedim>::InternalData::initialize_face,
- data,
- update_flags,
- std_cxx11::cref(q),
- quadrature.size())));
+ InternalData *data = new InternalData;
+
+ // build the Q1 and Qp internal data objects in parallel
+ Threads::Task<typename MappingQGeneric<dim,spacedim>::InternalData *>
+ do_get_data = Threads::new_task (&MappingQGeneric<dim,spacedim>::get_subface_data,
+ *qp_mapping,
+ update_flags,
+ quadrature);
+
if (!use_mapping_q_on_all_cells)
data->mapping_q1_data.reset (q1_mapping->get_subface_data (update_flags, quadrature));
- task.join ();
+ // wait for the task above to finish and use returned value
+ data->mapping_qp_data.reset (do_get_data.return_value());
return data;
}
*data.mapping_q1_data,
output_data);
else
- MappingQGeneric<dim,spacedim>::fill_fe_values(cell,
- updated_cell_similarity,
- quadrature,
- data,
- output_data);
+ qp_mapping->fill_fe_values(cell,
+ updated_cell_similarity,
+ quadrature,
+ *data.mapping_qp_data,
+ output_data);
return updated_cell_similarity;
}
*data.mapping_q1_data,
output_data);
else
- MappingQGeneric<dim,spacedim>::fill_fe_face_values(cell,
- face_no,
- quadrature,
- data,
- output_data);
+ qp_mapping->fill_fe_face_values(cell,
+ face_no,
+ quadrature,
+ *data.mapping_qp_data,
+ output_data);
}
*data.mapping_q1_data,
output_data);
else
- MappingQGeneric<dim,spacedim>::fill_fe_subface_values(cell,
- face_no,
- subface_no,
- quadrature,
- data,
- output_data);
+ qp_mapping->fill_fe_subface_values(cell,
+ face_no,
+ subface_no,
+ quadrature,
+ *data.mapping_qp_data,
+ output_data);
}
// check whether we should in fact work on the Q1 portion of it
if (data->use_mapping_q1_on_current_cell)
- return q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
+ q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
else
// otherwise use the full mapping
- MappingQGeneric<dim,spacedim>::transform(input, mapping_type, mapping_data, output);
+ qp_mapping->transform(input, mapping_type, *data->mapping_qp_data, output);
}
// check whether we should in fact work on the Q1 portion of it
if (data->use_mapping_q1_on_current_cell)
- return q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
+ q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
else
// otherwise use the full mapping
- MappingQGeneric<dim,spacedim>::transform(input, mapping_type, mapping_data, output);
+ qp_mapping->transform(input, mapping_type, *data->mapping_qp_data, output);
}
// check whether we should in fact work on the Q1 portion of it
if (data->use_mapping_q1_on_current_cell)
- return q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
+ q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
else
// otherwise use the full mapping
- MappingQGeneric<dim,spacedim>::transform(input, mapping_type, mapping_data, output);
+ qp_mapping->transform(input, mapping_type, *data->mapping_qp_data, output);
}
// check whether we should in fact work on the Q1 portion of it
if (data->use_mapping_q1_on_current_cell)
- return q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
+ q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
else
// otherwise use the full mapping
- MappingQGeneric<dim,spacedim>::transform(input, mapping_type, mapping_data, output);
+ qp_mapping->transform(input, mapping_type, *data->mapping_qp_data, output);
}
// check whether we should in fact work on the Q1 portion of it
if (data->use_mapping_q1_on_current_cell)
- return q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
+ q1_mapping->transform (input, mapping_type, *data->mapping_q1_data, output);
else
// otherwise use the full mapping
- MappingQGeneric<dim,spacedim>::transform(input, mapping_type, mapping_data, output);
+ qp_mapping->transform(input, mapping_type, *data->mapping_qp_data, output);
}
// mapping, then either use our own facilities or that of the Q1
// mapping we store
if (use_mapping_q_on_all_cells || cell->has_boundary_lines())
- return this->MappingQGeneric<dim,spacedim>::transform_unit_to_real_cell (cell, p);
+ return qp_mapping->transform_unit_to_real_cell (cell, p);
else
return q1_mapping->transform_unit_to_real_cell (cell, p);
}
use_mapping_q_on_all_cells
||
(dim!=spacedim) )
- return MappingQGeneric<dim,spacedim>::transform_real_to_unit_cell(cell, p);
+ return qp_mapping->transform_real_to_unit_cell(cell, p);
else
return q1_mapping->transform_real_to_unit_cell(cell, p);
}
// .... MAPPING Q EULERIAN CONSTRUCTOR
+template <int dim, class EulerVectorType, int spacedim>
+MappingQEulerian<dim, EulerVectorType, spacedim>::MappingQEulerianGeneric::
+MappingQEulerianGeneric (const unsigned int degree,
+ const MappingQEulerian<dim, EulerVectorType, spacedim> &mapping_q_eulerian)
+ :
+ MappingQGeneric<dim,spacedim>(degree),
+ mapping_q_eulerian (mapping_q_eulerian)
+{}
+
template <int dim, class EulerVectorType, int spacedim>
MappingQEulerian<dim, EulerVectorType, spacedim>::
MappingQEulerian (const unsigned int degree,
// current vector
this->q1_mapping.reset (new MappingQ1Eulerian<dim,EulerVectorType,spacedim>(euler_vector,
euler_dof_handler));
+
+ // also reset the qp mapping pointer with our own class
+ this->qp_mapping.reset (new MappingQEulerianGeneric(degree,*this));
}
// current vector
this->q1_mapping.reset (new MappingQ1Eulerian<dim,EulerVectorType,spacedim>(euler_vector,
euler_dof_handler));
+
+ // also reset the qp mapping pointer with our own class
+ this->qp_mapping.reset (new MappingQEulerianGeneric(degree,*this));
}
(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const
{
// get the vertices as the first 2^dim mapping support points
- const std::vector<Point<spacedim> > a = compute_mapping_support_points(cell);
+ const std::vector<Point<spacedim> > a
+ = dynamic_cast<const MappingQEulerianGeneric &>(*this->qp_mapping).compute_mapping_support_points(cell);
std_cxx11::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell> vertex_locations;
std::copy (a.begin(),
+template <int dim, class EulerVectorType, int spacedim>
+std_cxx11::array<Point<spacedim>, GeometryInfo<dim>::vertices_per_cell>
+MappingQEulerian<dim, EulerVectorType, spacedim>::MappingQEulerianGeneric::
+get_vertices
+(const typename Triangulation<dim,spacedim>::cell_iterator &cell) const
+{
+ return mapping_q_eulerian.get_vertices (cell);
+}
+
+
+
+
template <int dim, class EulerVectorType, int spacedim>
std::vector<Point<spacedim> >
-MappingQEulerian<dim, EulerVectorType, spacedim>::
+MappingQEulerian<dim, EulerVectorType, spacedim>::MappingQEulerianGeneric::
compute_mapping_support_points (const typename Triangulation<dim,spacedim>::cell_iterator &cell) const
{
// first, basic assertion with respect to vector size,
- const types::global_dof_index n_dofs = euler_dof_handler->n_dofs();
- const types::global_dof_index vector_size = euler_vector->size();
+ const types::global_dof_index n_dofs = mapping_q_eulerian.euler_dof_handler->n_dofs();
+ const types::global_dof_index vector_size = mapping_q_eulerian.euler_vector->size();
(void)n_dofs;
(void)vector_size;
// we then transform our tria iterator into a dof iterator so we can access
// data not associated with triangulations
- typename DoFHandler<dim,spacedim>::cell_iterator dof_cell(*cell, euler_dof_handler);
+ typename DoFHandler<dim,spacedim>::cell_iterator dof_cell(*cell,
+ mapping_q_eulerian.euler_dof_handler);
Assert (dof_cell->active() == true, ExcInactiveCell());
// that. this implies that the user should order components appropriately,
// or create a separate dof handler for the displacements.
- const unsigned int n_support_pts = support_quadrature.size();
- const unsigned int n_components = euler_dof_handler->get_fe().n_components();
+ const unsigned int n_support_pts = mapping_q_eulerian.support_quadrature.size();
+ const unsigned int n_components = mapping_q_eulerian.euler_dof_handler->get_fe().n_components();
Assert (n_components >= spacedim, ExcDimensionMismatch(n_components, spacedim) );
// fill shift vector for each support point using an fe_values object. make
// sure that the fe_values variable isn't used simultaneously from different
// threads
- Threads::Mutex::ScopedLock lock(fe_values_mutex);
- fe_values.reinit(dof_cell);
- fe_values.get_function_values(*euler_vector, shift_vector);
+ Threads::Mutex::ScopedLock lock(mapping_q_eulerian.fe_values_mutex);
+ mapping_q_eulerian.fe_values.reinit(dof_cell);
+ mapping_q_eulerian.fe_values.get_function_values(*mapping_q_eulerian.euler_vector, shift_vector);
// and finally compute the positions of the support points in the deformed
// configuration.
std::vector<Point<spacedim> > a(n_support_pts);
for (unsigned int q=0; q<n_support_pts; ++q)
{
- a[q] = fe_values.quadrature_point(q);
+ a[q] = mapping_q_eulerian.fe_values.quadrature_point(q);
for (unsigned int d=0; d<spacedim; ++d)
a[q](d) += shift_vector[q](d);
}
fe_q(dim == 3 ? new FE_Q<dim>(this->polynomial_degree) : 0),
support_point_weights_on_quad (compute_support_point_weights_on_quad<dim>(this->polynomial_degree)),
support_point_weights_on_hex (compute_support_point_weights_on_hex<dim>(this->polynomial_degree))
-{}
+{
+ Assert (p >= 1, ExcMessage ("It only makes sense to create polynomial mappings "
+ "with a polynomial degree greater or equal to one."));
+}
projects / dealii.git / commitdiff