* @author Timo Heister, 2012
*/
template <int dim>
-class TensorProductPolynomialsConst : public TensorProductPolynomials<dim>
+class TensorProductPolynomialsConst
{
public:
/**
template <class Pol>
TensorProductPolynomialsConst(const std::vector<Pol> &pols);
+ /**
+ * Print the list of <tt>tensor_polys</tt> indices to <tt>out</tt>.
+ */
+ void
+ output_indices(std::ostream &out) const;
+
+ /**
+ * Set the ordering of the polynomials. Requires
+ * <tt>renumber.size()==tensor_polys.n()</tt>. Stores a copy of
+ * <tt>renumber</tt>.
+ */
+ void
+ set_numbering(const std::vector<unsigned int> &renumber);
+
+ /**
+ * Give read access to the renumber vector.
+ */
+ const std::vector<unsigned int> &
+ get_numbering() const;
+
+ /**
+ * Give read access to the inverse renumber vector.
+ */
+ const std::vector<unsigned int> &
+ get_numbering_inverse() const;
+
/**
* Compute the value and the first and second derivatives of each tensor
* product polynomial at <tt>unit_point</tt>.
*/
unsigned int
n() const;
+
+private:
+ /**
+ * The TensorProductPolynomials object
+ */
+ TensorProductPolynomials<dim> tensor_polys;
+
+ /**
+ * Index map for reordering the polynomials.
+ */
+ std::vector<unsigned int> index_map;
+
+ /**
+ * Index map for reordering the polynomials.
+ */
+ std::vector<unsigned int> index_map_inverse;
};
/** @} */
template <class Pol>
inline TensorProductPolynomialsConst<dim>::TensorProductPolynomialsConst(
const std::vector<Pol> &pols)
- : TensorProductPolynomials<dim>(pols)
-{
- // append index for renumbering
- this->index_map.push_back(TensorProductPolynomials<dim>::n());
- this->index_map_inverse.push_back(TensorProductPolynomials<dim>::n());
-}
+ : tensor_polys(pols)
+ , index_map(tensor_polys.n() + 1)
+ , index_map_inverse(tensor_polys.n() + 1)
+{}
inline unsigned int
TensorProductPolynomialsConst<dim>::n() const
{
- return TensorProductPolynomials<dim>::n() + 1;
+ return tensor_polys.n() + 1;
+}
+
+
+
+template <int dim>
+inline const std::vector<unsigned int> &
+TensorProductPolynomialsConst<dim>::get_numbering() const
+{
+ return index_map;
+}
+
+
+template <int dim>
+inline const std::vector<unsigned int> &
+TensorProductPolynomialsConst<dim>::get_numbering_inverse() const
+{
+ return index_map_inverse;
}
const unsigned int i,
const Point<dim> & p) const
{
- const unsigned int max_indices = TensorProductPolynomials<dim>::n();
+ const unsigned int max_indices = tensor_polys.n();
Assert(i <= max_indices, ExcInternalError());
// treat the regular basis functions
if (i < max_indices)
- return this
- ->TensorProductPolynomials<dim>::template compute_derivative<order>(i, p);
+ return tensor_polys.template compute_derivative<order>(i, p);
else
// this is for the constant function
return Tensor<order, dim>();
/* ------------------- TensorProductPolynomialsConst -------------- */
+
+template <int dim>
+void
+TensorProductPolynomialsConst<dim>::output_indices(std::ostream &out) const
+{
+ unsigned int ix[dim];
+ for (unsigned int i = 0; i < tensor_polys.n(); ++i)
+ {
+ tensor_polys.compute_index(i, ix);
+ out << i << "\t";
+ for (unsigned int d = 0; d < dim; ++d)
+ out << ix[d] << " ";
+ out << std::endl;
+ }
+}
+
+
+
+template <int dim>
+void
+TensorProductPolynomialsConst<dim>::set_numbering(
+ const std::vector<unsigned int> &renumber)
+{
+ Assert(renumber.size() == index_map.size(),
+ ExcDimensionMismatch(renumber.size(), index_map.size()));
+
+ index_map = renumber;
+ for (unsigned int i = 0; i < index_map.size(); ++i)
+ index_map_inverse[index_map[i]] = i;
+
+ std::vector<unsigned int> renumber_base;
+ for (unsigned int i = 0; i < tensor_polys.n(); ++i)
+ renumber_base.push_back(renumber[i]);
+
+ tensor_polys.set_numbering(renumber_base);
+}
+
+
template <int dim>
double
TensorProductPolynomialsConst<dim>::compute_value(const unsigned int i,
const Point<dim> & p) const
{
- const unsigned int max_indices = TensorProductPolynomials<dim>::n();
+ const unsigned int max_indices = tensor_polys.n();
Assert(i <= max_indices, ExcInternalError());
// treat the regular basis functions
if (i < max_indices)
- return this->TensorProductPolynomials<dim>::compute_value(i, p);
+ return tensor_polys.compute_value(i, p);
else
// this is for the constant function
return 1.;
TensorProductPolynomialsConst<dim>::compute_grad(const unsigned int i,
const Point<dim> & p) const
{
- const unsigned int max_indices = TensorProductPolynomials<dim>::n();
+ const unsigned int max_indices = tensor_polys.n();
Assert(i <= max_indices, ExcInternalError());
// treat the regular basis functions
if (i < max_indices)
- return this->TensorProductPolynomials<dim>::compute_grad(i, p);
+ return tensor_polys.compute_grad(i, p);
else
// this is for the constant function
return Tensor<1, dim>();
TensorProductPolynomialsConst<dim>::compute_grad_grad(const unsigned int i,
const Point<dim> &p) const
{
- const unsigned int max_indices = TensorProductPolynomials<dim>::n();
+ const unsigned int max_indices = tensor_polys.n();
Assert(i <= max_indices, ExcInternalError());
// treat the regular basis functions
if (i < max_indices)
- return this->TensorProductPolynomials<dim>::compute_grad_grad(i, p);
+ return tensor_polys.compute_grad_grad(i, p);
else
// this is for the constant function
return Tensor<2, dim>();
std::vector<Tensor<3, dim>> &third_derivatives,
std::vector<Tensor<4, dim>> &fourth_derivatives) const
{
- Assert(values.size() == TensorProductPolynomials<dim>::n() + 1 ||
- values.size() == 0,
- ExcDimensionMismatch2(values.size(),
- TensorProductPolynomials<dim>::n() + 1,
- 0));
- Assert(grads.size() == TensorProductPolynomials<dim>::n() + 1 ||
- grads.size() == 0,
- ExcDimensionMismatch2(grads.size(),
- TensorProductPolynomials<dim>::n() + 1,
- 0));
- Assert(grad_grads.size() == TensorProductPolynomials<dim>::n() + 1 ||
- grad_grads.size() == 0,
- ExcDimensionMismatch2(grad_grads.size(),
- TensorProductPolynomials<dim>::n() + 1,
- 0));
- Assert(third_derivatives.size() == TensorProductPolynomials<dim>::n() + 1 ||
+ Assert(values.size() == tensor_polys.n() + 1 || values.size() == 0,
+ ExcDimensionMismatch2(values.size(), tensor_polys.n() + 1, 0));
+ Assert(grads.size() == tensor_polys.n() + 1 || grads.size() == 0,
+ ExcDimensionMismatch2(grads.size(), tensor_polys.n() + 1, 0));
+ Assert(grad_grads.size() == tensor_polys.n() + 1 || grad_grads.size() == 0,
+ ExcDimensionMismatch2(grad_grads.size(), tensor_polys.n() + 1, 0));
+ Assert(third_derivatives.size() == tensor_polys.n() + 1 ||
third_derivatives.size() == 0,
ExcDimensionMismatch2(third_derivatives.size(),
- TensorProductPolynomials<dim>::n() + 1,
+ tensor_polys.n() + 1,
0));
- Assert(fourth_derivatives.size() == TensorProductPolynomials<dim>::n() + 1 ||
+ Assert(fourth_derivatives.size() == tensor_polys.n() + 1 ||
fourth_derivatives.size() == 0,
ExcDimensionMismatch2(fourth_derivatives.size(),
- TensorProductPolynomials<dim>::n() + 1,
+ tensor_polys.n() + 1,
0));
// remove slot for const value, go into the base class compute method and
}
if (third_derivatives.empty() == false)
{
- third_derivatives.resize(TensorProductPolynomials<dim>::n());
+ third_derivatives.resize(tensor_polys.n());
do_3rd_derivatives = true;
}
if (fourth_derivatives.empty() == false)
{
- fourth_derivatives.resize(TensorProductPolynomials<dim>::n());
+ fourth_derivatives.resize(tensor_polys.n());
do_4th_derivatives = true;
}
- this->TensorProductPolynomials<dim>::evaluate(
+ tensor_polys.evaluate(
p, values, grads, grad_grads, third_derivatives, fourth_derivatives);
// for dgq node: values =1, grads=0, grads_grads=0, third_derivatives=0,