<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.
//---------------------------------------------------------------------------
// $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
* 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.
*
/**
- * 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
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.
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
// $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
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()));
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();
// 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;
}
}
// $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
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
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
}