From: Wolfgang Bangerth <bangerth@math.tamu.edu>
Date: Fri, 14 Aug 2015 01:52:45 +0000 (-0500)
Subject: Fix #1343.
X-Git-Tag: v8.4.0-rc2~610^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7ca15876dccb2bdd52b2d7c9f6d6460014863852;p=dealii.git

Fix #1343.

#1343 copied an *old* version of the code in .cc and .templates.h
files into the .h files, thereby undoing a significant portion of
the work of Hamed Maien in #1322. This restores this by using the
correct code and should thereby fix the ~700 testsuite failures that
#1343 introduced.
---

diff --git a/include/deal.II/fe/fe_face.h b/include/deal.II/fe/fe_face.h
index baeb46664b..352bc36399 100644
--- a/include/deal.II/fe/fe_face.h
+++ b/include/deal.II/fe/fe_face.h
@@ -248,8 +248,8 @@ protected:
                 const Quadrature<0> &quadrature) const
   {
     // generate a new data object and initialize some fields
-    typename FiniteElement<1, spacedim>::InternalDataBase *data =
-      new typename FiniteElement<1, spacedim>::InternalDataBase;
+    typename FiniteElement<1,spacedim>::InternalDataBase *data =
+      new typename FiniteElement<1,spacedim>::InternalDataBase;
 
     // check what needs to be initialized only once and what on every
     // cell/face/subface we visit
diff --git a/include/deal.II/fe/fe_poly.h b/include/deal.II/fe/fe_poly.h
index 4f90876b05..410327cfca 100644
--- a/include/deal.II/fe/fe_poly.h
+++ b/include/deal.II/fe/fe_poly.h
@@ -163,14 +163,16 @@ public:
                                                    const unsigned int component) const;
 
 protected:
-  /**
-  NOTE: The following function has its definition inlined into the class declaration
-    * because we otherwise run into a compiler error with MS Visual Studio.
-    */
+  /*
+   * NOTE: The following function has its definition inlined into the class declaration
+   * because we otherwise run into a compiler error with MS Visual Studio.
+   */
+
+
   virtual
   typename FiniteElement<dim,spacedim>::InternalDataBase *
   get_data(const UpdateFlags update_flags,
-           const Mapping<dim,spacedim> &mapping,
+           const Mapping<dim,spacedim> &/*mapping*/,
            const Quadrature<dim> &quadrature) const
   {
     // generate a new data object and
@@ -190,24 +192,24 @@ protected:
 
     // some scratch arrays
     std::vector<double> values(0);
-    std::vector<Tensor<1, dim> > grads(0);
-    std::vector<Tensor<2, dim> > grad_grads(0);
+    std::vector<Tensor<1,dim> > grads(0);
+    std::vector<Tensor<2,dim> > grad_grads(0);
 
     // initialize fields only if really
     // necessary. otherwise, don't
     // allocate memory
     if (flags & update_values)
       {
-        values.resize(this->dofs_per_cell);
-        data->shape_values.resize(this->dofs_per_cell,
-                                  std::vector<double>(n_q_points));
+        values.resize (this->dofs_per_cell);
+        data->shape_values.resize (this->dofs_per_cell,
+                                   std::vector<double> (n_q_points));
       }
 
     if (flags & update_gradients)
       {
-        grads.resize(this->dofs_per_cell);
-        data->shape_gradients.resize(this->dofs_per_cell,
-                                     std::vector<Tensor<1, dim> >(n_q_points));
+        grads.resize (this->dofs_per_cell);
+        data->shape_gradients.resize (this->dofs_per_cell,
+                                      std::vector<Tensor<1,dim> > (n_q_points));
       }
 
     // if second derivatives through
@@ -215,7 +217,12 @@ protected:
     // then initialize some objects for
     // that
     if (flags & update_hessians)
-      data->initialize_2nd(this, mapping, quadrature);
+      {
+        grad_grads.resize (this->dofs_per_cell);
+        data->shape_hessians.resize (this->dofs_per_cell,
+                                     std::vector<Tensor<2,dim> > (n_q_points));
+        data->untransformed_shape_hessians.resize (n_q_points);
+      }
 
     // next already fill those fields
     // of which we have information by
@@ -224,19 +231,23 @@ protected:
     // unit cell, and need to be
     // transformed when visiting an
     // actual cell
-    if (flags & (update_values | update_gradients))
-      for (unsigned int i = 0; i<n_q_points; ++i)
+    if (flags & (update_values | update_gradients | update_hessians))
+      for (unsigned int i=0; i<n_q_points; ++i)
         {
           poly_space.compute(quadrature.point(i),
                              values, grads, grad_grads);
 
           if (flags & update_values)
-            for (unsigned int k = 0; k<this->dofs_per_cell; ++k)
+            for (unsigned int k=0; k<this->dofs_per_cell; ++k)
               data->shape_values[k][i] = values[k];
 
           if (flags & update_gradients)
-            for (unsigned int k = 0; k<this->dofs_per_cell; ++k)
+            for (unsigned int k=0; k<this->dofs_per_cell; ++k)
               data->shape_gradients[k][i] = grads[k];
+
+          if (flags & update_hessians)
+            for (unsigned int k=0; k<this->dofs_per_cell; ++k)
+              data->shape_hessians[k][i] = grad_grads[k];
         }
     return data;
   }
diff --git a/include/deal.II/fe/fe_poly_face.h b/include/deal.II/fe/fe_poly_face.h
index c54b97f009..0dadbca149 100644
--- a/include/deal.II/fe/fe_poly_face.h
+++ b/include/deal.II/fe/fe_poly_face.h
@@ -113,23 +113,23 @@ protected:
 
     // some scratch arrays
     std::vector<double> values(0);
-    std::vector<Tensor<1, dim - 1> > grads(0);
-    std::vector<Tensor<2, dim - 1> > grad_grads(0);
+    std::vector<Tensor<1,dim-1> > grads(0);
+    std::vector<Tensor<2,dim-1> > grad_grads(0);
 
     // initialize fields only if really
     // necessary. otherwise, don't
     // allocate memory
     if (flags & update_values)
       {
-        values.resize(poly_space.n());
-        data->shape_values.resize(poly_space.n(),
-                                  std::vector<double>(n_q_points));
-        for (unsigned int i = 0; i<n_q_points; ++i)
+        values.resize (poly_space.n());
+        data->shape_values.resize (poly_space.n(),
+                                   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);
 
-            for (unsigned int k = 0; k<poly_space.n(); ++k)
+            for (unsigned int k=0; k<poly_space.n(); ++k)
               data->shape_values[k][i] = values[k];
           }
       }
diff --git a/include/deal.II/fe/fe_poly_tensor.h b/include/deal.II/fe/fe_poly_tensor.h
index 5c05271fb5..3252a07c01 100644
--- a/include/deal.II/fe/fe_poly_tensor.h
+++ b/include/deal.II/fe/fe_poly_tensor.h
@@ -200,27 +200,38 @@ protected:
     const unsigned int n_q_points = quadrature.size();
 
     // some scratch arrays
-    std::vector<Tensor<1, dim> > values(0);
-    std::vector<Tensor<2, dim> > grads(0);
-    std::vector<Tensor<3, dim> > grad_grads(0);
+    std::vector<Tensor<1,dim> > values(0);
+    std::vector<Tensor<2,dim> > grads(0);
+    std::vector<Tensor<3,dim> > grad_grads(0);
+
+    data->sign_change.resize (this->dofs_per_cell);
 
     // initialize fields only if really
     // necessary. otherwise, don't
     // allocate memory
     if (flags & update_values)
       {
-        values.resize(this->dofs_per_cell);
-        data->shape_values.resize(this->dofs_per_cell);
-        for (unsigned int i = 0; i<this->dofs_per_cell; ++i)
-          data->shape_values[i].resize(n_q_points);
+        values.resize (this->dofs_per_cell);
+        data->shape_values.resize (this->dofs_per_cell);
+        for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+          data->shape_values[i].resize (n_q_points);
+        if (mapping_type != mapping_none)
+          data->transformed_shape_values.resize (n_q_points);
       }
 
     if (flags & update_gradients)
       {
-        grads.resize(this->dofs_per_cell);
-        data->shape_grads.resize(this->dofs_per_cell);
-        for (unsigned int i = 0; i<this->dofs_per_cell; ++i)
-          data->shape_grads[i].resize(n_q_points);
+        grads.resize (this->dofs_per_cell);
+        data->shape_grads.resize (this->dofs_per_cell);
+        for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+          data->shape_grads[i].resize (n_q_points);
+        data->transformed_shape_grads.resize (n_q_points);
+        if ( (mapping_type == mapping_raviart_thomas)
+             ||
+             (mapping_type == mapping_piola)
+             ||
+             (mapping_type == mapping_nedelec))
+          data->untransformed_shape_grads.resize(n_q_points);
       }
 
     // if second derivatives through
@@ -229,8 +240,8 @@ protected:
     // that
     if (flags & update_hessians)
       {
-        //      grad_grads.resize (this->dofs_per_cell);
-        data->initialize_2nd(this, mapping, quadrature);
+//      grad_grads.resize (this->dofs_per_cell);
+        data->initialize_2nd (this, mapping, quadrature);
       }
 
     // Compute shape function values
@@ -240,7 +251,7 @@ protected:
     // N_i(v_j)=\delta_ij for all basis
     // functions v_j
     if (flags & (update_values | update_gradients))
-      for (unsigned int k = 0; k<n_q_points; ++k)
+      for (unsigned int k=0; k<n_q_points; ++k)
         {
           poly_space.compute(quadrature.point(k),
                              values, grads, grad_grads);
@@ -248,14 +259,14 @@ protected:
           if (flags & update_values)
             {
               if (inverse_node_matrix.n_cols() == 0)
-                for (unsigned int i = 0; i<this->dofs_per_cell; ++i)
+                for (unsigned int i=0; i<this->dofs_per_cell; ++i)
                   data->shape_values[i][k] = values[i];
               else
-                for (unsigned int i = 0; i<this->dofs_per_cell; ++i)
+                for (unsigned int i=0; i<this->dofs_per_cell; ++i)
                   {
-                    Tensor<1, dim> add_values;
-                    for (unsigned int j = 0; j<this->dofs_per_cell; ++j)
-                      add_values += inverse_node_matrix(j, i) * values[j];
+                    Tensor<1,dim> add_values;
+                    for (unsigned int j=0; j<this->dofs_per_cell; ++j)
+                      add_values += inverse_node_matrix(j,i) * values[j];
                     data->shape_values[i][k] = add_values;
                   }
             }
@@ -263,14 +274,14 @@ protected:
           if (flags & update_gradients)
             {
               if (inverse_node_matrix.n_cols() == 0)
-                for (unsigned int i = 0; i<this->dofs_per_cell; ++i)
+                for (unsigned int i=0; i<this->dofs_per_cell; ++i)
                   data->shape_grads[i][k] = grads[i];
               else
-                for (unsigned int i = 0; i<this->dofs_per_cell; ++i)
+                for (unsigned int i=0; i<this->dofs_per_cell; ++i)
                   {
-                    Tensor<2, dim> add_grads;
-                    for (unsigned int j = 0; j<this->dofs_per_cell; ++j)
-                      add_grads += inverse_node_matrix(j, i) * grads[j];
+                    Tensor<2,dim> add_grads;
+                    for (unsigned int j=0; j<this->dofs_per_cell; ++j)
+                      add_grads += inverse_node_matrix(j,i) * grads[j];
                     data->shape_grads[i][k] = add_grads;
                   }
             }