]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement the hp::DoFHandler cell dof indices cache in the same way as we already...
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 10 Sep 2013 18:13:31 +0000 (18:13 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 10 Sep 2013 18:13:31 +0000 (18:13 +0000)
git-svn-id: https://svn.dealii.org/trunk@30666 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/include/deal.II/dofs/dof_accessor.templates.h
deal.II/include/deal.II/hp/dof_level.h
deal.II/source/hp/dof_handler.cc

index a2a4c03f1088bd594ef7eb2faa977efeaa54ff8a..5fee314d9e6d206850048a96a974c54dda82bf29 100644 (file)
@@ -2556,15 +2556,44 @@ namespace internal
       template <int dim, int spacedim, bool lda>
       static
       void
-      update_cell_dof_indices_cache (const DoFCellAccessor<dealii::hp::DoFHandler<dim,spacedim>, lda> &)
+      update_cell_dof_indices_cache (const DoFCellAccessor<dealii::hp::DoFHandler<dim,spacedim>, lda> &accessor)
       {
-//TODO[WB]: should implement a dof indices cache for hp as well
+        // caches are only for cells with DoFs, i.e., for active ones
+        if (accessor.has_children())
+          return;
 
-        // not implemented, but should also
-        // not be called
-        Assert (false, ExcNotImplemented());
+        const unsigned int dofs_per_cell   = accessor.get_fe().dofs_per_cell;
+
+        // make sure the cache is at least
+        // as big as we need it when
+        // writing to the last element of
+        // this cell
+        Assert (static_cast<unsigned int>(accessor.present_index)
+                <
+                accessor.dof_handler->levels[accessor.present_level]
+                                             ->cell_cache_offsets.size(),
+                ExcInternalError());
+        Assert (accessor.dof_handler->levels[accessor.present_level]
+                ->cell_cache_offsets[accessor.present_index]
+                <=
+                accessor.dof_handler->levels[accessor.present_level]
+                ->cell_dof_indices_cache.size(),
+                ExcInternalError());
+
+        std::vector<types::global_dof_index> dof_indices (dofs_per_cell);
+        accessor.dealii::DoFAccessor<dim,dealii::hp::DoFHandler<dim,spacedim>,lda>::get_dof_indices
+        (dof_indices, accessor.active_fe_index());
+
+        types::global_dof_index *next_dof_index
+        = &accessor.dof_handler->levels[accessor.present_level]
+           ->cell_dof_indices_cache[accessor.dof_handler->levels[accessor.present_level]
+                                                                 ->cell_cache_offsets[accessor.present_index]];
+        for (unsigned int i=0; i<dofs_per_cell; ++i, ++next_dof_index)
+          *next_dof_index = dof_indices[i];
       }
 
+
+
       /**
        * Implement setting dof
        * indices on a
@@ -2763,14 +2792,15 @@ namespace internal
       get_dof_indices (const DoFCellAccessor<dealii::hp::DoFHandler<dim,spacedim>, lda> &accessor,
                        std::vector<types::global_dof_index>                            &dof_indices)
       {
-        // no caching for
-        // hp::DoFHandler
-        // implemented
-        typedef
-        dealii::DoFAccessor<dim,dealii::hp::DoFHandler<dim,spacedim>, lda>
-        DoFAccessor;
-        accessor.DoFAccessor::get_dof_indices (dof_indices,
-                                               accessor.active_fe_index());
+        AssertDimension (dof_indices.size(), accessor.get_fe().dofs_per_cell);
+
+        Assert (!accessor.has_children(),
+                ExcMessage ("Cell must be active."));
+
+        types::global_dof_index *cache = &accessor.dof_handler->levels[accessor.level()]
+                                         ->cell_dof_indices_cache[accessor.dof_handler->levels[accessor.level()]->cell_cache_offsets[accessor.present_index]];
+        for (unsigned int i=0; i<accessor.get_fe().dofs_per_cell; ++i, ++cache)
+          dof_indices[i] = *cache;
       }
 
       /**
index 58da6c6c1686378b27077d9040b469c75da54944..37f56ed13c2080b5db3efac37034f0c1c04d69d8 100644 (file)
@@ -42,6 +42,10 @@ namespace internal
       struct Implementation;
     }
   }
+  namespace DoFCellAccessor
+  {
+    struct Implementation;
+  }
 }
 
 
@@ -149,6 +153,18 @@ namespace internal
        */
       std::vector<types::global_dof_index> dof_indices;
 
+      /**
+       * The offsets for each cell of the cache that holds all DoF indices.
+       */
+      std::vector<offset_type> cell_cache_offsets;
+
+      /**
+       * Cache for the DoF indices
+       * on cells. The size of this
+       * array equals the sum over all cells of selected_fe[active_fe_index[cell]].dofs_per_cell.
+       */
+      std::vector<types::global_dof_index> cell_dof_indices_cache;
+
     public:
 
       /**
@@ -275,6 +291,7 @@ namespace internal
        */
       template <int, int> friend class dealii::hp::DoFHandler;
       friend struct dealii::internal::hp::DoFHandler::Implementation;
+      friend struct dealii::internal::DoFCellAccessor::Implementation;
     };
 
 
index eeec5b43e1b4b23b3e9f3e51be152c343db75bf5..d6afadf43c605ee894c1b7ee09275d482a0254c3 100644 (file)
@@ -580,8 +580,13 @@ namespace internal
                 = std::vector<DoFLevel::offset_type> (
                     dof_handler.tria->n_raw_lines(level),
                     (DoFLevel::offset_type)(-1));
+              dof_handler.levels[level]->cell_cache_offsets
+                = std::vector<DoFLevel::offset_type> (
+                    dof_handler.tria->n_raw_lines(level),
+                    (DoFLevel::offset_type)(-1));
 
               types::global_dof_index next_free_dof = 0;
+              types::global_dof_index cache_size = 0;
               for (typename DoFHandler<dim,spacedim>::active_cell_iterator
                    cell=dof_handler.begin_active(level);
                    cell!=dof_handler.end_active(level); ++cell)
@@ -589,11 +594,17 @@ namespace internal
                   {
                     dof_handler.levels[level]->dof_offsets[cell->index()] = next_free_dof;
                     next_free_dof += cell->get_fe().dofs_per_line;
+
+                    dof_handler.levels[level]->cell_cache_offsets[cell->index()] = cache_size;
+                    cache_size += cell->get_fe().dofs_per_cell;
                   }
 
               dof_handler.levels[level]->dof_indices
                 = std::vector<types::global_dof_index> (next_free_dof,
                                                         DoFHandler<dim,spacedim>::invalid_dof_index);
+              dof_handler.levels[level]->cell_dof_indices_cache
+                = std::vector<types::global_dof_index> (cache_size,
+                                                        DoFHandler<dim,spacedim>::invalid_dof_index);
             }
 
           // safety check: make sure that
@@ -697,8 +708,13 @@ namespace internal
                 = std::vector<DoFLevel::offset_type> (
                     dof_handler.tria->n_raw_quads(level),
                     (DoFLevel::offset_type)(-1));
+              dof_handler.levels[level]->cell_cache_offsets
+                = std::vector<DoFLevel::offset_type> (
+                    dof_handler.tria->n_raw_quads(level),
+                    (DoFLevel::offset_type)(-1));
 
               types::global_dof_index next_free_dof = 0;
+              types::global_dof_index cache_size = 0;
               for (typename DoFHandler<dim,spacedim>::active_cell_iterator
                    cell=dof_handler.begin_active(level);
                    cell!=dof_handler.end_active(level); ++cell)
@@ -706,11 +722,17 @@ namespace internal
                   {
                     dof_handler.levels[level]->dof_offsets[cell->index()] = next_free_dof;
                     next_free_dof += cell->get_fe().dofs_per_quad;
+
+                    dof_handler.levels[level]->cell_cache_offsets[cell->index()] = cache_size;
+                    cache_size += cell->get_fe().dofs_per_cell;
                   }
 
               dof_handler.levels[level]->dof_indices
                 = std::vector<types::global_dof_index> (next_free_dof,
                                                         DoFHandler<dim,spacedim>::invalid_dof_index);
+              dof_handler.levels[level]->cell_dof_indices_cache
+                = std::vector<types::global_dof_index> (cache_size,
+                                                        DoFHandler<dim,spacedim>::invalid_dof_index);
             }
 
           // safety check: make sure that
@@ -1059,11 +1081,16 @@ namespace internal
           for (unsigned int level=0; level<dof_handler.tria->n_levels(); ++level)
             {
               dof_handler.levels[level]->dof_offsets
-                = std::vector<DoFLevel::offset_type>
-                                                       (dof_handler.tria->n_raw_hexs(level),
-                                                       (DoFLevel::offset_type)(-1));
+                = std::vector<DoFLevel::offset_type> (
+                    dof_handler.tria->n_raw_hexs(level),
+                    (DoFLevel::offset_type)(-1));
+              dof_handler.levels[level]->cell_cache_offsets
+                = std::vector<DoFLevel::offset_type> (
+                    dof_handler.tria->n_raw_hexs(level),
+                    (DoFLevel::offset_type)(-1));
 
               types::global_dof_index next_free_dof = 0;
+              types::global_dof_index cache_size = 0;
               for (typename DoFHandler<dim,spacedim>::active_cell_iterator
                    cell=dof_handler.begin_active(level);
                    cell!=dof_handler.end_active(level); ++cell)
@@ -1071,11 +1098,17 @@ namespace internal
                   {
                     dof_handler.levels[level]->dof_offsets[cell->index()] = next_free_dof;
                     next_free_dof += cell->get_fe().dofs_per_hex;
+
+                    dof_handler.levels[level]->cell_cache_offsets[cell->index()] = cache_size;
+                    cache_size += cell->get_fe().dofs_per_cell;
                   }
 
               dof_handler.levels[level]->dof_indices
                 = std::vector<types::global_dof_index> (next_free_dof,
                                                         DoFHandler<dim,spacedim>::invalid_dof_index);
+              dof_handler.levels[level]->cell_dof_indices_cache
+                = std::vector<types::global_dof_index> (cache_size,
+                                                        DoFHandler<dim,spacedim>::invalid_dof_index);
             }
 
           // safety check: make sure that
@@ -2689,6 +2722,11 @@ namespace hp
       = std::vector<IndexSet> (1,
                                number_cache.locally_owned_dofs);
 
+    // update the cache used for cell dof indices and compress the data on the levels
+    for (active_cell_iterator cell = begin_active();
+         cell != end(); ++cell)
+      cell->update_cell_dof_indices_cache ();
+
     for (unsigned int level=0; level<levels.size(); ++level)
       levels[level]->compress_data (*finite_elements);
 
@@ -2736,6 +2774,11 @@ namespace hp
     // do the renumbering
     renumber_dofs_internal (new_numbers, dealii::internal::int2type<dim>());
 
+    // update the cache used for cell dof indices
+    for (active_cell_iterator cell = begin_active();
+         cell != end(); ++cell)
+      cell->update_cell_dof_indices_cache ();
+
     // now re-compress the dof indices
     for (unsigned int level=0; level<levels.size(); ++level)
       levels[level]->compress_data (*finite_elements);

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.