]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Reimplement ProductType do better deal with qualified number types. 4951/head
authorJean-Paul Pelteret <jppelteret@gmail.com>
Thu, 24 Aug 2017 15:21:51 +0000 (17:21 +0200)
committerJean-Paul Pelteret <jppelteret@gmail.com>
Thu, 24 Aug 2017 15:54:25 +0000 (17:54 +0200)
With reference to #4950.

doc/news/changes/incompatibilities/20170824Jean-PaulPelteretWolfgangBangerth [new file with mode: 0644]
include/deal.II/base/template_constraints.h
include/deal.II/differentiation/ad/sacado_product_types.h
include/deal.II/fe/fe_values.h
source/fe/fe_values.cc

diff --git a/doc/news/changes/incompatibilities/20170824Jean-PaulPelteretWolfgangBangerth b/doc/news/changes/incompatibilities/20170824Jean-PaulPelteretWolfgangBangerth
new file mode 100644 (file)
index 0000000..850f9c0
--- /dev/null
@@ -0,0 +1,6 @@
+Changed: Specialization of the ProductType class are now implemented through
+specialization of the class internal::ProductTypeImpl . This was done in order 
+to ensure that product operations performed with qualified number types do not
+result in the intended specializations being overlooked by the compiler.
+<br>
+(Jean-Paul Pelteret, Wolfgang Bangerth, 2017/08/24)
index 6171fd5ab3a9025feabaf8574e7d321a02c39422..5f2a01149dec598986f2489a80f52ca6c9df6413 100644 (file)
@@ -326,6 +326,31 @@ struct types_are_equal : std::is_same<T,U>
 
 
 
+namespace internal
+{
+
+  /**
+   * A struct that implements the default product type resulting from the
+   * multiplication of two types.
+   *
+   * @note Care should be taken when @p T or @p U have qualifiers (@p const or
+   * @p volatile) or are @p lvalue or @p rvalue references! It is recommended
+   * that specialization of this class is only made for unqualified (fully
+   * stripped) types and that the ProductType class be used to determine the
+   * result of operating with (potentially) qualified types.
+   *
+   * @author Wolfgang Bangerth, Jean-Paul Pelteret, 2017
+   */
+  template <typename T, typename U>
+  struct ProductTypeImpl
+  {
+    typedef decltype(std::declval<T>() * std::declval<U>()) type;
+  };
+
+}
+
+
+
 /**
  * A class with a local typedef that represents the type that results from the
  * product of two variables of type @p T and @p U. In other words, we would
@@ -372,54 +397,60 @@ struct types_are_equal : std::is_same<T,U>
  * used for the result of computing the product of unknowns and the values,
  * gradients, or other properties of shape functions.
  *
- * @author Wolfgang Bangerth, 2015
+ * @author Wolfgang Bangerth, 2015, 2017
  */
 template <typename T, typename U>
 struct ProductType
 {
-  typedef decltype(std::declval<T>() * std::declval<U>()) type;
+  typedef typename internal::ProductTypeImpl<
+  typename std::decay<T>::type, typename std::decay<U>::type>::type type;
 };
 
-// Annoyingly, there is no std::complex<T>::operator*(U) for scalars U
-// other than T (not even in C++11, or C++14). We provide our own overloads
-// in base/complex_overloads.h, but in order for them to work, we have to
-// manually specify all products we want to allow:
-
-template <typename T>
-struct ProductType<std::complex<T>,std::complex<T> >
+namespace internal
 {
-  typedef std::complex<T> type;
-};
 
-template <typename T, typename U>
-struct ProductType<std::complex<T>,std::complex<U> >
-{
-  typedef std::complex<typename ProductType<T,U>::type> type;
-};
+  // Annoyingly, there is no std::complex<T>::operator*(U) for scalars U
+  // other than T (not even in C++11, or C++14). We provide our own overloads
+  // in base/complex_overloads.h, but in order for them to work, we have to
+  // manually specify all products we want to allow:
 
-template <typename U>
-struct ProductType<double,std::complex<U> >
-{
-  typedef std::complex<typename ProductType<double,U>::type> type;
-};
+  template <typename T>
+  struct ProductTypeImpl<std::complex<T>,std::complex<T> >
+  {
+    typedef std::complex<T> type;
+  };
 
-template <typename T>
-struct ProductType<std::complex<T>,double>
-{
-  typedef std::complex<typename ProductType<T,double>::type> type;
-};
+  template <typename T, typename U>
+  struct ProductTypeImpl<std::complex<T>,std::complex<U> >
+  {
+    typedef std::complex<typename ProductType<T,U>::type> type;
+  };
 
-template <typename U>
-struct ProductType<float,std::complex<U> >
-{
-  typedef std::complex<typename ProductType<float,U>::type> type;
-};
+  template <typename U>
+  struct ProductTypeImpl<double,std::complex<U> >
+  {
+    typedef std::complex<typename ProductType<double,U>::type> type;
+  };
 
-template <typename T>
-struct ProductType<std::complex<T>,float>
-{
-  typedef std::complex<typename ProductType<T,float>::type> type;
-};
+  template <typename T>
+  struct ProductTypeImpl<std::complex<T>,double>
+  {
+    typedef std::complex<typename ProductType<T,double>::type> type;
+  };
+
+  template <typename U>
+  struct ProductTypeImpl<float,std::complex<U> >
+  {
+    typedef std::complex<typename ProductType<float,U>::type> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<std::complex<T>,float>
+  {
+    typedef std::complex<typename ProductType<T,float>::type> type;
+  };
+
+}
 
 
 
index 54de3762b19d2e603d701cefb78268586439256a..0ff6ce32fcbb7e4404088e9cd5aaeca9cabd6ccc 100644 (file)
@@ -29,47 +29,52 @@ DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
 
 DEAL_II_NAMESPACE_OPEN
 
-template <typename T>
-struct ProductType<Sacado::Fad::DFad<T>, float>
-{
-  typedef Sacado::Fad::DFad<T> type;
-};
-
-template <typename T>
-struct ProductType<float, Sacado::Fad::DFad<T> >
-{
-  typedef Sacado::Fad::DFad<T> type;
-};
-
-template <typename T>
-struct ProductType<Sacado::Fad::DFad<T>, double>
-{
-  typedef Sacado::Fad::DFad<T> type;
-};
-
-template <typename T>
-struct ProductType<double, Sacado::Fad::DFad<T> >
+namespace internal
 {
-  typedef Sacado::Fad::DFad<T> type;
-};
 
-template <typename T>
-struct ProductType<Sacado::Fad::DFad<T>, int>
-{
-  typedef Sacado::Fad::DFad<T> type;
-};
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Fad::DFad<T>, float>
+  {
+    typedef Sacado::Fad::DFad<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<float, Sacado::Fad::DFad<T> >
+  {
+    typedef Sacado::Fad::DFad<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Fad::DFad<T>, double>
+  {
+    typedef Sacado::Fad::DFad<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<double, Sacado::Fad::DFad<T> >
+  {
+    typedef Sacado::Fad::DFad<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<Sacado::Fad::DFad<T>, int>
+  {
+    typedef Sacado::Fad::DFad<T> type;
+  };
+
+  template <typename T>
+  struct ProductTypeImpl<int, Sacado::Fad::DFad<T> >
+  {
+    typedef Sacado::Fad::DFad<T> type;
+  };
+
+  template <typename T, typename U>
+  struct ProductTypeImpl<Sacado::Fad::DFad<T>, Sacado::Fad::DFad<U> >
+  {
+    typedef Sacado::Fad::DFad<typename ProductType<T,U>::type > type;
+  };
 
-template <typename T>
-struct ProductType<int, Sacado::Fad::DFad<T> >
-{
-  typedef Sacado::Fad::DFad<T> type;
-};
-
-template <typename T, typename U>
-struct ProductType<Sacado::Fad::DFad<T>, Sacado::Fad::DFad<U> >
-{
-  typedef Sacado::Fad::DFad<typename ProductType<T,U>::type > type;
-};
+}
 
 template <typename T>
 struct EnableIfScalar<Sacado::Fad::DFad<T> >
index 3839ea98d07d23e8c498da7aa5ee3ce99ebd9ac0..6faeb04ade25e3834a48819422cb9e37ca07f9cd 100644 (file)
@@ -178,31 +178,31 @@ namespace FEValuesViews
        * A typedef for the data type of the product of a @p Number and the
        * values of the view the Scalar class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Scalar<dim,spacedim>::value_type>::type value_type;
+      typedef typename ProductType<Number, typename Scalar<dim,spacedim>::value_type>::type value_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * gradients of the view the Scalar class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Scalar<dim,spacedim>::gradient_type>::type gradient_type;
+      typedef typename ProductType<Number, typename Scalar<dim,spacedim>::gradient_type>::type gradient_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * laplacians of the view the Scalar class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Scalar<dim,spacedim>::value_type>::type laplacian_type;
+      typedef typename ProductType<Number, typename Scalar<dim,spacedim>::value_type>::type laplacian_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * hessians of the view the Scalar class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Scalar<dim,spacedim>::hessian_type>::type hessian_type;
+      typedef typename ProductType<Number, typename Scalar<dim,spacedim>::hessian_type>::type hessian_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * third derivatives of the view the Scalar class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Scalar<dim,spacedim>::third_derivative_type>::type third_derivative_type;
+      typedef typename ProductType<Number, typename Scalar<dim,spacedim>::third_derivative_type>::type third_derivative_type;
     };
 
     /**
@@ -595,49 +595,49 @@ namespace FEValuesViews
        * A typedef for the data type of the product of a @p Number and the
        * values of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::value_type>::type value_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::value_type>::type value_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * gradients of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::gradient_type>::type gradient_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::gradient_type>::type gradient_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * symmetric gradients of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::symmetric_gradient_type>::type symmetric_gradient_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::symmetric_gradient_type>::type symmetric_gradient_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * divergences of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::divergence_type>::type divergence_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::divergence_type>::type divergence_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * laplacians of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::value_type>::type laplacian_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::value_type>::type laplacian_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * curls of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::curl_type>::type curl_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::curl_type>::type curl_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * hessians of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::hessian_type>::type hessian_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::hessian_type>::type hessian_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * third derivatives of the view the Vector class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Vector<dim,spacedim>::third_derivative_type>::type third_derivative_type;
+      typedef typename ProductType<Number, typename Vector<dim,spacedim>::third_derivative_type>::type third_derivative_type;
     };
 
     /**
@@ -1157,13 +1157,13 @@ namespace FEValuesViews
        * A typedef for the data type of the product of a @p Number and the
        * values of the view the SymmetricTensor class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename SymmetricTensor<2,dim,spacedim>::value_type>::type value_type;
+      typedef typename ProductType<Number, typename SymmetricTensor<2,dim,spacedim>::value_type>::type value_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * divergences of the view the SymmetricTensor class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename SymmetricTensor<2,dim,spacedim>::divergence_type>::type divergence_type;
+      typedef typename ProductType<Number, typename SymmetricTensor<2,dim,spacedim>::divergence_type>::type divergence_type;
     };
 
     /**
@@ -1413,13 +1413,13 @@ namespace FEValuesViews
        * A typedef for the data type of the product of a @p Number and the
        * values of the view the Tensor class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Tensor<2,dim,spacedim>::value_type>::type value_type;
+      typedef typename ProductType<Number, typename Tensor<2,dim,spacedim>::value_type>::type value_type;
 
       /**
        * A typedef for the data type of the product of a @p Number and the
        * divergences of the view the Tensor class.
        */
-      typedef typename ProductType<typename std::decay<Number>::type, typename Tensor<2,dim,spacedim>::divergence_type>::type divergence_type;
+      typedef typename ProductType<Number, typename Tensor<2,dim,spacedim>::divergence_type>::type divergence_type;
     };
 
     /**
index 0e2ac8ae7966c2af1558ca8079c1790562b642f5..f1ae74645b97f3756104eb102fd029dd99314bf0 100644 (file)
@@ -462,14 +462,15 @@ namespace FEValuesViews
     do_function_values (const ArrayView<Number> &dof_values,
                         const Table<2,double> &shape_values,
                         const std::vector<typename Scalar<dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                        std::vector<typename ProductType<typename std::decay<Number>::type,double>::type>            &values)
+                        std::vector<typename ProductType<Number,double>::type>            &values)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
                                                shape_values.n_cols() : values.size();
       AssertDimension (values.size(), n_quadrature_points);
 
-      std::fill (values.begin(), values.end(), dealii::internal::NumberType<typename std::decay<Number>::type>::value(0.0));
+      std::fill (values.begin(), values.end(),
+                 dealii::internal::NumberType<typename std::decay<Number>::type>::value(0.0));
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -495,7 +496,7 @@ namespace FEValuesViews
     do_function_derivatives (const ArrayView<Number> &dof_values,
                              const Table<2,dealii::Tensor<order,spacedim> > &shape_derivatives,
                              const std::vector<typename Scalar<dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                             std::vector<typename ProductType<typename std::decay<Number>::type,dealii::Tensor<order,spacedim> >::type> &derivatives)
+                             std::vector<typename ProductType<Number,dealii::Tensor<order,spacedim> >::type> &derivatives)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
@@ -503,7 +504,7 @@ namespace FEValuesViews
       AssertDimension (derivatives.size(), n_quadrature_points);
 
       std::fill (derivatives.begin(), derivatives.end(),
-                 typename ProductType<typename std::decay<Number>::type,dealii::Tensor<order,spacedim> >::type());
+                 typename ProductType<Number,dealii::Tensor<order,spacedim> >::type());
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -561,14 +562,15 @@ namespace FEValuesViews
     void do_function_values (const ArrayView<Number> &dof_values,
                              const Table<2,double>          &shape_values,
                              const std::vector<typename Vector<dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                             std::vector<typename ProductType<typename std::decay<Number>::type,dealii::Tensor<1,spacedim> >::type> &values)
+                             std::vector<typename ProductType<Number,dealii::Tensor<1,spacedim> >::type> &values)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
                                                shape_values.n_cols() : values.size();
       AssertDimension (values.size(), n_quadrature_points);
 
-      std::fill (values.begin(), values.end(), typename ProductType<typename std::decay<Number>::type,dealii::Tensor<1,spacedim> >::type());
+      std::fill (values.begin(), values.end(),
+                 typename ProductType<Number,dealii::Tensor<1,spacedim> >::type());
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -610,7 +612,7 @@ namespace FEValuesViews
     do_function_derivatives (const ArrayView<Number> &dof_values,
                              const Table<2,dealii::Tensor<order,spacedim> > &shape_derivatives,
                              const std::vector<typename Vector<dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                             std::vector<typename ProductType<typename std::decay<Number>::type,dealii::Tensor<order+1,spacedim> >::type> &derivatives)
+                             std::vector<typename ProductType<Number,dealii::Tensor<order+1,spacedim> >::type> &derivatives)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
@@ -618,7 +620,7 @@ namespace FEValuesViews
       AssertDimension (derivatives.size(), n_quadrature_points);
 
       std::fill (derivatives.begin(), derivatives.end(),
-                 typename ProductType<typename std::decay<Number>::type,dealii::Tensor<order+1,spacedim> >::type());
+                 typename ProductType<Number,dealii::Tensor<order+1,spacedim> >::type());
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -664,7 +666,7 @@ namespace FEValuesViews
     do_function_symmetric_gradients (const ArrayView<Number> &dof_values,
                                      const Table<2,dealii::Tensor<1,spacedim> > &shape_gradients,
                                      const std::vector<typename Vector<dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                                     std::vector<typename ProductType<typename std::decay<Number>::type,dealii::SymmetricTensor<2,spacedim> >::type> &symmetric_gradients)
+                                     std::vector<typename ProductType<Number,dealii::SymmetricTensor<2,spacedim> >::type> &symmetric_gradients)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
@@ -672,7 +674,7 @@ namespace FEValuesViews
       AssertDimension (symmetric_gradients.size(), n_quadrature_points);
 
       std::fill (symmetric_gradients.begin(), symmetric_gradients.end(),
-                 typename ProductType<typename std::decay<Number>::type,dealii::SymmetricTensor<2,spacedim> >::type());
+                 typename ProductType<Number,dealii::SymmetricTensor<2,spacedim> >::type());
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -768,7 +770,7 @@ namespace FEValuesViews
     do_function_curls (const ArrayView<Number> &dof_values,
                        const Table<2,dealii::Tensor<1,spacedim> > &shape_gradients,
                        const std::vector<typename Vector<dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                       std::vector<typename ProductType<typename std::decay<Number>::type,typename dealii::internal::CurlType<spacedim>::type>::type> &curls)
+                       std::vector<typename ProductType<Number,typename dealii::internal::CurlType<spacedim>::type>::type> &curls)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
@@ -776,7 +778,7 @@ namespace FEValuesViews
       AssertDimension (curls.size(), n_quadrature_points);
 
       std::fill (curls.begin(), curls.end(),
-                 typename ProductType<typename std::decay<Number>::type,typename dealii::internal::CurlType<spacedim>::type>::type());
+                 typename ProductType<Number,typename dealii::internal::CurlType<spacedim>::type>::type());
 
       switch (spacedim)
         {
@@ -1015,7 +1017,7 @@ namespace FEValuesViews
     do_function_values (const ArrayView<Number> &dof_values,
                         const dealii::Table<2,double>          &shape_values,
                         const std::vector<typename SymmetricTensor<2,dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                        std::vector<typename ProductType<typename std::decay<Number>::type,dealii::SymmetricTensor<2,spacedim> >::type> &values)
+                        std::vector<typename ProductType<Number,dealii::SymmetricTensor<2,spacedim> >::type> &values)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
@@ -1023,7 +1025,7 @@ namespace FEValuesViews
       AssertDimension (values.size(), n_quadrature_points);
 
       std::fill (values.begin(), values.end(),
-                 typename ProductType<typename std::decay<Number>::type,dealii::SymmetricTensor<2,spacedim> >::type());
+                 typename ProductType<Number,dealii::SymmetricTensor<2,spacedim> >::type());
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -1157,7 +1159,7 @@ namespace FEValuesViews
     do_function_values (const ArrayView<Number> &dof_values,
                         const dealii::Table<2,double>          &shape_values,
                         const std::vector<typename Tensor<2,dim,spacedim>::ShapeFunctionData> &shape_function_data,
-                        std::vector<typename ProductType<typename std::decay<Number>::type,dealii::Tensor<2,spacedim> >::type> &values)
+                        std::vector<typename ProductType<Number,dealii::Tensor<2,spacedim> >::type> &values)
     {
       const unsigned int dofs_per_cell = dof_values.size();
       const unsigned int n_quadrature_points = dofs_per_cell > 0 ?
@@ -1165,7 +1167,7 @@ namespace FEValuesViews
       AssertDimension (values.size(), n_quadrature_points);
 
       std::fill (values.begin(), values.end(),
-                 typename ProductType<typename std::decay<Number>::type,dealii::Tensor<2,spacedim> >::type());
+                 typename ProductType<Number,dealii::Tensor<2,spacedim> >::type());
 
       for (unsigned int shape_function=0;
            shape_function<dofs_per_cell; ++shape_function)
@@ -2674,7 +2676,8 @@ namespace internal
     AssertDimension(values.size(), n_quadrature_points);
 
     // initialize with zero
-    std::fill_n (values.begin(), n_quadrature_points, dealii::internal::NumberType<typename std::decay<Number>::type>::value(0.0));
+    std::fill_n (values.begin(), n_quadrature_points,
+                 dealii::internal::NumberType<typename std::decay<Number>::type>::value(0.0));
 
     // add up contributions of trial functions. note that here we deal with
     // scalar finite elements, so no need to check for non-primitivity of
@@ -2936,7 +2939,8 @@ namespace internal
     AssertDimension(laplacians.size(), n_quadrature_points);
 
     // initialize with zero
-    std::fill_n (laplacians.begin(), n_quadrature_points, dealii::internal::NumberType<typename std::decay<Number>::type>::value(0.0));
+    std::fill_n (laplacians.begin(), n_quadrature_points,
+                 dealii::internal::NumberType<typename std::decay<Number>::type>::value(0.0));
 
     // add up contributions of trial functions. note that here we deal with
     // scalar finite elements and also note that the Laplacian is

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.