* for a description of the orientation of the different faces.
*/
template <int dim>
-class QProjector : public Quadrature<dim>
+class QProjector
{
public:
-//TODO:[GK] remove this constructor again and return to old state with static-only functions. move this function to something in the finite element classes which is the only place where it is used.
- /**
- * Constructor for a quadrature rule on all (sub)faces.
- * The quadrature rule
- * @p{quadrature} is applied to
- * each face or subface, according
- * to the @{sub} flag.
- *
- * The weights of the new rule are
- * replications of the original
- * weights. This is not a proper
- * handling, but it is
- * consistent with later use. If
- * the (sub)face rule is applied to
- * the unity function, the result
- * is the number of (sub)faces.
- */
- QProjector (const Quadrature<dim-1>& quadrature,
- const bool sub);
-
+
/**
- * Compute the quadrature points on the
- * cell if the given quadrature formula
- * is used on face @p{face_no}. For further
- * details, see the general doc for
- * this class.
+ * Compute the quadrature points
+ * on the cell if the given
+ * quadrature formula is used on
+ * face @p{face_no}. For further
+ * details, see the general doc
+ * for this class.
*/
static void project_to_face (const Quadrature<dim-1> &quadrature,
const unsigned int face_no,
typename std::vector<Point<dim> > &q_points);
- /**
- * Projection to all faces.
- * Generate a formula that integrates
- * over all faces at the same time.
- */
- static void project_to_faces (const Quadrature<dim-1> &quadrature,
- typename std::vector<Point<dim> > &q_points);
-
/**
- * Compute the quadrature points on the
- * cell if the given quadrature formula
- * is used on face @p{face_no}, subface
- * number @p{subface_no}. For further
- * details, see the general doc for
- * this class.
+ * Compute the quadrature points
+ * on the cell if the given
+ * quadrature formula is used on
+ * face @p{face_no}, subface
+ * number @p{subface_no}. For
+ * further details, see the
+ * general doc for this class.
*/
static void project_to_subface (const Quadrature<dim-1> &quadrature,
const unsigned int face_no,
typename std::vector<Point<dim> > &q_points);
/**
- * Projection to all child faces.
- * Project to the children of all
- * faces at the same time. The
- * ordering is first by face,
- * then by subface
+ * Take a face quadrature formula
+ * and generate a cell quadrature
+ * formula from it where the
+ * quadrature points of the given
+ * argument are projected on all
+ * faces.
+ *
+ * The weights of the new rule
+ * are replications of the
+ * original weights. This is not
+ * a proper handling, in that the
+ * sum of weights does not equal
+ * one, but it is consistent with
+ * the use of this function,
+ * namely to generate sets of
+ * face quadrature points on a
+ * cell, one set of which will
+ * then be selected at each
+ * time. This is used in the
+ * @ref{FEFaceValues} class,
+ * where we initialize the
+ * values, derivatives, etc on
+ * all faces at once, while
+ * selecting the data of one
+ * particular face only happens
+ * later.
+ */
+ static Quadrature<dim>
+ project_to_all_faces (const Quadrature<dim-1> &quadrature);
+
+ /**
+ * This function is alike the
+ * previous one, but projects the
+ * given face quadrature formula
+ * to the subfaces of a cell,
+ * i.e. to the children of the
+ * faces of the unit cell.
*/
- static void project_to_subfaces (const Quadrature<dim-1> &quadrature,
- typename std::vector<Point<dim> > &q_points);
+ static Quadrature<dim>
+ project_to_all_subfaces (const Quadrature<dim-1> &quadrature);
};
void QProjector<3>::project_to_face (const Quadrature<2> &quadrature,
const unsigned int face_no,
std::vector<Point<3> > &q_points);
+
+template <>
+Quadrature<1> QProjector<1>::project_to_all_faces (const Quadrature<0> &quadrature);
+
+
template <>
void QProjector<1>::project_to_subface (const Quadrature<0> &,
const unsigned int,
const unsigned int face_no,
const unsigned int subface_no,
std::vector<Point<3> > &q_points);
+
+template <>
+Quadrature<1> QProjector<1>::project_to_all_subfaces (const Quadrature<0> &quadrature);
+
+
template <>
bool
QIterated<1>::uses_both_endpoints (const Quadrature<1> &base_quadrature);
//----------------------------------------------------------------------//
-template <int dim>
-QProjector<dim>::QProjector (const Quadrature<dim-1> &quadrature,
- const bool sub)
- :
- Quadrature<dim> (quadrature.n_quadrature_points
- * GeometryInfo<dim>::faces_per_cell
- * (sub ? GeometryInfo<dim>::subfaces_per_face : 1))
-{
- if (sub)
- project_to_subfaces (quadrature, quadrature_points);
- else
- project_to_faces (quadrature, quadrature_points);
-
- const unsigned int n = GeometryInfo<dim>::faces_per_cell
- * (sub ? GeometryInfo<dim>::subfaces_per_face : 1);
- unsigned int k=0;
- for (unsigned int i=0; i<n; ++i)
- for (unsigned int j=0; j<quadrature.n_quadrature_points; ++j)
- weights[k++] = quadrature.weight(j);
-}
-
template <>
void QProjector<1>::project_to_face (const Quadrature<0> &,
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_faces (const Quadrature<0> &)
+{
+ Assert (false, ExcImpossibleInDim(1));
+ return Quadrature<1>(0);
+};
+
+
+
template <int dim>
-void
-QProjector<dim>::project_to_faces (const Quadrature<dim-1> &quadrature,
- typename std::vector<Point<dim> > &q_points)
+Quadrature<dim>
+QProjector<dim>::project_to_all_faces (const Quadrature<dim-1> &quadrature)
{
- unsigned int npt = quadrature.n_quadrature_points;
- unsigned int nf = GeometryInfo<dim>::faces_per_cell;
-
- q_points.resize (npt*nf);
- std::vector <Point<dim> > help(npt);
+ const unsigned int n_points = quadrature.n_quadrature_points,
+ n_faces = GeometryInfo<dim>::faces_per_cell;
+
+ // first fix quadrature points
+ typename std::vector<Point<dim> > q_points (n_points * n_faces);
+ std::vector <Point<dim> > help(n_points);
- unsigned k=0;
- for (unsigned int i=0;i<nf;++i)
+ // project to each face and copy
+ // results
+ for (unsigned int face=0; face<n_faces; ++face)
{
- project_to_face(quadrature, i, help);
- for (unsigned int j=0;j<npt;++j)
- q_points[k++] = help[j];
+ project_to_face(quadrature, face, help);
+ std::copy (help.begin(), help.end(), q_points.begin()+n_points*face);
}
+
+ // next copy over weights
+ typename std::vector<double> weights (n_points * n_faces);
+ for (unsigned int face=0; face<n_faces; ++face)
+ std::copy (quadrature.get_weights().begin(),
+ quadrature.get_weights().end(),
+ weights.begin()+n_points*face);
+
+ return Quadrature<dim>(q_points, weights);
}
+template <>
+Quadrature<1>
+QProjector<1>::project_to_all_subfaces (const Quadrature<0> &)
+{
+ Assert (false, ExcImpossibleInDim(1));
+ return Quadrature<1>(0);
+};
+
+
template <int dim>
-void
-QProjector<dim>::project_to_subfaces (const Quadrature<dim-1> &quadrature,
- typename std::vector<Point<dim> > &q_points)
+Quadrature<dim>
+QProjector<dim>::project_to_all_subfaces (const Quadrature<dim-1> &quadrature)
{
- unsigned int npt = quadrature.n_quadrature_points;
- unsigned int nf = GeometryInfo<dim>::faces_per_cell;
- unsigned int nc = GeometryInfo<dim>::subfaces_per_face;
+ const unsigned int n_points = quadrature.n_quadrature_points,
+ n_faces = GeometryInfo<dim>::faces_per_cell,
+ subfaces_per_face = GeometryInfo<dim>::subfaces_per_face;
- q_points.resize (npt*nf*nc);
- std::vector <Point<dim> > help(npt);
+ // first fix quadrature points
+ typename std::vector<Point<dim> > q_points (n_points * n_faces * subfaces_per_face);
+ std::vector <Point<dim> > help(n_points);
- unsigned k=0;
- for (unsigned int i=0;i<nf;++i)
- for (unsigned int c=0;c<nc;++c)
+ // project to each face and copy
+ // results
+ for (unsigned int face=0; face<n_faces; ++face)
+ for (unsigned int subface=0; subface<subfaces_per_face; ++subface)
{
- project_to_subface(quadrature, i, c, help);
- for (unsigned int j=0;j<npt;++j)
- q_points[k++] = help[j];
- }
+ project_to_subface(quadrature, face, subface, help);
+ std::copy (help.begin(), help.end(),
+ q_points.begin()+n_points*(face*subfaces_per_face+subface));
+ };
+
+ // next copy over weights
+ typename std::vector<double> weights (n_points * n_faces * subfaces_per_face);
+ for (unsigned int face=0; face<n_faces; ++face)
+ for (unsigned int subface=0; subface<subfaces_per_face; ++subface)
+ std::copy (quadrature.get_weights().begin(),
+ quadrature.get_weights().end(),
+ weights.begin()+n_points*(face*subfaces_per_face+subface));
+
+ return Quadrature<dim>(q_points, weights);
}
const Mapping<dim> &mapping,
const Quadrature<dim-1> &quadrature) const
{
- QProjector<dim> q(quadrature, false);
- return get_data (flags, mapping, q);
+ return get_data (flags, mapping,
+ QProjector<dim>::project_to_all_faces(quadrature));
}
const Mapping<dim> &mapping,
const Quadrature<dim-1> &quadrature) const
{
- QProjector<dim> q(quadrature, true);
- return get_data (flags, mapping, q);
+ return get_data (flags, mapping,
+ QProjector<dim>::project_to_all_subfaces(quadrature));
}
MappingCartesian<dim>::get_face_data (const UpdateFlags update_flags,
const Quadrature<dim-1>& quadrature) const
{
- QProjector<dim> q (quadrature, false);
- InternalData* data = new InternalData (q);
+ InternalData* data
+ = new InternalData (QProjector<dim>::project_to_all_faces(quadrature));
data->update_once = update_once(update_flags);
data->update_each = update_each(update_flags);
MappingCartesian<dim>::get_subface_data (const UpdateFlags update_flags,
const Quadrature<dim-1> &quadrature) const
{
- QProjector<dim> q (quadrature, true);
- InternalData* data = new InternalData (q);
+ InternalData* data
+ = new InternalData (QProjector<dim>::project_to_all_subfaces(quadrature));
data->update_once = update_once(update_flags);
data->update_each = update_each(update_flags);
const Quadrature<dim-1>& quadrature) const
{
InternalData *data = new InternalData(n_shape_functions);
- QProjector<dim> q (quadrature, false);
+ Quadrature<dim> q (QProjector<dim>::project_to_all_faces(quadrature));
compute_face_data (update_flags, q,
quadrature.n_quadrature_points, *data);
if (!use_mapping_q_on_all_cells)
const Quadrature<dim-1>& quadrature) const
{
InternalData *data = new InternalData(n_shape_functions);
- QProjector<dim> q (quadrature, true);
+ Quadrature<dim> q (QProjector<dim>::project_to_all_subfaces(quadrature));
compute_face_data (update_flags, q,
quadrature.n_quadrature_points, *data);
if (!use_mapping_q_on_all_cells)
const Quadrature<dim-1> &quadrature) const
{
InternalData* data = new InternalData(n_shape_functions);
- QProjector<dim> q (quadrature, false);
- compute_face_data (update_flags, q, quadrature.n_quadrature_points, *data);
+ compute_face_data (update_flags,
+ QProjector<dim>::project_to_all_faces(quadrature),
+ quadrature.n_quadrature_points,
+ *data);
return data;
}
const Quadrature<dim-1>& quadrature) const
{
InternalData* data = new InternalData(n_shape_functions);
- QProjector<dim> q (quadrature, true);
- compute_face_data (update_flags, q, quadrature.n_quadrature_points, *data);
+ compute_face_data (update_flags,
+ QProjector<dim>::project_to_all_subfaces(quadrature),
+ quadrature.n_quadrature_points,
+ *data);
return data;
}
// quadrature formula which will be
// used to probe boundary points at
// curved faces
- QProjector<dim> *q_projector=0;
+ Quadrature<dim> *q_projector=0;
if (mapping!=0)
{
typename std::vector<Point<dim-1> > boundary_points(n_points);
boundary_points[i](0)= 1.*(i+1)/(n_points+1);
Quadrature<dim-1> quadrature(boundary_points);
- q_projector = new QProjector<dim> (quadrature, false);
+ q_projector = new Quadrature<dim> (QProjector<dim>::project_to_all_faces(quadrature));
}
for (; cell!=endc; ++cell)
boundary_points[i](0) = 1.*(i+1)/(n_points+1);
Quadrature<dim-1> quadrature (boundary_points);
- QProjector<dim> q_projector (quadrature, false);
+ Quadrature<dim> q_projector (QProjector<dim>::project_to_all_faces(quadrature));
// next loop over all
// boundary faces and
for (unsigned int n=0; n<quadratures.size(); ++n)
{
- QProjector<dim> quadrature(*quadratures[n], sub);
+ Quadrature<dim> quadrature (sub == false?
+ QProjector<dim>::project_to_all_faces(*quadratures[n]) :
+ QProjector<dim>::project_to_all_subfaces(*quadratures[n]));
const std::vector<Point<dim> > &points=quadrature.get_points();
const std::vector<double> &weights=quadrature.get_weights();