]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Patch by Christian Goll: Add an order parameter to GridRefinement::refine_and_coarsen...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 16 Dec 2011 14:06:08 +0000 (14:06 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 16 Dec 2011 14:06:08 +0000 (14:06 +0000)
git-svn-id: https://svn.dealii.org/trunk@24833 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 7ccff517e2eae7805f54368f87e06a14d01f0528..5a0bddcf522b012a3d682ff1b784dacec53837c7 100644 (file)
@@ -73,6 +73,14 @@ enabled due to a missing include file in file
 <h3>Specific improvements</h3>
 
 <ol>
+<li> Improved: The GridRefinement::refine_and_coarsen_optimize function
+assumed that the expected convergence order was 2. It has now gotten an
+argument by which the user can prescribe a different value. A bug has also
+been fixed in which the function incorrectly assumed in its algorithm that
+the mesh was two-dimensional.
+<br>
+(Christian Goll, 2011/12/16)
+
 <li> Fixed: Restriction and prolongation didn't work for elements of
 kind FE_Nothing. Consequently, many other parts of the library also
 didn't work, such as the SolutionTransfer class. This is now fixed.
index f67a38e72b80a50855d1da689b8cb0b238e63d79..88889d79538ab9f648da51df7347f01d12a22013 100644 (file)
@@ -1,7 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors
+//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -150,7 +150,7 @@ namespace GridRefinement
  * to sort the array of criteria.  Just like the other strategy
  * described above, this function only computes the threshold values
  * and then passes over to refine() and coarsen().
- *  
+ *
  * @arg @p criteria: the refinement criterion computed on each mesh
  * cell. Entries may not be negative.
  *
@@ -182,30 +182,34 @@ namespace GridRefinement
 
 
                                     /**
-                                     * Refine the triangulation by
-                                     * flagging certain cells to reach
-                                     * an optimal grid:
-                                     * We try to minimize the error
-                                     * multiplied with the number of
-                                     * cells in the new grid. All cells
-                                     * with large error indicator are
-                                     * refined to generate an optimal
-                                     * grid in the above sense.
-                                     * We assume that the error in one
-                                     * cell is reduced to a quarter
-                                     * after refinement.
-                                     * The new triangulation has three
-                                     * new cells for every flagged cell.
+                                     * Refine the triangulation by flagging
+                                     * certain cells to reach an optimal
+                                     * grid: We try to minimize the error
+                                     * multiplied with the number of cells in
+                                     * the new grid. All cells with large
+                                     * error indicator are refined to
+                                     * generate an optimal grid in the above
+                                     * sense.  We assume that the error in
+                                     * one cell is reduced to 1-2^{-order}
+                                     * after refinement, if 'order' is the
+                                     * expected order of convergence. This
+                                     * expected order of convergence must be
+                                     * passed as an argument but is defaulted
+                                     * to 2.  The new triangulation has
+                                     * ($2^d-1$) new cells for every flagged
+                                     * cell (the original cell is replaced by
+                                     * $2^d$ cells but it then made
+                                     * inactive).
                                      *
                                      * Refer to the general doc of
                                      * this class for more
                                      * information.
                                      */
-    
     template <int dim, class Vector, int spacedim>
     void
     refine_and_coarsen_optimize (Triangulation<dim,spacedim> &tria,
-                                 const Vector                &criteria);
+                                 const Vector                &criteria,
+                                 const unsigned int           order=2);
 
 /**
  * Flag all mesh cells for which the value in @p criteria exceeds @p
@@ -226,7 +230,7 @@ namespace GridRefinement
     void refine (Triangulation<dim,spacedim> &tria,
                 const Vector                &criteria,
                 const double                threshold);
-  
+
 /**
  * Flag all mesh cells for which the value in @p criteria
  * is less than @p threshold for coarsening.
@@ -246,14 +250,14 @@ namespace GridRefinement
     void coarsen (Triangulation<dim,spacedim> &tria,
                  const Vector                &criteria,
                  const double                threshold);
-  
+
                                   /**
                                    * An exception thrown if the
                                    * vector with cell criteria contains
                                    * negative values
                                    */
   DeclException0(ExcNegativeCriteria);
-  
+
                                   /**
                                    * One of the threshold parameters
                                    * causes trouble. Or the
index 6efcb42ac64fb18c315b1f2beae71f87010174b1..119c70b9f1e9a550e10f0debd55c7b0d117d1e44 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010 by the deal.II authors
+//    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -504,7 +504,8 @@ GridRefinement::refine_and_coarsen_fixed_fraction (Triangulation<dim,spacedim> &
 template <int dim, class Vector, int spacedim>
 void
 GridRefinement::refine_and_coarsen_optimize (Triangulation<dim,spacedim> &tria,
-                                            const Vector       &criteria)
+                                            const Vector       &criteria,
+                                            const unsigned int  order=2)
 {
   Assert (criteria.size() == tria.n_active_cells(),
          ExcDimensionMismatch(criteria.size(), tria.n_active_cells()));
@@ -518,7 +519,7 @@ GridRefinement::refine_and_coarsen_optimize (Triangulation<dim,spacedim> &tria,
 
   qsort_index(criteria,tmp,0,criteria.size()-1);
 
-  double s0 = 0.75 * criteria(tmp[0]);
+  double s0 = (1-std::pow(2,-order)) * criteria(tmp[0]);
   double E  = criteria.l1_norm();
 
   unsigned int N = criteria.size();
@@ -529,23 +530,24 @@ GridRefinement::refine_and_coarsen_optimize (Triangulation<dim,spacedim> &tria,
                                   // multiplied with the expected
                                   // number of cells.
                                   // We assume that the error is
-                                  // decreased by 3/4 a_K if the cell
+                                  // decreased by (1-2^(-order)) a_K if the cell
                                   // K with error indicator a_K is
-                                  // refined.
+                                  // refined and 'order' ist the expected
+                                   // order of convergence.
                                   // The expected number of cells is
-                                  // N+3*M (N is the current number
+                                  // N+(2^d-1)*M (N is the current number
                                   // of cells)
-  double min = (3.*(1.+M)+N) * (E-s0);
+  double min = ((std::pow(2.,dim)-1)*(1.+M)+N) * (E-s0);
 
   unsigned int minArg = N-1;
 
   for (M=1;M<criteria.size();++M)
     {
-      s0+= 0.75 * criteria(tmp[M]);
+      s0 += (1-std::pow(2.,-order)) * criteria(tmp[M]);
 
-      if ( (3.*(1+M)+N)*(E-s0) <= min)
+      if ( ((std::pow(2.,dim)-1)*(1+M)+N)*(E-s0) <= min)
        {
-         min = (3.*(1+M)+N)*(E-s0);
+         min = ((std::pow(2.,dim)-1)*(1+M)+N)*(E-s0);
          minArg = M;
        }
     }
index 61a4503643aad715ef2578d9cbbaf12c29e63b02..6d6b688c280685cd50f452105373bccf02a49692 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2010 by the deal.II authors
+//    Copyright (C) 2010, 2011 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -56,7 +56,8 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
   GridRefinement::
   refine_and_coarsen_optimize<deal_II_dimension,dealii::Vector<S>,deal_II_dimension>
   (Triangulation<deal_II_dimension> &,
-   const dealii::Vector<S> &);
+   const dealii::Vector<S> &,
+   const unsigned int);
 
 #if deal_II_dimension < 3
   template
@@ -100,6 +101,7 @@ for (S : REAL_SCALARS; deal_II_dimension : DIMENSIONS)
   GridRefinement::
   refine_and_coarsen_optimize<deal_II_dimension,dealii::Vector<S>,deal_II_dimension+1>
   (Triangulation<deal_II_dimension,deal_II_dimension+1> &,
-   const dealii::Vector<S> &);
+   const dealii::Vector<S> &,
+   const unsigned int);
 #endif
 }

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.