]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implementation of a complete polynomial space of order k based on monomials. Makes...
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 9 Mar 2004 14:45:44 +0000 (14:45 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 9 Mar 2004 14:45:44 +0000 (14:45 +0000)
git-svn-id: https://svn.dealii.org/trunk@8699 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/polynomials_p.h [new file with mode: 0644]
deal.II/base/source/polynomials_p.cc [new file with mode: 0644]

diff --git a/deal.II/base/include/base/polynomials_p.h b/deal.II/base/include/base/polynomials_p.h
new file mode 100644 (file)
index 0000000..ee11054
--- /dev/null
@@ -0,0 +1,113 @@
+//-------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-------------------------------------------------------------------
+#ifndef __deal2__polynomials_P_h
+#define __deal2__polynomials_P_h
+
+
+#include <base/config.h>
+#include <base/exceptions.h>
+#include <base/tensor.h>
+#include <base/point.h>
+#include <base/polynomial.h>
+#include <base/polynomial_space.h>
+#include <base/table.h>
+
+#include <vector>
+
+
+/**
+ * @brief The complete polynomial space of order <tt>k</tt> based on
+ * the monomials.
+ *
+ * This class implements the polynomial space of order <tt>k</tt>
+ * based on the monomials ${1,x,x^2,...}$. I.e. in <tt>d</tt>
+ * dimensions it constructs all polynomials of the form $\prod_{i=1}^d
+ * x_i^{n_i}$, where $\sum_i n_i\leq k$. The base polynomials are
+ * given a specific ordering, e.g. in 2 dimensions:
+ * ${1,x,y,xy,x^2,y^2,x^2y,xy^2,x^3,y^3,...}$. The ordering of the
+ * monomials in $P_k1$ matches the ordering of the monomials in $P_k2$
+ * for $k2>k1$.
+ *
+ * @author Ralf Hartmann, 2004
+ */
+template <int dim>
+class PolynomialsP: public PolynomialSpace<dim>
+{
+  public:
+                                    /**
+                                     * Constructor. Creates all basis
+                                     * functions of $P_k$.
+                                     * @arg k: the degree of the
+                                     * polynomial space
+                                     */
+    PolynomialsP (const unsigned int k);
+
+                                    /**
+                                     * Returns the degree <tt>k</tt>
+                                     * of the polynomial space
+                                     * <tt>P_k</tt>.
+                                     *
+                                     * Note, that this number is
+                                     * <tt>PolynomialSpace::degree()-1</tt>,
+                                     * compare definition in
+                                     * @ref{PolynomialSpace}.
+                                     */
+    unsigned int degree() const;
+
+                                    /**
+                                     * For the <tt>n</tt>th
+                                     * polynomial $p_n(x,y,z)=x^i y^j
+                                     * z^k$ this function gives the
+                                     * degrees i,j,k in the x,y,z
+                                     * directions.
+                                     */
+    void directional_degrees(unsigned int n,
+                            unsigned int (&degrees)[dim]) const;
+    
+  private:
+
+                                    /**
+                                     * Fills the <tt>index_map</tt>.
+                                     */
+    void create_polynomial_ordering(vector<unsigned int> &index_map) const;
+
+                                    /**
+                                     * Degree <tt>k<tt> of the
+                                     * polynomial space $P_k$,
+                                     * i.e. the number <tt>k<tt>
+                                     * which was given to the
+                                     * constructor.
+                                     */
+    const unsigned int k;
+};
+
+
+
+template <int dim>
+inline unsigned int
+PolynomialsP<dim>::degree() const
+{
+  return k;
+}
+
+
+template <int dim>
+inline void
+PolynomialsP<dim>::directional_degrees(unsigned int n,
+                                      unsigned int (&degrees)[dim]) const
+{
+  compute_index(n,degrees);
+}
+
+
+#endif
diff --git a/deal.II/base/source/polynomials_p.cc b/deal.II/base/source/polynomials_p.cc
new file mode 100644 (file)
index 0000000..e2ce678
--- /dev/null
@@ -0,0 +1,98 @@
+//-------------------------------------------------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//-------------------------------------------------------------------
+
+#include <base/polynomials_p.h>
+
+
+template <int dim>
+PolynomialsP<dim>::PolynomialsP (const unsigned int k)
+               :
+               PolynomialSpace<dim>(Polynomials::Monomial<double>::generate_complete_basis(k)),
+               k(k)
+{
+  vector<unsigned int> index_map(n());
+  create_polynomial_ordering(index_map);
+  set_polynomial_ordering(index_map);
+}
+
+
+template <>
+void PolynomialsP<1>::create_polynomial_ordering(
+  vector<unsigned int> &index_map) const
+{
+  Assert(index_map.size()==n(), ExcDimensionMismatch(index_map.size(), n()));
+
+                                  // identity
+  for (unsigned int i=0; i<n(); ++i)
+    index_map[i]=i;
+}
+
+
+
+const unsigned int imap2[6][21]=
+{{0},
+ {0,1,2},
+ {0,1,3,4,2,5},
+ {0,1,4,5,2,7,6,8,3,9},
+ {0,1,5,6,2,9,7,10,3,12,11,8,13,4,14},
+ {0,1,6,7,2,11,8,12,3,15,13,9,16,4,18,14,17,10,19,5,20}
+};
+
+
+template <>
+void PolynomialsP<2>::create_polynomial_ordering(
+  vector<unsigned int> &index_map) const
+{
+  Assert(index_map.size()==n(), ExcDimensionMismatch(index_map.size(), n()));
+  Assert(k<=5, ExcNotImplemented());
+  
+                                  // Given the number i of the
+                                  // polynomial in
+                                  // $1,x,y,xy,x2,y2,...$,
+                                  // index_map[i] gives the number of
+                                  // the polynomial in
+                                  // PolynomialSpace.
+  for (unsigned int i=0; i<n(); ++i)
+    index_map[i]=imap2[k][i];
+}
+
+
+const unsigned int imap3[4][20]=
+{{0},
+ {0,1,2,3},
+ {0,1,3,6,4,7,8,2,5,9},
+ {0,1,4,10,5,11,13,2,7,16,14,6,12,8,15,17,18,3,9,19}
+};
+
+template <>
+void PolynomialsP<3>::create_polynomial_ordering(
+  vector<unsigned int> &index_map) const
+{
+  Assert(index_map.size()==n(), ExcDimensionMismatch(index_map.size(), n()));
+  Assert(k<=3, ExcNotImplemented());
+  
+                                  // Given the number i of the
+                                  // polynomial in
+                                  // $1,x,y,xy,x2,y2,...$,
+                                  // index_map[i] gives the number of
+                                  // the polynomial in
+                                  // PolynomialSpace.
+  for (unsigned int i=0; i<n(); ++i)
+    index_map[i]=imap3[k][i];
+}
+
+
+
+template class PolynomialsP<1>;
+template class PolynomialsP<2>;
+template class PolynomialsP<3>;

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.