]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Polynomial is a template now
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 29 Jan 2001 06:17:53 +0000 (06:17 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 29 Jan 2001 06:17:53 +0000 (06:17 +0000)
git-svn-id: https://svn.dealii.org/trunk@3813 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/polynomial.h
deal.II/base/include/base/tensor_product_polynomials.h
deal.II/base/source/polynomial.cc
deal.II/base/source/tensor_product_polynomials.cc

index 445f31df5cb0f648ba35624277361efd0285c179..c0505d035db69c369b0127a465da70906d3852cf 100644 (file)
@@ -30,6 +30,7 @@
  *
  * @author Ralf Hartmann, Guido Kanschat, 2000
  */
+template <typename number>
 class Polynomial : public Subscriptor
 {
   public:
@@ -48,7 +49,7 @@ class Polynomial : public Subscriptor
                                      * the @p{coefficient} array
                                      * minus one.
                                      */
-    Polynomial (const std::vector<double> &coefficients);
+    Polynomial (const std::vector<number> &coefficients);
 
                                      /**
                                      * Default-Constructor.
@@ -63,7 +64,7 @@ class Polynomial : public Subscriptor
                                      * scheme for numerical stability
                                      * of the evaluation.
                                      */
-    double value (const double x) const;
+    number value (const number x) const;
     
                                     /**
                                      * Return the values and the
@@ -81,8 +82,8 @@ class Polynomial : public Subscriptor
                                      * 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
@@ -104,7 +105,7 @@ class Polynomial : public Subscriptor
                                      * passed down by derived
                                      * classes.
                                      */
-    std::vector<double> coefficients;
+    std::vector<number> coefficients;
 };
 
 
@@ -123,7 +124,7 @@ class Polynomial : public Subscriptor
  *
  * @author Ralf Hartmann, 2000
  */
-class LagrangeEquidistant: public Polynomial
+class LagrangeEquidistant: public Polynomial<double>
 {
   public:
                                     /**
@@ -141,15 +142,7 @@ class LagrangeEquidistant: public Polynomial
                                      /**
                                      * Default-constructor.
                                      */
-    LagrangeEquidistant ();
-    
-
-                                    /**
-                                     * Exception
-                                     */
-    DeclException1 (ExcInvalidSupportPoint,
-                   int,
-                   << "The support point " << arg1 << " is invalid.");
+    LagrangeEquidistant ();    
 
   private:
 
@@ -168,5 +161,46 @@ class LagrangeEquidistant: public Polynomial
 };
 
 
+/**
+ * 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
+
index 8827888568cdbd3cb4cf29a6a05fded714978ceb..2f7283f8e98da43dd3172a62f068f4b40fd2d809 100644 (file)
@@ -43,7 +43,7 @@ class TensorProductPolynomials
                                      * 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
@@ -81,7 +81,7 @@ class TensorProductPolynomials
                                      * Pointer to the @p{polynomials}
                                      * given to the constructor.
                                      */
-    std::vector<SmartPointer<Polynomial> > polynomials;
+    std::vector<SmartPointer<Polynomial<double> > > polynomials;
 
                                     /**
                                      * Number of tensor product
index 9dd7a5bd358c632ef28ac10709f464813d65046c..a85408df96408b8dcd7b033c68d24c021a7fc8e8 100644 (file)
 
 #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];
 
@@ -43,8 +46,10 @@ double Polynomial::value (const double x) const
 
 
 
-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());
@@ -68,7 +73,7 @@ void Polynomial::value (const double         x,
                                   // 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)
     {      
@@ -87,7 +92,7 @@ void Polynomial::value (const double         x,
 
 LagrangeEquidistant::LagrangeEquidistant (const unsigned int n,
                                          const unsigned int support_point):
-               Polynomial(compute_coefficients(n,support_point))
+               Polynomial<double>(compute_coefficients(n,support_point))
 {}
 
 
@@ -291,3 +296,7 @@ LagrangeEquidistant::compute_coefficients (const unsigned int n,
   
   return a;
 }
+
+template Polynomial<float>;
+template Polynomial<double>;
+template Polynomial<long double>;
index 80b390f3fa85eb54a8162282c07efc13b82bef50..4bb827db8e4e060d80bd11e919ce8164ba6708c0 100644 (file)
@@ -30,7 +30,7 @@ number power(number x, unsigned int y)
 
 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))
 {}

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.