]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Apply suggested changes
authorJulius Witte <julius.witte@iwr.uni-heidelberg.de>
Fri, 13 Mar 2020 08:02:34 +0000 (09:02 +0100)
committerJulius Witte <julius.witte@iwr.uni-heidelberg.de>
Fri, 13 Mar 2020 08:02:34 +0000 (09:02 +0100)
include/deal.II/matrix_free/shape_info.h
include/deal.II/matrix_free/shape_info.templates.h

index d143e4044177487cdf1eca36e160e756b0ae2485..26d22a0f7c806410c53e3924c2a10070114cfb26 100644 (file)
@@ -97,6 +97,8 @@ namespace internal
     template <typename Number>
     struct UnivariateShapeData
     {
+      UnivariateShapeData();
+
       /**
        * Return the memory consumption of this class in bytes.
        */
@@ -108,7 +110,7 @@ namespace internal
        * will select the most efficient algorithm based on the given element
        * type.
        */
-      ElementType element_type = ElementType::tensor_general;
+      ElementType element_type;
 
       /**
        * Stores the shape values of the 1D finite element evaluated on all 1D
@@ -241,18 +243,18 @@ namespace internal
       /**
        * Stores the degree of the element.
        */
-      unsigned int fe_degree = 0;
+      unsigned int fe_degree;
 
       /**
        * Stores the number of quadrature points per dimension.
        */
-      unsigned int n_q_points_1d = 0;
+      unsigned int n_q_points_1d;
 
       /**
        * Indicates whether the basis functions are nodal in 0 and 1, i.e., the
        * end points of the unit cell.
        */
-      bool nodal_at_cell_boundaries = false;
+      bool nodal_at_cell_boundaries;
     };
 
 
index da01d690184ab0d06af2249f1398393937c529f1..37c15930bbcd1fddac0facd3d39f5da9a629b4e2 100644 (file)
@@ -42,6 +42,18 @@ namespace internal
 {
   namespace MatrixFreeFunctions
   {
+    // ----------------- actual UnivariateShapeData functions
+    // --------------------
+
+    template <typename Number>
+    UnivariateShapeData<Number>::UnivariateShapeData()
+      : element_type(tensor_general)
+      , fe_degree(0)
+      , n_q_points_1d(0)
+      , nodal_at_cell_boundaries(false)
+    {}
+
+
     // ----------------- actual ShapeInfo functions --------------------
 
     template <typename Number>
@@ -87,31 +99,30 @@ namespace internal
 
       // assuming isotropy of dimensions and components
       data.resize(1);
-      UnivariateShapeData<Number> *univariate_shape_data = &data.front();
+      UnivariateShapeData<Number> &univariate_shape_data = data.front();
       data_access.reinit(n_dimensions, n_components);
-      data_access.fill(univariate_shape_data);
-      univariate_shape_data->quadrature    = quad;
-      univariate_shape_data->fe_degree     = fe->degree;
-      univariate_shape_data->n_q_points_1d = quad.size();
+      data_access.fill(&univariate_shape_data);
+      univariate_shape_data.quadrature    = quad;
+      univariate_shape_data.fe_degree     = fe->degree;
+      univariate_shape_data.n_q_points_1d = quad.size();
 
       // grant write access to common univariate shape data
-      auto &shape_values    = univariate_shape_data->shape_values;
-      auto &shape_gradients = univariate_shape_data->shape_gradients;
-      auto &shape_hessians  = univariate_shape_data->shape_hessians;
+      auto &shape_values    = univariate_shape_data.shape_values;
+      auto &shape_gradients = univariate_shape_data.shape_gradients;
+      auto &shape_hessians  = univariate_shape_data.shape_hessians;
       auto &shape_gradients_collocation =
-        univariate_shape_data->shape_gradients_collocation;
+        univariate_shape_data.shape_gradients_collocation;
       auto &shape_hessians_collocation =
-        univariate_shape_data->shape_hessians_collocation;
-      auto &inverse_shape_values = univariate_shape_data->inverse_shape_values;
-      auto &shape_data_on_face   = univariate_shape_data->shape_data_on_face;
-      auto &values_within_subface =
-        univariate_shape_data->values_within_subface;
+        univariate_shape_data.shape_hessians_collocation;
+      auto &inverse_shape_values  = univariate_shape_data.inverse_shape_values;
+      auto &shape_data_on_face    = univariate_shape_data.shape_data_on_face;
+      auto &values_within_subface = univariate_shape_data.values_within_subface;
       auto &gradients_within_subface =
-        univariate_shape_data->gradients_within_subface;
+        univariate_shape_data.gradients_within_subface;
       auto &hessians_within_subface =
-        univariate_shape_data->hessians_within_subface;
+        univariate_shape_data.hessians_within_subface;
       auto &nodal_at_cell_boundaries =
-        univariate_shape_data->nodal_at_cell_boundaries;
+        univariate_shape_data.nodal_at_cell_boundaries;
 
       const unsigned int fe_degree     = fe->degree;
       const unsigned int n_q_points_1d = quad.size();
@@ -418,9 +429,9 @@ namespace internal
         }
 
       if (element_type == tensor_general &&
-          check_1d_shapes_symmetric(*univariate_shape_data))
+          check_1d_shapes_symmetric(univariate_shape_data))
         {
-          if (check_1d_shapes_collocation(*univariate_shape_data))
+          if (check_1d_shapes_collocation(univariate_shape_data))
             element_type = tensor_symmetric_collocation;
           else
             element_type = tensor_symmetric;
@@ -439,7 +450,7 @@ namespace internal
             }
         }
       else if (element_type == tensor_symmetric_plus_dg0)
-        check_1d_shapes_symmetric(*univariate_shape_data);
+        check_1d_shapes_symmetric(univariate_shape_data);
 
       nodal_at_cell_boundaries = true;
       for (unsigned int i = 1; i < n_dofs_1d; ++i)
@@ -522,8 +533,7 @@ namespace internal
             }
         }
 
-      // TODO !!!
-      univariate_shape_data->element_type = this->element_type;
+      univariate_shape_data.element_type = this->element_type;
     }
 
 

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.