]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add the possibility of a maximal number of cells upon mesh refinement.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 24 Feb 2007 14:54:43 +0000 (14:54 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Sat, 24 Feb 2007 14:54:43 +0000 (14:54 +0000)
git-svn-id: https://svn.dealii.org/trunk@14498 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_refinement.h
deal.II/deal.II/source/grid/grid_refinement.cc
deal.II/doc/news/changes.html

index fd031eab82969a7f9b28edad1669acf629adaeed..6d17cbc12b45b5bcd6594e20815ea82979c958c9 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -67,6 +67,14 @@ template <class T> class Vector;
  *     want to chose a smaller value to avoid overrefinement in regions which
  *     do not contribute much to the error.
  *
+ *     The function takes an additional last argument that can be used to
+ *     specify a maximal number of cells. If this number is going to be
+ *     exceeded upon refinement, then refinement and coarsening fractions are
+ *     going to be adjusted in an attempt to reach the maximum number of
+ *     cells. In practice, it is complicated to reach this number exactly, so
+ *     the argument is only an indication. The default value of this argument
+ *     is to impose no limit on the number of cells.
+ *
  *   <li> @p refine_and_coarsen_fixed_fraction: this function computes the
  *     threshold such that the number of cells getting flagged for refinement
  *     makes up for a certain fraction of the total error. If this fraction is 50
@@ -90,6 +98,10 @@ template <class T> class Vector;
  *     functionals, but may lead to very slow convergence of the grid
  *     if only few cells are refined in each step.
  *
+ *     The function takes an additional parameter indicating the maximum
+ *     number of cells we want, just as the previous function. See there for
+ *     more information.
+ *
  *     From the point of view of implementation, this time we really need to
  *     sort the array of criteria.
  *     Just like the other strategy described above, this function only
@@ -193,6 +205,19 @@ class GridRefinement
                                      * @p fraction_of_cells shall be
                                      * a value between zero and one.
                                      *
+                                     * The last argument can be used to
+                                     * specify a maximal number of cells. If
+                                     * this number is going to be exceeded
+                                     * upon refinement, then refinement and
+                                     * coarsening fractions are going to be
+                                     * adjusted in an attempt to reach the
+                                     * maximum number of cells. In practice,
+                                     * it is complicated to reach this number
+                                     * exactly, so the argument is only an
+                                     * indication. The default value of this
+                                     * argument is to impose no limit on the
+                                     * number of cells.
+                                     *
                                      * Refer to the general doc of
                                      * this class for more
                                      * information.
@@ -203,7 +228,8 @@ class GridRefinement
     refine_and_coarsen_fixed_number (Triangulation<dim> &tria,
                                      const Vector       &criteria,
                                      const double        top_fraction_of_cells,
-                                     const double        bottom_fraction_of_cells);
+                                     const double        bottom_fraction_of_cells,
+                                    const unsigned int  max_n_cells = std::numeric_limits<unsigned int>::max());
     
                                     /**
                                      * Refine the triangulation by
@@ -220,6 +246,19 @@ class GridRefinement
                                      * <tt>*_fraction</tt> shall be a
                                      * values between zero and one.
                                      *
+                                     * The last argument can be used to
+                                     * specify a maximal number of cells. If
+                                     * this number is going to be exceeded
+                                     * upon refinement, then refinement and
+                                     * coarsening fractions are going to be
+                                     * adjusted in an attempt to reach the
+                                     * maximum number of cells. In practice,
+                                     * it is complicated to reach this number
+                                     * exactly, so the argument is only an
+                                     * indication. The default value of this
+                                     * argument is to impose no limit on the
+                                     * number of cells.
+                                     *
                                      * Refer to the general doc of
                                      * this class for more
                                      * information.
@@ -230,7 +269,8 @@ class GridRefinement
     refine_and_coarsen_fixed_fraction (Triangulation<dim> &tria,
                                        const Vector       &criteria,
                                        const double        top_fraction,
-                                       const double        bottom_fraction);
+                                       const double        bottom_fraction,
+                                      const unsigned int  max_n_cells = std::numeric_limits<unsigned int>::max());
 
 
 
index d624f626e990c0f4c22fb47d5272f1263e635f7b..136f42100733c79b58bcdec69441307ba7a96902 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
+//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -183,7 +183,8 @@ void
 GridRefinement::refine_and_coarsen_fixed_number (Triangulation<dim> &tria,
                                                 const Vector       &criteria,
                                                 const double        top_fraction,
-                                                const double        bottom_fraction)
+                                                const double        bottom_fraction,
+                                                const unsigned int  max_n_cells)
 {
                                   // correct number of cells is
                                   // checked in @p{refine}
@@ -192,9 +193,62 @@ GridRefinement::refine_and_coarsen_fixed_number (Triangulation<dim> &tria,
   Assert (top_fraction+bottom_fraction <= 1, ExcInvalidParameterValue());
   Assert (criteria.is_non_negative (), ExcInvalidParameterValue());
 
-  const int refine_cells=static_cast<int>(top_fraction*criteria.size());
-  const int coarsen_cells=static_cast<int>(bottom_fraction*criteria.size());
+  int refine_cells  = static_cast<int>(top_fraction*criteria.size());
+  int coarsen_cells = static_cast<int>(bottom_fraction*criteria.size());
 
+                                  // first we have to see whether we
+                                  // currently already exceed the target
+                                  // number of cells
+  if (tria.n_active_cells() >= max_n_cells)
+    {
+                                      // if yes, then we need to stop
+                                      // refining cells and instead try to
+                                      // only coarsen as many as it would
+                                      // take to get to the target
+      refine_cells  = 0;
+      coarsen_cells = (tria.n_active_cells() - max_n_cells) *
+                     GeometryInfo<dim>::children_per_cell /
+                     (GeometryInfo<dim>::children_per_cell - 1);
+    }
+                                  // otherwise, see if we would exceed the
+                                  // maximum desired number of cells with the
+                                  // number of cells that are likely going to
+                                  // result from refinement. here, each cell
+                                  // to be refined is replaced by
+                                  // C=GeometryInfo<dim>::children_per_cell
+                                  // new cells, i.e. there will be C-1 more
+                                  // cells than before. similarly, C cells
+                                  // will be replaced by 1
+  else if (static_cast<unsigned int>
+          (tria.n_active_cells()
+           + refine_cells * (GeometryInfo<dim>::children_per_cell - 1)
+           - (coarsen_cells *
+              (GeometryInfo<dim>::children_per_cell - 1) /
+              GeometryInfo<dim>::children_per_cell))
+          >
+          max_n_cells)
+    {
+                                      // we have to adjust the
+                                      // fractions. assume we want
+                                      // alpha*refine_fraction and
+                                      // alpha*coarsen_fraction as new
+                                      // fractions and the resulting number
+                                      // of cells to be equal to
+                                      // max_n_cells. this leads to the
+                                      // following equation for lambda
+      const double alpha
+       =
+       1. *
+       (max_n_cells - tria.n_active_cells())
+       /
+       (refine_cells * (GeometryInfo<dim>::children_per_cell - 1)
+        - (coarsen_cells *
+           (GeometryInfo<dim>::children_per_cell - 1) /
+           GeometryInfo<dim>::children_per_cell));
+      refine_cells  = static_cast<int> (refine_cells * alpha);
+      coarsen_cells = static_cast<int> (coarsen_cells * alpha);
+    }
+  
   if (refine_cells || coarsen_cells)
     {
       dealii::Vector<typename Vector::value_type> tmp(criteria);
@@ -204,16 +258,17 @@ GridRefinement::refine_and_coarsen_fixed_number (Triangulation<dim> &tria,
                            tmp.end(),
                            std::greater<double>());
          refine (tria, criteria, *(tmp.begin() + refine_cells));
-       };
+       }
 
       if (coarsen_cells)
        {
          std::nth_element (tmp.begin(), tmp.begin()+tmp.size()-coarsen_cells,
                            tmp.end(),
                            std::greater<double>());
-         coarsen (tria, criteria, *(tmp.begin() + tmp.size() - coarsen_cells));
-       };
-    };
+         coarsen (tria, criteria,
+                  *(tmp.begin() + tmp.size() - coarsen_cells));
+       }
+    }
 }
 
 
@@ -223,7 +278,8 @@ void
 GridRefinement::refine_and_coarsen_fixed_fraction (Triangulation<dim> &tria,
                                                   const Vector       &criteria,
                                                   const double        top_fraction,
-                                                  const double        bottom_fraction)
+                                                  const double        bottom_fraction,
+                                                  const unsigned int  max_n_cells)
 {
                                   // correct number of cells is
                                   // checked in @p{refine}
@@ -231,7 +287,7 @@ GridRefinement::refine_and_coarsen_fixed_fraction (Triangulation<dim> &tria,
   Assert ((bottom_fraction>=0) && (bottom_fraction<=1), ExcInvalidParameterValue());
   Assert (top_fraction+bottom_fraction <= 1, ExcInvalidParameterValue());
   Assert (criteria.is_non_negative (), ExcInvalidParameterValue());
-
+  
                                   // let tmp be the cellwise square of the
                                   // error, which is what we have to sum
                                   // up and compare with
@@ -244,19 +300,54 @@ GridRefinement::refine_and_coarsen_fixed_fraction (Triangulation<dim> &tria,
   std::sort (tmp.begin(), tmp.end(), std::greater<double>());
 
                                   // compute thresholds
-  typename dealii::Vector<typename Vector::value_type>::const_iterator pp=tmp.begin();
-  for (double sum=0; (sum<top_fraction*total_error) && (pp!=(tmp.end()-1)); ++pp)
+  typename dealii::Vector<typename Vector::value_type>::const_iterator
+    pp=tmp.begin();
+  for (double sum=0;
+       (sum<top_fraction*total_error) && (pp!=(tmp.end()-1));
+       ++pp)
     sum += *pp;
   double top_threshold = ( pp != tmp.begin () ?
                           (*pp+*(pp-1))/2 :
                           *pp );
-  typename dealii::Vector<typename Vector::value_type>::const_iterator qq=(tmp.end()-1);
-  for (double sum=0; (sum<bottom_fraction*total_error) && (qq!=tmp.begin()); --qq)
+  typename dealii::Vector<typename Vector::value_type>::const_iterator
+    qq=(tmp.end()-1);
+  for (double sum=0;
+       (sum<bottom_fraction*total_error) && (qq!=tmp.begin());
+       --qq)
     sum += *qq;
   double bottom_threshold = ( qq != (tmp.end()-1) ?
                              (*qq + *(qq+1))/2 :
                              0);
 
+                                  // we now have an idea how many cells we
+                                  // are going to refine and coarsen. we use
+                                  // this information to see whether we are
+                                  // over the limit and if so use a function
+                                  // that knows how to deal with this
+                                  // situation
+  {
+    const unsigned int refine_cells  = pp - tmp.begin(),
+                      coarsen_cells = tmp.end() - qq;
+
+    if (static_cast<unsigned int>
+       (tria.n_active_cells()
+        + refine_cells * (GeometryInfo<dim>::children_per_cell - 1)
+        - (coarsen_cells *
+           (GeometryInfo<dim>::children_per_cell - 1) /
+           GeometryInfo<dim>::children_per_cell))
+       >
+       max_n_cells)
+      {
+       refine_and_coarsen_fixed_number (tria,
+                                        criteria,
+                                        1.*refine_cells/criteria.size(),
+                                        1.*coarsen_cells/criteria.size(),
+                                        max_n_cells);
+       return;
+      }
+  }
+  
+
                                   // in some rare cases it may happen that
                                   // both thresholds are the same (e.g. if
                                   // there are many cells with the same
@@ -393,7 +484,8 @@ GridRefinement::
 refine_and_coarsen_fixed_number (Triangulation<deal_II_dimension> &,
                                  const Vector<double> &,
                                  const double,
-                                 const double);
+                                 const double,
+                                const unsigned int);
 
 template
 void
@@ -401,7 +493,8 @@ GridRefinement::
 refine_and_coarsen_fixed_number (Triangulation<deal_II_dimension> &,
                                  const Vector<float> &,
                                  const double,
-                                 const double);
+                                 const double,
+                                const unsigned int);
 
 template
 void
@@ -409,7 +502,8 @@ GridRefinement::
 refine_and_coarsen_fixed_fraction (Triangulation<deal_II_dimension> &,
                                    const Vector<double> &,
                                    const double,
-                                   const double);
+                                   const double,
+                                  const unsigned int);
 
 template
 void
@@ -417,7 +511,8 @@ GridRefinement::
 refine_and_coarsen_fixed_fraction (Triangulation<deal_II_dimension> &,
                                    const Vector<float> &,
                                    const double,
-                                   const double);
+                                   const double,
+                                  const unsigned int);
 
 template
 void
index 4ef647f001628f9ce15e5354505597017d090d0d..7a17e56e14fe02f3df17d7b5177a5b2d5263c517 100644 (file)
@@ -1016,23 +1016,35 @@ inconvenience this causes.
 <h3>deal.II</h3>
 
 <ol>
-   <li> <p>New: Added function <code class="class">GridGenerator</code>::<code 
-       class="member">hyper_cube_with_cylindrical_hole</code> that produces 
-          a square with a circular hole in the middle in 2d, and extrudes it along
-          the z-direction between 0 and L. 
+   <li> <p>Improved: The
+       <code>GridRefinement::refine_and_coarsen_fixed_number</code> and
+       <code>GridRefinement::refine_and_coarsen_fixed_fraction</code> functions
+       have gained an additional last argument that can be used to specify a
+       maximal number of cells that we would like to use in a
+       triangulation. Its default value is set to indicate that no limit is
+       desired, as is the previous behavior.
        <br>
-       (Luca Heltai 2007/02/15)
+       (WB 2007/02/20)
        </p>
 
-   <li> <p>WorkAround: The class <code class="class">GridOut</code>::<code 
-       class="member">write_msh</code> produces a mesh which can be visualized in 
-          the Gmsh reader. A bug in Gmsh prevented the boundary indicator to be 
-          properly visualized. The boundary indicator was added to the material flag of the
-          faces (which is ignored when reading back the mesh) in order for the Gmsh reader
-          to be able to display the boundary indicator in a proper way.
-       <br>
-       (Luca Heltai 2007/02/15)
-       </p>
+   <li> <p>New: Added function <code class="class">GridGenerator</code>::<code 
+        class="member">hyper_cube_with_cylindrical_hole</code> that produces 
+       a square with a circular hole in the middle in 2d, and extrudes it
+       along the z-direction between 0 and L. 
+        <br>
+        (Luca Heltai 2007/02/15)
+        </p>
+
+   <li> <p>Workaround: The class <code class="class">GridOut</code>::<code 
+        class="member">write_msh</code> produces a mesh which can be visualized
+       in the Gmsh reader. A bug in Gmsh prevented the boundary indicator to
+       be properly visualized. The boundary indicator was added to the
+       material flag of the faces (which is ignored when reading back the
+       mesh) in order for the Gmsh reader to be able to display the boundary
+       indicator in a proper way. 
+        <br>
+        (Luca Heltai 2007/02/15)
+        </p>
 
 
    <li> <p>Fixed: The function <code class="class">DoFTools</code>::<code 
@@ -1044,9 +1056,9 @@ inconvenience this causes.
        </p>
 
    <li> <p>Fixed: A local variable in
-       <code>TriaAccessor&lt;3,3&gt;::measure</code> and 
-          erroneously marked as static could lead to
-       wrong results and crashes in multithread mode. This is now fixed.
+       <code>TriaAccessor&lt;3,3&gt;::measure</code> was
+       erroneously marked as static could lead to
+       wrong results and crashes in multithreaded mode. This is now fixed.
        <br>
        (WB 2007/02/09)
        </p>

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.