]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Rename ShapeInfo::dofs_per_cell -> ShapeInfo::dofs_per_component_on_cell
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 28 Oct 2017 16:00:47 +0000 (18:00 +0200)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Sat, 28 Oct 2017 16:06:55 +0000 (18:06 +0200)
include/deal.II/matrix_free/evaluation_kernels.h
include/deal.II/matrix_free/fe_evaluation.h
include/deal.II/matrix_free/shape_info.h
include/deal.II/matrix_free/shape_info.templates.h

index 2c890e8cde29bbf34a1944c46090254aabc76389..f5bc0cd9e1ea1753f9cf230a82b1b9d893764bb1 100644 (file)
@@ -175,7 +175,7 @@ namespace internal
       {
         values_dofs = expanded_dof_values;
         for (unsigned int c=0; c<n_components; ++c)
-          expanded_dof_values[c] = scratch_data+2*(std::max(shape_info.dofs_per_cell,
+          expanded_dof_values[c] = scratch_data+2*(std::max(shape_info.dofs_per_component_on_cell,
                                                             shape_info.n_q_points)) +
                                    c*Utilities::fixed_power<dim>(shape_info.fe_degree+1);
         const int degree = fe_degree != -1 ? fe_degree : shape_info.fe_degree;
@@ -334,7 +334,7 @@ namespace internal
     if (type == MatrixFreeFunctions::tensor_symmetric_plus_dg0 && evaluate_values)
       for (unsigned int c=0; c<n_components; ++c)
         for (unsigned int q=0; q<shape_info.n_q_points; ++q)
-          values_quad[c][q] += values_dofs[c][shape_info.dofs_per_cell-1];
+          values_quad[c][q] += values_dofs[c][shape_info.dofs_per_component_on_cell-1];
   }
 
 
@@ -389,7 +389,7 @@ namespace internal
       {
         values_dofs = expanded_dof_values;
         for (unsigned int c=0; c<n_components; ++c)
-          expanded_dof_values[c] = scratch_data+2*(std::max(shape_info.dofs_per_cell,
+          expanded_dof_values[c] = scratch_data+2*(std::max(shape_info.dofs_per_component_on_cell,
                                                             shape_info.n_q_points)) +
                                    c*Utilities::fixed_power<dim>(shape_info.fe_degree+1);
       }
@@ -492,13 +492,13 @@ namespace internal
         if (integrate_values)
           for (unsigned int c=0; c<n_components; ++c)
             {
-              values_dofs[c][shape_info.dofs_per_cell-1] = values_quad[c][0];
+              values_dofs[c][shape_info.dofs_per_component_on_cell-1] = values_quad[c][0];
               for (unsigned int q=1; q<shape_info.n_q_points; ++q)
-                values_dofs[c][shape_info.dofs_per_cell-1] += values_quad[c][q];
+                values_dofs[c][shape_info.dofs_per_component_on_cell-1] += values_quad[c][q];
             }
         else
           for (unsigned int c=0; c<n_components; ++c)
-            values_dofs[c][shape_info.dofs_per_cell-1] = VectorizedArray<Number>();
+            values_dofs[c][shape_info.dofs_per_component_on_cell-1] = VectorizedArray<Number>();
       }
 
     if (type == MatrixFreeFunctions::truncated_tensor)
index f1121778efb26bc9dd8bcd9c63b52ce9d7532055..ca1eff12f45351641c3ff4b8ab4d8293ae0be92c 100644 (file)
@@ -1914,8 +1914,8 @@ public:
   static const unsigned int n_components  = n_components_;
   static const unsigned int static_n_q_points    = Utilities::fixed_int_power<n_q_points_1d,dim>::value;
   static const unsigned int static_dofs_per_component = Utilities::fixed_int_power<fe_degree+1,dim>::value;
-  static const unsigned int tensor_dofs_per_cell = Utilities::fixed_int_power<fe_degree+1,dim>::value;
-  static const unsigned int static_dofs_per_cell = n_components * Utilities::fixed_int_power<fe_degree+1,dim>::value;
+  static const unsigned int tensor_dofs_per_cell = n_components *static_dofs_per_component;
+  static const unsigned int static_dofs_per_cell = n_components *static_dofs_per_component;
 
   /**
    * Constructor. Takes all data stored in MatrixFree. If applied to problems
@@ -2147,7 +2147,7 @@ FEEvaluationBase<dim,n_components_,Number>
           ExcNotInitialized());
   AssertDimension (matrix_info->get_size_info().vectorization_length,
                    VectorizedArray<Number>::n_array_elements);
-  AssertDimension (data->dofs_per_cell*n_fe_components,
+  AssertDimension (data->dofs_per_component_on_cell*n_fe_components,
                    dof_info->dofs_per_cell[active_fe_index]);
   AssertDimension (data->n_q_points,
                    mapping_info->mapping_data_gen[quad_no].n_q_points[active_quad_index]);
@@ -2394,7 +2394,7 @@ FEEvaluationBase<dim,n_components_,Number>
 
   const unsigned int tensor_dofs_per_component =
     Utilities::fixed_power<dim>(this->data->fe_degree+1);
-  const unsigned int dofs_per_component = this->data->dofs_per_cell;
+  const unsigned int dofs_per_component = this->data->dofs_per_component_on_cell;
   const unsigned int n_quadrature_points = this->data->n_q_points;
 
   const unsigned int shift = std::max(tensor_dofs_per_component+1, dofs_per_component)*
@@ -3008,10 +3008,10 @@ FEEvaluationBase<dim,n_components_,Number>
     {
       Assert (!local_dof_indices.empty(), ExcNotInitialized());
 
-      unsigned int index = first_selected_component * this->data->dofs_per_cell;
+      unsigned int index = first_selected_component * this->data->dofs_per_component_on_cell;
       for (unsigned int comp = 0; comp<n_components; ++comp)
         {
-          for (unsigned int i=0; i<this->data->dofs_per_cell; ++i, ++index)
+          for (unsigned int i=0; i<this->data->dofs_per_component_on_cell; ++i, ++index)
             {
               operation.process_dof_global(local_dof_indices[this->data->lexicographic_numbering[index]],
                                            *src[0], values_dofs[comp][i][0]);
@@ -3036,7 +3036,7 @@ FEEvaluationBase<dim,n_components_,Number>
   const std::pair<unsigned short,unsigned short> *indicators_end =
     dof_info->end_indicators(cell);
   unsigned int ind_local = 0;
-  const unsigned int dofs_per_component = this->data->dofs_per_cell;
+  const unsigned int dofs_per_component = this->data->dofs_per_component_on_cell;
 
   const unsigned int n_irreg_components_filled = dof_info->row_starts[cell][2];
   const bool at_irregular_cell = n_irreg_components_filled > 0;
@@ -3462,7 +3462,7 @@ FEEvaluationBase<dim,n_components_,Number>
   // iterates over the elements of index_local_to_global and dof_indices
   // points to the global indices stored in index_local_to_global
   const unsigned int *dof_indices = dof_info->begin_indices_plain(cell);
-  const unsigned int dofs_per_component = this->data->dofs_per_cell;
+  const unsigned int dofs_per_component = this->data->dofs_per_component_on_cell;
 
   const unsigned int n_irreg_components_filled = dof_info->row_starts[cell][2];
   const bool at_irregular_cell = n_irreg_components_filled > 0;
@@ -3705,7 +3705,7 @@ Tensor<1,n_components_,VectorizedArray<Number> >
 FEEvaluationBase<dim,n_components_,Number>
 ::get_dof_value (const unsigned int dof) const
 {
-  AssertIndexRange (dof, this->data->dofs_per_cell);
+  AssertIndexRange (dof, this->data->dofs_per_component_on_cell);
   Tensor<1,n_components_,VectorizedArray<Number> > return_value;
   for (unsigned int comp=0; comp<n_components; comp++)
     return_value[comp] = this->values_dofs[comp][dof];
@@ -4067,7 +4067,7 @@ FEEvaluationBase<dim,n_components_,Number>
 #ifdef DEBUG
   this->dof_values_initialized = true;
 #endif
-  AssertIndexRange (dof, this->data->dofs_per_cell);
+  AssertIndexRange (dof, this->data->dofs_per_component_on_cell);
   for (unsigned int comp=0; comp<n_components; comp++)
     this->values_dofs[comp][dof] = val_in[comp];
 }
@@ -4287,7 +4287,7 @@ VectorizedArray<Number>
 FEEvaluationAccess<dim,1,Number>
 ::get_dof_value (const unsigned int dof) const
 {
-  AssertIndexRange (dof, this->data->dofs_per_cell);
+  AssertIndexRange (dof, this->data->dofs_per_component_on_cell);
   return this->values_dofs[0][dof];
 }
 
@@ -4389,7 +4389,7 @@ FEEvaluationAccess<dim,1,Number>
 {
 #ifdef DEBUG
   this->dof_values_initialized = true;
-  AssertIndexRange (dof, this->data->dofs_per_cell);
+  AssertIndexRange (dof, this->data->dofs_per_component_on_cell);
 #endif
   this->values_dofs[0][dof] = val_in;
 }
@@ -4910,7 +4910,7 @@ VectorizedArray<Number>
 FEEvaluationAccess<1,1,Number>
 ::get_dof_value (const unsigned int dof) const
 {
-  AssertIndexRange (dof, this->data->dofs_per_cell);
+  AssertIndexRange (dof, this->data->dofs_per_component_on_cell);
   return this->values_dofs[0][dof];
 }
 
@@ -5007,7 +5007,7 @@ FEEvaluationAccess<1,1,Number>
 {
 #ifdef DEBUG
   this->dof_values_initialized = true;
-  AssertIndexRange (dof, this->data->dofs_per_cell);
+  AssertIndexRange (dof, this->data->dofs_per_component_on_cell);
 #endif
   this->values_dofs[0][dof] = val_in;
 }
@@ -5099,8 +5099,8 @@ FEEvaluation<dim,fe_degree,n_q_points_1d,n_components_,Number>
                 const unsigned int quad_no)
   :
   BaseClass (data_in, fe_no, quad_no, fe_degree, static_n_q_points),
-  dofs_per_component (this->data->dofs_per_cell),
-  dofs_per_cell (this->data->dofs_per_cell*n_components_),
+  dofs_per_component (this->data->dofs_per_component_on_cell),
+  dofs_per_cell (this->data->dofs_per_component_on_cell *n_components_),
   n_q_points (this->data->n_q_points)
 {
   check_template_arguments(fe_no, 0);
@@ -5121,8 +5121,8 @@ FEEvaluation<dim,fe_degree,n_q_points_1d,n_components_,Number>
   BaseClass (mapping, fe, quadrature, update_flags,
              first_selected_component,
              static_cast<FEEvaluationBase<dim,1,Number>*>(nullptr)),
-  dofs_per_component (this->data->dofs_per_cell),
-  dofs_per_cell (this->data->dofs_per_cell*n_components_),
+  dofs_per_component (this->data->dofs_per_component_on_cell),
+  dofs_per_cell (this->data->dofs_per_component_on_cell *n_components_),
   n_q_points (this->data->n_q_points)
 {
   check_template_arguments(numbers::invalid_unsigned_int, 0);
@@ -5142,8 +5142,8 @@ FEEvaluation<dim,fe_degree,n_q_points_1d,n_components_,Number>
   BaseClass (StaticMappingQ1<dim>::mapping, fe, quadrature, update_flags,
              first_selected_component,
              static_cast<FEEvaluationBase<dim,1,Number>*>(nullptr)),
-  dofs_per_component (this->data->dofs_per_cell),
-  dofs_per_cell (this->data->dofs_per_cell*n_components_),
+  dofs_per_component (this->data->dofs_per_component_on_cell),
+  dofs_per_cell (this->data->dofs_per_component_on_cell *n_components_),
   n_q_points (this->data->n_q_points)
 {
   check_template_arguments(numbers::invalid_unsigned_int, 0);
@@ -5164,8 +5164,8 @@ FEEvaluation<dim,fe_degree,n_q_points_1d,n_components_,Number>
              fe, other.mapped_geometry->get_quadrature(),
              other.mapped_geometry->get_fe_values().get_update_flags(),
              first_selected_component, &other),
-  dofs_per_component (this->data->dofs_per_cell),
-  dofs_per_cell (this->data->dofs_per_cell*n_components_),
+  dofs_per_component (this->data->dofs_per_component_on_cell),
+  dofs_per_cell (this->data->dofs_per_component_on_cell *n_components_),
   n_q_points (this->data->n_q_points)
 {
   check_template_arguments(numbers::invalid_unsigned_int, 0);
@@ -5180,8 +5180,8 @@ FEEvaluation<dim,fe_degree,n_q_points_1d,n_components_,Number>
 ::FEEvaluation (const FEEvaluation &other)
   :
   BaseClass (other),
-  dofs_per_component (this->data->dofs_per_cell),
-  dofs_per_cell (this->data->dofs_per_cell*n_components_),
+  dofs_per_component (this->data->dofs_per_component_on_cell),
+  dofs_per_cell (this->data->dofs_per_component_on_cell *n_components_),
   n_q_points (this->data->n_q_points)
 {
   check_template_arguments(numbers::invalid_unsigned_int, 0);
@@ -5328,7 +5328,7 @@ FEEvaluation<dim,fe_degree,n_q_points_1d,n_components_,Number>
       AssertDimension (n_q_points,
                        this->mapping_info->mapping_data_gen[this->quad_no].
                        n_q_points[this->active_quad_index]);
-      AssertDimension (this->data->dofs_per_cell * this->n_fe_components,
+      AssertDimension (this->data->dofs_per_component_on_cell * this->n_fe_components,
                        this->dof_info->dofs_per_cell[this->active_fe_index]);
     }
 #endif
index acb90f0040d5fae09bfc8663f6726f770d7d3ba6..582617a8ac93345c652f2664bba78fcec3f0aaa9 100644 (file)
@@ -237,9 +237,10 @@ namespace internal
       unsigned int n_q_points;
 
       /**
-       * Stores the number of DoFs per cell in @p dim dimensions.
+       * Stores the number of DoFs per cell of the scalar element in @p dim
+       * dimensions.
        */
-      unsigned int dofs_per_cell;
+      unsigned int dofs_per_component_on_cell;
 
       /**
        * Stores the number of quadrature points per face in @p dim dimensions.
@@ -249,7 +250,7 @@ namespace internal
       /**
        * Stores the number of DoFs per face in @p dim dimensions.
        */
-      unsigned int dofs_per_face;
+      unsigned int dofs_per_component_on_face;
 
       /**
        * Indicates whether the basis functions are nodal in 0 and 1, i.e., the
@@ -364,9 +365,9 @@ namespace internal
       fe_degree (0),
       n_q_points_1d (0),
       n_q_points (0),
-      dofs_per_cell (0),
+      dofs_per_component_on_cell (0),
       n_q_points_face (0),
-      dofs_per_face (0),
+      dofs_per_component_on_face (0),
       nodal_at_cell_boundaries (false)
     {
       reinit (quad, fe_in, base_element_number);
index 489f51afdd4964096d5b41b5a573653de9fcebc5..8260f8f869d4cb236880d86ab01124306e2fa6ad 100644 (file)
@@ -62,9 +62,9 @@ namespace internal
       fe_degree (numbers::invalid_unsigned_int),
       n_q_points_1d(0),
       n_q_points (0),
-      dofs_per_cell (0),
+      dofs_per_component_on_cell (0),
       n_q_points_face (0),
-      dofs_per_face (0),
+      dofs_per_component_on_face (0),
       nodal_at_cell_boundaries (false)
     {}
 
@@ -181,9 +181,9 @@ namespace internal
       }
 
       n_q_points      = Utilities::fixed_power<dim>(n_q_points_1d);
-      dofs_per_cell   = fe->dofs_per_cell;
       n_q_points_face = dim>1?Utilities::fixed_power<dim-1>(n_q_points_1d):1;
-      dofs_per_face   = dim>1?Utilities::fixed_power<dim-1>(fe_degree+1):1;
+      dofs_per_component_on_cell = fe->dofs_per_cell;
+      dofs_per_component_on_face = dim>1?Utilities::fixed_power<dim-1>(fe_degree+1):1;
 
       const unsigned int array_size = n_dofs_1d*n_q_points_1d;
       this->shape_gradients.resize_fast (array_size);
@@ -302,7 +302,7 @@ namespace internal
       if (nodal_at_cell_boundaries == true)
         {
           face_to_cell_index_nodal.reinit(GeometryInfo<dim>::faces_per_cell,
-                                          dofs_per_face);
+                                          dofs_per_component_on_face);
           for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
             {
               const unsigned int direction = f/2;
@@ -313,7 +313,7 @@ namespace internal
               const unsigned int offset = (f%2)*fe_degree*shift;
 
               if (direction == 0 || direction == dim-1)
-                for (unsigned int i=0; i<dofs_per_face; ++i)
+                for (unsigned int i=0; i<dofs_per_component_on_face; ++i)
                   face_to_cell_index_nodal(f,i) = offset + i*stride;
               else
                 // local coordinate system on faces 2 and 3 is zx in
@@ -322,8 +322,8 @@ namespace internal
                 for (unsigned int j=0; j<=fe_degree; ++j)
                   for (unsigned int i=0; i<=fe_degree; ++i)
                     {
-                      const unsigned int ind = offset + j*dofs_per_face + i;
-                      AssertIndexRange(ind, dofs_per_cell);
+                      const unsigned int ind = offset + j*dofs_per_component_on_face + i;
+                      AssertIndexRange(ind, dofs_per_component_on_cell);
                       const unsigned int l = i*(fe_degree+1)+j;
                       face_to_cell_index_nodal(f,l) = ind;
                     }
@@ -333,7 +333,7 @@ namespace internal
       if (element_type == tensor_symmetric_hermite)
         {
           face_to_cell_index_hermite.reinit(GeometryInfo<dim>::faces_per_cell,
-                                            2*dofs_per_face);
+                                            2*dofs_per_component_on_face);
           for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
             {
               const unsigned int direction = f/2;
@@ -346,7 +346,7 @@ namespace internal
                 shift = -shift;
 
               if (direction == 0 || direction == dim-1)
-                for (unsigned int i=0; i<dofs_per_face; ++i)
+                for (unsigned int i=0; i<dofs_per_component_on_face; ++i)
                   {
                     face_to_cell_index_hermite(f,2*i) = offset + i*stride;
                     face_to_cell_index_hermite(f,2*i+1) = offset + i*stride + shift;
@@ -358,8 +358,8 @@ namespace internal
                 for (unsigned int j=0; j<=fe_degree; ++j)
                   for (unsigned int i=0; i<=fe_degree; ++i)
                     {
-                      const unsigned int ind = offset + j*dofs_per_face + i;
-                      AssertIndexRange(ind, dofs_per_cell);
+                      const unsigned int ind = offset + j*dofs_per_component_on_face + i;
+                      AssertIndexRange(ind, dofs_per_component_on_cell);
                       const unsigned int l = i*(fe_degree+1)+j;
                       face_to_cell_index_hermite(f,2*l) = ind;
                       face_to_cell_index_hermite(f,2*l+1) = ind+shift;
@@ -374,7 +374,7 @@ namespace internal
     bool
     ShapeInfo<Number>::check_1d_shapes_symmetric(const unsigned int n_q_points_1d)
     {
-      if (dofs_per_cell == 0)
+      if (dofs_per_component_on_cell == 0)
         return false;
 
       const double zero_tol =
@@ -476,7 +476,7 @@ namespace internal
     bool
     ShapeInfo<Number>::check_1d_shapes_collocation()
     {
-      if (dofs_per_cell != n_q_points)
+      if (dofs_per_component_on_cell != n_q_points)
         return false;
 
       const double zero_tol =

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.