]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove selected_fe from DoFHandler
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Mon, 21 Aug 2017 18:49:01 +0000 (20:49 +0200)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Tue, 22 Aug 2017 12:04:50 +0000 (14:04 +0200)
include/deal.II/dofs/dof_accessor.templates.h
include/deal.II/dofs/dof_handler.h
source/dofs/dof_handler.cc
source/dofs/dof_handler_policy.cc

index b1a118de8d32e5218c0222b274cca8f271b5beab..e8c442ccc9f9324776c295124ac8006356f43a83 100644 (file)
@@ -1140,15 +1140,12 @@ namespace internal
         (void)fe_index;
         Assert ((fe_index == dealii::DoFHandler<dim,spacedim>::default_fe_index),
                 ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-        Assert (dof_handler.selected_fe != nullptr,
-                ExcMessage ("No finite element collection is associated with "
-                            "this DoFHandler"));
-        Assert (local_index < dof_handler.selected_fe->dofs_per_vertex,
+        Assert (local_index < dof_handler.get_finite_element().dofs_per_vertex,
                 ExcIndexRange(local_index, 0,
-                              dof_handler.selected_fe->dofs_per_vertex));
+                              dof_handler.get_finite_element().dofs_per_vertex));
 
         dof_handler.vertex_dofs[vertex_index *
-                                dof_handler.selected_fe->dofs_per_vertex
+                                dof_handler.get_finite_element().dofs_per_vertex
                                 + local_index]
           = global_index;
       }
@@ -1229,16 +1226,13 @@ namespace internal
         (void)fe_index;
         Assert ((fe_index == dealii::DoFHandler<dim,spacedim>::default_fe_index),
                 ExcMessage ("Only the default FE index is allowed for non-hp DoFHandler objects"));
-        Assert (dof_handler.selected_fe != nullptr,
-                ExcMessage ("No finite element collection is associated with "
-                            "this DoFHandler"));
-        Assert (local_index < dof_handler.selected_fe->dofs_per_vertex,
+        Assert (local_index < dof_handler.get_finite_element().dofs_per_vertex,
                 ExcIndexRange(local_index, 0,
-                              dof_handler.selected_fe->dofs_per_vertex));
+                              dof_handler.get_finite_element().dofs_per_vertex));
 
         return
           dof_handler.vertex_dofs[vertex_index *
-                                  dof_handler.selected_fe->dofs_per_vertex
+                                  dof_handler.get_finite_element().dofs_per_vertex
                                   + local_index];
       }
 
index c3b2e206cea2247a55cf4c26d55d19257e1c42c1..5314a42c56843d3ed60fda31ea97597cfde4f4af 100644 (file)
@@ -354,11 +354,7 @@ public:
    * The purpose of this function is first discussed in the introduction to
    * the step-2 tutorial program.
    *
-   * @note A pointer of the finite element given as argument is stored.
-   * Therefore, the lifetime of the finite element object shall be longer than
-   * that of this object. If you don't want this behavior, you may want to
-   * call the @p clear member function which also releases the lock of this
-   * object to the finite element.
+   * @note A copy of the finite element given as argument is stored.
    */
   virtual void distribute_dofs (const FiniteElement<dim,spacedim> &fe);
 
@@ -402,9 +398,7 @@ public:
   void initialize_local_block_info();
 
   /**
-   * Clear all data of this object and especially delete the lock this object
-   * has to the finite element used the last time when @p distribute_dofs was
-   * called.
+   * Clear all data of this object.
    */
   virtual void clear ();
 
@@ -942,16 +936,6 @@ private:
   SmartPointer<const Triangulation<dim,spacedim>,DoFHandler<dim,spacedim> >
   tria;
 
-  /**
-   * Store a pointer to the finite element given latest for the distribution
-   * of dofs. In order to avoid destruction of the object before the lifetime
-   * of the DoF handler, we subscribe to the finite element object. To unlock
-   * the FE before the end of the lifetime of this DoF handler, use the
-   * <tt>clear()</tt> function (this clears all data of this object as well,
-   * though).
-   */
-  SmartPointer<const FiniteElement<dim,spacedim>,DoFHandler<dim,spacedim> >
-  selected_fe;
 
   /**
    * Store a pointer to a hp::FECollection object containing the (one)
@@ -1240,9 +1224,7 @@ inline
 const FiniteElement<dim,spacedim> &
 DoFHandler<dim,spacedim>::get_fe () const
 {
-  Assert(selected_fe!=nullptr,
-         ExcMessage("You are trying to access the DoFHandler's FiniteElement object before it has been initialized."));
-  return *selected_fe;
+  return this->get_finite_element();
 }
 
 
@@ -1255,9 +1237,7 @@ DoFHandler<dim,spacedim>::get_finite_element
 {
   (void) index;
   Assert(index == 0, ExcMessage("There is only one FiniteElement stored. The index must be zero!"));
-  Assert(selected_fe!=nullptr,
-         ExcMessage("You are trying to access the DoFHandler's FiniteElement object before it has been initialized."));
-  return *selected_fe;
+  return get_fe_collection()[0];
 }
 
 
@@ -1353,7 +1333,7 @@ void DoFHandler<dim,spacedim>::save (Archive &ar,
   // loading that this number is indeed correct; same with something that
   // identifies the FE and the policy
   unsigned int n_cells = tria->n_cells();
-  std::string  fe_name = selected_fe->get_name();
+  std::string  fe_name = this->get_finite_element().get_name();
   std::string  policy_name = internal::policy_to_string(*policy);
 
   ar &n_cells &fe_name &policy_name;
@@ -1403,7 +1383,7 @@ void DoFHandler<dim,spacedim>::load (Archive &ar,
   AssertThrow (n_cells == tria->n_cells(),
                ExcMessage ("The object being loaded into does not match the triangulation "
                            "that has been stored previously."));
-  AssertThrow (fe_name == selected_fe->get_name(),
+  AssertThrow (fe_name == this->get_finite_element().get_name(),
                ExcMessage ("The finite element associated with this DoFHandler does not match "
                            "the one that was associated with the DoFHandler previously stored."));
   AssertThrow (policy_name == internal::policy_to_string(*policy),
index 8210a0db24b43b0e65e4caf4eb9a065f9743c3cb..c129e3719393aaa9259d3cd12b9e714835296b30 100644 (file)
@@ -107,8 +107,8 @@ namespace internal
       unsigned int
       max_couplings_between_dofs (const DoFHandler<1,spacedim> &dof_handler)
       {
-        return std::min(static_cast<types::global_dof_index>(3*dof_handler.selected_fe->dofs_per_vertex +
-                                                             2*dof_handler.selected_fe->dofs_per_line),
+        return std::min(static_cast<types::global_dof_index>(3*dof_handler.get_finite_element().dofs_per_vertex +
+                                                             2*dof_handler.get_finite_element().dofs_per_line),
                         dof_handler.n_dofs());
       }
 
@@ -145,29 +145,29 @@ namespace internal
         switch (dof_handler.tria->max_adjacent_cells())
           {
           case 4:
-            max_couplings=19*dof_handler.selected_fe->dofs_per_vertex +
-                          28*dof_handler.selected_fe->dofs_per_line +
-                          8*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=19*dof_handler.get_finite_element().dofs_per_vertex +
+                          28*dof_handler.get_finite_element().dofs_per_line +
+                          8*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 5:
-            max_couplings=21*dof_handler.selected_fe->dofs_per_vertex +
-                          31*dof_handler.selected_fe->dofs_per_line +
-                          9*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=21*dof_handler.get_finite_element().dofs_per_vertex +
+                          31*dof_handler.get_finite_element().dofs_per_line +
+                          9*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 6:
-            max_couplings=28*dof_handler.selected_fe->dofs_per_vertex +
-                          42*dof_handler.selected_fe->dofs_per_line +
-                          12*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=28*dof_handler.get_finite_element().dofs_per_vertex +
+                          42*dof_handler.get_finite_element().dofs_per_line +
+                          12*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 7:
-            max_couplings=30*dof_handler.selected_fe->dofs_per_vertex +
-                          45*dof_handler.selected_fe->dofs_per_line +
-                          13*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=30*dof_handler.get_finite_element().dofs_per_vertex +
+                          45*dof_handler.get_finite_element().dofs_per_line +
+                          13*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 8:
-            max_couplings=37*dof_handler.selected_fe->dofs_per_vertex +
-                          56*dof_handler.selected_fe->dofs_per_line +
-                          16*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=37*dof_handler.get_finite_element().dofs_per_vertex +
+                          56*dof_handler.get_finite_element().dofs_per_line +
+                          16*dof_handler.get_finite_element().dofs_per_quad;
             break;
 
           // the following numbers are not based on actual counting but by
@@ -175,44 +175,44 @@ namespace internal
           // example, for dofs_per_vertex, the sequence above is 19, 21, 28,
           // 30, 37, and is continued as follows):
           case 9:
-            max_couplings=39*dof_handler.selected_fe->dofs_per_vertex +
-                          59*dof_handler.selected_fe->dofs_per_line +
-                          17*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=39*dof_handler.get_finite_element().dofs_per_vertex +
+                          59*dof_handler.get_finite_element().dofs_per_line +
+                          17*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 10:
-            max_couplings=46*dof_handler.selected_fe->dofs_per_vertex +
-                          70*dof_handler.selected_fe->dofs_per_line +
-                          20*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=46*dof_handler.get_finite_element().dofs_per_vertex +
+                          70*dof_handler.get_finite_element().dofs_per_line +
+                          20*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 11:
-            max_couplings=48*dof_handler.selected_fe->dofs_per_vertex +
-                          73*dof_handler.selected_fe->dofs_per_line +
-                          21*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=48*dof_handler.get_finite_element().dofs_per_vertex +
+                          73*dof_handler.get_finite_element().dofs_per_line +
+                          21*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 12:
-            max_couplings=55*dof_handler.selected_fe->dofs_per_vertex +
-                          84*dof_handler.selected_fe->dofs_per_line +
-                          24*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=55*dof_handler.get_finite_element().dofs_per_vertex +
+                          84*dof_handler.get_finite_element().dofs_per_line +
+                          24*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 13:
-            max_couplings=57*dof_handler.selected_fe->dofs_per_vertex +
-                          87*dof_handler.selected_fe->dofs_per_line +
-                          25*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=57*dof_handler.get_finite_element().dofs_per_vertex +
+                          87*dof_handler.get_finite_element().dofs_per_line +
+                          25*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 14:
-            max_couplings=63*dof_handler.selected_fe->dofs_per_vertex +
-                          98*dof_handler.selected_fe->dofs_per_line +
-                          28*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=63*dof_handler.get_finite_element().dofs_per_vertex +
+                          98*dof_handler.get_finite_element().dofs_per_line +
+                          28*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 15:
-            max_couplings=65*dof_handler.selected_fe->dofs_per_vertex +
-                          103*dof_handler.selected_fe->dofs_per_line +
-                          29*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=65*dof_handler.get_finite_element().dofs_per_vertex +
+                          103*dof_handler.get_finite_element().dofs_per_line +
+                          29*dof_handler.get_finite_element().dofs_per_quad;
             break;
           case 16:
-            max_couplings=72*dof_handler.selected_fe->dofs_per_vertex +
-                          114*dof_handler.selected_fe->dofs_per_line +
-                          32*dof_handler.selected_fe->dofs_per_quad;
+            max_couplings=72*dof_handler.get_finite_element().dofs_per_vertex +
+                          114*dof_handler.get_finite_element().dofs_per_line +
+                          32*dof_handler.get_finite_element().dofs_per_quad;
             break;
 
           default:
@@ -246,10 +246,10 @@ namespace internal
 
         types::global_dof_index max_couplings;
         if (max_adjacent_cells <= 8)
-          max_couplings=7*7*7*dof_handler.selected_fe->dofs_per_vertex +
-                        7*6*7*3*dof_handler.selected_fe->dofs_per_line +
-                        9*4*7*3*dof_handler.selected_fe->dofs_per_quad +
-                        27*dof_handler.selected_fe->dofs_per_hex;
+          max_couplings=7*7*7*dof_handler.get_finite_element().dofs_per_vertex +
+                        7*6*7*3*dof_handler.get_finite_element().dofs_per_line +
+                        9*4*7*3*dof_handler.get_finite_element().dofs_per_quad +
+                        27*dof_handler.get_finite_element().dofs_per_hex;
         else
           {
             Assert (false, ExcNotImplemented());
@@ -275,7 +275,7 @@ namespace internal
       {
         dof_handler.vertex_dofs
         .resize(dof_handler.tria->n_vertices() *
-                dof_handler.selected_fe->dofs_per_vertex,
+                dof_handler.get_finite_element().dofs_per_vertex,
                 numbers::invalid_dof_index);
 
         for (unsigned int i=0; i<dof_handler.tria->n_levels(); ++i)
@@ -284,12 +284,12 @@ namespace internal
 
             dof_handler.levels.back()->dof_object.dofs
             .resize (dof_handler.tria->n_raw_cells(i) *
-                     dof_handler.selected_fe->dofs_per_line,
+                     dof_handler.get_finite_element().dofs_per_line,
                      numbers::invalid_dof_index);
 
             dof_handler.levels.back()->cell_dof_indices_cache
             .resize (dof_handler.tria->n_raw_cells(i) *
-                     dof_handler.selected_fe->dofs_per_cell,
+                     dof_handler.get_finite_element().dofs_per_cell,
                      numbers::invalid_dof_index);
           }
       }
@@ -301,7 +301,7 @@ namespace internal
       {
         dof_handler.vertex_dofs
         .resize(dof_handler.tria->n_vertices() *
-                dof_handler.selected_fe->dofs_per_vertex,
+                dof_handler.get_finite_element().dofs_per_vertex,
                 numbers::invalid_dof_index);
 
         for (unsigned int i=0; i<dof_handler.tria->n_levels(); ++i)
@@ -310,12 +310,12 @@ namespace internal
 
             dof_handler.levels.back()->dof_object.dofs
             .resize (dof_handler.tria->n_raw_cells(i) *
-                     dof_handler.selected_fe->dofs_per_quad,
+                     dof_handler.get_finite_element().dofs_per_quad,
                      numbers::invalid_dof_index);
 
             dof_handler.levels.back()->cell_dof_indices_cache
             .resize (dof_handler.tria->n_raw_cells(i) *
-                     dof_handler.selected_fe->dofs_per_cell,
+                     dof_handler.get_finite_element().dofs_per_cell,
                      numbers::invalid_dof_index);
           }
 
@@ -325,7 +325,7 @@ namespace internal
           {
             dof_handler.faces->lines.dofs
             .resize (dof_handler.tria->n_raw_lines() *
-                     dof_handler.selected_fe->dofs_per_line,
+                     dof_handler.get_finite_element().dofs_per_line,
                      numbers::invalid_dof_index);
           }
       }
@@ -337,7 +337,7 @@ namespace internal
       {
         dof_handler.vertex_dofs
         .resize(dof_handler.tria->n_vertices() *
-                dof_handler.selected_fe->dofs_per_vertex,
+                dof_handler.get_finite_element().dofs_per_vertex,
                 numbers::invalid_dof_index);
 
         for (unsigned int i=0; i<dof_handler.tria->n_levels(); ++i)
@@ -346,12 +346,12 @@ namespace internal
 
             dof_handler.levels.back()->dof_object.dofs
             .resize (dof_handler.tria->n_raw_cells(i) *
-                     dof_handler.selected_fe->dofs_per_hex,
+                     dof_handler.get_finite_element().dofs_per_hex,
                      numbers::invalid_dof_index);
 
             dof_handler.levels.back()->cell_dof_indices_cache
             .resize (dof_handler.tria->n_raw_cells(i) *
-                     dof_handler.selected_fe->dofs_per_cell,
+                     dof_handler.get_finite_element().dofs_per_cell,
                      numbers::invalid_dof_index);
           }
         dof_handler.faces.reset (new internal::DoFHandler::DoFFaces<3>);
@@ -361,11 +361,11 @@ namespace internal
           {
             dof_handler.faces->lines.dofs
             .resize (dof_handler.tria->n_raw_lines() *
-                     dof_handler.selected_fe->dofs_per_line,
+                     dof_handler.get_finite_element().dofs_per_line,
                      numbers::invalid_dof_index);
             dof_handler.faces->quads.dofs
             .resize (dof_handler.tria->n_raw_quads() *
-                     dof_handler.selected_fe->dofs_per_quad,
+                     dof_handler.get_finite_element().dofs_per_quad,
                      numbers::invalid_dof_index);
           }
       }
@@ -653,7 +653,6 @@ template <int dim, int spacedim>
 DoFHandler<dim,spacedim>::DoFHandler (const Triangulation<dim,spacedim> &tria)
   :
   tria(&tria, typeid(*this).name()),
-  selected_fe(nullptr, typeid(*this).name()),
   fe_collection(nullptr),
   faces(nullptr),
   mg_faces (nullptr)
@@ -676,7 +675,6 @@ template <int dim, int spacedim>
 DoFHandler<dim,spacedim>::DoFHandler ()
   :
   tria(nullptr, typeid(*this).name()),
-  selected_fe(nullptr, typeid(*this).name()),
   fe_collection(nullptr)
 {}
 
@@ -965,7 +963,6 @@ std::size_t
 DoFHandler<dim,spacedim>::memory_consumption () const
 {
   std::size_t mem = (MemoryConsumption::memory_consumption (tria) +
-                     MemoryConsumption::memory_consumption (selected_fe) +
                      MemoryConsumption::memory_consumption (fe_collection) +
                      MemoryConsumption::memory_consumption (block_info_object) +
                      MemoryConsumption::memory_consumption (levels) +
@@ -994,7 +991,6 @@ DoFHandler<dim,spacedim>::memory_consumption () const
 template <int dim, int spacedim>
 void DoFHandler<dim,spacedim>::distribute_dofs (const FiniteElement<dim,spacedim> &ff)
 {
-  selected_fe = &ff;
   fe_collection = std_cxx14::make_unique<hp::FECollection<dim, spacedim>>(ff);
 
   // delete all levels and set them
@@ -1029,15 +1025,10 @@ void DoFHandler<dim,spacedim>::distribute_dofs (const FiniteElement<dim,spacedim
 
 
 template <int dim, int spacedim>
-void DoFHandler<dim, spacedim>::distribute_mg_dofs (const FiniteElement<dim, spacedim> &fe)
+void DoFHandler<dim, spacedim>::distribute_mg_dofs (const FiniteElement<dim, spacedim> &)
 {
-  (void)fe;
   Assert(levels.size()>0, ExcMessage("Distribute active DoFs using distribute_dofs() before calling distribute_mg_dofs()."));
 
-  const FiniteElement<dim, spacedim> *old_fe = selected_fe;
-  (void)old_fe;
-  Assert(old_fe == &fe, ExcMessage("You are required to use the same FE for level and active DoFs!") );
-
   Assert(((tria->get_mesh_smoothing() & Triangulation<dim, spacedim>::limit_level_difference_at_vertices)
           != Triangulation<dim, spacedim>::none),
          ExcMessage("The mesh smoothing requirement 'limit_level_difference_at_vertices' has to be set for using multigrid!"));
@@ -1091,7 +1082,6 @@ template <int dim, int spacedim>
 void DoFHandler<dim,spacedim>::clear ()
 {
   // release lock to old fe
-  selected_fe = nullptr;
   fe_collection.release();
 
   // release memory
index 4f3dcd98c254dac8334cb7a83660041e519d070a..7e4f82cf358477ae3133b4c4800c83a6b2d06b01 100644 (file)
@@ -1452,7 +1452,7 @@ namespace internal
               // really is unused
               Assert (dof_handler.get_triangulation()
                       .vertex_used((i-dof_handler.vertex_dofs.begin()) /
-                                   dof_handler.selected_fe->dofs_per_vertex)
+                                   dof_handler.get_finite_element().dofs_per_vertex)
                       == false,
                       ExcInternalError ());
         }

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.