From: wolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Thu, 23 Aug 2001 20:13:03 +0000 (+0000)
Subject: QProjector::project_to_child
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7f10011d4c32a7f39a9b37c4424524cc0e157481;p=dealii-svn.git

QProjector::project_to_child


git-svn-id: https://svn.dealii.org/trunk@4906 0785d39b-7218-0410-832d-ea1e28bc413d
---

diff --git a/deal.II/base/include/base/quadrature.h b/deal.II/base/include/base/quadrature.h
index 56bb445fbc..5a9d61f3f5 100644
--- a/deal.II/base/include/base/quadrature.h
+++ b/deal.II/base/include/base/quadrature.h
@@ -279,7 +279,6 @@ template <int dim>
 class QProjector
 {
   public:
-    
 				     /**
 				      * Compute the quadrature points
 				      * on the cell if the given
@@ -347,6 +346,30 @@ class QProjector
 				      */
     static Quadrature<dim>
     project_to_all_subfaces (const Quadrature<dim-1> &quadrature);
+
+				     /**
+				      * Project a give quadrature
+				      * formula to a child of a
+				      * cell. You may want to use this
+				      * function in case you want to
+				      * extend an integral only over
+				      * the area which a potential
+				      * child would occupy. The child
+				      * numbering is the same as the
+				      * children would be numbered
+				      * upon refinement of the cell.
+				      *
+				      * As integration using this
+				      * quadrature formula now only
+				      * extends over a fraction of the
+				      * cell, the weights of the
+				      * resulting object are scaled by
+				      * @p{1./GeometryInfo<dim>::children_per_cell}.
+				      */
+    static
+    Quadrature<dim>
+    project_to_child (const Quadrature<dim>  &quadrature,
+		      const unsigned int      child_no);
 };
 
 
@@ -369,6 +392,7 @@ template <>
 double Quadrature<0>::weight (const unsigned int) const;
 template <>
 const std::vector<double> & Quadrature<0>::get_weights () const;
+
 template <>
 void QProjector<1>::project_to_face (const Quadrature<0> &,
 				     const unsigned int,
diff --git a/deal.II/base/source/quadrature.cc b/deal.II/base/source/quadrature.cc
index 62ab8d2fd9..63e222cbc4 100644
--- a/deal.II/base/source/quadrature.cc
+++ b/deal.II/base/source/quadrature.cc
@@ -240,6 +240,8 @@ void QProjector<2>::project_to_face (const Quadrature<1>      &quadrature,
 {
   const unsigned int dim=2;
   Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
+  Assert (q_points.size() == quadrature.n_quadrature_points,
+	  ExcDimensionMismatch (q_points.size(), quadrature.n_quadrature_points));
   
   for (unsigned int p=0; p<quadrature.n_quadrature_points; ++p)
     switch (face_no)
@@ -270,6 +272,8 @@ void QProjector<3>::project_to_face (const Quadrature<2>    &quadrature,
 {
   const unsigned int dim=3;
   Assert (face_no<2*dim, ExcIndexRange (face_no, 0, 2*dim));
+  Assert (q_points.size() == quadrature.n_quadrature_points,
+	  ExcDimensionMismatch (q_points.size(), quadrature.n_quadrature_points));
   
   for (unsigned int p=0; p<quadrature.n_quadrature_points; ++p)
     switch (face_no)
@@ -332,6 +336,8 @@ void QProjector<2>::project_to_subface (const Quadrature<1>    &quadrature,
   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)));
+  Assert (q_points.size() == quadrature.n_quadrature_points,
+	  ExcDimensionMismatch (q_points.size(), quadrature.n_quadrature_points));
   
   for (unsigned int p=0; p<quadrature.n_quadrature_points; ++p)
     switch (face_no)
@@ -406,6 +412,8 @@ void QProjector<3>::project_to_subface (const Quadrature<2>    &quadrature,
   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)));
+  Assert (q_points.size() == quadrature.n_quadrature_points,
+	  ExcDimensionMismatch (q_points.size(), quadrature.n_quadrature_points));
 
 
 				   // for all faces and subfaces:
@@ -641,6 +649,60 @@ QProjector<dim>::project_to_all_subfaces (const Quadrature<dim-1> &quadrature)
 
 
 
+template <int dim>
+Quadrature<dim>
+QProjector<dim>::project_to_child (const Quadrature<dim>    &quadrature,
+				   const unsigned int        child_no)
+{
+  Assert (child_no < GeometryInfo<dim>::children_per_cell,
+	  ExcIndexRange (child_no, 0, GeometryInfo<dim>::children_per_cell));
+  
+  const unsigned int n_q_points = quadrature.n_quadrature_points;  
+    
+				   // projection is simple: copy over
+				   // old points, scale them, and if
+				   // necessary shift them
+  std::vector<Point<dim> > q_points = quadrature.get_points ();
+  for (unsigned int i=0; i<n_q_points; ++i)
+    q_points[i] /= 2;
+
+  switch (dim)
+    {
+      case 1:
+      {
+	if (child_no == 1)
+	  for (unsigned int i=0; i<n_q_points; ++i)
+	    q_points[i][0] += 0.5;
+	break;
+      };
+
+      case 2:
+      {
+	if ((child_no == 1) || (child_no == 2))
+	  for (unsigned int i=0; i<n_q_points; ++i)
+	    q_points[i][0] += 0.5;
+	if ((child_no == 2) || (child_no == 3))
+	  for (unsigned int i=0; i<n_q_points; ++i)
+	    q_points[i][1] += 0.5;
+	break;
+      };
+
+      default:
+	    Assert (false, ExcNotImplemented());
+    };
+
+				   // for the weights, things are
+				   // equally simple: copy them and
+				   // scale them
+  std::vector<double> weights = quadrature.get_weights ();
+  for (unsigned int i=0; i<n_q_points; ++i)
+    weights[i] *= (1./GeometryInfo<dim>::children_per_cell);
+
+  return Quadrature<dim> (q_points, weights);
+};
+
+
+
 // ------------------------------------------------------------ //
 
 
diff --git a/deal.II/doc/news/2001/c-3-1.html b/deal.II/doc/news/2001/c-3-1.html
index 77bb3a4f44..6da9f8736d 100644
--- a/deal.II/doc/news/2001/c-3-1.html
+++ b/deal.II/doc/news/2001/c-3-1.html
@@ -158,6 +158,13 @@ documentation, etc</a>.
 <h3>base</h3>
 
 <ol>
+  <li> <p> 
+       New: Function <code class="member">QProjector::project_to_child</code>
+       generates quadrature formulae which act on the area which a
+       child would occupy.
+       <br>
+       (WB 2001/08/23)
+
   <li> <p>
        Changed: The examples classes in the base directory are now
        moved into a namespace <code class="class">Functions</code> of