#include <deal.II/base/exceptions.h>
#include <deal.II/base/point.h>
#include <deal.II/base/polynomial.h>
+#include <deal.II/base/scalar_polynomials_base.h>
#include <deal.II/base/smartpointer.h>
#include <deal.II/base/tensor.h>
* 2005
*/
template <int dim>
-class PolynomialSpace
+class PolynomialSpace : public ScalarPolynomialsBase<dim>
{
public:
/**
* compute_grad_grad() functions, see below, in a loop over all polynomials.
*/
void
- compute(const Point<dim> & unit_point,
- std::vector<double> & values,
- std::vector<Tensor<1, dim>> &grads,
- std::vector<Tensor<2, dim>> &grad_grads,
- std::vector<Tensor<3, dim>> &third_derivatives,
- std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+ evaluate(const Point<dim> & unit_point,
+ std::vector<double> & values,
+ std::vector<Tensor<1, dim>> &grads,
+ std::vector<Tensor<2, dim>> &grad_grads,
+ std::vector<Tensor<3, dim>> &third_derivatives,
+ std::vector<Tensor<4, dim>> &fourth_derivatives) const override;
/**
* Compute the value of the <tt>i</tt>th polynomial at unit point
* given, then the result of this function is <i>N</i> in 1d,
* <i>N(N+1)/2</i> in 2d, and <i>N(N+1)(N+2)/6</i> in 3d.
*/
- unsigned int
- n() const;
+ static unsigned int
+ n_polynomials(const unsigned int n);
/**
- * Degree of the space. This is by definition the number of polynomials
- * given to the constructor, NOT the maximal degree of a polynomial in this
- * vector. The latter value is never checked and therefore left to the
- * application.
+ * Return the name of the space, which is <tt>PolynomialSpace</tt>.
*/
- unsigned int
- degree() const;
+ std::string
+ name() const override;
/**
- * Static function used in the constructor to compute the number of
- * polynomials.
- *
- * @warning The argument `n` is not the maximal degree, but the number of
- * onedimensional polynomials, thus the degree plus one.
+ * @copydoc ScalarPolynomialsBase<dim>::clone()
*/
- static unsigned int
- n_polynomials(const unsigned int n);
+ virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
+ clone() const override;
protected:
/**
*/
const std::vector<Polynomials::Polynomial<double>> polynomials;
- /**
- * Store the precomputed value which the <tt>n()</tt> function returns.
- */
- const unsigned int n_pols;
-
/**
* Index map for reordering the polynomials.
*/
template <int dim>
template <class Pol>
PolynomialSpace<dim>::PolynomialSpace(const std::vector<Pol> &pols)
- : polynomials(pols.begin(), pols.end())
- , n_pols(n_polynomials(polynomials.size()))
- , index_map(n_pols)
- , index_map_inverse(n_pols)
+ : ScalarPolynomialsBase<dim>(pols.size(), n_polynomials(pols.size()))
+ , polynomials(pols.begin(), pols.end())
+ , index_map(n_polynomials(pols.size()))
+ , index_map_inverse(n_polynomials(pols.size()))
{
// per default set this index map
// to identity. This map can be
// changed by the user through the
// set_numbering function
- for (unsigned int i = 0; i < n_pols; ++i)
+ for (unsigned int i = 0; i < this->n(); ++i)
{
index_map[i] = i;
index_map_inverse[i] = i;
}
-template <int dim>
-inline unsigned int
-PolynomialSpace<dim>::n() const
-{
- return n_pols;
-}
-
-
template <int dim>
-inline unsigned int
-PolynomialSpace<dim>::degree() const
+inline std::string
+PolynomialSpace<dim>::name() const
{
- return polynomials.size();
+ return "PolynomialSpace";
}
void
PolynomialSpace<dim>::output_indices(StreamType &out) const
{
- for (unsigned int i = 0; i < n_pols; ++i)
+ for (unsigned int i = 0; i < this->n(); ++i)
{
const std::array<unsigned int, dim> ix = compute_index(i);
out << i << "\t";
#define dealii_polynomials_rannacher_turek_h
#include <deal.II/base/point.h>
+#include <deal.II/base/scalar_polynomials_base.h>
#include <deal.II/base/tensor.h>
#include <vector>
* @date 2015
*/
template <int dim>
-class PolynomialsRannacherTurek
+class PolynomialsRannacherTurek : public ScalarPolynomialsBase<dim>
{
public:
/**
* zero. A size of zero means that we are not computing the vector entries.
*/
void
- compute(const Point<dim> & unit_point,
- std::vector<double> & values,
- std::vector<Tensor<1, dim>> &grads,
- std::vector<Tensor<2, dim>> &grad_grads,
- std::vector<Tensor<3, dim>> &third_derivatives,
- std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+ evaluate(const Point<dim> & unit_point,
+ std::vector<double> & values,
+ std::vector<Tensor<1, dim>> &grads,
+ std::vector<Tensor<2, dim>> &grad_grads,
+ std::vector<Tensor<3, dim>> &third_derivatives,
+ std::vector<Tensor<4, dim>> &fourth_derivatives) const override;
+
+ /**
+ * Return the name of the space, which is <tt>RannacherTurek</tt>.
+ */
+ std::string
+ name() const override;
+
+ /**
+ * @copydoc ScalarPolynomialsBase<dim>::clone()
+ */
+ virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
+ clone() const override;
};
}
+
+template <int dim>
+inline std::string
+PolynomialsRannacherTurek<dim>::name() const
+{
+ return "RannacherTurek";
+}
+
+
DEAL_II_NAMESPACE_CLOSE
#endif
* over all tensor product polynomials.
*/
void
- compute(const Point<dim> & unit_point,
- std::vector<double> & values,
- std::vector<Tensor<1, dim>> &grads,
- std::vector<Tensor<2, dim>> &grad_grads,
- std::vector<Tensor<3, dim>> &third_derivatives,
- std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+ evaluate(const Point<dim> & unit_point,
+ std::vector<double> & values,
+ std::vector<Tensor<1, dim>> &grads,
+ std::vector<Tensor<2, dim>> &grad_grads,
+ std::vector<Tensor<3, dim>> &third_derivatives,
+ std::vector<Tensor<4, dim>> &fourth_derivatives) const;
/**
* Compute the value of the <tt>i</tt>th tensor product polynomial at
* over all tensor product polynomials.
*/
void
- compute(const Point<dim> & unit_point,
- std::vector<double> & values,
- std::vector<Tensor<1, dim>> &grads,
- std::vector<Tensor<2, dim>> &grad_grads,
- std::vector<Tensor<3, dim>> &third_derivatives,
- std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+ evaluate(const Point<dim> & unit_point,
+ std::vector<double> & values,
+ std::vector<Tensor<1, dim>> &grads,
+ std::vector<Tensor<2, dim>> &grad_grads,
+ std::vector<Tensor<3, dim>> &third_derivatives,
+ std::vector<Tensor<4, dim>> &fourth_derivatives) const;
/**
* Compute the value of the <tt>i</tt>th tensor product polynomial at
* over all tensor product polynomials.
*/
void
- compute(const Point<dim> & unit_point,
- std::vector<double> & values,
- std::vector<Tensor<1, dim>> &grads,
- std::vector<Tensor<2, dim>> &grad_grads,
- std::vector<Tensor<3, dim>> &third_derivatives,
- std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+ evaluate(const Point<dim> & unit_point,
+ std::vector<double> & values,
+ std::vector<Tensor<1, dim>> &grads,
+ std::vector<Tensor<2, dim>> &grad_grads,
+ std::vector<Tensor<3, dim>> &third_derivatives,
+ std::vector<Tensor<4, dim>> &fourth_derivatives) const;
/**
* Compute the value of the <tt>i</tt>th tensor product polynomial at
* @code
* static const unsigned int dimension;
*
- * void compute (const Point<dim> &unit_point,
- * std::vector<double> &values,
- * std::vector<Tensor<1,dim> > &grads,
- * std::vector<Tensor<2,dim> > &grad_grads,
- * std::vector<Tensor<3,dim> > &third_derivatives,
- * std::vector<Tensor<4,dim> > &fourth_derivatives) const;
+ * void evaluate (const Point<dim> &unit_point,
+ * std::vector<double> &values,
+ * std::vector<Tensor<1,dim> > &grads,
+ * std::vector<Tensor<2,dim> > &grad_grads,
+ * std::vector<Tensor<3,dim> > &third_derivatives,
+ * std::vector<Tensor<4,dim> > &fourth_derivatives) const;
*
* double compute_value (const unsigned int i,
* const Point<dim> &p) const;
update_3rd_derivatives))
for (unsigned int i = 0; i < n_q_points; ++i)
{
- poly_space.compute(quadrature.point(i),
- values,
- grads,
- grad_grads,
- third_derivatives,
- fourth_derivatives);
+ poly_space.evaluate(quadrature.point(i),
+ values,
+ grads,
+ grad_grads,
+ third_derivatives,
+ fourth_derivatives);
// the values of shape functions at quadrature points don't change.
// consequently, write these values right into the output array if
std::vector<double>(n_q_points));
for (unsigned int i = 0; i < n_q_points; ++i)
{
- poly_space.compute(quadrature.point(i),
- values,
- grads,
- grad_grads,
- empty_vector_of_3rd_order_tensors,
- empty_vector_of_4th_order_tensors);
+ poly_space.evaluate(quadrature.point(i),
+ values,
+ grads,
+ grad_grads,
+ empty_vector_of_3rd_order_tensors,
+ empty_vector_of_4th_order_tensors);
for (unsigned int k = 0; k < poly_space.n(); ++k)
data.shape_values[k][i] = values[k];
#include <deal.II/base/exceptions.h>
#include <deal.II/base/polynomial_space.h>
+#include <deal.II/base/std_cxx14/memory.h>
#include <deal.II/base/table.h>
DEAL_II_NAMESPACE_OPEN
template <int dim>
void
-PolynomialSpace<dim>::compute(
+PolynomialSpace<dim>::evaluate(
const Point<dim> & p,
std::vector<double> & values,
std::vector<Tensor<1, dim>> &grads,
{
const unsigned int n_1d = polynomials.size();
- Assert(values.size() == n_pols || values.size() == 0,
- ExcDimensionMismatch2(values.size(), n_pols, 0));
- Assert(grads.size() == n_pols || grads.size() == 0,
- ExcDimensionMismatch2(grads.size(), n_pols, 0));
- Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
- ExcDimensionMismatch2(grad_grads.size(), n_pols, 0));
- Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
- ExcDimensionMismatch2(third_derivatives.size(), n_pols, 0));
- Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
- ExcDimensionMismatch2(fourth_derivatives.size(), n_pols, 0));
+ Assert(values.size() == this->n() || values.size() == 0,
+ ExcDimensionMismatch2(values.size(), this->n(), 0));
+ Assert(grads.size() == this->n() || grads.size() == 0,
+ ExcDimensionMismatch2(grads.size(), this->n(), 0));
+ Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+ ExcDimensionMismatch2(grad_grads.size(), this->n(), 0));
+ Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+ ExcDimensionMismatch2(third_derivatives.size(), this->n(), 0));
+ Assert(fourth_derivatives.size() == this->n() ||
+ fourth_derivatives.size() == 0,
+ ExcDimensionMismatch2(fourth_derivatives.size(), this->n(), 0));
unsigned int v_size = 0;
bool update_values = false, update_grads = false, update_grad_grads = false;
bool update_3rd_derivatives = false, update_4th_derivatives = false;
- if (values.size() == n_pols)
+ if (values.size() == this->n())
{
update_values = true;
v_size = 1;
}
- if (grads.size() == n_pols)
+ if (grads.size() == this->n())
{
update_grads = true;
v_size = 2;
}
- if (grad_grads.size() == n_pols)
+ if (grad_grads.size() == this->n())
{
update_grad_grads = true;
v_size = 3;
}
- if (third_derivatives.size() == n_pols)
+ if (third_derivatives.size() == this->n())
{
update_3rd_derivatives = true;
v_size = 4;
}
- if (fourth_derivatives.size() == n_pols)
+ if (fourth_derivatives.size() == this->n())
{
update_4th_derivatives = true;
v_size = 5;
}
+
+template <int dim>
+std::unique_ptr<ScalarPolynomialsBase<dim>>
+PolynomialSpace<dim>::clone() const
+{
+ return std_cxx14::make_unique<PolynomialSpace<dim>>(*this);
+}
+
+
template class PolynomialSpace<1>;
template class PolynomialSpace<2>;
template class PolynomialSpace<3>;
// will have first all polynomials
// in the x-component, then y and
// z.
- polynomial_space.compute(unit_point,
- p_values,
- p_grads,
- p_grad_grads,
- p_third_derivatives,
- p_fourth_derivatives);
+ polynomial_space.evaluate(unit_point,
+ p_values,
+ p_grads,
+ p_grad_grads,
+ p_third_derivatives,
+ p_fourth_derivatives);
std::fill(values.begin(), values.end(), Tensor<1, dim>());
for (unsigned int i = 0; i < p_values.size(); ++i)
#include <deal.II/base/geometry_info.h>
#include <deal.II/base/polynomials_rannacher_turek.h>
+#include <deal.II/base/std_cxx14/memory.h>
DEAL_II_NAMESPACE_OPEN
template <int dim>
PolynomialsRannacherTurek<dim>::PolynomialsRannacherTurek()
+ : ScalarPolynomialsBase<dim>(2, dealii::GeometryInfo<dim>::faces_per_cell)
{
Assert(dim == 2, ExcNotImplemented());
}
template <int dim>
void
-PolynomialsRannacherTurek<dim>::compute(
+PolynomialsRannacherTurek<dim>::evaluate(
const Point<dim> & unit_point,
std::vector<double> & values,
std::vector<Tensor<1, dim>> &grads,
std::vector<Tensor<3, dim>> &third_derivatives,
std::vector<Tensor<4, dim>> &fourth_derivatives) const
{
- const unsigned int n_pols = dealii::GeometryInfo<dim>::faces_per_cell;
+ const unsigned int n_pols = this->n();
Assert(values.size() == n_pols || values.size() == 0,
ExcDimensionMismatch(values.size(), n_pols));
Assert(grads.size() == n_pols || grads.size() == 0,
}
+
+template <int dim>
+std::unique_ptr<ScalarPolynomialsBase<dim>>
+PolynomialsRannacherTurek<dim>::clone() const
+{
+ return std_cxx14::make_unique<PolynomialsRannacherTurek<dim>>(*this);
+}
+
+
// explicit instantiations
#include "polynomials_rannacher_turek.inst"
template <int dim, typename PolynomialType>
void
-TensorProductPolynomials<dim, PolynomialType>::compute(
+TensorProductPolynomials<dim, PolynomialType>::evaluate(
const Point<dim> & p,
std::vector<double> & values,
std::vector<Tensor<1, dim>> &grads,
template <int dim>
void
-TensorProductPolynomialsBubbles<dim>::compute(
+TensorProductPolynomialsBubbles<dim>::evaluate(
const Point<dim> & p,
std::vector<double> & values,
std::vector<Tensor<1, dim>> &grads,
do_4th_derivatives = true;
}
- this->TensorProductPolynomials<dim>::compute(
+ this->TensorProductPolynomials<dim>::evaluate(
p, values, grads, grad_grads, third_derivatives, fourth_derivatives);
for (unsigned int i = this->n_tensor_pols;
template <int dim>
void
-TensorProductPolynomialsConst<dim>::compute(
+TensorProductPolynomialsConst<dim>::evaluate(
const Point<dim> & p,
std::vector<double> & values,
std::vector<Tensor<1, dim>> &grads,
do_4th_derivatives = true;
}
- this->TensorProductPolynomials<dim>::compute(
+ this->TensorProductPolynomials<dim>::evaluate(
p, values, grads, grad_grads, third_derivatives, fourth_derivatives);
// for dgq node: values =1, grads=0, grads_grads=0, third_derivatives=0,
for (unsigned int k = 0; k < quadrature.size(); ++k)
{
test_values[k].resize(poly.n());
- poly.compute(quadrature.point(k),
- test_values[k],
- dummy1,
- dummy2,
- dummy3,
- dummy4);
+ poly.evaluate(quadrature.point(k),
+ test_values[k],
+ dummy1,
+ dummy2,
+ dummy3,
+ dummy4);
for (unsigned int i = 0; i < poly.n(); ++i)
{
test_values[k][i] *= quadrature.weight(k);
if (fe_internal.update_each & (update_values | update_gradients))
for (unsigned int i = 0; i < n_q_points; ++i)
{
- polynomial_space.compute(mapping_data.quadrature_points[i],
- values,
- grads,
- grad_grads,
- empty_vector_of_3rd_order_tensors,
- empty_vector_of_4th_order_tensors);
+ polynomial_space.evaluate(mapping_data.quadrature_points[i],
+ values,
+ grads,
+ grad_grads,
+ empty_vector_of_3rd_order_tensors,
+ empty_vector_of_4th_order_tensors);
if (fe_internal.update_each & update_values)
for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
if (fe_internal.update_each & (update_values | update_gradients))
for (unsigned int i = 0; i < n_q_points; ++i)
{
- polynomial_space.compute(mapping_data.quadrature_points[i],
- values,
- grads,
- grad_grads,
- empty_vector_of_3rd_order_tensors,
- empty_vector_of_4th_order_tensors);
+ polynomial_space.evaluate(mapping_data.quadrature_points[i],
+ values,
+ grads,
+ grad_grads,
+ empty_vector_of_3rd_order_tensors,
+ empty_vector_of_4th_order_tensors);
if (fe_internal.update_each & update_values)
for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
if (fe_internal.update_each & (update_values | update_gradients))
for (unsigned int i = 0; i < n_q_points; ++i)
{
- polynomial_space.compute(mapping_data.quadrature_points[i],
- values,
- grads,
- grad_grads,
- empty_vector_of_3rd_order_tensors,
- empty_vector_of_4th_order_tensors);
+ polynomial_space.evaluate(mapping_data.quadrature_points[i],
+ values,
+ grads,
+ grad_grads,
+ empty_vector_of_3rd_order_tensors,
+ empty_vector_of_4th_order_tensors);
if (fe_internal.update_each & update_values)
for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
data.shape_fourth_derivatives.size() != 0)
for (unsigned int point = 0; point < n_points; ++point)
{
- tensor_pols.compute(
+ tensor_pols.evaluate(
unit_points[point], values, grads, grad2, grad3, grad4);
if (data.shape_values.size() != 0)
std::vector<Tensor<3, dim>> third1(n), third2(n);
std::vector<Tensor<4, dim>> fourth1(n), fourth2(n);
- p.compute(x, values1, gradients1, second1, third1, fourth1);
+ p.evaluate(x, values1, gradients1, second1, third1, fourth1);
q.compute(x, values2, gradients2, second2, third2, fourth2);
for (unsigned int k = 0; k < n; ++k)
std::vector<Tensor<3, dim>> third(n);
std::vector<Tensor<4, dim>> fourth(n);
- p.compute(x, values, gradients, second, third, fourth);
+ p.evaluate(x, values, gradients, second, third, fourth);
for (unsigned int k = 0; k < n; ++k)
{