Do not use a plain pointer in PolynomialsABF.

These objects are being copied by the FE_Poly* classes, so plain pointers without
dedicated copy constructors and operators are likely going to lead to memory
corruption. It's not clear to me how this ever worked, but it's easy to fix.

diff --git a/include/deal.II/base/polynomials_abf.h b/include/deal.II/base/polynomials_abf.h
index 3ce995f56ee1e25d8d4f5001b34703e136f9ef37..f183c406eebe1e096b89b3884d600f4fa22fb3ac 100644 (file)
--- a/include/deal.II/base/polynomials_abf.h
+++ b/include/deal.II/base/polynomials_abf.h
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2004 - 2015 by the deal.II authors
+// Copyright (C) 2004 - 2017 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -27,6 +27,7 @@
 #include <deal.II/base/table.h>
 #include <deal.II/base/thread_management.h>
 
+#include <deal.II/base/std_cxx11/shared_ptr.h>
 #include <vector>
 
 DEAL_II_NAMESPACE_OPEN
@@ -63,11 +64,6 @@ public:
    */
   PolynomialsABF (const unsigned int k);
 
-  /**
-   * Destructor deleting the polynomials.
-   */
-  ~PolynomialsABF ();
-
   /**
    * Compute the value and the first and second derivatives of each Raviart-
    * Thomas polynomial at @p unit_point.
@@ -118,9 +114,10 @@ private:
 
   /**
    * An object representing the polynomial space for a single component. We
-   * can re-use it by rotating the coordinates of the evaluation point.
+   * can re-use it for the other vector components by rotating the
+   * coordinates of the evaluation point.
    */
-  AnisotropicPolynomials<dim> *polynomial_space;
+  const AnisotropicPolynomials<dim> polynomial_space;
 
   /**
    * Number of Raviart-Thomas polynomials.

diff --git a/source/base/polynomials_abf.cc b/source/base/polynomials_abf.cc
index ae13827fa4cad300c80b61586829a343c264a45a..c420338a97d620288616c07791ad8bb0eb23e3cd 100644 (file)
--- a/source/base/polynomials_abf.cc
+++ b/source/base/polynomials_abf.cc
@@ -23,38 +23,42 @@
 DEAL_II_NAMESPACE_OPEN
 
 
+
+namespace
+{
+  template <int dim>
+  std::vector<std::vector< Polynomials::Polynomial< double > > >
+  get_abf_polynomials (const unsigned int k)
+  {
+    std::vector<std::vector< Polynomials::Polynomial< double > > > pols(dim);
+    pols[0] = Polynomials::LagrangeEquidistant::generate_complete_basis(k+2);
+
+    if (k == 0)
+      for (unsigned int d=1; d<dim; ++d)
+        pols[d] = Polynomials::Legendre::generate_complete_basis(0);
+    else
+      for (unsigned int d=1; d<dim; ++d)
+        pols[d] = Polynomials::LagrangeEquidistant::generate_complete_basis(k);
+
+    return pols;
+  }
+}
+
 template <int dim>
 PolynomialsABF<dim>::PolynomialsABF (const unsigned int k)
   :
   my_degree(k),
+  polynomial_space(get_abf_polynomials<dim>(k)),
   n_pols(compute_n_pols(k))
 {
-  std::vector<std::vector< Polynomials::Polynomial< double > > > pols(dim);
-  pols[0] = Polynomials::LagrangeEquidistant::generate_complete_basis(k+2);
-
-  if (k == 0)
-    for (unsigned int d=1; d<dim; ++d)
-      pols[d] = Polynomials::Legendre::generate_complete_basis(0);
-  else
-    for (unsigned int d=1; d<dim; ++d)
-      pols[d] = Polynomials::LagrangeEquidistant::generate_complete_basis(k);
-
-  polynomial_space = new AnisotropicPolynomials<dim>(pols);
-
   // check that the dimensions match. we only store one of the 'dim'
   // anisotropic polynomials that make up the vector-valued space, so
   // multiply by 'dim'
-  Assert (dim * polynomial_space->n() == compute_n_pols(k),
+  Assert (dim * polynomial_space.n() == compute_n_pols(k),
           ExcInternalError());
 }
 
 
-template <int dim>
-PolynomialsABF<dim>::~PolynomialsABF ()
-{
-  delete polynomial_space;
-}
-
 
 template <int dim>
 void
@@ -76,7 +80,7 @@ PolynomialsABF<dim>::compute (const Point<dim>            &unit_point,
   Assert(fourth_derivatives.size()==n_pols|| fourth_derivatives.size()==0,
          ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
 
-  const unsigned int n_sub = polynomial_space->n();
+  const unsigned int n_sub = polynomial_space.n();
   // guard access to the scratch
   // arrays in the following block
   // using a mutex to make sure they
@@ -107,8 +111,8 @@ PolynomialsABF<dim>::compute (const Point<dim>            &unit_point,
       for (unsigned int c=0; c<dim; ++c)
         p(c) = unit_point((c+d)%dim);
 
-      polynomial_space->compute (p, p_values, p_grads, p_grad_grads,
-                                 p_third_derivatives, p_fourth_derivatives);
+      polynomial_space.compute (p, p_values, p_grads, p_grad_grads,
+                                p_third_derivatives, p_fourth_derivatives);
 
       for (unsigned int i=0; i<p_values.size(); ++i)
         values[i+d*n_sub][d] = p_values[i];