From: Martin Kronbichler Date: Thu, 4 Feb 2010 15:04:13 +0000 (+0000) Subject: Use vector > instead of plain pointers in order to correctly... X-Git-Tag: v8.0.0~6553 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=380e3cf6ee9ff29475c552a38de85657f45bd5d2;p=dealii.git Use vector > instead of plain pointers in order to correctly free memory. git-svn-id: https://svn.dealii.org/trunk@20494 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/polynomial.h b/deal.II/base/include/base/polynomial.h index 8eda982788..81f8b79c47 100644 --- a/deal.II/base/include/base/polynomial.h +++ b/deal.II/base/include/base/polynomial.h @@ -18,6 +18,7 @@ #include #include #include +#include #include @@ -73,13 +74,13 @@ namespace Polynomials * polynomial of degree @p n. */ Polynomial (const unsigned int n); - + /** * Default constructor creating * an illegal object. */ Polynomial (); - + /** * Return the value of this * polynomial at the given point. @@ -89,7 +90,7 @@ namespace Polynomials * of the evaluation. */ number value (const number x) const; - + /** * Return the values and the * derivatives of the @@ -193,12 +194,12 @@ namespace Polynomials * Add a second polynomial. */ Polynomial& operator += (const Polynomial& p); - + /** * Subtract a second polynomial. */ Polynomial& operator -= (const Polynomial& p); - + /** * Print coefficients. */ @@ -226,7 +227,7 @@ namespace Polynomials */ static void multiply (std::vector& coefficients, const number factor); - + /** * Coefficients of the polynomial * $\sum_i a_i x^i$. This vector @@ -278,7 +279,7 @@ namespace Polynomials static std::vector > generate_complete_basis (const unsigned int degree); - + private: /** * Needed by constructor. @@ -286,7 +287,7 @@ namespace Polynomials static std::vector make_vector(unsigned int n, const double coefficient); }; - + /** * Lagrange polynomials with equidistant interpolation points in @@ -340,7 +341,7 @@ namespace Polynomials static std::vector > generate_complete_basis (const unsigned int degree); - + private: /** @@ -351,7 +352,7 @@ namespace Polynomials * called in the * constructor. */ - static + static void compute_coefficients (const unsigned int n, const unsigned int support_point, @@ -380,9 +381,9 @@ namespace Polynomials std::vector > generate_complete_basis (const std::vector >& points); }; - - - + + + /** * Legendre polynomials of arbitrary degree on [0,1]. * @@ -419,30 +420,33 @@ namespace Polynomials static std::vector > generate_complete_basis (const unsigned int degree); - + private: /** * Coefficients for the interval $[0,1]$. */ - static std::vector *> shifted_coefficients; - + static std::vector > > shifted_coefficients; + /** * Vector with already computed - * coefficients. For each degree - * of the polynomial, we keep one - * pointer to the list of - * coefficients; we do so rather - * than keeping a vector of + * coefficients. For each degree of the + * polynomial, we keep one pointer to + * the list of coefficients; we do so + * rather than keeping a vector of * vectors in order to simplify - * programming multithread-safe. + * programming multithread-safe. In + * order to avoid memory leak, we use a + * shared_ptr in order to correctly + * free the memory of the vectors when + * the global destructor is called. */ - static std::vector *> recursive_coefficients; - + static std::vector > > recursive_coefficients; + /** * Compute coefficients recursively. */ static void compute_coefficients (const unsigned int p); - + /** * Get coefficients for * constructor. This way, it can @@ -459,23 +463,23 @@ namespace Polynomials /** * Hierarchical polynomials of arbitrary degree on [0,1]. * - * When Constructing a Hierarchical polynomial of degree p, + * When Constructing a Hierarchical polynomial of degree p, * the coefficients will be computed by a recursion formula. The * coefficients are stored in a static data vector to be available * when needed next time. * - * These hierarchical polynomials are based on those of Demkowicz, Oden, + * These hierarchical polynomials are based on those of Demkowicz, Oden, * Rachowicz, and Hardy (CMAME 77 (1989) 79-112, Sec. 4). The first two - * polynomials are the standard linear shape functions given by + * polynomials are the standard linear shape functions given by * $\phi_{0}(x) = 1 - x$ and $\phi_{1}(x) = x$. For $l \geq 2$ * we use the definitions $\phi_{l}(x) = (2x-1)^l - 1, l = 2,4,6,...$ - * and $\phi_{l}(x) = (2x-1)^l - (2x-1), l = 3,5,7,...$. These satisfy the - * recursion relations $\phi_{l}(x) = (2x-1)\phi_{l-1}, l=3,5,7,...$ and - * $\phi_{l}(x) = (2x-1)\phi_{l-1} + \phi_{2}, l=4,6,8,...$. + * and $\phi_{l}(x) = (2x-1)^l - (2x-1), l = 3,5,7,...$. These satisfy the + * recursion relations $\phi_{l}(x) = (2x-1)\phi_{l-1}, l=3,5,7,...$ and + * $\phi_{l}(x) = (2x-1)\phi_{l-1} + \phi_{2}, l=4,6,8,...$. * - * The degrees of freedom are the values at the vertices and the + * The degrees of freedom are the values at the vertices and the * derivatives at the midpoint. Currently, we do not scale the - * polynomials in any way, although better conditioning of the + * polynomials in any way, although better conditioning of the * element stiffness matrix could possibly be achieved with scaling. * * Calling the constructor with a given index p will generate the @@ -529,7 +533,7 @@ namespace Polynomials static std::vector > generate_complete_basis (const unsigned int degree); - + private: /** * Compute coefficients recursively. @@ -545,22 +549,22 @@ namespace Polynomials */ static const std::vector & get_coefficients (const unsigned int p); - + static std::vector *> recursive_coefficients; - }; + }; } /** @} */ /* -------------------------- inline functions --------------------- */ -namespace Polynomials +namespace Polynomials { template inline - Polynomial::Polynomial () + Polynomial::Polynomial () {} - + template inline unsigned int diff --git a/deal.II/base/source/polynomial.cc b/deal.II/base/source/polynomial.cc index 2d864be150..4c4bf0880c 100644 --- a/deal.II/base/source/polynomial.cc +++ b/deal.II/base/source/polynomial.cc @@ -702,22 +702,12 @@ namespace Polynomials // ------------------ class Legendre --------------- // -//TODO:[?] This class leaks memory, but only at the very end of a program. -// Since it expands the Legendre::coefficients array, the elements -// of this static variable are not destroyed at the end of the program -// run. While this is not a problem (since the returned memory could -// not be used anyway then), it is a little confusing when looking at -// a memory checker such as "purify". Maybe, this should be handled somehow -// to avoid this confusion in future. - // Reserve space for polynomials up to degree 19. Should be sufficient // for the start. - std::vector *> - Legendre::recursive_coefficients(20, - static_cast*>(0)); - std::vector *> - Legendre::shifted_coefficients(20, - static_cast*>(0)); + std::vector > > + Legendre::recursive_coefficients(20); + std::vector > > + Legendre::shifted_coefficients(20); Legendre::Legendre (const unsigned int k) @@ -761,7 +751,7 @@ namespace Polynomials // no, then generate the // respective coefficients { - recursive_coefficients.resize (k+1, 0); + recursive_coefficients.resize (k+1); if (k<=1) { @@ -781,9 +771,15 @@ namespace Polynomials (*c1)[1] = 1.; // now make these arrays - // const - recursive_coefficients[0] = c0; - recursive_coefficients[1] = c1; + // const. use shared_ptr for + // recursive_coefficients because + // that avoids a memory leak that + // would appear if we used plain + // pointers. + recursive_coefficients[0] = + std_cxx1x::shared_ptr >(c0); + recursive_coefficients[1] = + std_cxx1x::shared_ptr >(c1); // Compute polynomials // orthogonal on [0,1] @@ -795,8 +791,8 @@ namespace Polynomials Polynomial::shift (*c1, -1.); Polynomial::scale(*c1, 2.); Polynomial::multiply(*c1, std::sqrt(3.)); - shifted_coefficients[0]=c0; - shifted_coefficients[1]=c1; + shifted_coefficients[0]=std_cxx1x::shared_ptr >(c0); + shifted_coefficients[1]=std_cxx1x::shared_ptr >(c1); } else { @@ -830,14 +826,16 @@ namespace Polynomials // created vector to the // const pointer in the // coefficients array - recursive_coefficients[k] = ck; + recursive_coefficients[k] = + std_cxx1x::shared_ptr >(ck); // and compute the // coefficients for [0,1] ck = new std::vector(*ck); Polynomial::shift (*ck, -1.); Polynomial::scale(*ck, 2.); Polynomial::multiply(*ck, std::sqrt(2.*k+1.)); - shifted_coefficients[k] = ck; + shifted_coefficients[k] = + std_cxx1x::shared_ptr >(ck); }; }; }