]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Let FETools::get_fe_by_name return a std::unique_ptr
authorDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Fri, 9 Feb 2018 22:32:41 +0000 (23:32 +0100)
committerDaniel Arndt <daniel.arndt@iwr.uni-heidelberg.de>
Sat, 10 Feb 2018 01:24:50 +0000 (02:24 +0100)
include/deal.II/fe/fe_tools.h
include/deal.II/fe/fe_tools.templates.h
source/fe/fe_tools.inst.in
tests/bits/fe_tools_09.cc
tests/bits/fe_tools_10.cc
tests/bits/fe_tools_11.cc
tests/fe/fe_tools_01.cc
tests/fe/get_fe_by_name_01.cc

index f4017d89cc23e41af21d53fc7cdb2247aef5a2b7..4c2aabfc86efe3d413294debdcff6f5be707f95c 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/base/geometry_info.h>
 #include <deal.II/base/tensor.h>
 #include <deal.II/base/symmetric_tensor.h>
+#include <deal.II/base/std_cxx14/memory.h>
 #include <deal.II/distributed/tria.h>
 #include <deal.II/fe/component_mask.h>
 #include <deal.II/lac/parallel_vector.h>
@@ -83,7 +84,7 @@ namespace FETools
     /**
      * Create a FiniteElement and return a pointer to it.
      */
-    virtual FiniteElement<dim,spacedim> *
+    virtual std::unique_ptr<FiniteElement<dim,spacedim> >
     get (const unsigned int degree) const = 0;
 
     /**
@@ -91,8 +92,9 @@ namespace FETools
      * implemented for FE_Q) and return a pointer to it.
      */
 
-    virtual FiniteElement<dim,spacedim> *
+    virtual std::unique_ptr<FiniteElement<dim,spacedim> >
     get (const Quadrature<1> &quad) const = 0;
+
     /**
      * Virtual destructor doing nothing but making the compiler happy.
      */
@@ -117,14 +119,14 @@ namespace FETools
     /**
      * Create a FiniteElement and return a pointer to it.
      */
-    virtual FiniteElement<FE::dimension,FE::space_dimension> *
+    virtual std::unique_ptr<FiniteElement<FE::dimension,FE::space_dimension> >
     get (const unsigned int degree) const;
 
     /**
      * Create a FiniteElement from a quadrature formula (currently only
      * implemented for FE_Q) and return a pointer to it.
      */
-    virtual FiniteElement<FE::dimension,FE::space_dimension> *
+    virtual std::unique_ptr<FiniteElement<FE::dimension,FE::space_dimension> >
     get (const Quadrature<1> &quad) const;
   };
 
@@ -1211,9 +1213,8 @@ namespace FETools
    * If no finite element can be reconstructed from this string, an exception
    * of type @p FETools::ExcInvalidFEName is thrown.
    *
-   * The function returns a pointer to a newly create finite element. It is in
-   * the caller's responsibility to destroy the object pointed to at an
-   * appropriate later time.
+   * The function returns a std::unique_ptr to a newly created finite element
+   * meaning the caller obtains ownership over the returned object.
    *
    * Since the value of the template argument can't be deduced from the
    * (string) argument given to this function, you have to explicitly specify
@@ -1226,7 +1227,7 @@ namespace FETools
    * one wants to get a codimension 1 finite element.
    */
   template <int dim, int spacedim = dim>
-  FiniteElement<dim, spacedim> *
+  std::unique_ptr<FiniteElement<dim, spacedim> >
   get_fe_by_name (const std::string &name);
 
 
@@ -1391,10 +1392,10 @@ namespace FETools
 namespace FETools
 {
   template <class FE>
-  FiniteElement<FE::dimension, FE::space_dimension> *
+  std::unique_ptr<FiniteElement<FE::dimension, FE::space_dimension> >
   FEFactory<FE>::get (const unsigned int degree) const
   {
-    return new FE(degree);
+    return std_cxx14::make_unique<FE>(degree);
   }
 
   namespace Compositing
index 41da17234b3a3234bb792cbbd0e7b5444335ea85..5079e80a5ce78882ccf122a50b149d9468c67cac 100644 (file)
@@ -21,6 +21,7 @@
 #include <deal.II/base/qprojector.h>
 #include <deal.II/base/thread_management.h>
 #include <deal.II/base/utilities.h>
+#include <deal.II/base/std_cxx14/memory.h>
 #include <deal.II/lac/full_matrix.h>
 #include <deal.II/lac/householder.h>
 #include <deal.II/lac/constraint_matrix.h>
@@ -938,7 +939,7 @@ namespace FETools
 
   // Not implemented in the general case.
   template <class FE>
-  FiniteElement<FE::dimension, FE::space_dimension> *
+  std::unique_ptr<FiniteElement<FE::dimension, FE::space_dimension> >
   FEFactory<FE>::get (const Quadrature<1> &) const
   {
     Assert(false, ExcNotImplemented());
@@ -947,111 +948,111 @@ namespace FETools
 
   // Specializations for FE_Q.
   template <>
-  FiniteElement<1, 1> *
+  std::unique_ptr<FiniteElement<1, 1> >
   FEFactory<FE_Q<1, 1> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q<1>(quad);
+    return std_cxx14::make_unique<FE_Q<1>>(quad);
   }
 
   template <>
-  FiniteElement<2, 2> *
+  std::unique_ptr<FiniteElement<2, 2> >
   FEFactory<FE_Q<2, 2> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q<2>(quad);
+    return std_cxx14::make_unique<FE_Q<2>>(quad);
   }
 
   template <>
-  FiniteElement<3, 3> *
+  std::unique_ptr<FiniteElement<3, 3> >
   FEFactory<FE_Q<3, 3> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q<3>(quad);
+    return std_cxx14::make_unique<FE_Q<3>>(quad);
   }
 
   // Specializations for FE_Q_DG0.
   template <>
-  FiniteElement<1, 1> *
+  std::unique_ptr<FiniteElement<1, 1> >
   FEFactory<FE_Q_DG0<1, 1> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q_DG0<1>(quad);
+    return std_cxx14::make_unique<FE_Q_DG0<1>>(quad);
   }
 
   template <>
-  FiniteElement<2, 2> *
+  std::unique_ptr<FiniteElement<2, 2> >
   FEFactory<FE_Q_DG0<2, 2> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q_DG0<2>(quad);
+    return std_cxx14::make_unique<FE_Q_DG0<2>>(quad);
   }
 
   template <>
-  FiniteElement<3, 3> *
+  std::unique_ptr<FiniteElement<3, 3> >
   FEFactory<FE_Q_DG0<3, 3> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q_DG0<3>(quad);
+    return std_cxx14::make_unique<FE_Q_DG0<3>>(quad);
   }
 
   // Specializations for FE_Q_Bubbles.
   template <>
-  FiniteElement<1, 1> *
+  std::unique_ptr<FiniteElement<1, 1> >
   FEFactory<FE_Q_Bubbles<1, 1> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q_Bubbles<1>(quad);
+    return std_cxx14::make_unique<FE_Q_Bubbles<1>>(quad);
   }
 
   template <>
-  FiniteElement<2, 2> *
+  std::unique_ptr<FiniteElement<2, 2> >
   FEFactory<FE_Q_Bubbles<2, 2> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q_Bubbles<2>(quad);
+    return std_cxx14::make_unique<FE_Q_Bubbles<2>>(quad);
   }
 
   template <>
-  FiniteElement<3, 3> *
+  std::unique_ptr<FiniteElement<3, 3> >
   FEFactory<FE_Q_Bubbles<3, 3> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_Q_Bubbles<3>(quad);
+    return std_cxx14::make_unique<FE_Q_Bubbles<3>>(quad);
   }
 
   // Specializations for FE_DGQArbitraryNodes.
   template <>
-  FiniteElement<1, 1> *
+  std::unique_ptr<FiniteElement<1, 1> >
   FEFactory<FE_DGQ<1> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_DGQArbitraryNodes<1>(quad);
+    return std_cxx14::make_unique<FE_DGQArbitraryNodes<1>>(quad);
   }
 
   template <>
-  FiniteElement<1, 2> *
+  std::unique_ptr<FiniteElement<1, 2> >
   FEFactory<FE_DGQ<1, 2> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_DGQArbitraryNodes<1, 2>(quad);
+    return std_cxx14::make_unique<FE_DGQArbitraryNodes<1, 2>>(quad);
   }
 
   template <>
-  FiniteElement<1, 3> *
+  std::unique_ptr<FiniteElement<1, 3> >
   FEFactory<FE_DGQ<1, 3> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_DGQArbitraryNodes<1, 3>(quad);
+    return std_cxx14::make_unique<FE_DGQArbitraryNodes<1, 3>>(quad);
   }
 
   template <>
-  FiniteElement<2, 2> *
+  std::unique_ptr<FiniteElement<2, 2> >
   FEFactory<FE_DGQ<2> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_DGQArbitraryNodes<2>(quad);
+    return std_cxx14::make_unique<FE_DGQArbitraryNodes<2>>(quad);
   }
 
   template <>
-  FiniteElement<2, 3> *
+  std::unique_ptr<FiniteElement<2, 3> >
   FEFactory<FE_DGQ<2, 3> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_DGQArbitraryNodes<2, 3>(quad);
+    return std_cxx14::make_unique<FE_DGQArbitraryNodes<2, 3>>(quad);
   }
 
   template <>
-  FiniteElement<3, 3> *
+  std::unique_ptr<FiniteElement<3, 3> >
   FEFactory<FE_DGQ<3> >::get (const Quadrature<1> &quad) const
   {
-    return new FE_DGQArbitraryNodes<3>(quad);
+    return std_cxx14::make_unique<FE_DGQArbitraryNodes<3>>(quad);
   }
 }
 
@@ -1065,9 +1066,9 @@ namespace
 
   template <int dim>
   void
-  fill_no_codim_fe_names (std::map<std::string,std::shared_ptr<const Subscriptor> > &result)
+  fill_no_codim_fe_names (std::map<std::string,std::unique_ptr<const Subscriptor> > &result)
   {
-    typedef std::shared_ptr<const Subscriptor> FEFactoryPointer;
+    typedef std::unique_ptr<const Subscriptor> FEFactoryPointer;
 
     result["FE_Q_Hierarchical"]
       = FEFactoryPointer(new FETools::FEFactory<FE_Q_Hierarchical<dim> >);
@@ -1130,9 +1131,9 @@ namespace
   // nonzero codimension.
   template <int dim, int spacedim>
   void
-  fill_codim_fe_names (std::map<std::string,std::shared_ptr<const Subscriptor> > &result)
+  fill_codim_fe_names (std::map<std::string,std::unique_ptr<const Subscriptor> > &result)
   {
-    typedef std::shared_ptr<const Subscriptor> FEFactoryPointer;
+    typedef std::unique_ptr<const Subscriptor> FEFactoryPointer;
 
     result["FE_Bernstein"]
       = FEFactoryPointer(new FETools::FEFactory<FE_Bernstein<dim,spacedim> >);
@@ -1166,12 +1167,12 @@ namespace
   // by the functions above.
   std::vector<std::vector<
   std::map<std::string,
-      std::shared_ptr<const Subscriptor> > > >
+      std::unique_ptr<const Subscriptor> > > >
       fill_default_map()
   {
     std::vector<std::vector<
     std::map<std::string,
-        std::shared_ptr<const Subscriptor> > > >
+        std::unique_ptr<const Subscriptor> > > >
         result(4);
 
     for (unsigned int d=0; d<4; ++d)
@@ -1211,7 +1212,7 @@ namespace
   // all entries of this vector are of different type, we store
   // pointers to generic objects and cast them when needed.
 
-  // We use a shared pointer to factory objects, to ensure that they
+  // We use a unique pointer to factory objects, to ensure that they
   // get deleted at the end of the program run and don't end up as
   // apparent memory leaks to programs like valgrind.
 
@@ -1223,7 +1224,7 @@ namespace
   static
   std::vector<std::vector<
   std::map<std::string,
-      std::shared_ptr<const Subscriptor> > > >
+      std::unique_ptr<const Subscriptor> > > >
       fe_name_map = fill_default_map();
 }
 
@@ -2293,7 +2294,7 @@ namespace FETools
     // Insert the normalized name into
     // the map
     fe_name_map[dim][spacedim][name] =
-      std::shared_ptr<const Subscriptor> (factory);
+      std::unique_ptr<const Subscriptor> (factory);
   }
 
 
@@ -2306,10 +2307,10 @@ namespace FETools
       // have a unique interface. could be done
       // smarter?
       template <int dim, int spacedim>
-      FiniteElement<dim,spacedim> *
+      std::unique_ptr<FiniteElement<dim,spacedim> >
       get_fe_by_name_ext (std::string &name,
                           const std::map<std::string,
-                          std::shared_ptr<const Subscriptor> >
+                          std::unique_ptr<const Subscriptor> >
                           &fe_name_map)
       {
         // Extract the name of the
@@ -2338,120 +2339,84 @@ namespace FETools
             // recursive calls if one of
             // these calls should throw
             // an exception
-            std::vector<const FiniteElement<dim,spacedim>*> base_fes;
+            std::vector<std::unique_ptr<const FiniteElement<dim,spacedim> > > base_fes;
             std::vector<unsigned int>        base_multiplicities;
-            try
+
+            // Now, just the [...]
+            // part should be left.
+            if (name.size() == 0 || name[0] != '[')
+              throw (std::string("Invalid first character in ") + name);
+            do
               {
-                // Now, just the [...]
-                // part should be left.
-                if (name.size() == 0 || name[0] != '[')
-                  throw (std::string("Invalid first character in ") + name);
-                do
+                // Erase the
+                // leading '[' or '-'
+                name.erase(0,1);
+                // Now, the name of the
+                // first base element is
+                // first... Let's get it
+                base_fes.push_back (get_fe_by_name_ext<dim,spacedim> (name,
+                                                                      fe_name_map));
+                // next check whether
+                // FESystem placed a
+                // multiplicity after
+                // the element name
+                if (name[0] == '^')
                   {
-                    // Erase the
-                    // leading '[' or '-'
+                    // yes. Delete the '^'
+                    // and read this
+                    // multiplicity
                     name.erase(0,1);
-                    // Now, the name of the
-                    // first base element is
-                    // first... Let's get it
-                    base_fes.push_back (get_fe_by_name_ext<dim,spacedim> (name,
-                                                                          fe_name_map));
-                    // next check whether
-                    // FESystem placed a
-                    // multiplicity after
-                    // the element name
-                    if (name[0] == '^')
-                      {
-                        // yes. Delete the '^'
-                        // and read this
-                        // multiplicity
-                        name.erase(0,1);
-
-                        const std::pair<int,unsigned int> tmp
-                          = Utilities::get_integer_at_position (name, 0);
-                        name.erase(0, tmp.second);
-                        // add to length,
-                        // including the '^'
-                        base_multiplicities.push_back (tmp.first);
-                      }
-                    else
-                      // no, so
-                      // multiplicity is
-                      // 1
-                      base_multiplicities.push_back (1);
-
-                    // so that's it for
-                    // this base
-                    // element. base
-                    // elements are
-                    // separated by '-',
-                    // and the list is
-                    // terminated by ']',
-                    // so loop while the
-                    // next character is
-                    // '-'
-                  }
-                while (name[0] == '-');
-
-                // so we got to the end
-                // of the '-' separated
-                // list. make sure that
-                // we actually had a ']'
-                // there
-                if (name.size() == 0 || name[0] != ']')
-                  throw (std::string("Invalid first character in ") + name);
-                name.erase(0,1);
-                // just one more sanity check
-                Assert ((base_fes.size() == base_multiplicities.size())
-                        &&
-                        (base_fes.size() > 0),
-                        ExcInternalError());
-
-                // ok, apparently
-                // everything went ok. so
-                // generate the composed
-                // element
-                FiniteElement<dim,spacedim> *system_element = nullptr;
-
-                // uses new FESystem constructor
-                // which is independent of
-                // the number of FEs in the system
-                system_element = new FESystem<dim,spacedim>(base_fes, base_multiplicities);
-
-                // now we don't need the
-                // list of base elements
-                // any more
-                for (unsigned int i=0; i<base_fes.size(); ++i)
-                  delete base_fes[i];
-
-                // finally return our
-                // findings
-                // Add the closing ']' to
-                // the length
-                return system_element;
 
+                    const std::pair<int,unsigned int> tmp
+                      = Utilities::get_integer_at_position (name, 0);
+                    name.erase(0, tmp.second);
+                    // add to length,
+                    // including the '^'
+                    base_multiplicities.push_back (tmp.first);
+                  }
+                else
+                  // no, so
+                  // multiplicity is
+                  // 1
+                  base_multiplicities.push_back (1);
+
+                // so that's it for
+                // this base
+                // element. base
+                // elements are
+                // separated by '-',
+                // and the list is
+                // terminated by ']',
+                // so loop while the
+                // next character is
+                // '-'
               }
-            catch (...)
-              {
-                // ups, some exception
-                // was thrown. prevent a
-                // memory leak, and then
-                // pass on the exception
-                // to the caller
-                for (unsigned int i=0; i<base_fes.size(); ++i)
-                  delete base_fes[i];
-                throw;
-              }
-
-            // this is a place where we
-            // should really never get,
-            // since above we have either
-            // returned from the
-            // try-clause, or have
-            // re-thrown in the catch
-            // clause. check that we
-            // never get here
-            Assert (false, ExcInternalError());
+            while (name[0] == '-');
+
+            // so we got to the end
+            // of the '-' separated
+            // list. make sure that
+            // we actually had a ']'
+            // there
+            if (name.size() == 0 || name[0] != ']')
+              throw (std::string("Invalid first character in ") + name);
+            name.erase(0,1);
+            // just one more sanity check
+            Assert ((base_fes.size() == base_multiplicities.size())
+                    &&
+                    (base_fes.size() > 0),
+                    ExcInternalError());
+
+            // this is a workaround since the constructor for FESystem
+            // only accepts raw pointers.
+            std::vector<const FiniteElement<dim,spacedim> *> raw_base_fes;
+            for (const auto &base_fe: base_fes)
+              raw_base_fes.push_back(base_fe.get());
+
+            // ok, apparently everything went ok. so generate the composed element.
+            // and return it.
+            return std_cxx14::make_unique<FESystem<dim,spacedim> >
+                   (raw_base_fes, base_multiplicities);
           }
         else if (name_part == "FE_Nothing")
           {
@@ -2556,7 +2521,7 @@ namespace FETools
 
 
       template <int dim,int spacedim>
-      FiniteElement<dim,spacedim> *get_fe_by_name (std::string &name)
+      std::unique_ptr<FiniteElement<dim,spacedim> >get_fe_by_name (std::string &name)
       {
         return get_fe_by_name_ext<dim,spacedim> (name, fe_name_map[dim][spacedim]);
       }
@@ -2569,13 +2534,13 @@ namespace FETools
   FiniteElement<dim> *
   get_fe_from_name (const std::string &parameter_name)
   {
-    return get_fe_by_name<dim,dim> (parameter_name);
+    return get_fe_by_name<dim,dim> (parameter_name).release();
   }
 
 
 
   template <int dim, int spacedim>
-  FiniteElement<dim, spacedim> *
+  std::unique_ptr<FiniteElement<dim, spacedim> >
   get_fe_by_name (const std::string &parameter_name)
   {
     std::string name = Utilities::trim(parameter_name);
@@ -2643,7 +2608,7 @@ namespace FETools
 
     try
       {
-        FiniteElement<dim,spacedim> *fe = internal::get_fe_by_name<dim,spacedim> (name);
+        auto fe = internal::get_fe_by_name<dim,spacedim> (name);
 
         // Make sure the auxiliary function
         // ate up all characters of the name.
index 30b1a8b698a4b6b2664f11f480e99018a45abc8b..e84477c400b3c753d39902f073aa817ea3e6229f 100644 (file)
@@ -130,7 +130,7 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension :  SPACE_DIMENSIONS
      const FiniteElement<deal_II_dimension,deal_II_space_dimension> &,
      FullMatrix<double> &);
 
-    template FiniteElement<deal_II_dimension,deal_II_space_dimension> *
+    template std::unique_ptr<FiniteElement<deal_II_dimension,deal_II_space_dimension> >
     get_fe_by_name<deal_II_dimension,deal_II_space_dimension> (const std::string &);
 
     template
index e598844729562a0b14fa0d4ce52f3080500168ee..54c54d4da7f58fe5eee80a4fe2d9bb4b7a128b38 100644 (file)
@@ -29,8 +29,6 @@ void
 check_this (const FiniteElement<dim> &fe1,
             const FiniteElement<dim> &fe2)
 {
-  FiniteElement<dim> *p1, *p2;
-
   // check that the name of the fe
   // and the name of the fe that we
   // re-create from this name are
@@ -38,18 +36,16 @@ check_this (const FiniteElement<dim> &fe1,
   // pretty good indication that the
   // two FEs are actually the same
   deallog << fe1.get_name();
-  p1 = FETools::get_fe_by_name<dim, dim> (fe1.get_name());
+  std::unique_ptr<FiniteElement<dim>> p1 = FETools::get_fe_by_name<dim, dim> (fe1.get_name());
   AssertThrow (fe1.get_name() == p1->get_name(),
                ExcInternalError());
   deallog << " ok" << std::endl;
-  delete p1;
 
   // same for fe2
   deallog << fe2.get_name();
-  p2 = FETools::get_fe_by_name<dim, dim> (fe2.get_name());
+  std::unique_ptr<FiniteElement<dim>> p2 = FETools::get_fe_by_name<dim, dim> (fe2.get_name());
   AssertThrow (fe2.get_name() == p2->get_name(),
                ExcInternalError());
   deallog << " ok" << std::endl;
-  delete p2;
 }
 
index e4f7cfa3bc04bd5d5c2b11937ce66ee0ff2dedc3..84d5d141765f989f9f373888ff99b588416389b0 100644 (file)
@@ -49,8 +49,6 @@ void
 check_this (const FiniteElement<dim> &fe1,
             const FiniteElement<dim> &fe2)
 {
-  FiniteElement<dim> *p1, *p2;
-
   // check that the name of the fe
   // and the name of the fe that we
   // re-create from this name are
@@ -58,18 +56,16 @@ check_this (const FiniteElement<dim> &fe1,
   // pretty good indication that the
   // two FEs are actually the same
   deallog << modify_name<dim> (fe1.get_name());
-  p1 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe1.get_name()));
+  std::unique_ptr<FiniteElement<dim> > p1 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe1.get_name()));
   AssertThrow (fe1.get_name() == p1->get_name(),
                ExcInternalError());
   deallog << " ok" << std::endl;
-  delete p1;
 
   // same for fe2
   deallog << modify_name<dim> (fe2.get_name());
-  p2 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe2.get_name()));
+  std::unique_ptr<FiniteElement<dim> > p2 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe2.get_name()));
   AssertThrow (fe2.get_name() == p2->get_name(),
                ExcInternalError());
   deallog << " ok" << std::endl;
-  delete p2;
 }
 
index 6298ef14508d77a9d853c2bca449abc2ca6793ac..4e2a6cab688c0ed4ab3c5fd0dabe7f0b5e8892fc 100644 (file)
@@ -48,8 +48,6 @@ void
 check_this (const FiniteElement<dim> &fe1,
             const FiniteElement<dim> &fe2)
 {
-  FiniteElement<dim> *p1, *p2;
-
   // check that the name of the fe
   // and the name of the fe that we
   // re-create from this name are
@@ -57,18 +55,16 @@ check_this (const FiniteElement<dim> &fe1,
   // pretty good indication that the
   // two FEs are actually the same
   deallog << modify_name<dim> (fe1.get_name());
-  p1 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe1.get_name()));
+  std::unique_ptr<FiniteElement<dim>> p1 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe1.get_name()));
   AssertThrow (fe1.get_name() == p1->get_name(),
                ExcInternalError());
   deallog << " ok" << std::endl;
-  delete p1;
 
   // same for fe2
   deallog << modify_name<dim> (fe2.get_name());
-  p2 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe2.get_name()));
+  std::unique_ptr<FiniteElement<dim>> p2 = FETools::get_fe_by_name<dim, dim> (modify_name<dim> (fe2.get_name()));
   AssertThrow (fe2.get_name() == p2->get_name(),
                ExcInternalError());
   deallog << " ok" << std::endl;
-  delete p2;
 }
 
index adc93b5ab6433f40719e5a58068fc0d5e16d3a66..26afa8595306d7f9d1309cebe38708c3ab0abfd3 100644 (file)
@@ -27,7 +27,7 @@
 template <int dim>
 void test_fe(const char *name)
 {
-  FiniteElement<dim> *fe = FETools::get_fe_by_name<dim, dim>(std::string(name));
+  std::unique_ptr<FiniteElement<dim> > fe = FETools::get_fe_by_name<dim, dim>(std::string(name));
 
   deallog << fe->get_name() << std::endl
           << '\t' << fe->dofs_per_cell
@@ -36,8 +36,6 @@ void test_fe(const char *name)
           << '\t' << fe->dofs_per_quad
           << '\t' << fe->dofs_per_hex
           << std::endl;
-
-  delete fe;
 }
 
 
index f61a491a4833de3d1587a43fe745857bca3342e7..25c80cdadb15d770953d49852c32bec81a732a1e 100644 (file)
@@ -35,13 +35,11 @@ public:
   {
     std::string name = myname;
 
-    FiniteElement<dim, spacedim> *fe = FETools::get_fe_by_name<dim, spacedim>(name);
+    std::unique_ptr<FiniteElement<dim, spacedim> >fe = FETools::get_fe_by_name<dim, spacedim>(name);
 
     deallog << "Read " << name << std::endl;
     deallog << "Generated :" << std::endl;
     deallog << fe->get_name()  << std::endl;
-
-    delete fe;
   }
 
 

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.