]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Change polynomial classes to derive from ScalarPolynomialsBase 8678/head
authorgrahambenharper <grahambenharper@gmail.com>
Sun, 1 Sep 2019 08:16:59 +0000 (02:16 -0600)
committergrahambenharper <grahambenharper@gmail.com>
Tue, 3 Sep 2019 00:11:17 +0000 (18:11 -0600)
18 files changed:
include/deal.II/base/polynomial_space.h
include/deal.II/base/polynomials_rannacher_turek.h
include/deal.II/base/tensor_product_polynomials.h
include/deal.II/base/tensor_product_polynomials_bubbles.h
include/deal.II/base/tensor_product_polynomials_const.h
include/deal.II/fe/fe_poly.h
include/deal.II/fe/fe_poly_face.h
source/base/polynomial_space.cc
source/base/polynomials_bdm.cc
source/base/polynomials_rannacher_turek.cc
source/base/tensor_product_polynomials.cc
source/base/tensor_product_polynomials_bubbles.cc
source/base/tensor_product_polynomials_const.cc
source/fe/fe_bdm.cc
source/fe/fe_dgp_nonparametric.cc
source/fe/mapping_q_generic.cc
tests/base/anisotropic_1.cc
tests/base/polynomial_test.cc

index 1e43d9a73a65744826c30e51fc69c05ba019bf37..368bda6d18b11b00ccb6d8cfd2acbfc75506fa2f 100644 (file)
@@ -22,6 +22,7 @@
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/point.h>
 #include <deal.II/base/polynomial.h>
+#include <deal.II/base/scalar_polynomials_base.h>
 #include <deal.II/base/smartpointer.h>
 #include <deal.II/base/tensor.h>
 
@@ -95,7 +96,7 @@ DEAL_II_NAMESPACE_OPEN
  * 2005
  */
 template <int dim>
-class PolynomialSpace
+class PolynomialSpace : public ScalarPolynomialsBase<dim>
 {
 public:
   /**
@@ -142,12 +143,12 @@ public:
    * compute_grad_grad() functions, see below, in a loop over all polynomials.
    */
   void
-  compute(const Point<dim> &           unit_point,
-          std::vector<double> &        values,
-          std::vector<Tensor<1, dim>> &grads,
-          std::vector<Tensor<2, dim>> &grad_grads,
-          std::vector<Tensor<3, dim>> &third_derivatives,
-          std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+  evaluate(const Point<dim> &           unit_point,
+           std::vector<double> &        values,
+           std::vector<Tensor<1, dim>> &grads,
+           std::vector<Tensor<2, dim>> &grad_grads,
+           std::vector<Tensor<3, dim>> &third_derivatives,
+           std::vector<Tensor<4, dim>> &fourth_derivatives) const override;
 
   /**
    * Compute the value of the <tt>i</tt>th polynomial at unit point
@@ -194,27 +195,20 @@ public:
    * given, then the result of this function is <i>N</i> in 1d,
    * <i>N(N+1)/2</i> in 2d, and <i>N(N+1)(N+2)/6</i> in 3d.
    */
-  unsigned int
-  n() const;
+  static unsigned int
+  n_polynomials(const unsigned int n);
 
   /**
-   * Degree of the space. This is by definition the number of polynomials
-   * given to the constructor, NOT the maximal degree of a polynomial in this
-   * vector. The latter value is never checked and therefore left to the
-   * application.
+   * Return the name of the space, which is <tt>PolynomialSpace</tt>.
    */
-  unsigned int
-  degree() const;
+  std::string
+  name() const override;
 
   /**
-   * Static function used in the constructor to compute the number of
-   * polynomials.
-   *
-   * @warning The argument `n` is not the maximal degree, but the number of
-   * onedimensional polynomials, thus the degree plus one.
+   * @copydoc ScalarPolynomialsBase<dim>::clone()
    */
-  static unsigned int
-  n_polynomials(const unsigned int n);
+  virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
+  clone() const override;
 
 protected:
   /**
@@ -234,11 +228,6 @@ private:
    */
   const std::vector<Polynomials::Polynomial<double>> polynomials;
 
-  /**
-   * Store the precomputed value which the <tt>n()</tt> function returns.
-   */
-  const unsigned int n_pols;
-
   /**
    * Index map for reordering the polynomials.
    */
@@ -270,16 +259,16 @@ PolynomialSpace<3>::compute_index(const unsigned int n) const;
 template <int dim>
 template <class Pol>
 PolynomialSpace<dim>::PolynomialSpace(const std::vector<Pol> &pols)
-  : polynomials(pols.begin(), pols.end())
-  , n_pols(n_polynomials(polynomials.size()))
-  , index_map(n_pols)
-  , index_map_inverse(n_pols)
+  : ScalarPolynomialsBase<dim>(pols.size(), n_polynomials(pols.size()))
+  , polynomials(pols.begin(), pols.end())
+  , index_map(n_polynomials(pols.size()))
+  , index_map_inverse(n_polynomials(pols.size()))
 {
   // per default set this index map
   // to identity. This map can be
   // changed by the user through the
   // set_numbering function
-  for (unsigned int i = 0; i < n_pols; ++i)
+  for (unsigned int i = 0; i < this->n(); ++i)
     {
       index_map[i]         = i;
       index_map_inverse[i] = i;
@@ -287,20 +276,12 @@ PolynomialSpace<dim>::PolynomialSpace(const std::vector<Pol> &pols)
 }
 
 
-template <int dim>
-inline unsigned int
-PolynomialSpace<dim>::n() const
-{
-  return n_pols;
-}
-
-
 
 template <int dim>
-inline unsigned int
-PolynomialSpace<dim>::degree() const
+inline std::string
+PolynomialSpace<dim>::name() const
 {
-  return polynomials.size();
+  return "PolynomialSpace";
 }
 
 
@@ -309,7 +290,7 @@ template <class StreamType>
 void
 PolynomialSpace<dim>::output_indices(StreamType &out) const
 {
-  for (unsigned int i = 0; i < n_pols; ++i)
+  for (unsigned int i = 0; i < this->n(); ++i)
     {
       const std::array<unsigned int, dim> ix = compute_index(i);
       out << i << "\t";
index 37b8b3ac813ad432814a98673c135f864f80f62b..bf3fa41c463c31a480b6445dd95a885e2a78355d 100644 (file)
@@ -18,6 +18,7 @@
 #define dealii_polynomials_rannacher_turek_h
 
 #include <deal.II/base/point.h>
+#include <deal.II/base/scalar_polynomials_base.h>
 #include <deal.II/base/tensor.h>
 
 #include <vector>
@@ -38,7 +39,7 @@ DEAL_II_NAMESPACE_OPEN
  * @date 2015
  */
 template <int dim>
-class PolynomialsRannacherTurek
+class PolynomialsRannacherTurek : public ScalarPolynomialsBase<dim>
 {
 public:
   /**
@@ -85,12 +86,24 @@ public:
    * zero. A size of zero means that we are not computing the vector entries.
    */
   void
-  compute(const Point<dim> &           unit_point,
-          std::vector<double> &        values,
-          std::vector<Tensor<1, dim>> &grads,
-          std::vector<Tensor<2, dim>> &grad_grads,
-          std::vector<Tensor<3, dim>> &third_derivatives,
-          std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+  evaluate(const Point<dim> &           unit_point,
+           std::vector<double> &        values,
+           std::vector<Tensor<1, dim>> &grads,
+           std::vector<Tensor<2, dim>> &grad_grads,
+           std::vector<Tensor<3, dim>> &third_derivatives,
+           std::vector<Tensor<4, dim>> &fourth_derivatives) const override;
+
+  /**
+   * Return the name of the space, which is <tt>RannacherTurek</tt>.
+   */
+  std::string
+  name() const override;
+
+  /**
+   * @copydoc ScalarPolynomialsBase<dim>::clone()
+   */
+  virtual std::unique_ptr<ScalarPolynomialsBase<dim>>
+  clone() const override;
 };
 
 
@@ -204,6 +217,15 @@ PolynomialsRannacherTurek<dim>::compute_derivative(const unsigned int i,
 }
 
 
+
+template <int dim>
+inline std::string
+PolynomialsRannacherTurek<dim>::name() const
+{
+  return "RannacherTurek";
+}
+
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 3a956a8b5ece5287e985368084c6f90ed5b61472..9c39dec516942ed5dafaa68046c37907f5290baf 100644 (file)
@@ -118,12 +118,12 @@ public:
    * over all tensor product polynomials.
    */
   void
-  compute(const Point<dim> &           unit_point,
-          std::vector<double> &        values,
-          std::vector<Tensor<1, dim>> &grads,
-          std::vector<Tensor<2, dim>> &grad_grads,
-          std::vector<Tensor<3, dim>> &third_derivatives,
-          std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+  evaluate(const Point<dim> &           unit_point,
+           std::vector<double> &        values,
+           std::vector<Tensor<1, dim>> &grads,
+           std::vector<Tensor<2, dim>> &grad_grads,
+           std::vector<Tensor<3, dim>> &third_derivatives,
+           std::vector<Tensor<4, dim>> &fourth_derivatives) const;
 
   /**
    * Compute the value of the <tt>i</tt>th tensor product polynomial at
index 80559955c36c5e357d305d3f925ea1c899403be5..7a0f3dd0187e5fec7f5b73f9be835952026e27eb 100644 (file)
@@ -69,12 +69,12 @@ public:
    * over all tensor product polynomials.
    */
   void
-  compute(const Point<dim> &           unit_point,
-          std::vector<double> &        values,
-          std::vector<Tensor<1, dim>> &grads,
-          std::vector<Tensor<2, dim>> &grad_grads,
-          std::vector<Tensor<3, dim>> &third_derivatives,
-          std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+  evaluate(const Point<dim> &           unit_point,
+           std::vector<double> &        values,
+           std::vector<Tensor<1, dim>> &grads,
+           std::vector<Tensor<2, dim>> &grad_grads,
+           std::vector<Tensor<3, dim>> &third_derivatives,
+           std::vector<Tensor<4, dim>> &fourth_derivatives) const;
 
   /**
    * Compute the value of the <tt>i</tt>th tensor product polynomial at
index ff2680bd45693ec30b04f92ef648110afd42a198..1671532506d65660c48060a01d5c2f3ed8ea6172 100644 (file)
@@ -75,12 +75,12 @@ public:
    * over all tensor product polynomials.
    */
   void
-  compute(const Point<dim> &           unit_point,
-          std::vector<double> &        values,
-          std::vector<Tensor<1, dim>> &grads,
-          std::vector<Tensor<2, dim>> &grad_grads,
-          std::vector<Tensor<3, dim>> &third_derivatives,
-          std::vector<Tensor<4, dim>> &fourth_derivatives) const;
+  evaluate(const Point<dim> &           unit_point,
+           std::vector<double> &        values,
+           std::vector<Tensor<1, dim>> &grads,
+           std::vector<Tensor<2, dim>> &grad_grads,
+           std::vector<Tensor<3, dim>> &third_derivatives,
+           std::vector<Tensor<4, dim>> &fourth_derivatives) const;
 
   /**
    * Compute the value of the <tt>i</tt>th tensor product polynomial at
index 57a4419522208dfecf3c9d21161afad8c7424228..d3917724d429cb580f30c2e3b4ec6d2a1b797881 100644 (file)
@@ -40,12 +40,12 @@ DEAL_II_NAMESPACE_OPEN
  * @code
  *  static const unsigned int dimension;
  *
- *  void compute (const Point<dim>            &unit_point,
- *                std::vector<double>         &values,
- *                std::vector<Tensor<1,dim> > &grads,
- *                std::vector<Tensor<2,dim> > &grad_grads,
- *                std::vector<Tensor<3,dim> > &third_derivatives,
- *                std::vector<Tensor<4,dim> > &fourth_derivatives) const;
+ *  void evaluate (const Point<dim>            &unit_point,
+ *                 std::vector<double>         &values,
+ *                 std::vector<Tensor<1,dim> > &grads,
+ *                 std::vector<Tensor<2,dim> > &grad_grads,
+ *                 std::vector<Tensor<3,dim> > &third_derivatives,
+ *                 std::vector<Tensor<4,dim> > &fourth_derivatives) const;
  *
  *  double compute_value (const unsigned int i,
  *                        const Point<dim> &p) const;
@@ -303,12 +303,12 @@ protected:
                         update_3rd_derivatives))
       for (unsigned int i = 0; i < n_q_points; ++i)
         {
-          poly_space.compute(quadrature.point(i),
-                             values,
-                             grads,
-                             grad_grads,
-                             third_derivatives,
-                             fourth_derivatives);
+          poly_space.evaluate(quadrature.point(i),
+                              values,
+                              grads,
+                              grad_grads,
+                              third_derivatives,
+                              fourth_derivatives);
 
           // the values of shape functions at quadrature points don't change.
           // consequently, write these values right into the output array if
index 4a79199f40dc6ba4c9db1cadbd89d863291f3a46..1c5780fcf9e7db598ff7837ce4cb1c83f0b1056c 100644 (file)
@@ -133,12 +133,12 @@ protected:
                                  std::vector<double>(n_q_points));
         for (unsigned int i = 0; i < n_q_points; ++i)
           {
-            poly_space.compute(quadrature.point(i),
-                               values,
-                               grads,
-                               grad_grads,
-                               empty_vector_of_3rd_order_tensors,
-                               empty_vector_of_4th_order_tensors);
+            poly_space.evaluate(quadrature.point(i),
+                                values,
+                                grads,
+                                grad_grads,
+                                empty_vector_of_3rd_order_tensors,
+                                empty_vector_of_4th_order_tensors);
 
             for (unsigned int k = 0; k < poly_space.n(); ++k)
               data.shape_values[k][i] = values[k];
index a87fcefa226de00c6dd37e30457d2b6a1bd940c5..95f87b2be32e77ea939be4ae03ae3560c1afe29e 100644 (file)
@@ -15,6 +15,7 @@
 
 #include <deal.II/base/exceptions.h>
 #include <deal.II/base/polynomial_space.h>
+#include <deal.II/base/std_cxx14/memory.h>
 #include <deal.II/base/table.h>
 
 DEAL_II_NAMESPACE_OPEN
@@ -200,7 +201,7 @@ PolynomialSpace<dim>::compute_grad_grad(const unsigned int i,
 
 template <int dim>
 void
-PolynomialSpace<dim>::compute(
+PolynomialSpace<dim>::evaluate(
   const Point<dim> &           p,
   std::vector<double> &        values,
   std::vector<Tensor<1, dim>> &grads,
@@ -210,41 +211,42 @@ PolynomialSpace<dim>::compute(
 {
   const unsigned int n_1d = polynomials.size();
 
-  Assert(values.size() == n_pols || values.size() == 0,
-         ExcDimensionMismatch2(values.size(), n_pols, 0));
-  Assert(grads.size() == n_pols || grads.size() == 0,
-         ExcDimensionMismatch2(grads.size(), n_pols, 0));
-  Assert(grad_grads.size() == n_pols || grad_grads.size() == 0,
-         ExcDimensionMismatch2(grad_grads.size(), n_pols, 0));
-  Assert(third_derivatives.size() == n_pols || third_derivatives.size() == 0,
-         ExcDimensionMismatch2(third_derivatives.size(), n_pols, 0));
-  Assert(fourth_derivatives.size() == n_pols || fourth_derivatives.size() == 0,
-         ExcDimensionMismatch2(fourth_derivatives.size(), n_pols, 0));
+  Assert(values.size() == this->n() || values.size() == 0,
+         ExcDimensionMismatch2(values.size(), this->n(), 0));
+  Assert(grads.size() == this->n() || grads.size() == 0,
+         ExcDimensionMismatch2(grads.size(), this->n(), 0));
+  Assert(grad_grads.size() == this->n() || grad_grads.size() == 0,
+         ExcDimensionMismatch2(grad_grads.size(), this->n(), 0));
+  Assert(third_derivatives.size() == this->n() || third_derivatives.size() == 0,
+         ExcDimensionMismatch2(third_derivatives.size(), this->n(), 0));
+  Assert(fourth_derivatives.size() == this->n() ||
+           fourth_derivatives.size() == 0,
+         ExcDimensionMismatch2(fourth_derivatives.size(), this->n(), 0));
 
   unsigned int v_size = 0;
   bool update_values = false, update_grads = false, update_grad_grads = false;
   bool update_3rd_derivatives = false, update_4th_derivatives = false;
-  if (values.size() == n_pols)
+  if (values.size() == this->n())
     {
       update_values = true;
       v_size        = 1;
     }
-  if (grads.size() == n_pols)
+  if (grads.size() == this->n())
     {
       update_grads = true;
       v_size       = 2;
     }
-  if (grad_grads.size() == n_pols)
+  if (grad_grads.size() == this->n())
     {
       update_grad_grads = true;
       v_size            = 3;
     }
-  if (third_derivatives.size() == n_pols)
+  if (third_derivatives.size() == this->n())
     {
       update_3rd_derivatives = true;
       v_size                 = 4;
     }
-  if (fourth_derivatives.size() == n_pols)
+  if (fourth_derivatives.size() == this->n())
     {
       update_4th_derivatives = true;
       v_size                 = 5;
@@ -396,6 +398,15 @@ PolynomialSpace<dim>::compute(
 }
 
 
+
+template <int dim>
+std::unique_ptr<ScalarPolynomialsBase<dim>>
+PolynomialSpace<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialSpace<dim>>(*this);
+}
+
+
 template class PolynomialSpace<1>;
 template class PolynomialSpace<2>;
 template class PolynomialSpace<3>;
index 80a2abe0c1d67925662eb6e047b545de201a5622..34005fcede4571796276a05ff6f64c2380be1941 100644 (file)
@@ -98,12 +98,12 @@ PolynomialsBDM<dim>::evaluate(
     // will have first all polynomials
     // in the x-component, then y and
     // z.
-    polynomial_space.compute(unit_point,
-                             p_values,
-                             p_grads,
-                             p_grad_grads,
-                             p_third_derivatives,
-                             p_fourth_derivatives);
+    polynomial_space.evaluate(unit_point,
+                              p_values,
+                              p_grads,
+                              p_grad_grads,
+                              p_third_derivatives,
+                              p_fourth_derivatives);
 
     std::fill(values.begin(), values.end(), Tensor<1, dim>());
     for (unsigned int i = 0; i < p_values.size(); ++i)
index a896143dab368fc48824d0c4e3f127f6a4390330..952c58e41a93e18e6d8a498fe31c471fb68fb1e8 100644 (file)
 
 #include <deal.II/base/geometry_info.h>
 #include <deal.II/base/polynomials_rannacher_turek.h>
+#include <deal.II/base/std_cxx14/memory.h>
 
 DEAL_II_NAMESPACE_OPEN
 
 
 template <int dim>
 PolynomialsRannacherTurek<dim>::PolynomialsRannacherTurek()
+  : ScalarPolynomialsBase<dim>(2, dealii::GeometryInfo<dim>::faces_per_cell)
 {
   Assert(dim == 2, ExcNotImplemented());
 }
@@ -141,7 +143,7 @@ PolynomialsRannacherTurek<dim>::compute_grad_grad(
 
 template <int dim>
 void
-PolynomialsRannacherTurek<dim>::compute(
+PolynomialsRannacherTurek<dim>::evaluate(
   const Point<dim> &           unit_point,
   std::vector<double> &        values,
   std::vector<Tensor<1, dim>> &grads,
@@ -149,7 +151,7 @@ PolynomialsRannacherTurek<dim>::compute(
   std::vector<Tensor<3, dim>> &third_derivatives,
   std::vector<Tensor<4, dim>> &fourth_derivatives) const
 {
-  const unsigned int n_pols = dealii::GeometryInfo<dim>::faces_per_cell;
+  const unsigned int n_pols = this->n();
   Assert(values.size() == n_pols || values.size() == 0,
          ExcDimensionMismatch(values.size(), n_pols));
   Assert(grads.size() == n_pols || grads.size() == 0,
@@ -187,6 +189,15 @@ PolynomialsRannacherTurek<dim>::compute(
 }
 
 
+
+template <int dim>
+std::unique_ptr<ScalarPolynomialsBase<dim>>
+PolynomialsRannacherTurek<dim>::clone() const
+{
+  return std_cxx14::make_unique<PolynomialsRannacherTurek<dim>>(*this);
+}
+
+
 // explicit instantiations
 #include "polynomials_rannacher_turek.inst"
 
index 71efbce20e260aba8a6842d097024dd3ce4b02f2..71d5f381861ff0ec68a063e80f35a443da7d3ecc 100644 (file)
@@ -241,7 +241,7 @@ TensorProductPolynomials<dim, PolynomialType>::compute_grad_grad(
 
 template <int dim, typename PolynomialType>
 void
-TensorProductPolynomials<dim, PolynomialType>::compute(
+TensorProductPolynomials<dim, PolynomialType>::evaluate(
   const Point<dim> &           p,
   std::vector<double> &        values,
   std::vector<Tensor<1, dim>> &grads,
index e927914d31412e53187892637fd74716474aaad3..afd88bfb827462945fc7ca94629d860374bc2f4f 100644 (file)
@@ -215,7 +215,7 @@ TensorProductPolynomialsBubbles<dim>::compute_grad_grad(
 
 template <int dim>
 void
-TensorProductPolynomialsBubbles<dim>::compute(
+TensorProductPolynomialsBubbles<dim>::evaluate(
   const Point<dim> &           p,
   std::vector<double> &        values,
   std::vector<Tensor<1, dim>> &grads,
@@ -273,7 +273,7 @@ TensorProductPolynomialsBubbles<dim>::compute(
       do_4th_derivatives = true;
     }
 
-  this->TensorProductPolynomials<dim>::compute(
+  this->TensorProductPolynomials<dim>::evaluate(
     p, values, grads, grad_grads, third_derivatives, fourth_derivatives);
 
   for (unsigned int i = this->n_tensor_pols;
index bbeef38690e4121d0396e68e7f57e876cae5f8a5..63bce45b98c7d18d460924b0871b09926b14c8da 100644 (file)
@@ -86,7 +86,7 @@ TensorProductPolynomialsConst<dim>::compute_grad_grad(const unsigned int i,
 
 template <int dim>
 void
-TensorProductPolynomialsConst<dim>::compute(
+TensorProductPolynomialsConst<dim>::evaluate(
   const Point<dim> &           p,
   std::vector<double> &        values,
   std::vector<Tensor<1, dim>> &grads,
@@ -141,7 +141,7 @@ TensorProductPolynomialsConst<dim>::compute(
       do_4th_derivatives = true;
     }
 
-  this->TensorProductPolynomials<dim>::compute(
+  this->TensorProductPolynomials<dim>::evaluate(
     p, values, grads, grad_grads, third_derivatives, fourth_derivatives);
 
   // for dgq node: values =1, grads=0, grads_grads=0, third_derivatives=0,
index 0a09f141f848d595f7303798e599e2e9f1c2803b..26b747877be79fae65adb66d103c1c51acdfeff2 100644 (file)
@@ -294,12 +294,12 @@ namespace internal
         for (unsigned int k = 0; k < quadrature.size(); ++k)
           {
             test_values[k].resize(poly.n());
-            poly.compute(quadrature.point(k),
-                         test_values[k],
-                         dummy1,
-                         dummy2,
-                         dummy3,
-                         dummy4);
+            poly.evaluate(quadrature.point(k),
+                          test_values[k],
+                          dummy1,
+                          dummy2,
+                          dummy3,
+                          dummy4);
             for (unsigned int i = 0; i < poly.n(); ++i)
               {
                 test_values[k][i] *= quadrature.weight(k);
index 9c3bc1b760ec3fe2df4cf6f185ac7f8dcbe667f0..f748a4ed0376f4229155a30ebe72b329ceaaab98 100644 (file)
@@ -332,12 +332,12 @@ FE_DGPNonparametric<dim, spacedim>::fill_fe_values(
   if (fe_internal.update_each & (update_values | update_gradients))
     for (unsigned int i = 0; i < n_q_points; ++i)
       {
-        polynomial_space.compute(mapping_data.quadrature_points[i],
-                                 values,
-                                 grads,
-                                 grad_grads,
-                                 empty_vector_of_3rd_order_tensors,
-                                 empty_vector_of_4th_order_tensors);
+        polynomial_space.evaluate(mapping_data.quadrature_points[i],
+                                  values,
+                                  grads,
+                                  grad_grads,
+                                  empty_vector_of_3rd_order_tensors,
+                                  empty_vector_of_4th_order_tensors);
 
         if (fe_internal.update_each & update_values)
           for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
@@ -388,12 +388,12 @@ FE_DGPNonparametric<dim, spacedim>::fill_fe_face_values(
   if (fe_internal.update_each & (update_values | update_gradients))
     for (unsigned int i = 0; i < n_q_points; ++i)
       {
-        polynomial_space.compute(mapping_data.quadrature_points[i],
-                                 values,
-                                 grads,
-                                 grad_grads,
-                                 empty_vector_of_3rd_order_tensors,
-                                 empty_vector_of_4th_order_tensors);
+        polynomial_space.evaluate(mapping_data.quadrature_points[i],
+                                  values,
+                                  grads,
+                                  grad_grads,
+                                  empty_vector_of_3rd_order_tensors,
+                                  empty_vector_of_4th_order_tensors);
 
         if (fe_internal.update_each & update_values)
           for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
@@ -445,12 +445,12 @@ FE_DGPNonparametric<dim, spacedim>::fill_fe_subface_values(
   if (fe_internal.update_each & (update_values | update_gradients))
     for (unsigned int i = 0; i < n_q_points; ++i)
       {
-        polynomial_space.compute(mapping_data.quadrature_points[i],
-                                 values,
-                                 grads,
-                                 grad_grads,
-                                 empty_vector_of_3rd_order_tensors,
-                                 empty_vector_of_4th_order_tensors);
+        polynomial_space.evaluate(mapping_data.quadrature_points[i],
+                                  values,
+                                  grads,
+                                  grad_grads,
+                                  empty_vector_of_3rd_order_tensors,
+                                  empty_vector_of_4th_order_tensors);
 
         if (fe_internal.update_each & update_values)
           for (unsigned int k = 0; k < this->dofs_per_cell; ++k)
index 4c118116d39857bcfdf5038ff53f768a5f671745..b4e5d22f8b8307918251289332c4d916e4be7317 100644 (file)
@@ -296,7 +296,7 @@ namespace internal
             data.shape_fourth_derivatives.size() != 0)
           for (unsigned int point = 0; point < n_points; ++point)
             {
-              tensor_pols.compute(
+              tensor_pols.evaluate(
                 unit_points[point], values, grads, grad2, grad3, grad4);
 
               if (data.shape_values.size() != 0)
index 66cce181ef464fcfd9ea2d43aef0147ab4b6e41a..e6f3327a91333c93b8fb5785ef896292e3ae8859 100644 (file)
@@ -42,7 +42,7 @@ check_poly(const Point<dim> &     x,
   std::vector<Tensor<3, dim>> third1(n), third2(n);
   std::vector<Tensor<4, dim>> fourth1(n), fourth2(n);
 
-  p.compute(x, values1, gradients1, second1, third1, fourth1);
+  p.evaluate(x, values1, gradients1, second1, third1, fourth1);
   q.compute(x, values2, gradients2, second2, third2, fourth2);
 
   for (unsigned int k = 0; k < n; ++k)
index 2e18d9fe37595d76b574315adf165f902905b7a3..7229f64a37ea819dabfc23c1a34a4217ea0295c9 100644 (file)
@@ -41,7 +41,7 @@ check_poly(const Point<dim> &x, const PolynomialType &p)
   std::vector<Tensor<3, dim>> third(n);
   std::vector<Tensor<4, dim>> fourth(n);
 
-  p.compute(x, values, gradients, second, third, fourth);
+  p.evaluate(x, values, gradients, second, third, fourth);
 
   for (unsigned int k = 0; k < n; ++k)
     {

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.