*
* @author Ralf Hartmann, Guido Kanschat, 2000
*/
+template <typename number>
class Polynomial : public Subscriptor
{
public:
* the @p{coefficient} array
* minus one.
*/
- Polynomial (const std::vector<double> &coefficients);
+ Polynomial (const std::vector<number> &coefficients);
/**
* Default-Constructor.
* scheme for numerical stability
* of the evaluation.
*/
- double value (const double x) const;
+ number value (const number x) const;
/**
* Return the values and the
* scheme for numerical stability
* of the evaluation.
*/
- void value (const double x,
- std::vector<double> &values) const;
+ void value (const number x,
+ std::vector<number> &values) const;
/**
* Exception
* passed down by derived
* classes.
*/
- std::vector<double> coefficients;
+ std::vector<number> coefficients;
};
*
* @author Ralf Hartmann, 2000
*/
-class LagrangeEquidistant: public Polynomial
+class LagrangeEquidistant: public Polynomial<double>
{
public:
/**
/**
* Default-constructor.
*/
- LagrangeEquidistant ();
-
-
- /**
- * Exception
- */
- DeclException1 (ExcInvalidSupportPoint,
- int,
- << "The support point " << arg1 << " is invalid.");
+ LagrangeEquidistant ();
private:
};
+/**
+ * Legendre polynomials of arbitrary order
+ *
+ * Constructing a Legendre polynomial of order @p{k}, the coefficients
+ * will be computed by the three-term recursion formula. The
+ * coefficients are stored in a static data vector to be available
+ * when needed next time.
+ *
+ * @author Guido Kanschat, 2000
+ */
+template <typename number>
+class Legendre : public Polynomial<number>
+{
+public:
+ /**
+ * Constructor for polynomial of
+ * order @p{k}.
+ */
+ Legendre (unsigned int k);
+private:
+ /**
+ * Vector with already computed
+ * coefficients.
+ */
+ static std::vector<vector<number> > coefficients
+
+ /**
+ * Compute coefficients recursively.
+ */
+ static void compute_coeficients (unsigned int k);
+
+ /**
+ * Get coefficients for
+ * constructor. This way, it can
+ * use the non-standard constructor
+ * of @ref{Polynomial}.
+ */
+ static const std::vector<number>& get_coefficients (unsigned int k);
+}
+
#endif
+
* and will be copied into the
* member variable @p{polynomials}.
*/
- TensorProductPolynomials(const std::vector<SmartPointer<Polynomial> > &pols);
+ TensorProductPolynomials(const std::vector<SmartPointer<Polynomial<double> > > &pols);
/**
* Calculates the polynomials
* Pointer to the @p{polynomials}
* given to the constructor.
*/
- std::vector<SmartPointer<Polynomial> > polynomials;
+ std::vector<SmartPointer<Polynomial<double> > > polynomials;
/**
* Number of tensor product
#include <base/polynomial.h>
-
-Polynomial::Polynomial (const std::vector<double> &a):
+template <typename number>
+Polynomial<number>::Polynomial (const std::vector<number> &a):
coefficients(a)
{}
-Polynomial::Polynomial ()
+template <typename number>
+Polynomial<number>::Polynomial ()
:
coefficients(0)
{}
-double Polynomial::value (const double x) const
+template <typename number>
+number
+Polynomial<number>::value (const number x) const
{
Assert (coefficients.size() > 0, ExcVoidPolynomial());
const unsigned int m=coefficients.size();
// Horner scheme
- double value = coefficients.back();
+ number value = coefficients.back();
for (int k=m-2; k>=0; --k)
value = value*x + coefficients[k];
-void Polynomial::value (const double x,
- std::vector<double> &values) const
+template <typename number>
+void
+Polynomial<number>::value (const number x,
+ std::vector<number> &values) const
{
Assert (coefficients.size() > 0, ExcVoidPolynomial());
Assert (values.size() > 0, ExcEmptyArray());
// then do it properly by the
// full Horner scheme
const unsigned int m=coefficients.size();
- std::vector<double> a(coefficients);
+ std::vector<number> a(coefficients);
unsigned int j_faculty=1;
for (unsigned int j=0; j<values_size; ++j)
{
LagrangeEquidistant::LagrangeEquidistant (const unsigned int n,
const unsigned int support_point):
- Polynomial(compute_coefficients(n,support_point))
+ Polynomial<double>(compute_coefficients(n,support_point))
{}
return a;
}
+
+template Polynomial<float>;
+template Polynomial<double>;
+template Polynomial<long double>;
template <int dim>
TensorProductPolynomials<dim>::TensorProductPolynomials(
- const std::vector<SmartPointer<Polynomial> > &pols):
+ const std::vector<SmartPointer<Polynomial<double> > > &pols):
polynomials(pols),
n_tensor_pols(power(polynomials.size(), dim))
{}