* minus one.
*/
Polynomial (const typename std::vector<number> &coefficients);
-
- /**
- * Default-Constructor.
- */
- Polynomial ();
/**
* Return the value of this
* of this class and may be
* passed down by derived
* classes.
+ *
+ * This vector cannot be constant
+ * since we want to allow copying
+ * of polynomials.
*/
typename std::vector<number> coefficients;
};
/**
- * Lagrange polynomials with equidistant interpolation
- * points in [0,1]. The polynomial of degree @p{n} has got @p{n+1} interpolation
+ * Lagrange polynomials with equidistant interpolation points in
+ * [0,1]. The polynomial of degree @p{n} has got @p{n+1} interpolation
* points. The interpolation points are sorted in ascending
* order. This order gives an index to each interpolation point. A
- * Lagrangian polynomial equals to 1 at its `support point',
- * and 0 at all other interpolation
- * points. For example, if the degree is 3, and the support point is 1,
- * then the polynomial represented by this object is cubic and its
- * value is 1 at the point @p{x=1/3}, and zero at the point @p{x=0},
- * @p{x=2/3}, and @p{x=1}.
+ * Lagrangian polynomial equals to 1 at its `support point', and 0 at
+ * all other interpolation points. For example, if the degree is 3,
+ * and the support point is 1, then the polynomial represented by this
+ * object is cubic and its value is 1 at the point @p{x=1/3}, and zero
+ * at the point @p{x=0}, @p{x=2/3}, and @p{x=1}. All the polynomials
+ * have polynomial order equal to @p{degree}, but together they span
+ * the entire space of polynomials of degree less than or equal
+ * @p{degree}.
*
* The Lagrange polynomials are implemented up to degree 10.
*
LagrangeEquidistant (const unsigned int n,
const unsigned int support_point);
- /**
- * Default-constructor.
+ /**
+ * Return a vector of polynomial
+ * objects of order @p{degree},
+ * which then spans the full
+ * space of polynomials up to the
+ * given degree. The polynomials
+ * are generated by calling the
+ * destructor of this class with
+ * the same degree but support
+ * point running from zero to
+ * @p{degree}. This function may
+ * be used to initialize the
+ * @ref{TensorProductPolynomials}
+ * and @ref{PolynomialSpace}
+ * classes.
*/
- LagrangeEquidistant ();
-
+ static
+ std::vector<Polynomial<double> >
+ generate_complete_basis (const unsigned int degree);
+
private:
/**
*/
Legendre (const unsigned int k);
+ /**
+ * Return a vector of Legendre
+ * polynomial objects of orders
+ * zero through @p{degree}, which
+ * then spans the full space of
+ * polynomials up to the given
+ * degree. This function may be
+ * used to initialize the
+ * @ref{TensorProductPolynomials}
+ * and @ref{PolynomialSpace}
+ * classes.
+ */
+ static
+ std::vector<Polynomial<number> >
+ generate_complete_basis (const unsigned int degree);
+
private:
/**
* Vector with already computed
};
+/* -------------------------- inline functions --------------------- */
+
template <typename number>
inline
unsigned int
/**
- * Polynomial space of degree at most n in higher dimensions.
+ * Representation of the space of polynomials of degree at most n in
+ * higher dimensions.
*
* Given a vector of @{n} one-dimensional polynomials @{P0} to @{Pn},
* where @{Pi} has degree @p{i}, this class generates all polynomials
- * the form @p{ Pijk(x,y,z) = Pi(x)Pj(y)Pk(z)}, where the sum of
- * @p{i}, @p{j} and @p{k} is below/equal @p{n}.
+ * of the form @p{ Pijk(x,y,z) = Pi(x)Pj(y)Pk(z)}, where the sum of
+ * @p{i}, @p{j} and @p{k} is less than or equal @p{n}.
*
* @author Guido Kanschat, 2002
*/
* vector of pointers to
* one-dimensional polynomials
* and will be copied into the
- * member variable @p{polynomials}.
+ * member variable
+ * @p{polynomials}. The static
+ * type of the template argument
+ * @p{pols} needs to be
+ * convertible to
+ * @p{Polynomial<double>},
+ * i.e. should usually be a
+ * derived class of
+ * @p{Polynomial<double>}.
*/
template <class Pol>
PolynomialSpace(const typename std::vector<Pol> &pols);
* functions, see below, in a
* loop over all polynomials.
*/
- void compute(const Point<dim> &unit_point,
- std::vector<double> &values,
- typename std::vector<Tensor<1,dim> > &grads,
- typename std::vector<Tensor<2,dim> > &grad_grads) const;
+ void compute (const Point<dim> &unit_point,
+ std::vector<double> &values,
+ typename std::vector<Tensor<1,dim> > &grads,
+ typename std::vector<Tensor<2,dim> > &grad_grads) const;
/**
* Computes the value of the
const Point<dim> &p) const;
/**
- * Returns the number of tensor
- * product polynomials. For $n$
- * 1d polynomials this is $n^dim$.
+ * Return the number of
+ * polynomials spanning the space
+ * represented by this
+ * class. Here, if @p{N} is the
+ * number of one-dimensional
+ * polynomials given, then the
+ * result of this function is
+ * @p{N} in 1d, @p{N(N+1)/2} in
+ * 2d, and @p{N(N+1)(N+2)/6 in
+ * 3d.
*/
unsigned int n() const;
* polynomials given to the
* constructor.
*/
- std::vector<Polynomial<double> > polynomials;
+ const std::vector<Polynomial<double> > polynomials;
/**
- * Number of tensor product
- * polynomials. For $n$ 1d
- * polynomials this is $n^dim$.
+ * Store the precomputed value
+ * which the @p{n()} function
+ * returns.
*/
- unsigned int n_pols;
-
+ const unsigned int n_pols;
+
/**
- * Computes @p{x} to the power of
- * @p{y} for unsigned int @p{x}
- * and @p{y}. It is a private
- * function as it is only used in
- * this class.
+ * Static function used in the
+ * constructor to compute the
+ * number of polynomials.
*/
- static unsigned int power(const unsigned int x, const unsigned int y);
+ static unsigned int compute_n_pols (const unsigned int n);
};
template <int dim>
template <class Pol>
-PolynomialSpace<dim>::PolynomialSpace(
- const typename std::vector<Pol> &pols):
- polynomials (pols.begin(), pols.end())
-{
- const unsigned int n=polynomials.size();
-
- n_pols = n;
- for (unsigned int i=1;i<dim;++i)
- {
- n_pols *= (n+i);
- n_pols /= (i+1);
- }
-}
+PolynomialSpace<dim>::
+PolynomialSpace (const typename std::vector<Pol> &pols)
+ :
+ polynomials (pols.begin(), pols.end()),
+ n_pols (compute_n_pols(polynomials.size()))
+{}
* polynomials given to the
* constructor.
*/
- std::vector<Polynomial<double> > polynomials;
+ const std::vector<Polynomial<double> > polynomials;
/**
* Number of tensor product
template <int dim>
template <class Pol>
-TensorProductPolynomials<dim>::TensorProductPolynomials(
- const typename std::vector<Pol> &pols):
+TensorProductPolynomials<dim>::
+TensorProductPolynomials(const typename std::vector<Pol> &pols)
+ :
polynomials (pols.begin(), pols.end()),
n_tensor_pols(power(pols.size(), dim)),
- n_pols_to(dim+1)
+ n_pols_to(dim+1, 0)
{
const unsigned int n_pols=polynomials.size();
+template <typename number>
+std::vector<Polynomial<number> >
+Legendre<number>::generate_complete_basis (const unsigned int degree)
+{
+ std::vector<Polynomial<double> > v;
+ v.reserve(degree+1);
+ for (unsigned int i=0; i<=degree; ++i)
+ v.push_back (Legendre<double>(i));
+ return v;
+};
+
+
+
// explicit instantiations
template class Legendre<double>;
-template <typename number>
-Polynomial<number>::Polynomial ()
- :
- coefficients(0)
-{}
-
-
-
template <typename number>
number
Polynomial<number>::value (const number x) const
-LagrangeEquidistant::LagrangeEquidistant ()
-{}
-
-
-
std::vector<double>
LagrangeEquidistant::compute_coefficients (const unsigned int n,
const unsigned int support_point)
return a;
}
+
+std::vector<Polynomial<double> >
+LagrangeEquidistant::
+generate_complete_basis (const unsigned int degree)
+{
+ if (degree==0)
+ // create constant polynomial
+ return std::vector<Polynomial<double> >
+ (1, Polynomial<double> (std::vector<double> (1,1.)));
+ else
+ {
+ // create array of Lagrange
+ // polynomials
+ std::vector<Polynomial<double> > v;
+ for (unsigned int i=0; i<=degree; ++i)
+ v.push_back(LagrangeEquidistant(degree,i));
+ return v;
+ };
+};
+
+
+
template class Polynomial<float>;
template class Polynomial<double>;
template class Polynomial<long double>;
#include <base/polynomial_space.h>
-
template <int dim>
-unsigned int PolynomialSpace<dim>::power(const unsigned int x,
- const unsigned int y)
+unsigned int
+PolynomialSpace<dim>::compute_n_pols (const unsigned int n)
{
- unsigned int value=1;
- for (unsigned int i=0; i<y; ++i)
- value*=x;
- return value;
-}
+ unsigned int n_pols = n;
+ for (unsigned int i=1;i<dim;++i)
+ {
+ n_pols *= (n+i);
+ n_pols /= (i+1);
+ }
+ return n_pols;
+};
<h3>base</h3>
<ol>
+ <li> <p>
+ New: The <code class="class">Legendre</code> and
+ <code class="class">LagrangeEquidistant</code> classes now have
+ static member functions <code
+ class="member">generate_complete_basis<code> which returns an
+ array of polynomial objects spanning the complete space up to a
+ specified order in 1d. This may be used to generate the
+ respective polynomial spaces in higher space dimensions.
+ <br>
+ (WB 2002/05/27)
+ </p>
+
+ <li> <p>
+ Changed: The <code class="class">Polynomial</code> and
+ <code class="class">LagrangeEquidistant</code> classes have lost
+ their default constructor, as that did not make much sense
+ anyway.
+ <br>
+ (WB 2002/05/27)
+ </p>
+
<li> <p>
Fixed: When forward declaring the <code
class="class">Tensor</code> class, we now also forward declare
// $Id$
// Version: $Name$
//
-// Copyright (C) 2000, 2001 by the deal.II authors
+// Copyright (C) 2000, 2001, 2002 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
deallog.attach(logfile);
deallog.depth_console(0);
- std::vector<Polynomial<double> > p (15);
+ std::vector<Polynomial<double> > p;
deallog << "Legendre" << endl;
- for (unsigned int i=0;i<p.size();++i)
- p[i] = Legendre<double>(i);
+ for (unsigned int i=0;i<15;++i)
+ p.push_back (Legendre<double>(i));
for (unsigned int i=0;i<p.size();++i)
for (unsigned int j=0;j<=i;++j)
deallog << "LagrangeEquidistant" << endl;
- p.resize(6);
- for (unsigned int i=0;i<p.size();++i)
- p[i] = LagrangeEquidistant(p.size(), i);
+ p.clear();
+ for (unsigned int i=0;i<6;++i)
+ p.push_back(LagrangeEquidistant(6, i));
// We add 1.0001 bacuse of bugs in
// the ostream classes
deallog.attach(logfile);
deallog.depth_console(0);
- vector<Polynomial<double> > p(3);
- for (unsigned int i=0;i<p.size();++i)
- p[i] = LagrangeEquidistant(p.size(), i);
+ vector<Polynomial<double> > p;
+ for (unsigned int i=0;i<3;++i)
+ p.push_back (LagrangeEquidistant(3, i));
check_dimensions(p);
- for (unsigned int i=0;i<p.size();++i)
- p[i] = Legendre<double>(i);
+ p.clear ();
+ for (unsigned int i=0;i<3;++i)
+ p.push_back (Legendre<double>(i));
check_dimensions(p);
}