]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Prepare structure for coarsening of grid (not yet implemented, not documented).
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 28 May 1998 13:58:50 +0000 (13:58 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Thu, 28 May 1998 13:58:50 +0000 (13:58 +0000)
git-svn-id: https://svn.dealii.org/trunk@371 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/Todo
deal.II/deal.II/include/grid/tria.h
deal.II/deal.II/include/grid/tria_accessor.h
deal.II/deal.II/source/grid/tria.cc
deal.II/deal.II/source/grid/tria_accessor.cc

index 79ea4e7166f9212909e57076eb9700e99cf4247f..114b8ecd12f8757adaa7a8dd89ff737bd4e2078b 100644 (file)
@@ -106,6 +106,9 @@ Remove the this-> coding in tria_iterator.templates.h. These were
 Re-enable printing of a preamble to ucd files in data_io.cc.
 
 
+Implement coarsening of grids and update docs for that.
+
+
 
 
 DEAL:
index ca89c5e6e117306860eacbbcd1a55e95264d7887..781f87189a2b018e39d280c5b95cf816204244ed 100644 (file)
@@ -83,6 +83,13 @@ class TriangulationLevel<0> {
                                      *  that of the #quads# vector, etc.
                                      */
     vector<bool> refine_flags;
+
+                                    /**
+                                     * Same meaning as the one above, but
+                                     * specifies whether a cell must be
+                                     * coarsened.
+                                     */
+    vector<bool> coarsen_flags;
     
                                     /**
                                      *  Levels and indices of the neighbors
@@ -1075,7 +1082,11 @@ class TriaDimensionInfo<2> {
  *
  *   You may write other information to the output file between different sets
  *   of refinement information, as long as you read it upon re-creation of the
- *   grid.
+ *   grid. You should make sure that the other information in the new
+ *   triangulation which is to be created from the saves flags, matches that of
+ *   the old triangulation, for example the smoothing level; if not, the
+ *   cells actually created from the flags may be other ones, since smoothing
+ *   adds additional cells, but the number depending on the smoothing level.
  *
  *
  *   \subsection{User flags}
@@ -1507,6 +1518,21 @@ class Triangulation : public TriaDimensionInfo<dim> {
                                      *  This function is dimension specific.
                                      */ 
     void execute_refinement ();
+
+                                    /**
+                                     * Coarsen all cells which were flagged for
+                                     * coarsening, or rather: delete all
+                                     * children of those cells of which all
+                                     * child cells are flagged for coarsening
+                                     * and several other constraints hold (see
+                                     * the general doc of this class).
+                                     *
+                                     * The function resets all refinement
+                                     * flags to false.
+                                     *
+                                     * This function is dimension specific.
+                                     */
+    void execute_coarsening ();
                                     /*@}*/
 
                                     /**
index 28566d98c60cd950a4b9f917a1064a742d27f71a..042ea3a921c4abb3c27eb6daa5908945dbf9081f 100644 (file)
@@ -884,8 +884,9 @@ class CellAccessor :  public TriaSubstructAccessor<dim> {
     bool refine_flag_set () const;
 
                                     /**
-                                     *  Flag the cell pointed to
-                                     *  for refinement.
+                                     *  Flag the cell pointed to fo
+                                     *  refinement. This function is only
+                                     *  allowed for active cells.
                                      */
     void set_refine_flag () const;
 
@@ -894,6 +895,25 @@ class CellAccessor :  public TriaSubstructAccessor<dim> {
                                      */
     void clear_refine_flag () const;
 
+
+                                    /**
+                                     *  Return whether the coarsen flag
+                                     *  is set or not.
+                                     */
+    bool coarsen_flag_set () const;
+
+                                    /**
+                                     *  Flag the cell pointed to for
+                                     *  coarsening. This function is only
+                                     *  allowed for active cells.
+                                     */
+    void set_coarsen_flag () const;
+
+                                    /**
+                                     *  Clear the coarsen flag.
+                                     */
+    void clear_coarsen_flag () const;
+
                                     /**
                                      *  Return a pointer to the #i#th
                                      *  child. Overloaded version which returns
index 0bd077b5d43971565e6fd35fec8081553dfaff36..c133ffef7790db4313b6dfe24d2c0193acbc7f9b 100644 (file)
@@ -2137,6 +2137,17 @@ void Triangulation<1>::execute_refinement () {
 #ifdef DEBUG
   for (unsigned int level=0; level<levels.size(); ++level) 
     levels[level]->monitor_memory (1);
+
+                                  // check whether really all refinement flags
+                                  // are reset (also of previously non-active
+                                  // cells which we may not have touched. If
+                                  // the refinement flag of a non-active cell
+                                  // is set, something went wrong since the
+                                  // cell-accessors should have caught this)
+  cell_iterator cell = begin(),
+               endc = end();
+  while (cell != endc)
+    Assert (!(cell++)->refine_flag_set(), ExcInternalError ());  
 #endif
 };
 
@@ -2710,10 +2721,28 @@ void Triangulation<2>::execute_refinement () {
 #ifdef DEBUG
   for (unsigned int level=0; level<levels.size(); ++level) 
     levels[level]->monitor_memory (2);
+
+                                  // check whether really all refinement flags
+                                  // are reset (also of previously non-active
+                                  // cells which we may not have touched. If
+                                  // the refinement flag of a non-active cell
+                                  // is set, something went wrong since the
+                                  // cell-accessors should have caught this)
+  cell_iterator cell = begin(),
+               endc = end();
+  while (cell != endc)
+    Assert (!(cell++)->refine_flag_set(), ExcInternalError ());
 #endif
 };
 
-  
+
+
+
+template <int dim>
+void Triangulation<dim>::execute_coarsening () {
+  Assert (false, ExcInternalError());
+};
+
 
 
 template <int dim>
@@ -2815,6 +2844,11 @@ void TriangulationLevel<0>::reserve_space (const unsigned int total_cells,
   refine_flags.insert (refine_flags.end(),
                       total_cells - refine_flags.size(),
                       false);
+  
+  coarsen_flags.reserve (total_cells);
+  coarsen_flags.insert (coarsen_flags.end(),
+                       total_cells - coarsen_flags.size(),
+                       false);
 
   neighbors.reserve (total_cells*(2*dimension));
   neighbors.insert (neighbors.end(),
@@ -2830,12 +2864,18 @@ void TriangulationLevel<0>::monitor_memory (const unsigned int true_dimension) c
 //       refine_flags.size()<256,
 //       ExcMemoryWasted ("refine_flags",
 //                        refine_flags.size(), refine_flags.capacity()));
+//  Assert (coarsen_flags.size() == coarsen_flags.capacity() ||
+//       coarsen_flags.size()<256,
+//       ExcMemoryWasted ("coarsen_flags",
+//                        coarsen_flags.size(), coarsen_flags.capacity()));
 //  Assert (neighbors.size() ==  neighbors.capacity() ||
 //       neighbors.size()<256,
 //       ExcMemoryWasted ("neighbors",
 //                        neighbors.size(), neighbors.capacity()));
   Assert (2*true_dimension*refine_flags.size() == neighbors.size(),
          ExcMemoryInexact (refine_flags.size(), neighbors.size()));
+  Assert (2*true_dimension*coarsen_flags.size() == neighbors.size(),
+         ExcMemoryInexact (coarsen_flags.size(), neighbors.size()));
 };
 
 
index 9f23965d26dfd67b974ddf1206c629c84f8df7d4..75d49dbf646a55716765919a4b58139be316f549 100644 (file)
@@ -651,7 +651,13 @@ bool CellAccessor<dim>::at_boundary (const unsigned int i) const {
 
 template <int dim>
 bool CellAccessor<dim>::refine_flag_set () const {
-  Assert (used() && active(),
+  Assert (used(), typename TriaSubstructAccessor<dim>::ExcCellNotUsed());
+                                  // cells flagged for refinement must be active
+                                  // (the #set_refine_flag# function checks this,
+                                  // but activity may change when refinement is
+                                  // executed and for some reason the refine
+                                  // flag is not cleared).
+  Assert (active() ||  !tria->levels[present_level]->refine_flags[present_index],
          typename TriaSubstructAccessor<dim>::ExcRefineCellNotActive());
   return tria->levels[present_level]->refine_flags[present_index];
 };
@@ -676,6 +682,39 @@ void CellAccessor<dim>::clear_refine_flag () const {
 
 
 
+template <int dim>
+bool CellAccessor<dim>::coarsen_flag_set () const {
+  Assert (used(), typename TriaSubstructAccessor<dim>::ExcCellNotUsed());
+                                  // cells flagged for coarsening must be active
+                                  // (the #set_refine_flag# function checks this,
+                                  // but activity may change when refinement is
+                                  // executed and for some reason the refine
+                                  // flag is not cleared).
+  Assert (active() ||  !tria->levels[present_level]->coarsen_flags[present_index],
+         typename TriaSubstructAccessor<dim>::ExcRefineCellNotActive());
+  return tria->levels[present_level]->coarsen_flags[present_index];
+};
+
+
+
+template <int dim>
+void CellAccessor<dim>::set_coarsen_flag () const {
+  Assert (used() && active(),
+         typename TriaSubstructAccessor<dim>::ExcRefineCellNotActive());
+  tria->levels[present_level]->coarsen_flags[present_index] = true;
+};
+
+
+
+template <int dim>
+void CellAccessor<dim>::clear_coarsen_flag () const {
+  Assert (used() && active(),
+         typename TriaSubstructAccessor<dim>::ExcRefineCellNotActive());
+  tria->levels[present_level]->coarsen_flags[present_index] = false;
+};
+
+
+
 template <int dim>
 TriaIterator<dim,CellAccessor<dim> >
 CellAccessor<dim>::neighbor (const unsigned int i) const {

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.