]> https://gitweb.dealii.org/ - dealii.git/commitdiff
use auto iterator
authornivesh <nidhinivesh.d@gmail.com>
Mon, 25 Jun 2018 14:18:05 +0000 (16:18 +0200)
committernivesh <nidhinivesh.d@gmail.com>
Mon, 25 Jun 2018 14:18:05 +0000 (16:18 +0200)
include/deal.II/fe/fe_enriched.h
source/fe/fe_enriched.cc
source/fe/fe_enriched.inst.in

index ecd651be019d8a50bb595641dccc7c66fba446d2..e5a801b7be1a1db0b7a351f1ddd4a8668b535a25 100644 (file)
@@ -767,7 +767,7 @@ namespace ColorEnriched
      * operator()) defining the subdomain 1. The function takes in a cell and
      * returns a boolean.
      * @param[in] predicate_2 Same as @p predicate_1 but defines subdomain 2.
-     * @return A boolean "true" if the subdomains share atleast a vertex.
+     * @return A boolean "true" if the subdomains share at least a vertex.
      */
     template <int dim, int spacedim>
     bool
@@ -779,7 +779,8 @@ namespace ColorEnriched
     /**
      * Assign colors to subdomains using Graph coloring algorithm where each
      * subdomain is considered as a graph node. Subdomains which are
-     * connected i.e share atleast a vertex have different color. Each subdomain
+     * connected i.e share at least a vertex have different color. Each
+     * subdomain
      * is defined using a predicate function of @p predicates.
      *
      * @param[in] dof_handler a hp::DoFHandler object
@@ -897,10 +898,10 @@ namespace ColorEnriched
         &fe_sets, // total list of color sets possible
       const std::vector<std::function<const Function<spacedim> *(
         const typename Triangulation<dim, spacedim>::cell_iterator &)>>
-        &color_enrichments,               // color wise enrichment functions
-      const FE_Q<dim, spacedim> &fe_base, // basic FE element
-      const FE_Q<dim, spacedim>
-        &fe_enriched, // FE element multiplied by enrichment function
+        &color_enrichments, // color wise enrichment functions
+      const FiniteElement<dim, spacedim> &fe_base, // basic FE element
+      const FiniteElement<dim, spacedim>
+        &fe_enriched, // FE multiplied by enrichment function
       const FE_Nothing<dim, spacedim> &fe_nothing,
       hp::FECollection<dim, spacedim> &fe_collection);
   } // namespace internal
@@ -1067,6 +1068,11 @@ namespace ColorEnriched
    * const hp::FECollection<dim>&
    * fe_collection(FE_helper.build_fe_collection(dof_handler));
    * @endcode
+   *
+   * @warning The current implementation relies on assigning each cell a
+   * material id, which shall not be modified after the setup
+   * and h-adaptive refinement. For a given cell, the material id is used
+   * to define color predicate map, which doesn't change with refinement.
    */
   template <int dim, int spacedim = dim>
   struct Helper
@@ -1081,8 +1087,8 @@ namespace ColorEnriched
      * belongs to a sub-domain with index (i).
      * @param enrichments std::vector of enrichment functions
      */
-    Helper(const FE_Q<dim, spacedim> &                             fe_base,
-           const FE_Q<dim, spacedim> &                             fe_enriched,
+    Helper(const FiniteElement<dim, spacedim> &                    fe_base,
+           const FiniteElement<dim, spacedim> &                    fe_enriched,
            const std::vector<predicate_function<dim, spacedim>> &  predicates,
            const std::vector<std::shared_ptr<Function<spacedim>>> &enrichments);
 
@@ -1091,8 +1097,6 @@ namespace ColorEnriched
      * mesh cells are initialized to work with
      * ColorEnriched::Helper<dim,spacedim>::fe_collection.
      *
-     * Returns an hp::FECollection object.
-     *
      * @param dof_handler an hp::DoFHandler object
      * @return hp::FECollection, a collection of
      * finite elements needed by @p dof_handler.
@@ -1111,13 +1115,13 @@ namespace ColorEnriched
      * A base FiniteElement used for constructing FE_Enriched
      * object required by ColorEnriched::Helper<dim,spacedim>::fe_collection.
      */
-    const FE_Q<dim, spacedim> &fe_base;
+    const FiniteElement<dim, spacedim> &fe_base;
 
     /**
      * An enriched FiniteElement used for constructing FE_Enriched
      * object required by ColorEnriched::Helper<dim,spacedim>::fe_collection.
      */
-    const FE_Q<dim, spacedim> &fe_enriched;
+    const FiniteElement<dim, spacedim> &fe_enriched;
 
     /**
      * A finite element with zero degrees of freedom used for
@@ -1157,10 +1161,9 @@ namespace ColorEnriched
      * and return a function pointer. These are needed while constructing
      * fe_collection.
      *
-     * color_enrichments[i](cell_iterator) calls the correct enrichment function
-     * (i.e. whose corresponding predicate has the color i) for the cell.
-     * Note that this call to cell_iterator returns the enrichment function
-     * which is a pointer to class derived from Function<dim>.
+     * color_enrichments[i](cell_iterator) returns a pointer to
+     * the correct enrichment function (i.e. whose corresponding
+     * predicate has the color i) for the cell.
      */
     std::vector<cell_iterator_function> color_enrichments;
 
index adf4ea21d4f90ad405353a572886ce7a4612b45f..2a920d9c5d8f0059830eb888fd57a6a5d0d8c34f 100644 (file)
@@ -1165,10 +1165,8 @@ namespace ColorEnriched
        * function ensures. The number of pairs is equal to the number
        * of predicates active in the given cell.
        */
-      unsigned int map_index(0);
-      auto         cell = dof_handler.begin_active();
-      auto         endc = dof_handler.end();
-      for (; cell != endc; ++cell)
+      unsigned int map_index = 0;
+      for (auto cell : dof_handler.active_cell_iterators())
         {
           // set default FE index ==> no enrichment and no active predicates
           cell->set_active_fe_index(0);
@@ -1189,6 +1187,8 @@ namespace ColorEnriched
           // pairs (1, 4) and (2, 5) at key 100 (i.e unique id of cell is
           // mapped with a map which associates color with predicate id)
           // Note that color list for the cell would be {1,2}.
+          std::map<unsigned int, unsigned int> &cell_map =
+            cellwise_color_predicate_map[map_index];
           for (unsigned int i = 0; i < predicates.size(); ++i)
             {
               if (predicates[i](cell))
@@ -1197,7 +1197,7 @@ namespace ColorEnriched
                    * create a pair predicate color and predicate id and add it
                    * to a map associated with each enriched cell
                    */
-                  auto ret = cellwise_color_predicate_map[map_index].insert(
+                  auto ret = cell_map.insert(
                     std::pair<unsigned int, unsigned int>(predicate_colors[i],
                                                           i));
 
@@ -1255,7 +1255,7 @@ namespace ColorEnriched
        * If we don't take further actions, we may find a dominating FE that
        * is too restrictive, i.e. enriched FE consisting of only FE_Nothing.
        * New elements needs to be added to FECollection object to help
-       * find the correct enriched FE element underlying the spaces in the
+       * find the correct enriched FE underlying the spaces in the
        * adjacent cells. This is done by creating an appropriate
        * set in fe_sets and a call to the function
        * make_fe_collection_from_colored_enrichments at a later stage.
@@ -1276,48 +1276,48 @@ namespace ColorEnriched
        */
 
       // loop through faces
-      for (auto cell = dof_handler.begin_active(); cell != dof_handler.end();
-           ++cell)
-        for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
-             ++face)
-          {
-            unsigned int fe_index = cell->active_fe_index();
-            auto         fe_set   = fe_sets[fe_index];
-
-            // find active indices of neighboring cells
-            if (!cell->at_boundary(face))
-              {
-                unsigned int nbr_fe_index =
-                  cell->neighbor(face)->active_fe_index();
-                // find corresponding fe set
-                auto nbr_fe_set = fe_sets[nbr_fe_index];
-
-                // find intersection of the fe_sets
-                std::set<unsigned int> intersection_set;
-                for (auto s : fe_set)
-                  if (nbr_fe_set.find(s) != nbr_fe_set.end())
-                    intersection_set.insert(s);
-
-                // add only the new sets
-                bool found = false;
-                if (!intersection_set.empty())
-                  {
-                    for (unsigned int j = 0; j < fe_sets.size(); ++j)
-                      {
-                        if (fe_sets[j] == intersection_set)
-                          {
-                            found = true;
-                            break;
-                          }
-                      }
-
-                    if (!found)
-                      {
-                        fe_sets.push_back(intersection_set);
-                      }
-                  }
-              }
-          }
+      for (auto cell : dof_handler.active_cell_iterators())
+        {
+          const unsigned int           fe_index = cell->active_fe_index();
+          const std::set<unsigned int> fe_set   = fe_sets.at(fe_index);
+          for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell;
+               ++face)
+            {
+              // find active indices of neighboring cells
+              if (!cell->at_boundary(face))
+                {
+                  const unsigned int nbr_fe_index =
+                    cell->neighbor(face)->active_fe_index();
+                  // find corresponding fe set
+                  auto nbr_fe_set = fe_sets.at(nbr_fe_index);
+
+                  // find intersection of the fe_sets
+                  std::set<unsigned int> intersection_set;
+                  for (auto s : fe_set)
+                    if (nbr_fe_set.find(s) != nbr_fe_set.end())
+                      intersection_set.insert(s);
+
+                  // add only the new sets
+                  bool found = false;
+                  if (!intersection_set.empty())
+                    {
+                      for (unsigned int j = 0; j < fe_sets.size(); ++j)
+                        {
+                          if (fe_sets[j] == intersection_set)
+                            {
+                              found = true;
+                              break;
+                            }
+                        }
+
+                      if (!found)
+                        {
+                          fe_sets.push_back(intersection_set);
+                        }
+                    }
+                }
+            }
+        }
     }
 
 
@@ -1373,11 +1373,11 @@ namespace ColorEnriched
       const std::vector<std::set<unsigned int>> &fe_sets,
       const std::vector<std::function<const Function<spacedim> *(
         const typename Triangulation<dim, spacedim>::cell_iterator &)>>
-        &                              color_enrichments,
-      const FE_Q<dim, spacedim> &      fe_base,
-      const FE_Q<dim, spacedim> &      fe_enriched,
-      const FE_Nothing<dim, spacedim> &fe_nothing,
-      hp::FECollection<dim, spacedim> &fe_collection)
+        &                                 color_enrichments,
+      const FiniteElement<dim, spacedim> &fe_base,
+      const FiniteElement<dim, spacedim> &fe_enriched,
+      const FE_Nothing<dim, spacedim> &   fe_nothing,
+      hp::FECollection<dim, spacedim> &   fe_collection)
     {
       // define dummy function which is associated with FE_Nothing
       using cell_function = std::function<const Function<spacedim> *(
@@ -1443,8 +1443,8 @@ namespace ColorEnriched
 
   template <int dim, int spacedim>
   Helper<dim, spacedim>::Helper(
-    const FE_Q<dim, spacedim> &                             fe_base,
-    const FE_Q<dim, spacedim> &                             fe_enriched,
+    const FiniteElement<dim, spacedim> &                    fe_base,
+    const FiniteElement<dim, spacedim> &                    fe_enriched,
     const std::vector<predicate_function<dim, spacedim>> &  predicates,
     const std::vector<std::shared_ptr<Function<spacedim>>> &enrichments)
     : fe_base(fe_base)
index 9371ca024b416bf190655f3e53f77caee41d51d1..cf527af200dc48ae3df71e0177731d91c58eeeea 100644 (file)
@@ -78,8 +78,10 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
               const typename Triangulation<
                 deal_II_dimension,
                 deal_II_space_dimension>::cell_iterator &)>> &color_enrichments,
-          const FE_Q<deal_II_dimension, deal_II_space_dimension> &fe_base,
-          const FE_Q<deal_II_dimension, deal_II_space_dimension> &fe_enriched,
+          const FiniteElement<deal_II_dimension, deal_II_space_dimension>
+            &fe_base,
+          const FiniteElement<deal_II_dimension, deal_II_space_dimension>
+            &fe_enriched,
           const FE_Nothing<deal_II_dimension, deal_II_space_dimension>
             &fe_nothing,
           hp::FECollection<deal_II_dimension, deal_II_space_dimension>

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.