]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Change Polynomials classes to derive from TensorPolynomialsBase 8588/head
authorgrahambenharper <grahambenharper@gmail.com>
Fri, 16 Aug 2019 07:37:23 +0000 (01:37 -0600)
committergrahambenharper <grahambenharper@gmail.com>
Sat, 17 Aug 2019 21:09:41 +0000 (15:09 -0600)
12 files changed:
include/deal.II/base/polynomials_abf.h
include/deal.II/base/polynomials_bdm.h
include/deal.II/base/polynomials_bernardi_raugel.h
include/deal.II/base/polynomials_nedelec.h
include/deal.II/base/polynomials_raviart_thomas.h
include/deal.II/base/polynomials_rt_bubbles.h
source/base/polynomials_abf.cc
source/base/polynomials_bdm.cc
source/base/polynomials_bernardi_raugel.cc
source/base/polynomials_nedelec.cc
source/base/polynomials_raviart_thomas.cc
source/base/polynomials_rt_bubbles.cc

index 1be8c69decfeac53a9dcbec27c1d5079009df29f..44ef70688b208841c1b5c82f6a5b12ee6ab75953 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/polynomial_space.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_polynomials_base.h>
 #include <deal.II/base/tensor_product_polynomials.h>
 #include <deal.II/base/thread_management.h>
 
@@ -51,7 +52,7 @@ DEAL_II_NAMESPACE_OPEN
  * @date 2006
  */
 template <int dim>
-class PolynomialsABF
+class PolynomialsABF : public TensorPolynomialsBase<dim>
 {
 public:
   /**
@@ -82,26 +83,13 @@ public:
           std::vector<Tensor<2, dim>> &grads,
           std::vector<Tensor<3, dim>> &grad_grads,
           std::vector<Tensor<4, dim>> &third_derivatives,
-          std::vector<Tensor<5, dim>> &fourth_derivatives) const;
-
-  /**
-   * Return the number of ABF polynomials.
-   */
-  unsigned int
-  n() const;
-
-  /**
-   * Return the degree of the ABF space, which is two less than the highest
-   * polynomial degree.
-   */
-  unsigned int
-  degree() const;
+          std::vector<Tensor<5, dim>> &fourth_derivatives) const override;
 
   /**
    * Return the name of the space, which is <tt>ABF</tt>.
    */
   std::string
-  name() const;
+  name() const override;
 
   /**
    * Return the number of polynomials in the space <tt>RT(degree)</tt> without
@@ -111,12 +99,13 @@ public:
   static unsigned int
   compute_n_pols(unsigned int degree);
 
-private:
   /**
-   * The degree of this object as given to the constructor.
+   * @copydoc TensorPolynomialsBase<dim>::clone()
    */
-  const unsigned int my_degree;
+  virtual std::unique_ptr<TensorPolynomialsBase<dim>>
+  clone() const override;
 
+private:
   /**
    * An object representing the polynomial space for a single component. We
    * can re-use it for the other vector components by rotating the
@@ -124,11 +113,6 @@ private:
    */
   const AnisotropicPolynomials<dim> polynomial_space;
 
-  /**
-   * Number of Raviart-Thomas polynomials.
-   */
-  unsigned int n_pols;
-
   /**
    * A mutex that guards the following scratch arrays.
    */
@@ -161,22 +145,6 @@ private:
 };
 
 
-template <int dim>
-inline unsigned int
-PolynomialsABF<dim>::n() const
-{
-  return n_pols;
-}
-
-
-template <int dim>
-inline unsigned int
-PolynomialsABF<dim>::degree() const
-{
-  return my_degree;
-}
-
-
 template <int dim>
 inline std::string
 PolynomialsABF<dim>::name() const
index 45d962ae5b2b6dbe3c524adf05cd9c36aa9b7a20..be37ecc584e44e9115c9a46f97c038054abeb7a6 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/polynomial_space.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_polynomials_base.h>
 #include <deal.II/base/thread_management.h>
 
 #include <vector>
@@ -97,7 +98,7 @@ DEAL_II_NAMESPACE_OPEN
  * @date 2003, 2005, 2009
  */
 template <int dim>
-class PolynomialsBDM
+class PolynomialsBDM : public TensorPolynomialsBase<dim>
 {
 public:
   /**
@@ -128,13 +129,7 @@ public:
           std::vector<Tensor<2, dim>> &grads,
           std::vector<Tensor<3, dim>> &grad_grads,
           std::vector<Tensor<4, dim>> &third_derivatives,
-          std::vector<Tensor<5, dim>> &fourth_derivatives) const;
-
-  /**
-   * Return the number of BDM polynomials.
-   */
-  unsigned int
-  n() const;
+          std::vector<Tensor<5, dim>> &fourth_derivatives) const override;
 
   /**
    * Return the degree of the BDM space, which is one less than the highest
@@ -147,7 +142,7 @@ public:
    * Return the name of the space, which is <tt>BDM</tt>.
    */
   std::string
-  name() const;
+  name() const override;
 
   /**
    * Return the number of polynomials in the space <tt>BDM(degree)</tt>
@@ -157,6 +152,12 @@ public:
   static unsigned int
   compute_n_pols(unsigned int degree);
 
+  /**
+   * @copydoc TensorPolynomialsBase<dim>::clone()
+   */
+  virtual std::unique_ptr<TensorPolynomialsBase<dim>>
+  clone() const override;
+
 private:
   /**
    * An object representing the polynomial space used here. The constructor
@@ -170,11 +171,6 @@ private:
    */
   std::vector<Polynomials::Polynomial<double>> monomials;
 
-  /**
-   * Number of BDM polynomials.
-   */
-  unsigned int n_pols;
-
   /**
    * A mutex that guards the following scratch arrays.
    */
@@ -207,14 +203,6 @@ private:
 };
 
 
-template <int dim>
-inline unsigned int
-PolynomialsBDM<dim>::n() const
-{
-  return n_pols;
-}
-
-
 template <int dim>
 inline unsigned int
 PolynomialsBDM<dim>::degree() const
index eb37c351061faed576e460a0c7e444de38955744..3c4c605256a32a720a2100210f8ca5e7b7609738 100644 (file)
@@ -22,6 +22,7 @@
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/polynomial_space.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_polynomials_base.h>
 #include <deal.II/base/tensor_product_polynomials.h>
 
 #include <vector>
@@ -81,7 +82,7 @@ DEAL_II_NAMESPACE_OPEN
  * @date 2018
  */
 template <int dim>
-class PolynomialsBernardiRaugel
+class PolynomialsBernardiRaugel : public TensorPolynomialsBase<dim>
 {
 public:
   /**
@@ -93,26 +94,11 @@ public:
    */
   PolynomialsBernardiRaugel(const unsigned int k);
 
-  /**
-   * Return the number of Bernardi-Raugel polynomials.
-   */
-  unsigned int
-  n() const;
-
-
-  /**
-   * Return the degree of Bernardi-Raugel polynomials.
-   * Since the bubble functions are quadratic in at least one variable,
-   * the degree of the Bernardi-Raugel polynomials is two.
-   */
-  unsigned int
-  degree() const;
-
   /**
    * Return the name of the space, which is <tt>BernardiRaugel</tt>.
    */
   std::string
-  name() const;
+  name() const override;
 
   /**
    * Compute the value and derivatives of each Bernardi-Raugel
@@ -132,7 +118,7 @@ public:
           std::vector<Tensor<2, dim>> &grads,
           std::vector<Tensor<3, dim>> &grad_grads,
           std::vector<Tensor<4, dim>> &third_derivatives,
-          std::vector<Tensor<5, dim>> &fourth_derivatives) const;
+          std::vector<Tensor<5, dim>> &fourth_derivatives) const override;
 
   /**
    * Return the number of polynomials in the space <tt>BR(degree)</tt> without
@@ -142,17 +128,13 @@ public:
   static unsigned int
   compute_n_pols(const unsigned int k);
 
-private:
-  /**
-   * The degree of this object given to the constructor (must be 1).
-   */
-  const unsigned int my_degree;
-
   /**
-   * The number of Bernardi-Raugel polynomials.
+   * @copydoc TensorPolynomialsBase<dim>::clone()
    */
-  const unsigned int n_pols;
+  virtual std::unique_ptr<TensorPolynomialsBase<dim>>
+  clone() const override;
 
+private:
   /**
    * An object representing the polynomial space of Q
    * functions which forms the <tt>BR</tt> polynomials through
@@ -184,24 +166,6 @@ private:
 };
 
 
-template <int dim>
-inline unsigned int
-PolynomialsBernardiRaugel<dim>::n() const
-{
-  return n_pols;
-}
-
-
-
-template <int dim>
-inline unsigned int
-PolynomialsBernardiRaugel<dim>::degree() const
-{
-  return 2;
-}
-
-
-
 template <int dim>
 inline std::string
 PolynomialsBernardiRaugel<dim>::name() const
index d082478075555d2f37339a4799e5ed351c379c8b..ef1e061a261348422e76d9d24c1a38582a765b7d 100644 (file)
@@ -25,6 +25,7 @@
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/polynomial_space.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_polynomials_base.h>
 #include <deal.II/base/tensor_product_polynomials.h>
 
 #include <vector>
@@ -49,7 +50,7 @@ DEAL_II_NAMESPACE_OPEN
  * @date 2009, 2010
  */
 template <int dim>
-class PolynomialsNedelec
+class PolynomialsNedelec : public TensorPolynomialsBase<dim>
 {
 public:
   /**
@@ -79,26 +80,13 @@ public:
           std::vector<Tensor<2, dim>> &grads,
           std::vector<Tensor<3, dim>> &grad_grads,
           std::vector<Tensor<4, dim>> &third_derivatives,
-          std::vector<Tensor<5, dim>> &fourth_derivatives) const;
-
-  /**
-   * Return the number of Nédélec polynomials.
-   */
-  unsigned int
-  n() const;
-
-  /**
-   * Return the degree of the Nédélec space, which is one less than the
-   * highest polynomial degree.
-   */
-  unsigned int
-  degree() const;
+          std::vector<Tensor<5, dim>> &fourth_derivatives) const override;
 
   /**
    * Return the name of the space, which is <tt>Nedelec</tt>.
    */
   std::string
-  name() const;
+  name() const override;
 
   /**
    * Return the number of polynomials in the space <tt>N(degree)</tt> without
@@ -108,23 +96,19 @@ public:
   static unsigned int
   compute_n_pols(unsigned int degree);
 
-private:
   /**
-   * The degree of this object as given to the constructor.
+   * @copydoc TensorPolynomialsBase<dim>::clone()
    */
-  const unsigned int my_degree;
+  virtual std::unique_ptr<TensorPolynomialsBase<dim>>
+  clone() const override;
 
+private:
   /**
    * An object representing the polynomial space for a single component. We
    * can re-use it by rotating the coordinates of the evaluation point.
    */
   const AnisotropicPolynomials<dim> polynomial_space;
 
-  /**
-   * Number of Nédélec polynomials.
-   */
-  const unsigned int n_pols;
-
   /**
    * A static member function that creates the polynomial space we use to
    * initialize the #polynomial_space member variable.
@@ -134,22 +118,6 @@ private:
 };
 
 
-template <int dim>
-inline unsigned int
-PolynomialsNedelec<dim>::n() const
-{
-  return n_pols;
-}
-
-
-template <int dim>
-inline unsigned int
-PolynomialsNedelec<dim>::degree() const
-{
-  return my_degree;
-}
-
-
 template <int dim>
 inline std::string
 PolynomialsNedelec<dim>::name() const
index 15501e71b9a4ca746d3328124a33d7d484e31ee8..c0cd7c7bc43fc73745d051d2b5a3769d00bf0e7a 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/polynomial_space.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_polynomials_base.h>
 #include <deal.II/base/tensor_product_polynomials.h>
 
 #include <vector>
@@ -46,7 +47,7 @@ DEAL_II_NAMESPACE_OPEN
  * @date 2005
  */
 template <int dim>
-class PolynomialsRaviartThomas
+class PolynomialsRaviartThomas : public TensorPolynomialsBase<dim>
 {
 public:
   /**
@@ -77,26 +78,13 @@ public:
           std::vector<Tensor<2, dim>> &grads,
           std::vector<Tensor<3, dim>> &grad_grads,
           std::vector<Tensor<4, dim>> &third_derivatives,
-          std::vector<Tensor<5, dim>> &fourth_derivatives) const;
-
-  /**
-   * Return the number of Raviart-Thomas polynomials.
-   */
-  unsigned int
-  n() const;
-
-  /**
-   * Return the degree of the Raviart-Thomas space, which is one less than
-   * the highest polynomial degree.
-   */
-  unsigned int
-  degree() const;
+          std::vector<Tensor<5, dim>> &fourth_derivatives) const override;
 
   /**
    * Return the name of the space, which is <tt>RaviartThomas</tt>.
    */
   std::string
-  name() const;
+  name() const override;
 
   /**
    * Return the number of polynomials in the space <tt>RT(degree)</tt> without
@@ -106,23 +94,19 @@ public:
   static unsigned int
   compute_n_pols(unsigned int degree);
 
-private:
   /**
-   * The degree of this object as given to the constructor.
+   * @copydoc TensorPolynomialsBase<dim>::clone()
    */
-  const unsigned int my_degree;
+  virtual std::unique_ptr<TensorPolynomialsBase<dim>>
+  clone() const override;
 
+private:
   /**
    * An object representing the polynomial space for a single component. We
    * can re-use it by rotating the coordinates of the evaluation point.
    */
   const AnisotropicPolynomials<dim> polynomial_space;
 
-  /**
-   * Number of Raviart-Thomas polynomials.
-   */
-  const unsigned int n_pols;
-
   /**
    * A static member function that creates the polynomial space we use to
    * initialize the #polynomial_space member variable.
@@ -132,25 +116,6 @@ private:
 };
 
 
-
-template <int dim>
-inline unsigned int
-PolynomialsRaviartThomas<dim>::n() const
-{
-  return n_pols;
-}
-
-
-
-template <int dim>
-inline unsigned int
-PolynomialsRaviartThomas<dim>::degree() const
-{
-  return my_degree;
-}
-
-
-
 template <int dim>
 inline std::string
 PolynomialsRaviartThomas<dim>::name() const
index 6643d789a98394e678f3e58f80366fdf39c78f19..a77f1af3607b4c768df9d3818a75f7e04830e40e 100644 (file)
@@ -24,6 +24,7 @@
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/polynomials_raviart_thomas.h>
 #include <deal.II/base/tensor.h>
+#include <deal.II/base/tensor_polynomials_base.h>
 
 #include <vector>
 
@@ -84,7 +85,7 @@ DEAL_II_NAMESPACE_OPEN
  */
 
 template <int dim>
-class PolynomialsRT_Bubbles
+class PolynomialsRT_Bubbles : public TensorPolynomialsBase<dim>
 {
 public:
   /**
@@ -111,26 +112,13 @@ public:
           std::vector<Tensor<2, dim>> &grads,
           std::vector<Tensor<3, dim>> &grad_grads,
           std::vector<Tensor<4, dim>> &third_derivatives,
-          std::vector<Tensor<5, dim>> &fourth_derivatives) const;
-
-  /**
-   * Returns the number of enhanced Raviart-Thomas polynomials.
-   */
-  unsigned int
-  n() const;
-
-  /**
-   * Returns the degree of the RT_bubble space, which is one less than the
-   * highest polynomial degree.
-   */
-  unsigned int
-  degree() const;
+          std::vector<Tensor<5, dim>> &fourth_derivatives) const override;
 
   /**
    * Return the name of the space, which is <tt>RT_Bubbles</tt>.
    */
   std::string
-  name() const;
+  name() const override;
 
   /**
    * Return the number of polynomials in the space <tt>RT_Bubbles(degree)</tt>
@@ -140,12 +128,13 @@ public:
   static unsigned int
   compute_n_pols(const unsigned int degree);
 
-private:
   /**
-   * The degree of this object as given to the constructor.
+   * @copydoc TensorPolynomialsBase<dim>::clone()
    */
-  const unsigned int my_degree;
+  virtual std::unique_ptr<TensorPolynomialsBase<dim>>
+  clone() const override;
 
+private:
   /**
    * An object representing the Raviart-Thomas part of the space
    */
@@ -156,33 +145,9 @@ private:
    * to <i>k+1</i>.
    */
   std::vector<Polynomials::Polynomial<double>> monomials;
-
-  /**
-   * Number of RT_Bubbles polynomials.
-   */
-  unsigned int n_pols;
 };
 
 
-
-template <int dim>
-inline unsigned int
-PolynomialsRT_Bubbles<dim>::n() const
-{
-  return n_pols;
-}
-
-
-
-template <int dim>
-inline unsigned int
-PolynomialsRT_Bubbles<dim>::degree() const
-{
-  return my_degree;
-}
-
-
-
 template <int dim>
 inline std::string
 PolynomialsRT_Bubbles<dim>::name() const
index 4615c881568eede89c37d4e798c0aff39939ff35..6bd98f574717c41f801f9fdfe4374867433ec98e 100644 (file)
@@ -16,6 +16,7 @@
 
 #include <deal.II/base/polynomials_abf.h>
 #include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/std_cxx14/memory.h>
 
 #include <iomanip>
 #include <iostream>
@@ -47,9 +48,8 @@ namespace
 
 template <int dim>
 PolynomialsABF<dim>::PolynomialsABF(const unsigned int k)
-  : my_degree(k)
+  : TensorPolynomialsBase<dim>(k, compute_n_pols(k))
   , polynomial_space(get_abf_polynomials<dim>(k))
-  , n_pols(compute_n_pols(k))
 {
   // check that the dimensions match. we only store one of the 'dim'
   // anisotropic polynomials that make up the vector-valued space, so
@@ -69,16 +69,17 @@ PolynomialsABF<dim>::compute(
   std::vector<Tensor<4, dim>> &third_derivatives,
   std::vector<Tensor<5, dim>> &fourth_derivatives) const
 {
-  Assert(values.size() == n_pols || values.size() == 0,
-         ExcDimensionMismatch(values.size(), n_pols));
-  Assert(grads.size() == n_pols || grads.size() == 0,
-         ExcDimensionMismatch(grads.size(), n_pols));
-  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
-         ExcDimensionMismatch(grad_grads.size(), n_pols));
-  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
-         ExcDimensionMismatch(third_derivatives.size(), n_pols));
-  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
-         ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
+  Assert(values.size() == this->n() || values.size() == 0,
+         ExcDimensionMismatch(values.size(), this->n()));
+  Assert(grads.size() == this->n() || grads.size() == 0,
+         ExcDimensionMismatch(grads.size(), this->n()));
+  Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+         ExcDimensionMismatch(grad_grads.size(), this->n()));
+  Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+         ExcDimensionMismatch(third_derivatives.size(), this->n()));
+  Assert(fourth_derivatives.size() == this->n() ||
+           fourth_derivatives.size() == 0,
+         ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
 
   const unsigned int n_sub = polynomial_space.n();
   // guard access to the scratch
@@ -180,6 +181,14 @@ PolynomialsABF<dim>::compute_n_pols(const unsigned int k)
 }
 
 
+template <int dim>
+std::unique_ptr<TensorPolynomialsBase<dim>>
+PolynomialsABF<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialsABF<dim>>(*this);
+}
+
+
 template class PolynomialsABF<1>;
 template class PolynomialsABF<2>;
 template class PolynomialsABF<3>;
index 51338da71dad0a665205d2dff0fc1daecd47a040..4669b0ddbb10ab760f93c73627bc3c0c432d5f9c 100644 (file)
@@ -18,6 +18,7 @@
 #include <deal.II/base/polynomial_space.h>
 #include <deal.II/base/polynomials_bdm.h>
 #include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/std_cxx14/memory.h>
 
 #include <iomanip>
 #include <iostream>
@@ -27,9 +28,9 @@ DEAL_II_NAMESPACE_OPEN
 
 template <int dim>
 PolynomialsBDM<dim>::PolynomialsBDM(const unsigned int k)
-  : polynomial_space(Polynomials::Legendre::generate_complete_basis(k))
+  : TensorPolynomialsBase<dim>(k, compute_n_pols(k))
+  , polynomial_space(Polynomials::Legendre::generate_complete_basis(k))
   , monomials((dim == 2) ? (1) : (k + 2))
-  , n_pols(compute_n_pols(k))
   , p_values(polynomial_space.n())
   , p_grads(polynomial_space.n())
   , p_grad_grads(polynomial_space.n())
@@ -60,16 +61,17 @@ PolynomialsBDM<dim>::compute(
   std::vector<Tensor<4, dim>> &third_derivatives,
   std::vector<Tensor<5, dim>> &fourth_derivatives) const
 {
-  Assert(values.size() == n_pols || values.size() == 0,
-         ExcDimensionMismatch(values.size(), n_pols));
-  Assert(grads.size() == n_pols || grads.size() == 0,
-         ExcDimensionMismatch(grads.size(), n_pols));
-  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
-         ExcDimensionMismatch(grad_grads.size(), n_pols));
-  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
-         ExcDimensionMismatch(third_derivatives.size(), n_pols));
-  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
-         ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
+  Assert(values.size() == this->n() || values.size() == 0,
+         ExcDimensionMismatch(values.size(), this->n()));
+  Assert(grads.size() == this->n() || grads.size() == 0,
+         ExcDimensionMismatch(grads.size(), this->n()));
+  Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+         ExcDimensionMismatch(grad_grads.size(), this->n()));
+  Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+         ExcDimensionMismatch(third_derivatives.size(), this->n()));
+  Assert(fourth_derivatives.size() == this->n() ||
+           fourth_derivatives.size() == 0,
+         ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
 
   // third and fourth derivatives not implemented
   (void)third_derivatives;
@@ -353,7 +355,7 @@ PolynomialsBDM<dim>::compute(
               grad_grads[start + 2][2][2][2] = 0.;
             }
         }
-      Assert(start == n_pols, ExcInternalError());
+      Assert(start == this->n(), ExcInternalError());
     }
 }
 
@@ -440,6 +442,14 @@ PolynomialsBDM<dim>::compute_n_pols(unsigned int k)
 }
 
 
+template <int dim>
+std::unique_ptr<TensorPolynomialsBase<dim>>
+PolynomialsBDM<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialsBDM<dim>>(*this);
+}
+
+
 template class PolynomialsBDM<1>;
 template class PolynomialsBDM<2>;
 template class PolynomialsBDM<3>;
index af1570d297422fad8b51bd1c3d40e3899ec0ed62..62e643274eb787ea3f5ce7bc36f9ac45e359b565 100644 (file)
 
 
 #include <deal.II/base/polynomials_bernardi_raugel.h>
+#include <deal.II/base/std_cxx14/memory.h>
 
 DEAL_II_NAMESPACE_OPEN
 
 
 template <int dim>
 PolynomialsBernardiRaugel<dim>::PolynomialsBernardiRaugel(const unsigned int k)
-  : my_degree(k)
-  , n_pols(compute_n_pols(k))
+  : TensorPolynomialsBase<dim>(k + 1, compute_n_pols(k))
   , polynomial_space_Q(create_polynomials_Q())
   , polynomial_space_bubble(create_polynomials_bubble())
 {}
@@ -75,16 +75,17 @@ PolynomialsBernardiRaugel<dim>::compute(
   std::vector<Tensor<4, dim>> &third_derivatives,
   std::vector<Tensor<5, dim>> &fourth_derivatives) const
 {
-  Assert(values.size() == n_pols || values.size() == 0,
-         ExcDimensionMismatch(values.size(), n_pols));
-  Assert(grads.size() == n_pols || grads.size() == 0,
-         ExcDimensionMismatch(grads.size(), n_pols));
-  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
-         ExcDimensionMismatch(grad_grads.size(), n_pols));
-  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
-         ExcDimensionMismatch(third_derivatives.size(), n_pols));
-  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
-         ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
+  Assert(values.size() == this->n() || values.size() == 0,
+         ExcDimensionMismatch(values.size(), this->n()));
+  Assert(grads.size() == this->n() || grads.size() == 0,
+         ExcDimensionMismatch(grads.size(), this->n()));
+  Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+         ExcDimensionMismatch(grad_grads.size(), this->n()));
+  Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+         ExcDimensionMismatch(third_derivatives.size(), this->n()));
+  Assert(fourth_derivatives.size() == this->n() ||
+           fourth_derivatives.size() == 0,
+         ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
 
   std::vector<double>         Q_values;
   std::vector<Tensor<1, dim>> Q_grads;
@@ -248,6 +249,14 @@ PolynomialsBernardiRaugel<dim>::compute_n_pols(const unsigned int k)
   return 0;
 }
 
+
+template <int dim>
+std::unique_ptr<TensorPolynomialsBase<dim>>
+PolynomialsBernardiRaugel<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialsBernardiRaugel<dim>>(*this);
+}
+
 template class PolynomialsBernardiRaugel<1>; // to prevent errors
 template class PolynomialsBernardiRaugel<2>;
 template class PolynomialsBernardiRaugel<3>;
index 1b48d076af6eebd7408a433755d5c752b93d0482..de531fada5e69b6a2784e1e6a275d6cb7f51bc40 100644 (file)
@@ -17,6 +17,7 @@
 #include <deal.II/base/polynomial.h>
 #include <deal.II/base/polynomials_nedelec.h>
 #include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/std_cxx14/memory.h>
 
 #include <iomanip>
 #include <iostream>
@@ -26,9 +27,8 @@ DEAL_II_NAMESPACE_OPEN
 
 template <int dim>
 PolynomialsNedelec<dim>::PolynomialsNedelec(const unsigned int k)
-  : my_degree(k)
+  : TensorPolynomialsBase<dim>(k, compute_n_pols(k))
   , polynomial_space(create_polynomials(k))
-  , n_pols(compute_n_pols(k))
 {}
 
 template <int dim>
@@ -59,16 +59,17 @@ PolynomialsNedelec<dim>::compute(
   std::vector<Tensor<4, dim>> &third_derivatives,
   std::vector<Tensor<5, dim>> &fourth_derivatives) const
 {
-  Assert(values.size() == n_pols || values.size() == 0,
-         ExcDimensionMismatch(values.size(), n_pols));
-  Assert(grads.size() == n_pols || grads.size() == 0,
-         ExcDimensionMismatch(grads.size(), n_pols));
-  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
-         ExcDimensionMismatch(grad_grads.size(), n_pols));
-  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
-         ExcDimensionMismatch(third_derivatives.size(), n_pols));
-  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
-         ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
+  Assert(values.size() == this->n() || values.size() == 0,
+         ExcDimensionMismatch(values.size(), this->n()));
+  Assert(grads.size() == this->n() || grads.size() == 0,
+         ExcDimensionMismatch(grads.size(), this->n()));
+  Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+         ExcDimensionMismatch(grad_grads.size(), this->n()));
+  Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+         ExcDimensionMismatch(third_derivatives.size(), this->n()));
+  Assert(fourth_derivatives.size() == this->n() ||
+           fourth_derivatives.size() == 0,
+         ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
 
   // third and fourth derivatives not implemented
   (void)third_derivatives;
@@ -80,7 +81,8 @@ PolynomialsNedelec<dim>::compute(
   // and second derivatives vectors of
   // <tt>polynomial_space</tt> at
   // <tt>unit_point</tt>
-  const unsigned int  n_basis = polynomial_space.n();
+  const unsigned int  n_basis   = polynomial_space.n();
+  const unsigned int  my_degree = this->degree();
   std::vector<double> unit_point_values((values.size() == 0) ? 0 : n_basis);
   std::vector<Tensor<1, dim>> unit_point_grads((grads.size() == 0) ? 0 :
                                                                      n_basis);
@@ -1506,6 +1508,14 @@ PolynomialsNedelec<dim>::compute_n_pols(unsigned int k)
 }
 
 
+template <int dim>
+std::unique_ptr<TensorPolynomialsBase<dim>>
+PolynomialsNedelec<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialsNedelec<dim>>(*this);
+}
+
+
 template class PolynomialsNedelec<1>;
 template class PolynomialsNedelec<2>;
 template class PolynomialsNedelec<3>;
index ed4f3525b7a7e2ab4ab4290bfb0e83a687e711c2..d15e210ca8a855194d1b7e0f1a979d20b8e5ed19 100644 (file)
@@ -16,6 +16,7 @@
 
 #include <deal.II/base/polynomials_raviart_thomas.h>
 #include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/std_cxx14/memory.h>
 #include <deal.II/base/thread_management.h>
 
 #include <iomanip>
@@ -35,9 +36,8 @@ DEAL_II_NAMESPACE_OPEN
 
 template <int dim>
 PolynomialsRaviartThomas<dim>::PolynomialsRaviartThomas(const unsigned int k)
-  : my_degree(k)
+  : TensorPolynomialsBase<dim>(k, compute_n_pols(k))
   , polynomial_space(create_polynomials(k))
-  , n_pols(compute_n_pols(k))
 {}
 
 
@@ -69,16 +69,17 @@ PolynomialsRaviartThomas<dim>::compute(
   std::vector<Tensor<4, dim>> &third_derivatives,
   std::vector<Tensor<5, dim>> &fourth_derivatives) const
 {
-  Assert(values.size() == n_pols || values.size() == 0,
-         ExcDimensionMismatch(values.size(), n_pols));
-  Assert(grads.size() == n_pols || grads.size() == 0,
-         ExcDimensionMismatch(grads.size(), n_pols));
-  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
-         ExcDimensionMismatch(grad_grads.size(), n_pols));
-  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
-         ExcDimensionMismatch(third_derivatives.size(), n_pols));
-  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
-         ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
+  Assert(values.size() == this->n() || values.size() == 0,
+         ExcDimensionMismatch(values.size(), this->n()));
+  Assert(grads.size() == this->n() || grads.size() == 0,
+         ExcDimensionMismatch(grads.size(), this->n()));
+  Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+         ExcDimensionMismatch(grad_grads.size(), this->n()));
+  Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+         ExcDimensionMismatch(third_derivatives.size(), this->n()));
+  Assert(fourth_derivatives.size() == this->n() ||
+           fourth_derivatives.size() == 0,
+         ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
 
   // have a few scratch
   // arrays. because we don't want to
@@ -183,6 +184,14 @@ PolynomialsRaviartThomas<dim>::compute_n_pols(unsigned int k)
 }
 
 
+template <int dim>
+std::unique_ptr<TensorPolynomialsBase<dim>>
+PolynomialsRaviartThomas<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialsRaviartThomas<dim>>(*this);
+}
+
+
 template class PolynomialsRaviartThomas<1>;
 template class PolynomialsRaviartThomas<2>;
 template class PolynomialsRaviartThomas<3>;
index 117d244723462c7738e2a9d7a50b17cc9b8c49ac..c4412c930c2c14fcee5bc90ed9242f1b4c95b199 100644 (file)
@@ -17,6 +17,7 @@
 #include "deal.II/base/polynomials_rt_bubbles.h"
 
 #include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/std_cxx14/memory.h>
 #include <deal.II/base/thread_management.h>
 
 #include <iomanip>
@@ -27,10 +28,9 @@ DEAL_II_NAMESPACE_OPEN
 
 template <int dim>
 PolynomialsRT_Bubbles<dim>::PolynomialsRT_Bubbles(const unsigned int k)
-  : my_degree(k)
+  : TensorPolynomialsBase<dim>(k, compute_n_pols(k))
   , raviart_thomas_space(k - 1)
   , monomials(k + 2)
-  , n_pols(compute_n_pols(k))
 {
   Assert(dim >= 2, ExcImpossibleInDim(dim));
 
@@ -50,16 +50,17 @@ PolynomialsRT_Bubbles<dim>::compute(
   std::vector<Tensor<4, dim>> &third_derivatives,
   std::vector<Tensor<5, dim>> &fourth_derivatives) const
 {
-  Assert(values.size() == n_pols || values.size() == 0,
-         ExcDimensionMismatch(values.size(), n_pols));
-  Assert(grads.size() == n_pols || grads.size() == 0,
-         ExcDimensionMismatch(grads.size(), n_pols));
-  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
-         ExcDimensionMismatch(grad_grads.size(), n_pols));
-  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
-         ExcDimensionMismatch(third_derivatives.size(), n_pols));
-  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
-         ExcDimensionMismatch(fourth_derivatives.size(), n_pols));
+  Assert(values.size() == this->n() || values.size() == 0,
+         ExcDimensionMismatch(values.size(), this->n()));
+  Assert(grads.size() == this->n() || grads.size() == 0,
+         ExcDimensionMismatch(grads.size(), this->n()));
+  Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+         ExcDimensionMismatch(grad_grads.size(), this->n()));
+  Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+         ExcDimensionMismatch(third_derivatives.size(), this->n()));
+  Assert(fourth_derivatives.size() == this->n() ||
+           fourth_derivatives.size() == 0,
+         ExcDimensionMismatch(fourth_derivatives.size(), this->n()));
 
   // Third and fourth derivatives are not implemented
   (void)third_derivatives;
@@ -67,7 +68,8 @@ PolynomialsRT_Bubbles<dim>::compute(
   (void)fourth_derivatives;
   Assert(fourth_derivatives.size() == 0, ExcNotImplemented());
 
-  const unsigned int n_sub = raviart_thomas_space.n();
+  const unsigned int n_sub     = raviart_thomas_space.n();
+  const unsigned int my_degree = this->degree();
 
   // Guard access to the scratch arrays in the following block
   // using a mutex to make sure they are not used by multiple threads
@@ -194,7 +196,7 @@ PolynomialsRT_Bubbles<dim>::compute(
                 monoval_plus[0][1] * monoval_i[1][2];
             }
         }
-      Assert(start == n_pols - my_degree - 1, ExcInternalError());
+      Assert(start == this->n() - my_degree - 1, ExcInternalError());
     }
   else if (dim == 3)
     {
@@ -825,7 +827,7 @@ PolynomialsRT_Bubbles<dim>::compute(
                 start += 2;
             }
         }
-      Assert(start == n_pols - 2 * n_curls, ExcInternalError());
+      Assert(start == this->n() - 2 * n_curls, ExcInternalError());
     }
 }
 
@@ -843,6 +845,14 @@ PolynomialsRT_Bubbles<dim>::compute_n_pols(const unsigned int k)
 }
 
 
+template <int dim>
+std::unique_ptr<TensorPolynomialsBase<dim>>
+PolynomialsRT_Bubbles<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialsRT_Bubbles<dim>>(*this);
+}
+
+
 template class PolynomialsRT_Bubbles<1>;
 template class PolynomialsRT_Bubbles<2>;
 template class PolynomialsRT_Bubbles<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.