]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add MappingCartesian::maybe_update_jacobians()
authorSimon Sticko <simon@sticko.se>
Tue, 16 Nov 2021 13:32:22 +0000 (14:32 +0100)
committerSimon Sticko <simon@sticko.se>
Tue, 16 Nov 2021 15:14:40 +0000 (16:14 +0100)
To remove some duplication in the class.

include/deal.II/fe/mapping_cartesian.h
source/fe/mapping_cartesian.cc

index d62039f4137f9b8b1b55aef035f297edd1d4b066..71e0aa385b140981abf4bf2909c14a1aae1fc1d6 100644 (file)
@@ -376,6 +376,17 @@ private:
       &output_data) const;
 
 
+  /**
+   * Compute the Jacobians if the UpdateFlags of the incoming
+   * InternalData object say that they should be updated.
+   */
+  void
+  maybe_update_jacobians(
+    const InternalData &             data,
+    const CellSimilarity::Similarity cell_similarity,
+    internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
+      &output_data) const;
+
   /**
    * Compute the inverse Jacobians if the UpdateFlags of the incoming
    * InternalData object say that they should be updated.
index f5b44264e5ad3f00b50f343619c2ddcccfc65d88..8963c8465b906ddb2f68f64899bcc7f87791a006 100644 (file)
@@ -385,6 +385,28 @@ MappingCartesian<dim, spacedim>::maybe_update_jacobian_derivatives(
 
 
 
+template <int dim, int spacedim>
+void
+MappingCartesian<dim, spacedim>::maybe_update_jacobians(
+  const InternalData &             data,
+  const CellSimilarity::Similarity cell_similarity,
+  internal::FEValuesImplementation::MappingRelatedData<dim, spacedim>
+    &output_data) const
+{
+  // "compute" Jacobian at the quadrature points, which are all the
+  // same
+  if (data.update_each & update_jacobians)
+    if (cell_similarity != CellSimilarity::translation)
+      for (unsigned int i = 0; i < output_data.jacobians.size(); ++i)
+        {
+          output_data.jacobians[i] = DerivativeForm<1, dim, spacedim>();
+          for (unsigned int j = 0; j < dim; ++j)
+            output_data.jacobians[i][j][j] = data.cell_extents[j];
+        }
+}
+
+
+
 template <int dim, int spacedim>
 void
 MappingCartesian<dim, spacedim>::maybe_update_inverse_jacobians(
@@ -443,17 +465,9 @@ MappingCartesian<dim, spacedim>::fill_fe_values(
           for (unsigned int i = 0; i < output_data.JxW_values.size(); ++i)
             output_data.JxW_values[i] = J * quadrature.weight(i);
       }
-  // "compute" Jacobian at the quadrature points, which are all the
-  // same
-  if (data.update_each & update_jacobians)
-    if (cell_similarity != CellSimilarity::translation)
-      for (unsigned int i = 0; i < output_data.jacobians.size(); ++i)
-        {
-          output_data.jacobians[i] = DerivativeForm<1, dim, spacedim>();
-          for (unsigned int j = 0; j < dim; ++j)
-            output_data.jacobians[i][j][j] = data.cell_extents[j];
-        }
 
+
+  maybe_update_jacobians(data, cell_similarity, output_data);
   maybe_update_jacobian_derivatives(data, cell_similarity, output_data);
   maybe_update_inverse_jacobians(data, cell_similarity, output_data);
 
@@ -480,7 +494,6 @@ MappingCartesian<dim, spacedim>::fill_mapping_data_for_generic_points(
          ExcNotImplemented());
 
   output_data.initialize(unit_points.size(), update_flags);
-  const unsigned int n_points = unit_points.size();
 
   InternalData data;
   data.update_each = update_flags;
@@ -493,16 +506,7 @@ MappingCartesian<dim, spacedim>::fill_mapping_data_for_generic_points(
                                       data,
                                       output_data.quadrature_points);
 
-  for (unsigned int i = 0; i < n_points; ++i)
-    {
-      if (update_flags & update_jacobians)
-        {
-          output_data.jacobians[i] = DerivativeForm<1, dim, spacedim>();
-          for (unsigned int j = 0; j < dim; ++j)
-            output_data.jacobians[i][j][j] = data.cell_extents[j];
-        }
-    }
-
+  maybe_update_jacobians(data, CellSimilarity::none, output_data);
   maybe_update_inverse_jacobians(data, CellSimilarity::none, output_data);
 }
 
@@ -560,14 +564,7 @@ MappingCartesian<dim, spacedim>::fill_fe_face_values(
       data.volume_element = J;
     }
 
-  if (data.update_each & update_jacobians)
-    for (unsigned int i = 0; i < output_data.jacobians.size(); ++i)
-      {
-        output_data.jacobians[i] = DerivativeForm<1, dim, spacedim>();
-        for (unsigned int d = 0; d < dim; ++d)
-          output_data.jacobians[i][d][d] = data.cell_extents[d];
-      }
-
+  maybe_update_jacobians(data, CellSimilarity::none, output_data);
   maybe_update_jacobian_derivatives(data, CellSimilarity::none, output_data);
   maybe_update_inverse_jacobians(data, CellSimilarity::none, output_data);
 }
@@ -631,14 +628,7 @@ MappingCartesian<dim, spacedim>::fill_fe_subface_values(
       data.volume_element = J;
     }
 
-  if (data.update_each & update_jacobians)
-    for (unsigned int i = 0; i < output_data.jacobians.size(); ++i)
-      {
-        output_data.jacobians[i] = DerivativeForm<1, dim, spacedim>();
-        for (unsigned int d = 0; d < dim; ++d)
-          output_data.jacobians[i][d][d] = data.cell_extents[d];
-      }
-
+  maybe_update_jacobians(data, CellSimilarity::none, output_data);
   maybe_update_jacobian_derivatives(data, CellSimilarity::none, output_data);
   maybe_update_inverse_jacobians(data, CellSimilarity::none, output_data);
 }

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.