// $Id$
// Version: $Name$
//
-// Copyright (C) 2005, 2006, 2007, 2008 by the deal.II authors
+// Copyright (C) 2005, 2006, 2007, 2008, 2009 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
*
* 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
- * 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.
+ * original weights. Thus, the
+ * sum of the weights is not one,
+ * but the number of faces, which
+ * is the surface of the
+ * reference cell.
+ *
+ * This in particular allows us
+ * to extract a subset of points
+ * corresponding to a single face
+ * and use it as a quadrature on
+ * this face, as is done in
+ * FEFaceValues.
*
* @note In 3D, this function
* produces eight sets of
project_to_all_faces (const SubQuadrature &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.
+ * 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
+ * subfaces.
+ *
+ * Like in project_to_all_faces(),
+ * the weights of the new rule
+ * sum up to the number of faces
+ * (not subfaces), which
+ * is the surface of the
+ * reference cell.
+ *
+ * This in particular allows us
+ * to extract a subset of points
+ * corresponding to a single subface
+ * and use it as a quadrature on
+ * this face, as is done in
+ * FESubfaceValues.
*/
static Quadrature<dim>
project_to_all_subfaces (const SubQuadrature &quadrature);
/**
- * Project a give quadrature
+ * Project a given quadrature
* formula to a child of a
* cell. You may want to use this
* function in case you want to
project_to_child (const Quadrature<dim> &quadrature,
const unsigned int child_no);
+ /**
+ * Project a quadrature rule to
+ * all children of a
+ * cell. Similarly to
+ * project_to_all_subfaces(),
+ * this function replicates the
+ * formula generated by
+ * project_to_child() for all
+ * children, such that the
+ * weights sum up to one, the
+ * volume of the total cell
+ * again.
+ *
+ * The child
+ * numbering is the same as the
+ * children would be numbered
+ * upon refinement of the cell.
+ */
+ static
+ Quadrature<dim>
+ project_to_all_children (const Quadrature<dim> &quadrature);
+
/**
* Project the onedimensional
* rule <tt>quadrature</tt> to
std::vector<Point<1> > &q_points)
{
const unsigned int dim=1;
- Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
+ AssertIndexRange (face_no, GeometryInfo<dim>::faces_per_cell);
AssertDimension (q_points.size(), 1);
q_points[0] = Point<dim>((double) face_no);
std::vector<Point<2> > &q_points)
{
const unsigned int dim=2;
- Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
+ AssertIndexRange (face_no, GeometryInfo<dim>::faces_per_cell);
Assert (q_points.size() == quadrature.size(),
ExcDimensionMismatch (q_points.size(), quadrature.size()));
std::vector<Point<3> > &q_points)
{
const unsigned int dim=3;
- Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
+ AssertIndexRange (face_no, GeometryInfo<dim>::faces_per_cell);
Assert (q_points.size() == quadrature.size(),
ExcDimensionMismatch (q_points.size(), quadrature.size()));
const RefinementCase<0> &)
{
const unsigned int dim=1;
- Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
+ AssertIndexRange (face_no, GeometryInfo<dim>::faces_per_cell);
AssertDimension (q_points.size(), 1);
q_points[0] = Point<dim>((double) face_no);
const RefinementCase<1> &)
{
const unsigned int dim=2;
- Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
- Assert (subface_no<(1<<(dim-1)), ExcIndexRange (face_no, 0, 1<<(dim-1)));
+ AssertIndexRange (face_no, GeometryInfo<dim>::faces_per_cell);
+ AssertIndexRange (subface_no, GeometryInfo<dim>::max_children_per_face);
+
Assert (q_points.size() == quadrature.size(),
ExcDimensionMismatch (q_points.size(), quadrature.size()));
const RefinementCase<2> &ref_case)
{
const unsigned int dim=3;
- Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
- Assert (subface_no<(1<<(dim-1)), ExcIndexRange (face_no, 0, 1<<(dim-1)));
+ AssertIndexRange (face_no, GeometryInfo<dim>::faces_per_cell);
+ AssertIndexRange (subface_no, GeometryInfo<dim>::max_children_per_face);
Assert (q_points.size() == quadrature.size(),
ExcDimensionMismatch (q_points.size(), quadrature.size()));
}
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_all_children (const Quadrature<dim> &quadrature)
+{
+ const unsigned int n_points = quadrature.size(),
+ n_children = GeometryInfo<dim>::max_children_per_cell;
+
+ std::vector<Point<dim> > q_points(n_points * n_children);
+ std::vector<double> weights(n_points * n_children);
+
+ // project to each child and copy
+ // results
+ for (unsigned int child=0; child<n_children; ++child)
+ {
+ Quadrature<dim> help = project_to_child(quadrature, child);
+ for (unsigned int i=0;i<n_points;++i)
+ {
+ q_points[child*n_points+i] = help.point(i);
+ weights[child*n_points+i] = help.weight(i);
+ }
+ }
+ return Quadrature<dim>(q_points, weights);
+}
+
+
+
template <int dim>
Quadrature<dim>