]> https://gitweb.dealii.org/ - dealii.git/commitdiff
std::unique_ptr<const CellIteratorBase> -> CellIteratorContainer 12946/head
authorPeter Munch <peterrmuench@gmail.com>
Sat, 13 Nov 2021 21:05:19 +0000 (22:05 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Sat, 13 Nov 2021 22:41:21 +0000 (23:41 +0100)
include/deal.II/fe/fe_values.h
source/fe/fe_values.cc

index 3b6977c482de102056b30bcf791ecd34d929422d..07c9cbfa9d6fed465bc9749a606b907b97eaaa0b 100644 (file)
@@ -3759,6 +3759,14 @@ public:
     << "> flag to be set, but it was apparently not specified "
     << "upon construction.");
 
+  /**
+   * FEValues::reinit() has not been called for any cell.
+   *
+   * @ingroup Exceptions
+   */
+  DeclExceptionMsg(ExcNotReinited,
+                   "FEValues object is not reinit'ed to any cell");
+
   /**
    * Mismatch between the FEValues FiniteElement and
    * cell->get_dof_handler().get_fe()
@@ -3795,77 +3803,48 @@ public:
     "only works for those. See FiniteElement::is_primitive() for more information.");
 
 protected:
-  /**
-   * Reset a unique_ptr. If we can, do not de-allocate the previously
-   * held memory but re-use it for the next item to avoid the repeated
-   * memory allocation. We do this because FEValues objects are heavily
-   * used in multithreaded contexts where memory allocations are evil.
-   */
-  template <typename Type, typename Pointer, typename Iterator>
-  static void
-  reset_pointer_in_place_if_possible(std::unique_ptr<Pointer> &present_cell,
-                                     const Iterator &          new_cell);
-
   /**
    * Objects of the FEValues class need to store an iterator
    * to the present cell in order to be able to extract the values of the
    * degrees of freedom on this cell in the get_function_values() and assorted
-   * functions. On the other hand, this class should also work for different
-   * iterators, as long as they have the same interface to extract the DoF
-   * values (i.e., for example, they need to have a @p
-   * get_interpolated_dof_values function).
-   *
-   * This calls for a common base class of iterator classes, and making the
-   * functions we need here @p virtual. On the other hand, this is the only
-   * place in the library where we need this, and introducing a base class of
-   * iterators and making a function virtual penalizes <em>all</em> users of
-   * the iterators, which are basically intended as very fast accessor
-   * functions. So we do not want to do this. Rather, what we do here is
-   * making the functions we need virtual only for use with <em>this
-   * class</em>. The idea is the following: have a common base class which
-   * declares some pure virtual functions, and for each possible iterator
-   * type, we have a derived class which stores the iterator to the cell and
-   * implements these functions. Since the iterator classes have the same
-   * interface, we can make the derived classes a template, templatized on the
-   * iterator type.
-   *
-   * This way, the use of virtual functions is restricted to only this class,
-   * and other users of iterators do not have to bear the negative effects.
-   *
-   * @note This class is an example of the
-   * <a href="https://www.artima.com/cppsource/type_erasure.html">type
-   * erasure</a> design pattern.
-   */
-  class CellIteratorBase
+   * functions.
+   */
+  class CellIteratorContainer
   {
   public:
     DeclExceptionMsg(
       ExcNeedsDoFHandler,
-      "You have previously called the FEValues::reinit function with a "
+      "You have previously called the FEValues::reinit() function with a "
       "cell iterator of type Triangulation<dim,spacedim>::cell_iterator. However, "
       "when you do this, you cannot call some functions in the FEValues "
       "class, such as the get_function_values/gradients/hessians/third_derivatives "
       "functions. If you need these functions, then you need to call "
-      "FEValues::reinit with an iterator type that allows to extract "
+      "FEValues::reinit() with an iterator type that allows to extract "
       "degrees of freedom, such as DoFHandler<dim,spacedim>::cell_iterator.");
 
+    /**
+     * Constructor.
+     */
+    CellIteratorContainer();
+
     /**
      * Constructor.
      */
     template <bool lda>
-    CellIteratorBase(
+    CellIteratorContainer(
       const TriaIterator<DoFCellAccessor<dim, spacedim, lda>> &cell);
 
     /**
      * Constructor.
      */
-    CellIteratorBase(
+    CellIteratorContainer(
       const typename Triangulation<dim, spacedim>::cell_iterator &cell);
 
     /**
-     * Destructor.
+     * Indicate whether FEValues::reinit() was called.
      */
-    virtual ~CellIteratorBase() = default;
+    bool
+    is_initialized() const;
 
     /**
      * Conversion operator to an iterator for triangulations. This
@@ -3901,9 +3880,10 @@ protected:
                                 Vector<IndexSet::value_type> &out) const;
 
   private:
-    const typename Triangulation<dim, spacedim>::cell_iterator cell;
-    const DoFHandler<dim, spacedim> *                          dof_handler;
-    const bool                                                 level_dof_access;
+    bool                                                 initialized;
+    typename Triangulation<dim, spacedim>::cell_iterator cell;
+    const DoFHandler<dim, spacedim> *                    dof_handler;
+    bool                                                 level_dof_access;
   };
 
   /**
@@ -3911,7 +3891,7 @@ protected:
    * is necessary for the <tt>get_function_*</tt> functions as well as the
    * functions of same name in the extractor classes.
    */
-  std::unique_ptr<const CellIteratorBase> present_cell;
+  CellIteratorContainer present_cell;
 
   /**
    * A signal connection we use to ensure we get informed whenever the
@@ -5582,8 +5562,7 @@ FEValuesBase<dim, spacedim>::shape_value(const unsigned int i,
   Assert(this->update_flags & update_values,
          ExcAccessToUninitializedField("update_values"));
   Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   // if the entire FE is primitive,
   // then we can take a short-cut:
   if (fe->is_primitive())
@@ -5619,8 +5598,7 @@ FEValuesBase<dim, spacedim>::shape_value_component(
   Assert(this->update_flags & update_values,
          ExcAccessToUninitializedField("update_values"));
   AssertIndexRange(component, fe->n_components());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   // check whether the shape function
   // is non-zero at all within
@@ -5648,8 +5626,7 @@ FEValuesBase<dim, spacedim>::shape_grad(const unsigned int i,
   Assert(this->update_flags & update_gradients,
          ExcAccessToUninitializedField("update_gradients"));
   Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   // if the entire FE is primitive,
   // then we can take a short-cut:
   if (fe->is_primitive())
@@ -5685,8 +5662,7 @@ FEValuesBase<dim, spacedim>::shape_grad_component(
   Assert(this->update_flags & update_gradients,
          ExcAccessToUninitializedField("update_gradients"));
   AssertIndexRange(component, fe->n_components());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   // check whether the shape function
   // is non-zero at all within
   // this component:
@@ -5713,8 +5689,7 @@ FEValuesBase<dim, spacedim>::shape_hessian(const unsigned int i,
   Assert(this->update_flags & update_hessians,
          ExcAccessToUninitializedField("update_hessians"));
   Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   // if the entire FE is primitive,
   // then we can take a short-cut:
   if (fe->is_primitive())
@@ -5750,8 +5725,7 @@ FEValuesBase<dim, spacedim>::shape_hessian_component(
   Assert(this->update_flags & update_hessians,
          ExcAccessToUninitializedField("update_hessians"));
   AssertIndexRange(component, fe->n_components());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   // check whether the shape function
   // is non-zero at all within
   // this component:
@@ -5778,8 +5752,7 @@ FEValuesBase<dim, spacedim>::shape_3rd_derivative(const unsigned int i,
   Assert(this->update_flags & update_3rd_derivatives,
          ExcAccessToUninitializedField("update_3rd_derivatives"));
   Assert(fe->is_primitive(i), ExcShapeFunctionNotPrimitive(i));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   // if the entire FE is primitive,
   // then we can take a short-cut:
   if (fe->is_primitive())
@@ -5815,8 +5788,7 @@ FEValuesBase<dim, spacedim>::shape_3rd_derivative_component(
   Assert(this->update_flags & update_3rd_derivatives,
          ExcAccessToUninitializedField("update_3rd_derivatives"));
   AssertIndexRange(component, fe->n_components());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   // check whether the shape function
   // is non-zero at all within
   // this component:
@@ -5867,8 +5839,7 @@ FEValuesBase<dim, spacedim>::get_quadrature_points() const
 {
   Assert(this->update_flags & update_quadrature_points,
          ExcAccessToUninitializedField("update_quadrature_points"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.quadrature_points;
 }
 
@@ -5880,8 +5851,7 @@ FEValuesBase<dim, spacedim>::get_JxW_values() const
 {
   Assert(this->update_flags & update_JxW_values,
          ExcAccessToUninitializedField("update_JxW_values"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.JxW_values;
 }
 
@@ -5893,8 +5863,7 @@ FEValuesBase<dim, spacedim>::get_jacobians() const
 {
   Assert(this->update_flags & update_jacobians,
          ExcAccessToUninitializedField("update_jacobians"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobians;
 }
 
@@ -5906,8 +5875,7 @@ FEValuesBase<dim, spacedim>::get_jacobian_grads() const
 {
   Assert(this->update_flags & update_jacobian_grads,
          ExcAccessToUninitializedField("update_jacobians_grads"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_grads;
 }
 
@@ -5920,8 +5888,7 @@ FEValuesBase<dim, spacedim>::jacobian_pushed_forward_grad(
 {
   Assert(this->update_flags & update_jacobian_pushed_forward_grads,
          ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_pushed_forward_grads[i];
 }
 
@@ -5933,8 +5900,7 @@ FEValuesBase<dim, spacedim>::get_jacobian_pushed_forward_grads() const
 {
   Assert(this->update_flags & update_jacobian_pushed_forward_grads,
          ExcAccessToUninitializedField("update_jacobian_pushed_forward_grads"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_pushed_forward_grads;
 }
 
@@ -5946,8 +5912,7 @@ FEValuesBase<dim, spacedim>::jacobian_2nd_derivative(const unsigned int i) const
 {
   Assert(this->update_flags & update_jacobian_2nd_derivatives,
          ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_2nd_derivatives[i];
 }
 
@@ -5959,8 +5924,7 @@ FEValuesBase<dim, spacedim>::get_jacobian_2nd_derivatives() const
 {
   Assert(this->update_flags & update_jacobian_2nd_derivatives,
          ExcAccessToUninitializedField("update_jacobian_2nd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_2nd_derivatives;
 }
 
@@ -5974,8 +5938,7 @@ FEValuesBase<dim, spacedim>::jacobian_pushed_forward_2nd_derivative(
   Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
          ExcAccessToUninitializedField(
            "update_jacobian_pushed_forward_2nd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_pushed_forward_2nd_derivatives[i];
 }
 
@@ -5988,8 +5951,7 @@ FEValuesBase<dim, spacedim>::get_jacobian_pushed_forward_2nd_derivatives() const
   Assert(this->update_flags & update_jacobian_pushed_forward_2nd_derivatives,
          ExcAccessToUninitializedField(
            "update_jacobian_pushed_forward_2nd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_pushed_forward_2nd_derivatives;
 }
 
@@ -6001,8 +5963,7 @@ FEValuesBase<dim, spacedim>::jacobian_3rd_derivative(const unsigned int i) const
 {
   Assert(this->update_flags & update_jacobian_3rd_derivatives,
          ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_3rd_derivatives[i];
 }
 
@@ -6014,8 +5975,7 @@ FEValuesBase<dim, spacedim>::get_jacobian_3rd_derivatives() const
 {
   Assert(this->update_flags & update_jacobian_3rd_derivatives,
          ExcAccessToUninitializedField("update_jacobian_3rd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_3rd_derivatives;
 }
 
@@ -6029,8 +5989,7 @@ FEValuesBase<dim, spacedim>::jacobian_pushed_forward_3rd_derivative(
   Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
          ExcAccessToUninitializedField(
            "update_jacobian_pushed_forward_3rd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_pushed_forward_3rd_derivatives[i];
 }
 
@@ -6043,8 +6002,7 @@ FEValuesBase<dim, spacedim>::get_jacobian_pushed_forward_3rd_derivatives() const
   Assert(this->update_flags & update_jacobian_pushed_forward_3rd_derivatives,
          ExcAccessToUninitializedField(
            "update_jacobian_pushed_forward_3rd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.jacobian_pushed_forward_3rd_derivatives;
 }
 
@@ -6056,8 +6014,7 @@ FEValuesBase<dim, spacedim>::get_inverse_jacobians() const
 {
   Assert(this->update_flags & update_inverse_jacobians,
          ExcAccessToUninitializedField("update_inverse_jacobians"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   return this->mapping_output.inverse_jacobians;
 }
 
@@ -6112,8 +6069,7 @@ FEValuesBase<dim, spacedim>::quadrature_point(const unsigned int i) const
   Assert(this->update_flags & update_quadrature_points,
          ExcAccessToUninitializedField("update_quadrature_points"));
   AssertIndexRange(i, this->mapping_output.quadrature_points.size());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   return this->mapping_output.quadrature_points[i];
 }
@@ -6127,8 +6083,7 @@ FEValuesBase<dim, spacedim>::JxW(const unsigned int i) const
   Assert(this->update_flags & update_JxW_values,
          ExcAccessToUninitializedField("update_JxW_values"));
   AssertIndexRange(i, this->mapping_output.JxW_values.size());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   return this->mapping_output.JxW_values[i];
 }
@@ -6142,8 +6097,7 @@ FEValuesBase<dim, spacedim>::jacobian(const unsigned int i) const
   Assert(this->update_flags & update_jacobians,
          ExcAccessToUninitializedField("update_jacobians"));
   AssertIndexRange(i, this->mapping_output.jacobians.size());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   return this->mapping_output.jacobians[i];
 }
@@ -6157,8 +6111,7 @@ FEValuesBase<dim, spacedim>::jacobian_grad(const unsigned int i) const
   Assert(this->update_flags & update_jacobian_grads,
          ExcAccessToUninitializedField("update_jacobians_grads"));
   AssertIndexRange(i, this->mapping_output.jacobian_grads.size());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   return this->mapping_output.jacobian_grads[i];
 }
@@ -6172,8 +6125,7 @@ FEValuesBase<dim, spacedim>::inverse_jacobian(const unsigned int i) const
   Assert(this->update_flags & update_inverse_jacobians,
          ExcAccessToUninitializedField("update_inverse_jacobians"));
   AssertIndexRange(i, this->mapping_output.inverse_jacobians.size());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   return this->mapping_output.inverse_jacobians[i];
 }
@@ -6188,8 +6140,7 @@ FEValuesBase<dim, spacedim>::normal_vector(const unsigned int i) const
          (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
            "update_normal_vectors")));
   AssertIndexRange(i, this->mapping_output.normal_vectors.size());
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   return this->mapping_output.normal_vectors[i];
 }
index 4a6b14bc7c48e79b9f90e5d9b22ac32b1d9786e7..af32d666f74396e98dd00bf768479b1d6d65eb62 100644 (file)
@@ -1552,17 +1552,17 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell and call internal worker
     // function
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_values<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_values,
@@ -1583,8 +1583,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_values<dim, spacedim>(
@@ -1607,16 +1607,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_derivatives<1, dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -1637,8 +1637,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_derivatives<1, dim, spacedim>(
@@ -1661,16 +1661,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_hessians,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_hessians")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_derivatives<2, dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_hessians,
@@ -1691,8 +1691,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_hessians,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_hessians")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_derivatives<2, dim, spacedim>(
@@ -1715,16 +1715,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_hessians,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_hessians")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_laplacians<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_hessians,
@@ -1745,8 +1745,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_hessians,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_hessians")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_laplacians<dim, spacedim>(
@@ -1770,16 +1770,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_3rd_derivatives,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_3rd_derivatives")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_derivatives<3, dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_3rd_derivatives,
@@ -1801,8 +1801,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_3rd_derivatives,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_3rd_derivatives")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_derivatives<3, dim, spacedim>(
@@ -1825,16 +1825,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_values<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_values,
@@ -1855,8 +1855,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_values<dim, spacedim>(
@@ -1879,16 +1879,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_derivatives<1, dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -1909,8 +1909,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_derivatives<1, dim, spacedim>(
@@ -1934,16 +1934,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_symmetric_gradients<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -1965,8 +1965,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_symmetric_gradients<dim, spacedim>(
@@ -1989,17 +1989,17 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs
     // on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_divergences<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -2020,8 +2020,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_divergences<dim, spacedim>(
@@ -2044,16 +2044,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
+    Assert(fe_values->present_cell.is_initialized(),
            ExcMessage("FEValues object is not reinited to any cell"));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_curls<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -2074,7 +2074,7 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
+    Assert(fe_values->present_cell.is_initialized(),
            ExcMessage("FEValues object is not reinited to any cell"));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
@@ -2098,16 +2098,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_hessians,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_hessians")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_derivatives<2, dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_hessians,
@@ -2128,8 +2128,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_hessians,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_hessians")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_derivatives<2, dim, spacedim>(
@@ -2155,18 +2155,18 @@ namespace FEValuesViews
     Assert(laplacians.size() == fe_values->n_quadrature_points,
            ExcDimensionMismatch(laplacians.size(),
                                 fe_values->n_quadrature_points));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     Assert(
-      fe_function.size() == fe_values->present_cell->n_dofs_for_dof_handler(),
+      fe_function.size() == fe_values->present_cell.n_dofs_for_dof_handler(),
       ExcDimensionMismatch(fe_function.size(),
-                           fe_values->present_cell->n_dofs_for_dof_handler()));
+                           fe_values->present_cell.n_dofs_for_dof_handler()));
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_laplacians<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_hessians,
@@ -2190,8 +2190,8 @@ namespace FEValuesViews
     Assert(laplacians.size() == fe_values->n_quadrature_points,
            ExcDimensionMismatch(laplacians.size(),
                                 fe_values->n_quadrature_points));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_laplacians<dim, spacedim>(
@@ -2215,16 +2215,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_3rd_derivatives,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_3rd_derivatives")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_derivatives<3, dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_3rd_derivatives,
@@ -2246,8 +2246,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_3rd_derivatives,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_3rd_derivatives")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_derivatives<3, dim, spacedim>(
@@ -2270,16 +2270,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_values<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_values,
@@ -2300,8 +2300,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_values<dim, spacedim>(
@@ -2324,17 +2324,17 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs
     // on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_divergences<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -2356,8 +2356,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_divergences<dim, spacedim>(
@@ -2380,16 +2380,16 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_values<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_values,
@@ -2410,8 +2410,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_values,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_values")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_values<dim, spacedim>(
@@ -2434,17 +2434,17 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs
     // on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_divergences<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -2465,8 +2465,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_divergences<dim, spacedim>(
@@ -2489,17 +2489,17 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(fe_function.size(),
-                    fe_values->present_cell->n_dofs_for_dof_handler());
+                    fe_values->present_cell.n_dofs_for_dof_handler());
 
     // get function values of dofs
     // on this cell
     dealii::Vector<typename InputVector::value_type> dof_values(
       fe_values->dofs_per_cell);
-    fe_values->present_cell->get_interpolated_dof_values(fe_function,
-                                                         dof_values);
+    fe_values->present_cell.get_interpolated_dof_values(fe_function,
+                                                        dof_values);
     internal::do_function_gradients<dim, spacedim>(
       make_array_view(dof_values.begin(), dof_values.end()),
       fe_values->finite_element_output.shape_gradients,
@@ -2520,8 +2520,8 @@ namespace FEValuesViews
     Assert(fe_values->update_flags & update_gradients,
            (typename FEValuesBase<dim, spacedim>::ExcAccessToUninitializedField(
              "update_gradients")));
-    Assert(fe_values->present_cell.get() != nullptr,
-           ExcMessage("FEValues object is not reinit'ed to any cell"));
+    Assert(fe_values->present_cell.is_initialized(),
+           (typename FEValuesBase<dim, spacedim>::ExcNotReinited()));
     AssertDimension(dof_values.size(), fe_values->dofs_per_cell);
 
     internal::do_function_gradients<dim, spacedim>(
@@ -2589,13 +2589,25 @@ namespace internal
 } // namespace internal
 
 
-/* ---------------- FEValuesBase<dim,spacedim>::CellIteratorBase --------- */
+/* ---------------- FEValuesBase<dim,spacedim>::CellIteratorContainer ---------
+ */
+
+template <int dim, int spacedim>
+FEValuesBase<dim, spacedim>::CellIteratorContainer::CellIteratorContainer()
+  : initialized(false)
+  , cell(typename Triangulation<dim, spacedim>::cell_iterator(nullptr, -1, -1))
+  , dof_handler(nullptr)
+  , level_dof_access(false)
+{}
+
+
 
 template <int dim, int spacedim>
 template <bool lda>
-FEValuesBase<dim, spacedim>::CellIteratorBase::CellIteratorBase(
+FEValuesBase<dim, spacedim>::CellIteratorContainer::CellIteratorContainer(
   const TriaIterator<DoFCellAccessor<dim, spacedim, lda>> &cell)
-  : cell(cell)
+  : initialized(true)
+  , cell(cell)
   , dof_handler(&cell->get_dof_handler())
   , level_dof_access(lda)
 {}
@@ -2603,9 +2615,10 @@ FEValuesBase<dim, spacedim>::CellIteratorBase::CellIteratorBase(
 
 
 template <int dim, int spacedim>
-FEValuesBase<dim, spacedim>::CellIteratorBase::CellIteratorBase(
+FEValuesBase<dim, spacedim>::CellIteratorContainer::CellIteratorContainer(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell)
-  : cell(cell)
+  : initialized(true)
+  , cell(cell)
   , dof_handler(nullptr)
   , level_dof_access(false)
 {}
@@ -2613,9 +2626,20 @@ FEValuesBase<dim, spacedim>::CellIteratorBase::CellIteratorBase(
 
 
 template <int dim, int spacedim>
-FEValuesBase<dim, spacedim>::CellIteratorBase::
+bool
+FEValuesBase<dim, spacedim>::CellIteratorContainer::is_initialized() const
+{
+  return initialized;
+}
+
+
+
+template <int dim, int spacedim>
+FEValuesBase<dim, spacedim>::CellIteratorContainer::
 operator typename Triangulation<dim, spacedim>::cell_iterator() const
 {
+  Assert(is_initialized(), ExcNotReinited());
+
   return cell;
 }
 
@@ -2623,9 +2647,12 @@ operator typename Triangulation<dim, spacedim>::cell_iterator() const
 
 template <int dim, int spacedim>
 types::global_dof_index
-FEValuesBase<dim, spacedim>::CellIteratorBase::n_dofs_for_dof_handler() const
+FEValuesBase<dim, spacedim>::CellIteratorContainer::n_dofs_for_dof_handler()
+  const
 {
+  Assert(is_initialized(), ExcNotReinited());
   Assert(dof_handler != nullptr, ExcNeedsDoFHandler());
+
   return dof_handler->n_dofs();
 }
 
@@ -2634,10 +2661,11 @@ FEValuesBase<dim, spacedim>::CellIteratorBase::n_dofs_for_dof_handler() const
 template <int dim, int spacedim>
 template <typename VectorType>
 void
-FEValuesBase<dim, spacedim>::CellIteratorBase::get_interpolated_dof_values(
+FEValuesBase<dim, spacedim>::CellIteratorContainer::get_interpolated_dof_values(
   const VectorType &                       in,
   Vector<typename VectorType::value_type> &out) const
 {
+  Assert(is_initialized(), ExcNotReinited());
   Assert(dof_handler != nullptr, ExcNeedsDoFHandler());
 
   if (level_dof_access)
@@ -2658,10 +2686,11 @@ FEValuesBase<dim, spacedim>::CellIteratorBase::get_interpolated_dof_values(
 
 template <int dim, int spacedim>
 void
-FEValuesBase<dim, spacedim>::CellIteratorBase::get_interpolated_dof_values(
+FEValuesBase<dim, spacedim>::CellIteratorContainer::get_interpolated_dof_values(
   const IndexSet &              in,
   Vector<IndexSet::value_type> &out) const
 {
+  Assert(is_initialized(), ExcNotReinited());
   Assert(dof_handler != nullptr, ExcNeedsDoFHandler());
   Assert(level_dof_access == false, ExcNotImplemented());
 
@@ -3344,13 +3373,12 @@ FEValuesBase<dim, spacedim>::get_function_values(
   Assert(this->update_flags & update_values,
          ExcAccessToUninitializedField("update_values"));
   AssertDimension(fe->n_components(), 1);
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_values(dof_values.begin(),
                                this->finite_element_output.shape_values,
                                values);
@@ -3390,16 +3418,15 @@ FEValuesBase<dim, spacedim>::get_function_values(
   std::vector<Vector<typename InputVector::value_type>> &values) const
 {
   using Number = typename InputVector::value_type;
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
 
   Assert(this->update_flags & update_values,
          ExcAccessToUninitializedField("update_values"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_values(
     dof_values.begin(),
     this->finite_element_output.shape_values,
@@ -3486,13 +3513,12 @@ FEValuesBase<dim, spacedim>::get_function_gradients(
   Assert(this->update_flags & update_gradients,
          ExcAccessToUninitializedField("update_gradients"));
   AssertDimension(fe->n_components(), 1);
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_derivatives(dof_values.begin(),
                                     this->finite_element_output.shape_gradients,
                                     gradients);
@@ -3537,13 +3563,12 @@ FEValuesBase<dim, spacedim>::get_function_gradients(
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_gradients,
          ExcAccessToUninitializedField("update_gradients"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_derivatives(
     dof_values.begin(),
     this->finite_element_output.shape_gradients,
@@ -3599,13 +3624,12 @@ FEValuesBase<dim, spacedim>::get_function_hessians(
   AssertDimension(fe->n_components(), 1);
   Assert(this->update_flags & update_hessians,
          ExcAccessToUninitializedField("update_hessians"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_derivatives(dof_values.begin(),
                                     this->finite_element_output.shape_hessians,
                                     hessians);
@@ -3625,7 +3649,7 @@ FEValuesBase<dim, spacedim>::get_function_hessians(
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_hessians,
          ExcAccessToUninitializedField("update_hessians"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
   AssertDimension(indices.size(), dofs_per_cell);
 
   boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
@@ -3651,13 +3675,12 @@ FEValuesBase<dim, spacedim>::get_function_hessians(
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_hessians,
          ExcAccessToUninitializedField("update_hessians"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_derivatives(
     dof_values.begin(),
     this->finite_element_output.shape_hessians,
@@ -3711,13 +3734,12 @@ FEValuesBase<dim, spacedim>::get_function_laplacians(
   Assert(this->update_flags & update_hessians,
          ExcAccessToUninitializedField("update_hessians"));
   AssertDimension(fe->n_components(), 1);
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_laplacians(dof_values.begin(),
                                    this->finite_element_output.shape_hessians,
                                    laplacians);
@@ -3757,15 +3779,14 @@ FEValuesBase<dim, spacedim>::get_function_laplacians(
   std::vector<Vector<typename InputVector::value_type>> &laplacians) const
 {
   using Number = typename InputVector::value_type;
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
+  Assert(present_cell.is_initialized(), ExcNotReinited());
   Assert(this->update_flags & update_hessians,
          ExcAccessToUninitializedField("update_hessians"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_laplacians(
     dof_values.begin(),
     this->finite_element_output.shape_hessians,
@@ -3849,13 +3870,12 @@ FEValuesBase<dim, spacedim>::get_function_third_derivatives(
   AssertDimension(fe->n_components(), 1);
   Assert(this->update_flags & update_3rd_derivatives,
          ExcAccessToUninitializedField("update_3rd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_derivatives(
     dof_values.begin(),
     this->finite_element_output.shape_3rd_derivatives,
@@ -3876,7 +3896,7 @@ FEValuesBase<dim, spacedim>::get_function_third_derivatives(
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_3rd_derivatives,
          ExcAccessToUninitializedField("update_3rd_derivatives"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
   AssertDimension(indices.size(), dofs_per_cell);
 
   boost::container::small_vector<Number, 200> dof_values(dofs_per_cell);
@@ -3903,13 +3923,12 @@ FEValuesBase<dim, spacedim>::get_function_third_derivatives(
   using Number = typename InputVector::value_type;
   Assert(this->update_flags & update_3rd_derivatives,
          ExcAccessToUninitializedField("update_3rd_derivatives"));
-  Assert(present_cell.get() != nullptr,
-         ExcMessage("FEValues object is not reinit'ed to any cell"));
-  AssertDimension(fe_function.size(), present_cell->n_dofs_for_dof_handler());
+  Assert(present_cell.is_initialized(), ExcNotReinited());
+  AssertDimension(fe_function.size(), present_cell.n_dofs_for_dof_handler());
 
   // get function values of dofs on this cell
   Vector<Number> dof_values(dofs_per_cell);
-  present_cell->get_interpolated_dof_values(fe_function, dof_values);
+  present_cell.get_interpolated_dof_values(fe_function, dof_values);
   internal::do_function_derivatives(
     dof_values.begin(),
     this->finite_element_output.shape_3rd_derivatives,
@@ -3956,7 +3975,7 @@ template <int dim, int spacedim>
 const typename Triangulation<dim, spacedim>::cell_iterator
 FEValuesBase<dim, spacedim>::get_cell() const
 {
-  return *present_cell;
+  return present_cell;
 }
 
 
@@ -4020,14 +4039,14 @@ FEValuesBase<dim, spacedim>::invalidate_present_cell()
 {
   // if there is no present cell, then we shouldn't be
   // connected via a signal to a triangulation
-  Assert(present_cell.get() != nullptr, ExcInternalError());
+  Assert(present_cell.is_initialized(), ExcInternalError());
 
   // so delete the present cell and
   // disconnect from the signal we have with
   // it
   tria_listener_refinement.disconnect();
   tria_listener_mesh_transform.disconnect();
-  present_cell.reset();
+  present_cell = {};
 }
 
 
@@ -4037,11 +4056,11 @@ void
 FEValuesBase<dim, spacedim>::maybe_invalidate_previous_present_cell(
   const typename Triangulation<dim, spacedim>::cell_iterator &cell)
 {
-  if (present_cell.get() != nullptr)
+  if (present_cell.is_initialized())
     {
       if (&cell->get_triangulation() !=
           &present_cell
-             ->
+             .
              operator typename Triangulation<dim, spacedim>::cell_iterator()
              ->get_triangulation())
         {
@@ -4102,7 +4121,7 @@ FEValuesBase<dim, spacedim>::check_cell_similarity(
     }
 
   // case that there has not been any cell before
-  if (this->present_cell.get() == nullptr)
+  if (this->present_cell.is_initialized() == false)
     cell_similarity = CellSimilarity::none;
   else
     // in MappingQ, data can have been modified during the last call. Then, we
@@ -4113,14 +4132,14 @@ FEValuesBase<dim, spacedim>::check_cell_similarity(
     cell_similarity =
       (cell->is_translation_of(
          static_cast<const typename Triangulation<dim, spacedim>::cell_iterator
-                       &>(*this->present_cell)) ?
+                       &>(this->present_cell)) ?
          CellSimilarity::translation :
          CellSimilarity::none);
 
   if ((dim < spacedim) && (cell_similarity == CellSimilarity::translation))
     {
       if (static_cast<const typename Triangulation<dim, spacedim>::cell_iterator
-                        &>(*this->present_cell)
+                        &>(this->present_cell)
             ->direction_flag() != cell->direction_flag())
         cell_similarity = CellSimilarity::inverted_translation;
     }
@@ -4265,32 +4284,6 @@ FEValues<dim, spacedim>::initialize(const UpdateFlags update_flags)
 
 
 
-template <int dim, int spacedim>
-template <typename Type, typename Pointer, typename Iterator>
-void
-FEValuesBase<dim, spacedim>::reset_pointer_in_place_if_possible(
-  std::unique_ptr<Pointer> &present_cell,
-  const Iterator &          new_cell)
-{
-  // see if the existing pointer is non-null and if the type of
-  // the old object pointed to matches that of the one we'd
-  // like to create
-  if (present_cell.get() && (typeid(*present_cell.get()) == typeid(Type)))
-    {
-      // call destructor of the old object
-      static_cast<const Type *>(present_cell.get())->~Type();
-
-      // then construct a new object in-place
-      new (const_cast<void *>(static_cast<const void *>(present_cell.get())))
-        Type(new_cell);
-    }
-  else
-    // if the types don't match, there is nothing we can do here
-    present_cell = std::make_unique<Type>(new_cell);
-}
-
-
-
 template <int dim, int spacedim>
 void
 FEValues<dim, spacedim>::reinit(
@@ -4308,9 +4301,7 @@ FEValues<dim, spacedim>::reinit(
   this->maybe_invalidate_previous_present_cell(cell);
   this->check_cell_similarity(cell);
 
-  this->template reset_pointer_in_place_if_possible<
-    typename FEValuesBase<dim, spacedim>::CellIteratorBase>(this->present_cell,
-                                                            cell);
+  this->present_cell = {cell};
 
   // this was the part of the work that is dependent on the actual
   // data type of the iterator. now pass on to the function doing
@@ -4342,9 +4333,7 @@ FEValues<dim, spacedim>::reinit(
   this->maybe_invalidate_previous_present_cell(cell);
   this->check_cell_similarity(cell);
 
-  this->template reset_pointer_in_place_if_possible<
-    typename FEValuesBase<dim, spacedim>::CellIteratorBase>(this->present_cell,
-                                                            cell);
+  this->present_cell = {cell};
 
   // this was the part of the work that is dependent on the actual
   // data type of the iterator. now pass on to the function doing
@@ -4365,7 +4354,7 @@ FEValues<dim, spacedim>::do_reinit()
   if (this->update_flags & update_mapping)
     {
       this->cell_similarity =
-        this->get_mapping().fill_fe_values(*this->present_cell,
+        this->get_mapping().fill_fe_values(this->present_cell,
                                            this->cell_similarity,
                                            quadrature,
                                            *this->mapping_data,
@@ -4376,7 +4365,7 @@ FEValues<dim, spacedim>::do_reinit()
   // already filled by the mapping, let it compute the
   // data for the mapped shape function values, gradients,
   // etc.
-  this->get_fe().fill_fe_values(*this->present_cell,
+  this->get_fe().fill_fe_values(this->present_cell,
                                 this->cell_similarity,
                                 this->quadrature,
                                 this->get_mapping(),
@@ -4609,9 +4598,7 @@ FEFaceValues<dim, spacedim>::reinit(
   AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
 
   this->maybe_invalidate_previous_present_cell(cell);
-  this->template reset_pointer_in_place_if_possible<
-    typename FEValuesBase<dim, spacedim>::CellIteratorBase>(this->present_cell,
-                                                            cell);
+  this->present_cell = {cell};
 
   // this was the part of the work that is dependent on the actual
   // data type of the iterator. now pass on to the function doing
@@ -4643,9 +4630,7 @@ FEFaceValues<dim, spacedim>::reinit(
   AssertIndexRange(face_no, GeometryInfo<dim>::faces_per_cell);
 
   this->maybe_invalidate_previous_present_cell(cell);
-  this->template reset_pointer_in_place_if_possible<
-    typename FEValuesBase<dim, spacedim>::CellIteratorBase>(this->present_cell,
-                                                            cell);
+  this->present_cell = {cell};
 
   // this was the part of the work that is dependent on the actual
   // data type of the iterator. now pass on to the function doing
@@ -4675,19 +4660,19 @@ FEFaceValues<dim, spacedim>::do_reinit(const unsigned int face_no)
 
   // first of all, set the present_face_index (if available)
   const typename Triangulation<dim, spacedim>::cell_iterator cell =
-    *this->present_cell;
+    this->present_cell;
   this->present_face_index = cell->face_index(face_no);
 
   if (this->update_flags & update_mapping)
     {
-      this->get_mapping().fill_fe_face_values(*this->present_cell,
+      this->get_mapping().fill_fe_face_values(this->present_cell,
                                               face_no,
                                               this->quadrature,
                                               *this->mapping_data,
                                               this->mapping_output);
     }
 
-  this->get_fe().fill_fe_face_values(*this->present_cell,
+  this->get_fe().fill_fe_face_values(this->present_cell,
                                      face_no,
                                      this->quadrature,
                                      this->get_mapping(),
@@ -4856,9 +4841,7 @@ FESubfaceValues<dim, spacedim>::reinit(
                     "instead in these cases."));
 
   this->maybe_invalidate_previous_present_cell(cell);
-  this->template reset_pointer_in_place_if_possible<
-    typename FEValuesBase<dim, spacedim>::CellIteratorBase>(this->present_cell,
-                                                            cell);
+  this->present_cell = {cell};
 
   // this was the part of the work that is dependent on the actual
   // data type of the iterator. now pass on to the function doing
@@ -4903,9 +4886,7 @@ FESubfaceValues<dim, spacedim>::reinit(
                       cell->face(face_no)->n_children()));
 
   this->maybe_invalidate_previous_present_cell(cell);
-  this->template reset_pointer_in_place_if_possible<
-    typename FEValuesBase<dim, spacedim>::CellIteratorBase>(this->present_cell,
-                                                            cell);
+  this->present_cell = {cell};
 
   // this was the part of the work that is dependent on the actual
   // data type of the iterator. now pass on to the function doing
@@ -4938,7 +4919,7 @@ FESubfaceValues<dim, spacedim>::do_reinit(const unsigned int face_no,
 
   // first of all, set the present_face_index (if available)
   const typename Triangulation<dim, spacedim>::cell_iterator cell =
-    *this->present_cell;
+    this->present_cell;
 
   if (!cell->face(face_no)->has_children())
     // no subfaces at all, so set present_face_index to this face rather
@@ -5008,7 +4989,7 @@ FESubfaceValues<dim, spacedim>::do_reinit(const unsigned int face_no,
   // now ask the mapping and the finite element to do the actual work
   if (this->update_flags & update_mapping)
     {
-      this->get_mapping().fill_fe_subface_values(*this->present_cell,
+      this->get_mapping().fill_fe_subface_values(this->present_cell,
                                                  face_no,
                                                  subface_no,
                                                  this->quadrature[0],
@@ -5016,7 +4997,7 @@ FESubfaceValues<dim, spacedim>::do_reinit(const unsigned int face_no,
                                                  this->mapping_output);
     }
 
-  this->get_fe().fill_fe_subface_values(*this->present_cell,
+  this->get_fe().fill_fe_subface_values(this->present_cell,
                                         face_no,
                                         subface_no,
                                         this->quadrature[0],

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.