From: bangerth Date: Sat, 24 Feb 2007 14:54:43 +0000 (+0000) Subject: Add the possibility of a maximal number of cells upon mesh refinement. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=5ac186178666a34f44e54953d8e64333ff9da420;p=dealii-svn.git Add the possibility of a maximal number of cells upon mesh refinement. git-svn-id: https://svn.dealii.org/trunk@14498 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/grid/grid_refinement.h b/deal.II/deal.II/include/grid/grid_refinement.h index fd031eab82..6d17cbc12b 100644 --- a/deal.II/deal.II/include/grid/grid_refinement.h +++ b/deal.II/deal.II/include/grid/grid_refinement.h @@ -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 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. + * *
  • @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 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 &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::max()); /** * Refine the triangulation by @@ -220,6 +246,19 @@ class GridRefinement * *_fraction 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 &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::max()); diff --git a/deal.II/deal.II/source/grid/grid_refinement.cc b/deal.II/deal.II/source/grid/grid_refinement.cc index d624f626e9..136f421007 100644 --- a/deal.II/deal.II/source/grid/grid_refinement.cc +++ b/deal.II/deal.II/source/grid/grid_refinement.cc @@ -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 &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 &tria, Assert (top_fraction+bottom_fraction <= 1, ExcInvalidParameterValue()); Assert (criteria.is_non_negative (), ExcInvalidParameterValue()); - const int refine_cells=static_cast(top_fraction*criteria.size()); - const int coarsen_cells=static_cast(bottom_fraction*criteria.size()); + int refine_cells = static_cast(top_fraction*criteria.size()); + int coarsen_cells = static_cast(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::children_per_cell / + (GeometryInfo::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::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 + (tria.n_active_cells() + + refine_cells * (GeometryInfo::children_per_cell - 1) + - (coarsen_cells * + (GeometryInfo::children_per_cell - 1) / + GeometryInfo::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::children_per_cell - 1) + - (coarsen_cells * + (GeometryInfo::children_per_cell - 1) / + GeometryInfo::children_per_cell)); + refine_cells = static_cast (refine_cells * alpha); + coarsen_cells = static_cast (coarsen_cells * alpha); + } + if (refine_cells || coarsen_cells) { dealii::Vector tmp(criteria); @@ -204,16 +258,17 @@ GridRefinement::refine_and_coarsen_fixed_number (Triangulation &tria, tmp.end(), std::greater()); 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()); - 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 &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 &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 &tria, std::sort (tmp.begin(), tmp.end(), std::greater()); // compute thresholds - typename dealii::Vector::const_iterator pp=tmp.begin(); - for (double sum=0; (sum::const_iterator + pp=tmp.begin(); + for (double sum=0; + (sum::const_iterator qq=(tmp.end()-1); - for (double sum=0; (sum::const_iterator + qq=(tmp.end()-1); + for (double sum=0; + (sum + (tria.n_active_cells() + + refine_cells * (GeometryInfo::children_per_cell - 1) + - (coarsen_cells * + (GeometryInfo::children_per_cell - 1) / + GeometryInfo::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 &, const Vector &, const double, - const double); + const double, + const unsigned int); template void @@ -401,7 +493,8 @@ GridRefinement:: refine_and_coarsen_fixed_number (Triangulation &, const Vector &, const double, - const double); + const double, + const unsigned int); template void @@ -409,7 +502,8 @@ GridRefinement:: refine_and_coarsen_fixed_fraction (Triangulation &, const Vector &, const double, - const double); + const double, + const unsigned int); template void @@ -417,7 +511,8 @@ GridRefinement:: refine_and_coarsen_fixed_fraction (Triangulation &, const Vector &, const double, - const double); + const double, + const unsigned int); template void diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index 4ef647f001..7a17e56e14 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -1016,23 +1016,35 @@ inconvenience this causes.

    deal.II

      -
    1. New: Added function GridGenerator::hyper_cube_with_cylindrical_hole that produces - a square with a circular hole in the middle in 2d, and extrudes it along - the z-direction between 0 and L. +

    2. Improved: The + GridRefinement::refine_and_coarsen_fixed_number and + GridRefinement::refine_and_coarsen_fixed_fraction 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.
      - (Luca Heltai 2007/02/15) + (WB 2007/02/20)

      -
    3. WorkAround: The class GridOut::write_msh 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. -
      - (Luca Heltai 2007/02/15) -

      +
    4. New: Added function GridGenerator::hyper_cube_with_cylindrical_hole that produces + a square with a circular hole in the middle in 2d, and extrudes it + along the z-direction between 0 and L. +
      + (Luca Heltai 2007/02/15) +

      + +
    5. Workaround: The class GridOut::write_msh 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. +
      + (Luca Heltai 2007/02/15) +

    6. Fixed: The function DoFTools::

    7. Fixed: A local variable in - TriaAccessor<3,3>::measure and - erroneously marked as static could lead to - wrong results and crashes in multithread mode. This is now fixed. + TriaAccessor<3,3>::measure was + erroneously marked as static could lead to + wrong results and crashes in multithreaded mode. This is now fixed.
      (WB 2007/02/09)