* called for the same cell first!
*/
virtual void
- fill_fe_values (const Mapping<dim,spacedim> &mapping,
+ fill_fe_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
- typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
- typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
- FEValuesData<dim,spacedim> &data) const = 0;
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
+ typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
+ FEValuesData<dim,spacedim> &data) const = 0;
/**
* Fill the fields of
virtual void
fill_fe_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
FEValuesData<dim,spacedim>& data) const;
* FiniteElement.
*/
virtual void
- fill_fe_values (const Mapping<dim,spacedim> &mapping,
+ fill_fe_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
- typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
- typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
+ typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
FEValuesData<dim,spacedim>& data) const;
/**
const Quadrature<dim>& quadrature) const ;
virtual void
- fill_fe_values (const Mapping<dim,spacedim> &mapping,
+ fill_fe_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
- FEValuesData<dim,spacedim>& data) const;
+ FEValuesData<dim,spacedim>& data) const;
virtual void
fill_fe_face_values (const Mapping<dim,spacedim> &mapping,
template <class POLY, int dim, int spacedim>
typename Mapping<dim,spacedim>::InternalDataBase *
FE_Poly<POLY,dim,spacedim>::get_data (const UpdateFlags update_flags,
- const Mapping<dim,spacedim> &mapping,
- const Quadrature<dim> &quadrature) const
+ const Mapping<dim,spacedim> &mapping,
+ const Quadrature<dim> &quadrature) const
{
// generate a new data object and
// initialize some fields
template <>
inline
void
-FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_values (const Mapping<1,2> &mapping,
- const Triangulation<1,2>::cell_iterator &/*cell*/,
- const Quadrature<1> &quadrature,
- Mapping<1,2>::InternalDataBase &mapping_data,
- Mapping<1,2>::InternalDataBase &fedata,
- FEValuesData<1,2> &data) const
+FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_values
+ (const Mapping<1,2> &mapping,
+ const Triangulation<1,2>::cell_iterator &/*cell*/,
+ const Quadrature<1> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ Mapping<1,2>::InternalDataBase &mapping_data,
+ Mapping<1,2>::InternalDataBase &fedata,
+ FEValuesData<1,2> &data) const
{
// convert data object to internal
// data for this class. fails with
for (unsigned int i=0; i<quadrature.size(); ++i)
data.shape_values(k,i) = fe_data.shape_values[k][i];
- if (flags & update_gradients)
+ if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
mapping_data, mapping_covariant);
}
// const typename QProjector<1>::DataSetDescriptor dsd;
- if (flags & update_hessians)
+ if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
{
AssertThrow(false, ExcNotImplemented());
/* this->compute_2nd (mapping, cell, dsd.cell(), */
template <>
inline
void
-FE_Poly<TensorProductPolynomials<2>,2,3>::fill_fe_values (const Mapping<2,3> &mapping,
- const Triangulation<2,3>::cell_iterator &/* cell */,
- const Quadrature<2> &quadrature,
- Mapping<2,3>::InternalDataBase & mapping_data,
- Mapping<2,3>::InternalDataBase &fedata,
- FEValuesData<2,3> &data) const
+FE_Poly<TensorProductPolynomials<2>,2,3>::fill_fe_values
+ (const Mapping<2,3> &mapping,
+ const Triangulation<2,3>::cell_iterator &/* cell */,
+ const Quadrature<2> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ Mapping<2,3>::InternalDataBase & mapping_data,
+ Mapping<2,3>::InternalDataBase &fedata,
+ FEValuesData<2,3> &data) const
{
// assert that the following dynamics
for (unsigned int i=0; i<quadrature.size(); ++i)
data.shape_values(k,i) = fe_data.shape_values[k][i];
- if (flags & update_gradients)
+ if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
mapping_data, mapping_covariant);
}
/* const QProjector<2>::DataSetDescriptor dsd; */
- if (flags & update_hessians)
+ if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
{
AssertThrow(false, ExcNotImplemented());
/* this->compute_2nd (mapping, cell, dsd.cell(), */
const Mapping<1,2>& mapping,
const Triangulation<1,2>::cell_iterator& /*cell*/,
const Quadrature<1>& quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
Mapping<1,2>::InternalDataBase& mapping_data,
Mapping<1,2>::InternalDataBase& fedata,
FEValuesData<1,2>& data) const
// TODO: I would think this should work. Guido
- if (flags & update_gradients)
+ if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
mapping_data, mapping_covariant);
}
// const typename QProjector<1>::DataSetDescriptor dsd;
- if (flags & update_hessians)
+ if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
{
AssertThrow(false, ExcNotImplemented());
/* this->compute_2nd (mapping, cell, dsd.cell(), */
template <>
inline
void
-FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_values (const Mapping<2,3>& mapping,
- const Triangulation<2,3>::cell_iterator &/* cell */,
- const Quadrature<2>& quadrature,
- Mapping<2,3>::InternalDataBase& mapping_data,
- Mapping<2,3>::InternalDataBase &fedata,
- FEValuesData<2,3> &data) const
+FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_values
+ (const Mapping<2,3>& mapping,
+ const Triangulation<2,3>::cell_iterator &/* cell */,
+ const Quadrature<2>& quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ Mapping<2,3>::InternalDataBase& mapping_data,
+ Mapping<2,3>::InternalDataBase &fedata,
+ FEValuesData<2,3> &data) const
{
Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
InternalData &fe_data = static_cast<InternalData &> (fedata);
data.shape_values(k,i) = fe_data.shape_values[k][i];
- if (flags & update_gradients)
+ if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
mapping_data, mapping_covariant);
}
/* const QProjector<2>::DataSetDescriptor dsd; */
- if (flags & update_hessians)
+ if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
{
AssertThrow(false, ExcNotImplemented());
/* this->compute_2nd (mapping, cell, dsd.cell(), */
template <class POLY, int dim, int spacedim>
void
-FE_Poly<POLY,dim,spacedim>::fill_fe_values (const Mapping<dim,spacedim> &mapping,
- const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
- typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
- typename Mapping<dim,spacedim>::InternalDataBase &fedata,
- FEValuesData<dim,spacedim> &data) const
+FE_Poly<POLY,dim,spacedim>::fill_fe_values
+ (const Mapping<dim,spacedim> &mapping,
+ const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ typename Mapping<dim,spacedim>::InternalDataBase &fedata,
+ FEValuesData<dim,spacedim> &data) const
{
// convert data object to internal
// data for this class. fails with
// an exception if that is not
- // possible
-
+ // possible
Assert (dynamic_cast<InternalData *> (&fedata) != 0, ExcInternalError());
InternalData &fe_data = static_cast<InternalData &> (fedata);
for (unsigned int i=0; i<quadrature.size(); ++i)
data.shape_values(k,i) = fe_data.shape_values[k][i];
- if (flags & update_gradients)
+ if (flags & update_gradients && cell_similarity != CellSimilarity::translation)
mapping.transform(fe_data.shape_gradients[k], data.shape_gradients[k],
mapping_data, mapping_covariant);
}
const typename QProjector<dim>::DataSetDescriptor dsd;
- if (flags & update_hessians)
+ if (flags & update_hessians && cell_similarity != CellSimilarity::translation)
this->compute_2nd (mapping, cell, dsd.cell(),
mapping_data, fe_data, data);
}
template <>
inline
void
-FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_face_values (const Mapping<1,2> &,
- const Triangulation<1,2>::cell_iterator &,
- const unsigned int,
- const Quadrature<0> &,
- Mapping<1,2>::InternalDataBase &,
- Mapping<1,2>::InternalDataBase &,
- FEValuesData<1,2> &) const
+FE_Poly<TensorProductPolynomials<1>,1,2>::fill_fe_face_values
+ (const Mapping<1,2> &,
+ const Triangulation<1,2>::cell_iterator &,
+ const unsigned int,
+ const Quadrature<0> &,
+ Mapping<1,2>::InternalDataBase &,
+ Mapping<1,2>::InternalDataBase &,
+ FEValuesData<1,2> &) const
{
AssertThrow(false, ExcNotImplemented());
}
inline
void
FE_Poly<PolynomialSpace<1>,1,2>::fill_fe_face_values (const Mapping<1,2> &,
- const Triangulation<1,2>::cell_iterator &,
- const unsigned int,
- const Quadrature<0> &,
- Mapping<1,2>::InternalDataBase &,
- Mapping<1,2>::InternalDataBase &,
- FEValuesData<1,2> &) const
+ const Triangulation<1,2>::cell_iterator &,
+ const unsigned int,
+ const Quadrature<0> &,
+ Mapping<1,2>::InternalDataBase &,
+ Mapping<1,2>::InternalDataBase &,
+ FEValuesData<1,2> &) const
{
AssertThrow(false, ExcNotImplemented());
}
inline
void
FE_Poly<PolynomialSpace<2>,2,3>::fill_fe_face_values (const Mapping<2,3> &,
- const Triangulation<2,3>::cell_iterator &,
- const unsigned int,
- const Quadrature<1> &,
- Mapping<2,3>::InternalDataBase &,
- Mapping<2,3>::InternalDataBase &,
- FEValuesData<2,3> &) const
+ const Triangulation<2,3>::cell_iterator &,
+ const unsigned int,
+ const Quadrature<1> &,
+ Mapping<2,3>::InternalDataBase &,
+ Mapping<2,3>::InternalDataBase &,
+ FEValuesData<2,3> &) const
{
AssertThrow(false, ExcNotImplemented());
}
template <class POLY, int dim, int spacedim>
void
-FE_Poly<POLY,dim,spacedim>::fill_fe_face_values (const Mapping<dim,spacedim> &mapping,
- const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const unsigned int face,
- const Quadrature<dim-1> &quadrature,
- typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
- typename Mapping<dim,spacedim>::InternalDataBase &fedata,
- FEValuesData<dim,spacedim> &data) const
+FE_Poly<POLY,dim,spacedim>::
+fill_fe_face_values (const Mapping<dim,spacedim> &mapping,
+ const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+ const unsigned int face,
+ const Quadrature<dim-1> &quadrature,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ typename Mapping<dim,spacedim>::InternalDataBase &fedata,
+ FEValuesData<dim,spacedim> &data) const
{
// convert data object to internal
// data for this class. fails with
const Quadrature<dim>& quadrature) const ;
virtual void
- fill_fe_values (const Mapping<dim,spacedim> &mapping,
+ fill_fe_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_internal,
typename Mapping<dim,spacedim>::InternalDataBase &fe_internal,
FEValuesData<dim,spacedim>& data) const;
virtual void
fill_fe_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
FEValuesData<dim,spacedim> &data) const;
* <tt>sub_no!=invalid_face_no</tt>.
*/
template <int dim_1>
- void compute_fill (const Mapping<dim,spacedim> &mapping,
+ void compute_fill (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const unsigned int face_no,
- const unsigned int sub_no,
- const Quadrature<dim_1> &quadrature,
+ const unsigned int face_no,
+ const unsigned int sub_no,
+ const Quadrature<dim_1> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
- FEValuesData<dim,spacedim> &data) const ;
+ FEValuesData<dim,spacedim> &data) const ;
private:
*/
const FiniteElement<dim,spacedim> & get_fe () const;
+ /**
+ * Return the update flags set
+ * for this object.
+ */
+ UpdateFlags get_update_flags () const;
+
/**
* Return a triangulation
* iterator to the current cell.
*/
void initialize (const UpdateFlags update_flags);
+ /**
+ * An enum variable that can store
+ * different states of the current cell
+ * in comparison to the previously
+ * visited cell. If wanted, additional
+ * states can be checked here and used
+ * in one of the methods used during
+ * reinit.
+ */
+ enum CellSimilarity::Similarity cell_similarity;
+
+ /**
+ * A function that checks whether the
+ * new cell is similar to the one
+ * previously used. Then, a significant
+ * amount of the data can be reused,
+ * e.g. the derivatives of the basis
+ * functions in real space, shape_grad.
+ */
+ void
+ check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell);
+
/**
* The reinit() functions do
* only that part of the work
template <int dim, int spacedim>
inline
+UpdateFlags
+FEValuesBase<dim,spacedim>::get_update_flags () const
+{
+ return this->update_flags;
+}
+
+
+
+template <int dim, int spacedim>
const typename Triangulation<dim,spacedim>::cell_iterator
FEValuesBase<dim,spacedim>::get_cell () const
{
template <int dim, int spacedim> class FEFaceValues;
template <int dim, int spacedim> class FESubfaceValues;
+ /**
+ * An enum variable that can store
+ * different states of the current cell
+ * in comparison to the previously
+ * visited cell. If wanted, additional
+ * states can be checked here and used
+ * in one of the methods used during
+ * reinit.
+ */
+namespace CellSimilarity
+{
+ enum Similarity
+ {
+ no_similarity,
+ translation
+ };
+}
+
/**
* The transformation type used
* for the Mapping::transform() functions.
*/
virtual void
fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
- InternalDataBase &internal,
- std::vector<Point<spacedim> > &quadrature_points,
- std::vector<double> &JxW_values,
- std::vector<Tensor<2,spacedim> > &jacobians,
- std::vector<Tensor<3,spacedim> > &jacobian_grads,
- std::vector<Tensor<2,spacedim> > &inverse_jacobians,
- std::vector<Point<spacedim> > &cell_normal_vectors
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ InternalDataBase &internal,
+ std::vector<Point<spacedim> > &quadrature_points,
+ std::vector<double> &JxW_values,
+ std::vector<Tensor<2,spacedim> > &jacobians,
+ std::vector<Tensor<3,spacedim> > &jacobian_grads,
+ std::vector<Tensor<2,spacedim> > &inverse_jacobians,
+ std::vector<Point<spacedim> > &cell_normal_vectors
) const=0;
virtual void
fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim>& quadrature,
- typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
- std::vector<Point<spacedim> > &quadrature_points,
- std::vector<double> &JxW_values,
- std::vector<Tensor<2,spacedim> > &jacobians,
- std::vector<Tensor<3,spacedim> > &jacobian_grads,
- std::vector<Tensor<2,spacedim> > &inverse_jacobians,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim, spacedim>::InternalDataBase &mapping_data,
+ std::vector<Point<spacedim> > &quadrature_points,
+ std::vector<double> &JxW_values,
+ std::vector<Tensor<2,spacedim> > &jacobians,
+ std::vector<Tensor<3,spacedim> > &jacobian_grads,
+ std::vector<Tensor<2,spacedim> > &inverse_jacobians,
std::vector<Point<spacedim> > &) const ;
void compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const unsigned int face_no,
const unsigned int sub_no,
+ const enum CellSimilarity::Similarity cell_similarity,
InternalData& data,
std::vector<Point<dim> > &quadrature_points,
std::vector<Point<dim> >& normal_vectors) const;
virtual void
fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const Quadrature<dim> &quadrature,
+ enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename std::vector<Point<spacedim> > &quadrature_points,
std::vector<double> &JxW_values,
*/
virtual void
fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> & quadrature,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename std::vector<Point<spacedim> > &quadrature_points,
std::vector<double> &JxW_values,
void compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const unsigned int npts,
const DataSetDescriptor data_set,
+ const enum CellSimilarity::Similarity cell_similarity,
InternalData &data,
std::vector<Point<spacedim> > &quadrature_points) const;
// $Id$
// Version: $Name$
//
-// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 2008 by the deal.II authors
+// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006, 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
protected:
+ /**
+ * Implementation of the interface in
+ * MappingQ1. Overrides the function in
+ * the base class, since we cannot use
+ * any cell similarity for this class.
+ */
+ virtual void
+ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+ const Quadrature<dim> &quadrature,
+ enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ typename std::vector<Point<spacedim> > &quadrature_points,
+ std::vector<double> &JxW_values,
+ std::vector<Tensor<2,spacedim> > &jacobians,
+ std::vector<Tensor<3,spacedim> > &jacobian_grads,
+ std::vector<Tensor<2,spacedim> > &inverse_jacobians,
+ std::vector<Point<spacedim> > &cell_normal_vectors) const;
/**
* Reference to the vector of
<< "-- expected size " << arg2);
protected:
+ /**
+ * Implementation of the interface in
+ * MappingQ. Overrides the function in
+ * the base class, since we cannot use
+ * any cell similarity for this class.
+ */
+ virtual void
+ fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+ const Quadrature<dim> &quadrature,
+ enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ typename std::vector<Point<spacedim> > &quadrature_points,
+ std::vector<double> &JxW_values,
+ std::vector<Tensor<2,spacedim> > &jacobians,
+ std::vector<Tensor<3,spacedim> > &jacobian_grads,
+ std::vector<Tensor<2,spacedim> > &inverse_jacobians,
+ std::vector<Point<spacedim> > &cell_normal_vectors) const;
+
/**
* Reference to the vector of
* shifts.
const Mapping<dim,spacedim>&,
const typename Triangulation<dim,spacedim>::cell_iterator&,
const Quadrature<dim>&,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase&,
typename Mapping<dim,spacedim>::InternalDataBase& fedata,
FEValuesData<dim,spacedim>& data) const
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
-
+ Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+ ExcInternalError());
+ InternalData &fe_data = static_cast<InternalData &> (fedata);
+
const UpdateFlags flags(fe_data.current_update_flags());
Assert (flags & update_quadrature_points, ExcInternalError());
-
+
const unsigned int n_q_points = data.quadrature_points.size();
-
+
if (flags & (update_values | update_gradients))
for (unsigned int i=0; i<n_q_points; ++i)
{
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+ Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+ ExcInternalError());
+ InternalData &fe_data = static_cast<InternalData &> (fedata);
const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
Assert (flags & update_quadrature_points, ExcInternalError());
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+ Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+ ExcInternalError());
+ InternalData &fe_data = static_cast<InternalData &> (fedata);
const UpdateFlags flags(fe_data.update_once | fe_data.update_each);
Assert (flags & update_quadrature_points, ExcInternalError());
template <int dim, int spacedim>
void
-FE_Nedelec<dim,spacedim>::fill_fe_values (const Mapping<dim,spacedim> &mapping,
- const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
- typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
- typename Mapping<dim,spacedim>::InternalDataBase &fedata,
- FEValuesData<dim,spacedim> &data) const
+FE_Nedelec<dim,spacedim>::fill_fe_values
+ (const Mapping<dim,spacedim> &mapping,
+ const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ typename Mapping<dim,spacedim>::InternalDataBase &fedata,
+ FEValuesData<dim,spacedim> &data) const
{
// convert data object to internal
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+ Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+ ExcInternalError());
+ InternalData &fe_data = static_cast<InternalData &> (fedata);
// get the flags indicating the
// fields that have to be filled
template <class POLY, int dim, int spacedim>
void
FE_PolyTensor<POLY,dim,spacedim>::fill_fe_values (
- const Mapping<dim,spacedim> &mapping,
+ const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename Mapping<dim,spacedim>::InternalDataBase &fedata,
- FEValuesData<dim,spacedim> &data) const
+ FEValuesData<dim,spacedim> &data) const
{
// convert data object to internal
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+ Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+ ExcInternalError());
+ InternalData &fe_data = static_cast<InternalData &> (fedata);
const unsigned int n_q_points = quadrature.size();
const UpdateFlags flags(fe_data.current_update_flags());
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+ Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+ ExcInternalError());
+ InternalData &fe_data = static_cast<InternalData &> (fedata);
const unsigned int n_q_points = quadrature.size();
// offset determines which data set
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &fe_data = dynamic_cast<InternalData &> (fedata);
+ Assert (dynamic_cast<InternalData *> (&fedata) != 0,
+ ExcInternalError());
+ InternalData &fe_data = static_cast<InternalData &> (fedata);
const unsigned int n_q_points = quadrature.size();
// offset determines which data set
template <int dim, int spacedim>
FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe,
- const unsigned int n_elements) :
+ const unsigned int n_elements) :
FiniteElement<dim,spacedim> (multiply_dof_numbers(fe, n_elements),
compute_restriction_is_additive_flags (fe, n_elements),
compute_nonzero_components(fe, n_elements)),
template <int dim, int spacedim>
FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe1,
- const unsigned int n1,
- const FiniteElement<dim,spacedim> &fe2,
- const unsigned int n2) :
+ const unsigned int n1,
+ const FiniteElement<dim,spacedim> &fe2,
+ const unsigned int n2) :
FiniteElement<dim,spacedim> (multiply_dof_numbers(fe1, n1, fe2, n2),
compute_restriction_is_additive_flags (fe1, n1,
fe2, n2),
template <int dim, int spacedim>
FESystem<dim,spacedim>::FESystem (const FiniteElement<dim,spacedim> &fe1,
- const unsigned int n1,
- const FiniteElement<dim,spacedim> &fe2,
- const unsigned int n2,
- const FiniteElement<dim,spacedim> &fe3,
- const unsigned int n3) :
- FiniteElement<dim,spacedim> (multiply_dof_numbers(fe1, n1,
+ const unsigned int n1,
+ const FiniteElement<dim,spacedim> &fe2,
+ const unsigned int n2,
+ const FiniteElement<dim,spacedim> &fe3,
+ const unsigned int n3) :
+ FiniteElement<dim,spacedim> (multiply_dof_numbers(fe1, n1,
fe2, n2,
fe3, n3),
compute_restriction_is_additive_flags (fe1, n1,
template <int dim, int spacedim>
double
FESystem<dim,spacedim>::shape_value (const unsigned int i,
- const Point<dim> &p) const
+ const Point<dim> &p) const
{
Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
Assert (this->is_primitive(i),
template <int dim, int spacedim>
double
FESystem<dim,spacedim>::shape_value_component (const unsigned int i,
- const Point<dim> &p,
- const unsigned int component) const
+ const Point<dim> &p,
+ const unsigned int component) const
{
Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
Assert (component < this->n_components(),
template <int dim, int spacedim>
Tensor<1,dim>
FESystem<dim,spacedim>::shape_grad (const unsigned int i,
- const Point<dim> &p) const
+ const Point<dim> &p) const
{
Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
Assert (this->is_primitive(i),
template <int dim, int spacedim>
Tensor<2,dim>
FESystem<dim,spacedim>::shape_grad_grad (const unsigned int i,
- const Point<dim> &p) const
+ const Point<dim> &p) const
{
Assert (i<this->dofs_per_cell, ExcIndexRange(i, 0, this->dofs_per_cell));
Assert (this->is_primitive(i),
template <int dim, int spacedim>
typename Mapping<dim,spacedim>::InternalDataBase *
FESystem<dim,spacedim>::get_data (const UpdateFlags flags_,
- const Mapping<dim,spacedim> &mapping,
- const Quadrature<dim> &quadrature) const
+ const Mapping<dim,spacedim> &mapping,
+ const Quadrature<dim> &quadrature) const
{
UpdateFlags flags = flags_;
template <int dim, int spacedim>
void
FESystem<dim,spacedim>::
-fill_fe_values (const Mapping<dim,spacedim> &mapping,
+fill_fe_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &quadrature,
+ const Quadrature<dim> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
- FEValuesData<dim,spacedim> &data) const
+ FEValuesData<dim,spacedim> &data) const
{
compute_fill(mapping, cell, invalid_face_number, invalid_face_number,
- quadrature, mapping_data, fe_data, data);
+ quadrature, cell_similarity, mapping_data, fe_data, data);
}
typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
FEValuesData<dim,spacedim> &data) const
{
- compute_fill (mapping, cell, face_no, invalid_face_number,
- quadrature, mapping_data, fe_data, data);
+ compute_fill (mapping, cell, face_no, invalid_face_number, quadrature,
+ CellSimilarity::no_similarity, mapping_data, fe_data, data);
}
template <int dim, int spacedim>
void
FESystem<dim,spacedim>::
-fill_fe_subface_values (const Mapping<dim,spacedim> &mapping,
+fill_fe_subface_values (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const unsigned int face_no,
- const unsigned int sub_no,
- const Quadrature<dim-1> &quadrature,
+ const unsigned int face_no,
+ const unsigned int sub_no,
+ const Quadrature<dim-1> &quadrature,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename Mapping<dim,spacedim>::InternalDataBase &fe_data,
- FEValuesData<dim,spacedim> &data) const
+ FEValuesData<dim,spacedim> &data) const
{
- compute_fill (mapping, cell, face_no, sub_no,
- quadrature, mapping_data, fe_data, data);
+ compute_fill (mapping, cell, face_no, sub_no, quadrature,
+ CellSimilarity::no_similarity, mapping_data, fe_data, data);
}
template <int dim_1>
void
FESystem<dim,spacedim>::
-compute_fill (const Mapping<dim,spacedim> &mapping,
+compute_fill (const Mapping<dim,spacedim> &mapping,
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const unsigned int face_no,
- const unsigned int sub_no,
- const Quadrature<dim_1> &quadrature,
+ const unsigned int face_no,
+ const unsigned int sub_no,
+ const Quadrature<dim_1> &quadrature,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
typename Mapping<dim,spacedim>::InternalDataBase &fedata,
- FEValuesData<dim,spacedim> &data) const
+ FEValuesData<dim,spacedim> &data) const
{
const unsigned int n_q_points = quadrature.size();
// data for this class. fails with
// an exception if that is not
// possible
- InternalData & fe_data = dynamic_cast<InternalData&> (fedata);
+ Assert (dynamic_cast<InternalData*> (&fedata) != 0, ExcInternalError());
+ InternalData & fe_data = static_cast<InternalData&> (fedata);
// Either dim_1==dim
// (fill_fe_values) or dim_1==dim-1
if (face_no==invalid_face_number)
{
Assert(dim_1==dim, ExcDimensionMismatch(dim_1,dim));
+ Assert (dynamic_cast<const Quadrature<dim> *>(quadrature_base_pointer) != 0,
+ ExcInternalError());
cell_quadrature
- = dynamic_cast<const Quadrature<dim> *>(quadrature_base_pointer);
- Assert (cell_quadrature != 0, ExcInternalError());
+ = static_cast<const Quadrature<dim> *>(quadrature_base_pointer);
}
else
{
Assert(dim_1==dim-1, ExcDimensionMismatch(dim_1,dim-1));
+ Assert (dynamic_cast<const Quadrature<dim-1> *>(quadrature_base_pointer) != 0,
+ ExcInternalError());
face_quadrature
- = dynamic_cast<const Quadrature<dim-1> *>(quadrature_base_pointer);
- Assert (face_quadrature != 0, ExcInternalError());
+ = static_cast<const Quadrature<dim-1> *>(quadrature_base_pointer);
}
// let base elements update the
// therefore use
// base_fe_data.update_flags.
if (face_no==invalid_face_number)
- base_fe.fill_fe_values(mapping, cell,
- *cell_quadrature, mapping_data, base_fe_data, base_data);
+ base_fe.fill_fe_values(mapping, cell, *cell_quadrature, cell_similarity,
+ mapping_data, base_fe_data, base_data);
else if (sub_no==invalid_face_number)
base_fe.fill_fe_face_values(mapping, cell, face_no,
*face_quadrature, mapping_data, base_fe_data, base_data);
// base element
const UpdateFlags base_flags(dim_1==dim ?
base_fe_data.current_update_flags() :
- base_fe_data.update_flags);
- for (unsigned int system_index=0; system_index<this->dofs_per_cell;
- ++system_index)
- if (this->system_to_base_table[system_index].first.first == base_no)
- {
- const unsigned int
- base_index = this->system_to_base_table[system_index].second;
- Assert (base_index<base_fe.dofs_per_cell, ExcInternalError());
+ base_fe_data.update_flags);
+
+ // if the current cell is just a
+ // translation of the previous
+ // one, the underlying data has
+ // not changed, and we don't even
+ // need to enter this section
+ if (cell_similarity != CellSimilarity::translation)
+ for (unsigned int system_index=0; system_index<this->dofs_per_cell;
+ ++system_index)
+ if (this->system_to_base_table[system_index].first.first == base_no)
+ {
+ const unsigned int
+ base_index = this->system_to_base_table[system_index].second;
+ Assert (base_index<base_fe.dofs_per_cell, ExcInternalError());
// now copy. if the
// shape function is
// this one value, and
// to which index to
// put
- unsigned int out_index = 0;
- for (unsigned int i=0; i<system_index; ++i)
- out_index += this->n_nonzero_components(i);
- unsigned int in_index = 0;
- for (unsigned int i=0; i<base_index; ++i)
- in_index += base_fe.n_nonzero_components(i);
+ unsigned int out_index = 0;
+ for (unsigned int i=0; i<system_index; ++i)
+ out_index += this->n_nonzero_components(i);
+ unsigned int in_index = 0;
+ for (unsigned int i=0; i<base_index; ++i)
+ in_index += base_fe.n_nonzero_components(i);
// then loop over the
// number of components
// to be copied
- Assert (this->n_nonzero_components(system_index) ==
- base_fe.n_nonzero_components(base_index),
- ExcInternalError());
- for (unsigned int s=0; s<this->n_nonzero_components(system_index); ++s)
- {
- if (base_flags & update_values)
- for (unsigned int q=0; q<n_q_points; ++q)
- data.shape_values[out_index+s][q] =
- base_data.shape_values(in_index+s,q);
+ Assert (this->n_nonzero_components(system_index) ==
+ base_fe.n_nonzero_components(base_index),
+ ExcInternalError());
+ for (unsigned int s=0; s<this->n_nonzero_components(system_index); ++s)
+ {
+ if (base_flags & update_values)
+ for (unsigned int q=0; q<n_q_points; ++q)
+ data.shape_values[out_index+s][q] =
+ base_data.shape_values(in_index+s,q);
- if (base_flags & update_gradients)
- for (unsigned int q=0; q<n_q_points; ++q)
- data.shape_gradients[out_index+s][q]=
- base_data.shape_gradients[in_index+s][q];
+ if (base_flags & update_gradients)
+ for (unsigned int q=0; q<n_q_points; ++q)
+ data.shape_gradients[out_index+s][q]=
+ base_data.shape_gradients[in_index+s][q];
// _we_ handle
// computation of
// should not
// have computed
// them!
- Assert (!(base_flags & update_hessians),
- ExcInternalError());
- };
- };
+ Assert (!(base_flags & update_hessians),
+ ExcInternalError());
+ };
+ };
};
if (fe_data.is_first_cell())
}
}
- if (fe_data.compute_hessians)
+ if (fe_data.compute_hessians && cell_similarity != CellSimilarity::translation)
{
unsigned int offset = 0;
if (face_no != invalid_face_number)
}
+
/*------------------------------- FEValuesBase ---------------------------*/
}
-
/*------------------------------- FEValues -------------------------------*/
template <int dim, int spacedim>
FEValues<dim,spacedim>::FEValues (const Mapping<dim,spacedim> &mapping,
- const FiniteElement<dim,spacedim> &fe,
- const Quadrature<dim> &q,
- const UpdateFlags update_flags)
+ const FiniteElement<dim,spacedim> &fe,
+ const Quadrature<dim> &q,
+ const UpdateFlags update_flags)
:
FEValuesBase<dim,spacedim> (q.size(),
fe.dofs_per_cell,
template <int dim, int spacedim>
FEValues<dim,spacedim>::FEValues (const FiniteElement<dim,spacedim> &fe,
- const Quadrature<dim> &q,
- const UpdateFlags update_flags)
+ const Quadrature<dim> &q,
+ const UpdateFlags update_flags)
:
FEValuesBase<dim,spacedim> (q.size(),
- fe.dofs_per_cell,
- update_default,
+ fe.dofs_per_cell,
+ update_default,
StaticMappingQ1<dim,spacedim>::mapping,
- fe),
- quadrature (q)
+ fe),
+ quadrature (q)
{
Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping"));
initialize (update_flags);
typedef FEValuesBase<dim,spacedim> FEVB;
Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
- static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
+ static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
typename FEVB::ExcFEDontMatch());
+ check_cell_similarity(cell);
+
// set new cell. auto_ptr will take
// care that old object gets
// destroyed and also that this
// passed to the constructor and
// used by the DoFHandler used by
// this cell, are the same
- typedef FEValuesBase<dim,spacedim> FEVB;
+ typedef FEValuesBase<dim,spacedim> FEVB;
Assert (static_cast<const FiniteElementData<dim>&>(*this->fe) ==
static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
typename FEVB::ExcFEDontMatch());
+ check_cell_similarity(cell);
+
// set new cell. auto_ptr will take
// care that old object gets
// destroyed and also that this
// static_cast<const FiniteElementData<dim>&>(cell->get_fe()),
// typename FEValuesBase<dim,spacedim>::ExcFEDontMatch());
+ check_cell_similarity(cell);
+
// set new cell. auto_ptr will take
// care that old object gets
// destroyed and also that this
template <int dim, int spacedim>
void FEValues<dim,spacedim>::reinit (const typename Triangulation<dim,spacedim>::cell_iterator &cell)
{
- // no FE in this cell, so no check
+ // no FE in this cell, so no assertion
// necessary here
+ check_cell_similarity(cell);
// set new cell. auto_ptr will take
// care that old object gets
+template <int dim, int spacedim>
+inline
+void
+FEValues<dim,spacedim>::check_cell_similarity (const typename Triangulation<dim,spacedim>::cell_iterator &cell)
+{
+ // case that there has not been any cell
+ // before
+ if (&*this->present_cell == 0)
+ {
+ cell_similarity = CellSimilarity::no_similarity;
+ return;
+ }
+
+ const typename Triangulation<dim,spacedim>::cell_iterator & present_cell =
+ *this->present_cell;
+
+ // test for translation
+ {
+ // otherwise, go through the vertices and
+ // check...
+ bool is_translation = true;
+ const Point<spacedim> dist = cell->vertex(0) - present_cell->vertex(0);
+ 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() > 1e-28)
+ {
+ is_translation = false;
+ break;
+ }
+ }
+
+ cell_similarity = (is_translation == true
+ ?
+ CellSimilarity::translation
+ :
+ CellSimilarity::no_similarity);
+ return;
+ }
+
+ // TODO: here, one could implement other
+ // checks for similarity, e.g. for
+ // children of a parallelogram.
+}
+
+
+
template <int dim, int spacedim>
void FEValues<dim,spacedim>::do_reinit ()
{
this->get_mapping().fill_fe_values(*this->present_cell,
quadrature,
+ this->cell_similarity,
*this->mapping_data,
this->quadrature_points,
this->JxW_values,
this->get_fe().fill_fe_values(this->get_mapping(),
*this->present_cell,
quadrature,
+ this->cell_similarity,
*this->mapping_data,
*this->fe_data,
*this);
template <int dim, int spacedim>
FEFaceValuesBase<dim,spacedim>::FEFaceValuesBase (const unsigned int n_q_points,
- const unsigned int dofs_per_cell,
- const UpdateFlags,
- const Mapping<dim,spacedim> &mapping,
- const FiniteElement<dim,spacedim> &fe,
- const Quadrature<dim-1>& quadrature)
+ const unsigned int dofs_per_cell,
+ const UpdateFlags,
+ const Mapping<dim,spacedim> &mapping,
+ const FiniteElement<dim,spacedim> &fe,
+ const Quadrature<dim-1>& quadrature)
:
FEValuesBase<dim,spacedim> (n_q_points,
dofs_per_cell,
template<int dim, int spacedim>
typename Mapping<dim, spacedim>::InternalDataBase *
MappingCartesian<dim, spacedim>::get_data (const UpdateFlags update_flags,
- const Quadrature<dim> &q) const
+ const Quadrature<dim> &q) const
{
InternalData* data = new InternalData (q);
template<int dim, int spacedim>
typename Mapping<dim, spacedim>::InternalDataBase *
MappingCartesian<dim, spacedim>::get_face_data (const UpdateFlags update_flags,
- const Quadrature<dim-1>& quadrature) const
+ const Quadrature<dim-1>& quadrature) const
{
InternalData* data
= new InternalData (QProjector<dim>::project_to_all_faces(quadrature));
template<int dim, int spacedim>
typename Mapping<dim, spacedim>::InternalDataBase *
MappingCartesian<dim, spacedim>::get_subface_data (const UpdateFlags update_flags,
- const Quadrature<dim-1> &quadrature) const
+ const Quadrature<dim-1> &quadrature) const
{
InternalData* data
= new InternalData (QProjector<dim>::project_to_all_subfaces(quadrature));
template<int dim, int spacedim>
void
MappingCartesian<dim, spacedim>::compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const unsigned int face_no,
- const unsigned int sub_no,
- InternalData &data,
- std::vector<Point<dim> > &quadrature_points,
- std::vector<Point<dim> > &normal_vectors) const
+ const unsigned int face_no,
+ const unsigned int sub_no,
+ const enum CellSimilarity::Similarity cell_similarity,
+ InternalData &data,
+ std::vector<Point<dim> > &quadrature_points,
+ std::vector<Point<dim> > &normal_vectors) const
{
const UpdateFlags update_flags(data.current_update_flags());
MappingCartesian<dim, spacedim>::
fill_fe_values (const typename Triangulation<dim,spacedim>::cell_iterator& cell,
const Quadrature<dim>& q,
+ const enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase& mapping_data,
std::vector<Point<spacedim> >& quadrature_points,
std::vector<double>& JxW_values,
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+ Assert (dynamic_cast<InternalData*> (&mapping_data) != 0, ExcInternalError());
+ InternalData &data = static_cast<InternalData&> (mapping_data);
std::vector<Point<dim> > dummy;
- compute_fill (cell, invalid_face_number, invalid_face_number,
+ compute_fill (cell, invalid_face_number, invalid_face_number, cell_similarity,
data,
quadrature_points,
dummy);
// local lengths in each coordinate
// direction
if (data.current_update_flags() & (update_JxW_values | update_volume_elements))
- {
- double J = data.length[0];
- for (unsigned int d=1;d<dim;++d)
- J *= data.length[d];
- data.volume_element = J;
- if (data.current_update_flags() & update_JxW_values)
- for (unsigned int i=0; i<JxW_values.size();++i)
- JxW_values[i] = J * q.weight(i);
- }
+ if (cell_similarity != CellSimilarity::translation)
+ {
+ double J = data.length[0];
+ for (unsigned int d=1;d<dim;++d)
+ J *= data.length[d];
+ data.volume_element = J;
+ if (data.current_update_flags() & update_JxW_values)
+ for (unsigned int i=0; i<JxW_values.size();++i)
+ JxW_values[i] = J * q.weight(i);
+ }
// "compute" Jacobian at the quadrature
// points, which are all the same
if (data.current_update_flags() & update_jacobians)
- for (unsigned int i=0; i<jacobians.size();++i)
- {
- jacobians[i]=Tensor<2,dim>();
- for (unsigned int j=0; j<dim; ++j)
- jacobians[j][j]=data.length[j];
- }
+ if (cell_similarity != CellSimilarity::translation)
+ for (unsigned int i=0; i<jacobians.size();++i)
+ {
+ jacobians[i]=Tensor<2,dim>();
+ for (unsigned int j=0; j<dim; ++j)
+ jacobians[j][j]=data.length[j];
+ }
// "compute" the derivative of the Jacobian
// at the quadrature points, which are all
// zero of course
if (data.current_update_flags() & update_jacobian_grads)
- for (unsigned int i=0; i<jacobian_grads.size();++i)
- jacobian_grads[i]=Tensor<3,dim>();
+ if (cell_similarity != CellSimilarity::translation)
+ for (unsigned int i=0; i<jacobian_grads.size();++i)
+ jacobian_grads[i]=Tensor<3,dim>();
// "compute" inverse Jacobian at the
// quadrature points, which are all
// the same
if (data.current_update_flags() & update_inverse_jacobians)
- for (unsigned int i=0; i<inverse_jacobians.size();++i)
- {
- inverse_jacobians[i]=Tensor<2,dim>();
- for (unsigned int j=0; j<dim; ++j)
- inverse_jacobians[j][j]=1./data.length[j];
- }
+ if (cell_similarity != CellSimilarity::translation)
+ for (unsigned int i=0; i<inverse_jacobians.size();++i)
+ {
+ inverse_jacobians[i]=Tensor<2,dim>();
+ for (unsigned int j=0; j<dim; ++j)
+ inverse_jacobians[j][j]=1./data.length[j];
+ }
}
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+ Assert (dynamic_cast<InternalData*> (&mapping_data) != 0,
+ ExcInternalError());
+ InternalData &data = static_cast<InternalData&> (mapping_data);
- compute_fill (cell, face_no, invalid_face_number,
+ compute_fill (cell, face_no, invalid_face_number,
+ CellSimilarity::no_similarity,
data,
quadrature_points,
normal_vectors);
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+ Assert (dynamic_cast<InternalData*> (&mapping_data) != 0, ExcInternalError());
+ InternalData &data = static_cast<InternalData&> (mapping_data);
- compute_fill (cell, face_no, sub_no,
+ compute_fill (cell, face_no, sub_no, CellSimilarity::no_similarity,
data,
quadrature_points,
normal_vectors);
AssertDimension (input.size(), output.size());
Assert (dynamic_cast<const InternalData *>(&mapping_data) != 0,
ExcInternalError());
- const InternalData &data = dynamic_cast<const InternalData&> (mapping_data);
+ const InternalData &data = static_cast<const InternalData&> (mapping_data);
switch (mapping_type)
{
template<int dim, int spacedim>
void
MappingQ<dim,spacedim>::compute_shapes_virtual (const std::vector<Point<dim> > &unit_points,
- typename MappingQ1<dim,spacedim>::InternalData &data) const
+ typename MappingQ1<dim,spacedim>::InternalData &data) const
{
const unsigned int n_points=unit_points.size();
std::vector<double> values;
template<int dim, int spacedim>
typename Mapping<dim,spacedim>::InternalDataBase *
MappingQ<dim,spacedim>::get_data (const UpdateFlags update_flags,
- const Quadrature<dim> &quadrature) const
+ const Quadrature<dim> &quadrature) const
{
InternalData *data = new InternalData(n_shape_functions);
this->compute_data (update_flags, quadrature,
template<int dim, int spacedim>
typename Mapping<dim,spacedim>::InternalDataBase *
MappingQ<dim,spacedim>::get_face_data (const UpdateFlags update_flags,
- const Quadrature<dim-1>& quadrature) const
+ const Quadrature<dim-1>& quadrature) const
{
InternalData *data = new InternalData(n_shape_functions);
const Quadrature<dim> q (QProjector<dim>::project_to_all_faces(quadrature));
template<int dim, int spacedim>
typename Mapping<dim,spacedim>::InternalDataBase *
MappingQ<dim,spacedim>::get_subface_data (const UpdateFlags update_flags,
- const Quadrature<dim-1>& quadrature) const
+ const Quadrature<dim-1>& quadrature) const
{
InternalData *data = new InternalData(n_shape_functions);
const Quadrature<dim> q (QProjector<dim>::project_to_all_subfaces(quadrature));
}
-
+ // Note that the CellSimilarity flag is
+ // modifyable, since MappingQ can need to
+ // recalculate data even when cells are
+ // similar.
template<int dim, int spacedim>
void
MappingQ<dim,spacedim>::fill_fe_values (
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const Quadrature<dim> &q,
+ enum CellSimilarity::Similarity cell_similarity,
typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
std::vector<Point<spacedim> > &quadrature_points,
std::vector<double> &JxW_values,
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+ Assert (dynamic_cast<InternalData*> (&mapping_data) != 0, ExcInternalError());
+ InternalData &data = static_cast<InternalData&> (mapping_data);
// check whether this cell needs
// the full mapping or can be
data.use_mapping_q1_on_current_cell = !(use_mapping_q_on_all_cells
|| cell->has_boundary_lines());
- // depending on this result, use
- // this or the other data object
- // for the mapping
+ // depending on this result, use this or
+ // the other data object for the
+ // mapping. furthermore, we need to
+ // ensure that the flag indicating
+ // whether we can use some similarity has
+ // to be modified - for a general
+ // MappingQ, the data needs to be
+ // recomputed anyway since then the
+ // mapping changes the data. this needs
+ // to be known also for later operations,
+ // so modify the variable here. TODO:
+ // still need to check how this will work
+ // when the previous cell disabled this
+ // flag.
typename MappingQ1<dim,spacedim>::InternalData *p_data=0;
if (data.use_mapping_q1_on_current_cell)
p_data=&data.mapping_q1_data;
else
- p_data=&data;
+ {
+ p_data=&data;
+ cell_similarity = CellSimilarity::no_similarity;
+ }
- MappingQ1<dim,spacedim>::fill_fe_values(cell, q, *p_data,
+ MappingQ1<dim,spacedim>::fill_fe_values(cell, q, cell_similarity, *p_data,
quadrature_points, JxW_values,
jacobians, jacobian_grads, inverse_jacobians,
cell_normal_vectors);
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+ Assert (dynamic_cast<InternalData*> (&mapping_data) != 0,
+ ExcInternalError());
+ InternalData &data = static_cast<InternalData&> (mapping_data);
// check whether this cell needs
// the full mapping or can be
// data for this class. fails with
// an exception if that is not
// possible
- InternalData &data = dynamic_cast<InternalData&> (mapping_data);
+ Assert (dynamic_cast<InternalData*> (&mapping_data) != 0,
+ ExcInternalError());
+ InternalData &data = static_cast<InternalData&> (mapping_data);
// check whether this cell needs
// the full mapping or can be
template<int dim, int spacedim>
void
MappingQ1<dim,spacedim>::compute_shapes (const std::vector<Point<dim> > &unit_points,
- InternalData &data) const
+ InternalData &data) const
{
// choose either the function implemented
// in this class, or whatever a virtual
template<int dim, int spacedim>
void
MappingQ1<dim,spacedim>::compute_data (const UpdateFlags update_flags,
- const Quadrature<dim> &q,
- const unsigned int n_original_q_points,
- InternalData &data) const
+ const Quadrature<dim> &q,
+ const unsigned int n_original_q_points,
+ InternalData &data) const
{
const unsigned int n_q_points = q.size();
template<int dim, int spacedim>
void
MappingQ1<dim,spacedim>::compute_face_data (const UpdateFlags update_flags,
- const Quadrature<dim>& q,
- const unsigned int n_original_q_points,
- InternalData& data) const
+ const Quadrature<dim>& q,
+ const unsigned int n_original_q_points,
+ InternalData& data) const
{
compute_data (update_flags, q, n_original_q_points, data);
template<int dim, int spacedim>
typename Mapping<dim,spacedim>::InternalDataBase *
MappingQ1<dim,spacedim>::get_face_data (const UpdateFlags update_flags,
- const Quadrature<dim-1> &quadrature) const
+ const Quadrature<dim-1> &quadrature) const
{
InternalData* data = new InternalData(n_shape_functions);
compute_face_data (update_flags,
MappingQ1<dim,spacedim>::compute_fill (const typename Triangulation<dim,spacedim>::cell_iterator &cell,
const unsigned int n_q_points,
const DataSetDescriptor data_set,
+ const enum CellSimilarity::Similarity cell_similarity,
InternalData &data,
std::vector<Point<spacedim> > &quadrature_points) const
{
const UpdateFlags update_flags(data.current_update_flags());
- if (update_flags & update_quadrature_points)
- {
- Assert (quadrature_points.size() == n_q_points,
- ExcDimensionMismatch(quadrature_points.size(), n_q_points));
- std::fill(quadrature_points.begin(), quadrature_points.end(),
- Point<spacedim>());
- }
-
- if (update_flags & update_contravariant_transformation)
- {
- Assert (data.contravariant.size() == n_q_points,
- ExcDimensionMismatch(data.contravariant.size(), n_q_points));
- std::fill(data.contravariant.begin(), data.contravariant.end(),
- Tensor<2,spacedim>());
- }
-
- if (update_flags & update_covariant_transformation)
- {
- Assert (data.covariant.size() == n_q_points,
- ExcDimensionMismatch(data.covariant.size(), n_q_points));
- }
// if necessary, recompute the
// support points of the
// transformation of this cell
// first compute quadrature points
if (update_flags & update_quadrature_points)
- for (unsigned int point=0; point<n_q_points; ++point)
- for (unsigned int k=0; k<data.n_shape_functions; ++k)
- quadrature_points[point]
- += data.shape(point+data_set,k) * data.mapping_support_points[k];
-
+ {
+ Assert (quadrature_points.size() == n_q_points,
+ ExcDimensionMismatch(quadrature_points.size(), n_q_points));
+ std::fill(quadrature_points.begin(), quadrature_points.end(),
+ Point<spacedim>());
+
+ for (unsigned int point=0; point<n_q_points; ++point)
+ for (unsigned int k=0; k<data.n_shape_functions; ++k)
+ quadrature_points[point]
+ += data.shape(point+data_set,k) * data.mapping_support_points[k];
+ }
+
// then Jacobians
if (update_flags & update_contravariant_transformation)
- for (unsigned int point=0; point<n_q_points; ++point)
- for (unsigned int k=0; k<data.n_shape_functions; ++k)
+ {
+ Assert (data.contravariant.size() == n_q_points,
+ ExcDimensionMismatch(data.contravariant.size(), n_q_points));
+
+ // if the current cell is just a
+ // translation of the previous one, no
+ // need to recompute jacobians...
+ if (cell_similarity != CellSimilarity::translation)
{
+ std::fill(data.contravariant.begin(), data.contravariant.end(),
+ Tensor<2,spacedim>());
+ for (unsigned int point=0; point<n_q_points; ++point)
+ for (unsigned int k=0; k<data.n_shape_functions; ++k)
+ {
// some compilers seem to have problems
// to detected that the two innermost
// loops just use the data of the same
// tensors, so get a reference to them by
// hand
- const Tensor<1,dim> &data_derv = data.derivative(point+data_set, k);
- const Tensor<1,spacedim> &supp_pts = data.mapping_support_points[k];
+ const Tensor<1,dim> &data_derv = data.derivative(point+data_set, k);
+ const Tensor<1,spacedim> &supp_pts = data.mapping_support_points[k];
- for (unsigned int i=0; i<spacedim; ++i)
- for (unsigned int j=0; j<dim; ++j)
- data.contravariant[point][i][j] += data_derv[j] * supp_pts[i];
+ for (unsigned int i=0; i<spacedim; ++i)
+ for (unsigned int j=0; j<dim; ++j)
+ data.contravariant[point][i][j] += data_derv[j] * supp_pts[i];
+ }
}
+ }
- if (update_flags & update_volume_elements)
- for (unsigned int point=0; point<n_q_points; ++point)
- data.volume_elements[point] = determinant(data.contravariant[point]);
-
+ if (update_flags & update_covariant_transformation)
+ {
+ Assert (data.covariant.size() == n_q_points,
+ ExcDimensionMismatch(data.covariant.size(), n_q_points));
+ if (cell_similarity != CellSimilarity::translation)
+ {
+ if (dim == spacedim)
// invert contravariant for
// covariant transformation
// matrices
- if ((update_flags & update_covariant_transformation)&&(dim==spacedim))
- for (unsigned int point=0; point<n_q_points; ++point)
- data.covariant[point] = invert(data.contravariant[point]);
+ for (unsigned int point=0; point<n_q_points; ++point)
+ data.covariant[point] = invert(data.contravariant[point]);
+ else if (dim == spacedim - 1)
+ {
// CODIMENSION 1
// calculate left-inversion of the
// contravariant matrices to obtain
// covariant ones (auxiliary
// rectangular fullmatrices are used)
- else if ((update_flags & update_covariant_transformation)&&(dim==spacedim-1))
- {
- FullMatrix<double> contravariant_matrix(spacedim,dim);
- FullMatrix<double> covariant_matrix(dim,spacedim);
+ FullMatrix<double> contravariant_matrix(spacedim,dim);
+ FullMatrix<double> covariant_matrix(dim,spacedim);
- for (unsigned int point=0; point<n_q_points; ++point)
- {
- contravariant_matrix.
- copy_from(data.contravariant[point],0,spacedim-1,0,dim-1);
- covariant_matrix.
- left_invert(contravariant_matrix);
- covariant_matrix.
- copy_to(data.covariant[point],0,dim-1,0,spacedim-1);
+ for (unsigned int point=0; point<n_q_points; ++point)
+ {
+ contravariant_matrix.
+ copy_from(data.contravariant[point],0,spacedim-1,0,dim-1);
+ covariant_matrix.
+ left_invert(contravariant_matrix);
+ covariant_matrix.
+ copy_to(data.covariant[point],0,dim-1,0,spacedim-1);
+ }
+ }
}
}
+ if (update_flags & update_volume_elements)
+ if (cell_similarity != CellSimilarity::translation)
+ for (unsigned int point=0; point<n_q_points; ++point)
+ data.volume_elements[point] = determinant(data.contravariant[point]);
}
void
MappingQ1<dim,spacedim>::fill_fe_values (
const typename Triangulation<dim,spacedim>::cell_iterator &cell,
- const Quadrature<dim> &q,
- typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
- std::vector<Point<spacedim> > &quadrature_points,
- std::vector<double> &JxW_values,
- std::vector<Tensor<2,spacedim> > &jacobians,
- std::vector<Tensor<3,spacedim> > &jacobian_grads,
- std::vector<Tensor<2,spacedim> > &inverse_jacobians,
- std::vector<Point<spacedim> > &cell_normal_vectors) const
+ const Quadrature<dim> &q,
+ const enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ std::vector<Point<spacedim> > &quadrature_points,
+ std::vector<double> &JxW_values,
+ std::vector<Tensor<2,spacedim> > &jacobians,
+ std::vector<Tensor<3,spacedim> > &jacobian_grads,
+ std::vector<Tensor<2,spacedim> > &inverse_jacobians,
+ std::vector<Point<spacedim> > &cell_normal_vectors) const
{
// ensure that the following cast is really correct:
Assert (dynamic_cast<InternalData *>(&mapping_data) != 0,
const unsigned int n_q_points=q.size();
- compute_fill (cell, n_q_points, DataSetDescriptor::cell (),
+ compute_fill (cell, n_q_points, DataSetDescriptor::cell (), cell_similarity,
data, quadrature_points);
const UpdateFlags update_flags(data.current_update_flags());
const std::vector<double> &weights=q.get_weights();
-
// Multiply quadrature weights by
// Jacobian determinants or by norm
// of the cross product of the
Assert( !(update_flags & update_cell_normal_vectors ) ||
(cell_normal_vectors.size() == n_q_points),
ExcDimensionMismatch(cell_normal_vectors.size(), n_q_points));
+ if (cell_similarity != CellSimilarity::translation)
for (unsigned int point=0; point<n_q_points; ++point)
{
{
Assert (jacobians.size() == n_q_points,
ExcDimensionMismatch(jacobians.size(), n_q_points));
+ if (cell_similarity != CellSimilarity::translation)
for (unsigned int point=0; point<n_q_points; ++point)
jacobians[point]
= data.contravariant[point];
{
Assert (jacobian_grads.size() == n_q_points,
ExcDimensionMismatch(jacobian_grads.size(), n_q_points));
- std::fill(jacobian_grads.begin(),
- jacobian_grads.end(),
- Tensor<3,spacedim>());
- for (unsigned int point=0; point<n_q_points; ++point)
- for (unsigned int k=0; k<data.n_shape_functions; ++k)
- for (unsigned int i=0; i<dim; ++i)
- for (unsigned int j=0; j<dim; ++j)
- for (unsigned int l=0; l<dim; ++l)
- jacobian_grads[point][i][j][l]
- += (data.second_derivative(point+DataSetDescriptor::cell (), k)[j][l]
- *
- data.mapping_support_points[k][i]);
+ if (cell_similarity != CellSimilarity::translation)
+ {
+ std::fill(jacobian_grads.begin(),
+ jacobian_grads.end(),
+ Tensor<3,spacedim>());
+
+ for (unsigned int point=0; point<n_q_points; ++point)
+ for (unsigned int k=0; k<data.n_shape_functions; ++k)
+ for (unsigned int i=0; i<dim; ++i)
+ for (unsigned int j=0; j<dim; ++j)
+ for (unsigned int l=0; l<dim; ++l)
+ jacobian_grads[point][i][j][l]
+ += (data.second_derivative(point+DataSetDescriptor::cell (), k)[j][l]
+ *
+ data.mapping_support_points[k][i]);
+ }
}
// copy values from InternalData to vector
// given by reference
{
Assert (inverse_jacobians.size() == n_q_points,
ExcDimensionMismatch(inverse_jacobians.size(), n_q_points));
- for (unsigned int point=0; point<n_q_points; ++point)
- inverse_jacobians[point]
- = transpose(data.covariant[point]);
+ if (cell_similarity != CellSimilarity::translation)
+ for (unsigned int point=0; point<n_q_points; ++point)
+ inverse_jacobians[point]
+ = transpose(data.covariant[point]);
}
}
std::vector<Tensor<1,dim> > &boundary_forms,
std::vector<Point<spacedim> > &normal_vectors) const
{
- compute_fill (cell, n_q_points, data_set, data, quadrature_points);
+ compute_fill (cell, n_q_points, data_set, CellSimilarity::no_similarity,
+ data, quadrature_points);
const UpdateFlags update_flags(data.current_update_flags());
MappingQ1Eulerian<dim, EulerVectorType, spacedim>::
MappingQ1Eulerian (const EulerVectorType &euler_transform_vectors,
const DoFHandler<dim> &shiftmap_dof_handler)
- :
- euler_transform_vectors(euler_transform_vectors),
- shiftmap_dof_handler(&shiftmap_dof_handler)
+ :
+ euler_transform_vectors(euler_transform_vectors),
+ shiftmap_dof_handler(&shiftmap_dof_handler)
{}
}
+
+template<int dim, class EulerVectorType, int spacedim>
+void
+MappingQ1Eulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
+ const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+ const Quadrature<dim> &q,
+ enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ std::vector<Point<spacedim> > &quadrature_points,
+ std::vector<double> &JxW_values,
+ std::vector<Tensor<2,spacedim> > &jacobians,
+ std::vector<Tensor<3,spacedim> > &jacobian_grads,
+ std::vector<Tensor<2,spacedim> > &inverse_jacobians,
+ std::vector<Point<spacedim> > &cell_normal_vectors) const
+{
+ // disable any previously detected
+ // similarity and hand on to the
+ // respective function of the base class.
+ cell_similarity = CellSimilarity::no_similarity;
+ MappingQ1<dim,spacedim>::fill_fe_values (cell, q, cell_similarity, mapping_data,
+ quadrature_points, JxW_values, jacobians,
+ jacobian_grads, inverse_jacobians,
+ cell_normal_vectors);
+}
+
+
+
// explicit instantiation
template class MappingQ1Eulerian<deal_II_dimension, Vector<double> >;
#ifdef DEAL_II_USE_PETSC
+template<int dim, class EulerVectorType, int spacedim>
+void
+MappingQEulerian<dim,EulerVectorType,spacedim>::fill_fe_values (
+ const typename Triangulation<dim,spacedim>::cell_iterator &cell,
+ const Quadrature<dim> &q,
+ enum CellSimilarity::Similarity cell_similarity,
+ typename Mapping<dim,spacedim>::InternalDataBase &mapping_data,
+ std::vector<Point<spacedim> > &quadrature_points,
+ std::vector<double> &JxW_values,
+ std::vector<Tensor<2,spacedim> > &jacobians,
+ std::vector<Tensor<3,spacedim> > &jacobian_grads,
+ std::vector<Tensor<2,spacedim> > &inverse_jacobians,
+ std::vector<Point<spacedim> > &cell_normal_vectors) const
+{
+ // disable any previously detected
+ // similarity and hand on to the respective
+ // function of the base class.
+ cell_similarity = CellSimilarity::no_similarity;
+ MappingQ<dim,spacedim>::fill_fe_values (cell, q, cell_similarity, mapping_data,
+ quadrature_points, JxW_values, jacobians,
+ jacobian_grads, inverse_jacobians,
+ cell_normal_vectors);
+}
+
+
+
// explicit instantiation
template class MappingQEulerian<deal_II_dimension, Vector<double> >;
#ifdef DEAL_II_USE_PETSC