]> https://gitweb.dealii.org/ - dealii.git/commitdiff
add QProjector::project_to_all_children
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 25 Nov 2009 04:25:16 +0000 (04:25 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Wed, 25 Nov 2009 04:25:16 +0000 (04:25 +0000)
git-svn-id: https://svn.dealii.org/trunk@20163 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/qprojector.h
deal.II/base/source/quadrature.cc

index 60e85bac46a5e36038a6c33a968917faf5ea3ec1..706f28d4ebbc2d5faa6bdda5a7ff9260154ab5ae 100644 (file)
@@ -2,7 +2,7 @@
 //    $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
@@ -158,23 +158,18 @@ class QProjector
                                      *
                                      * 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
@@ -187,18 +182,32 @@ class QProjector
     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
@@ -221,6 +230,28 @@ class QProjector
     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
index afa20c62e642ff3c6eda3032efd52b7d82015525..795b9484d7685c7b11f7fa0f6ef449e855ccd47a 100644 (file)
@@ -416,7 +416,7 @@ QProjector<1>::project_to_face (const Quadrature<0> &,
                                 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);
@@ -431,7 +431,7 @@ QProjector<2>::project_to_face (const Quadrature<1>      &quadrature,
                                 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()));
   
@@ -464,7 +464,7 @@ QProjector<3>::project_to_face (const Quadrature<2>    &quadrature,
                                 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()));
   
@@ -518,7 +518,7 @@ QProjector<1>::project_to_subface (const Quadrature<0> &,
                                   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);
@@ -535,8 +535,9 @@ QProjector<2>::project_to_subface (const Quadrature<1>    &quadrature,
                                   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()));
   
@@ -614,8 +615,8 @@ QProjector<3>::project_to_subface (const Quadrature<2>    &quadrature,
                                   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()));
 
@@ -1021,6 +1022,32 @@ QProjector<dim>::project_to_child (const Quadrature<dim>    &quadrature,
 }
 
 
+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>

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.