From: Wolfgang Bangerth Date: Fri, 12 Jan 2018 18:12:08 +0000 (-0700) Subject: Use automatic memory management, rather than doing it by hand. X-Git-Tag: v9.0.0-rc1~576^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c5e05d3381eda6d04def136fa2c52fd1aaa414a4;p=dealii.git Use automatic memory management, rather than doing it by hand. --- diff --git a/source/fe/fe_raviart_thomas.cc b/source/fe/fe_raviart_thomas.cc index 4e7a0ba180..95c16ec86a 100644 --- a/source/fe/fe_raviart_thomas.cc +++ b/source/fe/fe_raviart_thomas.cc @@ -193,9 +193,8 @@ FE_RaviartThomas::initialize_support_points (const unsigned int deg) if (deg==0) return; - // Create Legendre basis for the - // space D_xi Q_k - std::vector* > polynomials(dim); + // Create Legendre basis for the space D_xi Q_k + std::unique_ptr> polynomials[dim]; for (unsigned int dd=0; dd > > poly(dim); @@ -203,7 +202,7 @@ FE_RaviartThomas::initialize_support_points (const unsigned int deg) poly[d] = Polynomials::Legendre::generate_complete_basis(deg); poly[dd] = Polynomials::Legendre::generate_complete_basis(deg-1); - polynomials[dd] = new AnisotropicPolynomials(poly); + polynomials[dd] = std_cxx14::make_unique>(poly); } interior_weights.reinit(TableIndices<3>(n_interior_points, polynomials[0]->n(), dim)); @@ -217,9 +216,6 @@ FE_RaviartThomas::initialize_support_points (const unsigned int deg) * polynomials[d]->compute_value(i,cell_quadrature.point(k)); } - for (unsigned int d=0; dgeneralized_support_points.size(), ExcInternalError()); } @@ -323,10 +319,9 @@ FE_RaviartThomas::initialize_restriction() if (this->degree == 1) return; - // Create Legendre basis for the - // space D_xi Q_k. Here, we cannot + // Create Legendre basis for the space D_xi Q_k. Here, we cannot // use the shape functions - std::vector* > polynomials(dim); + std::unique_ptr> polynomials[dim]; for (unsigned int dd=0; dd > > poly(dim); @@ -334,7 +329,7 @@ FE_RaviartThomas::initialize_restriction() poly[d] = Polynomials::Legendre::generate_complete_basis(this->degree-1); poly[dd] = Polynomials::Legendre::generate_complete_basis(this->degree-2); - polynomials[dd] = new AnisotropicPolynomials(poly); + polynomials[dd] = std_cxx14::make_unique>(poly); } QGauss q_cell(this->degree); @@ -366,9 +361,6 @@ FE_RaviartThomas::initialize_restriction() * polynomials[d]->compute_value(i_weight, q_sub.point(k)); } } - - for (unsigned int d=0; d::initialize_support_points (const unsigned int deg) for (unsigned int d=0; d *quadrature; - if (dim == 1) quadrature = new QAnisotropic(high); - if (dim == 2) quadrature = new QAnisotropic(((d==0) ? low : high), - ((d==1) ? low : high)); - if (dim == 3) quadrature = new QAnisotropic(((d==0) ? low : high), - ((d==1) ? low : high), - ((d==2) ? low : high)); - Assert(dim<=3, ExcNotImplemented()); + std::unique_ptr> quadrature; + switch (dim) + { + case 1: + quadrature = std_cxx14::make_unique>(high); + break; + case 2: + quadrature = std_cxx14::make_unique>(((d==0) ? low : high), + ((d==1) ? low : high)); + break; + case 3: + quadrature = std_cxx14::make_unique>(((d==0) ? low : high), + ((d==1) ? low : high), + ((d==2) ? low : high)); + break; + default: + Assert(false, ExcNotImplemented()); + } for (unsigned int k=0; ksize(); ++k) this->generalized_support_points[current++] = quadrature->point(k); - delete quadrature; } Assert (current == this->dofs_per_cell, ExcInternalError()); }