]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Mappings: refactor the covariant gradient. 18577/head
authorDavid Wells <drwells@email.unc.edu>
Thu, 19 Jun 2025 14:30:02 +0000 (10:30 -0400)
committerDavid Wells <drwells@email.unc.edu>
Thu, 19 Jun 2025 14:36:00 +0000 (10:36 -0400)
include/deal.II/fe/mapping_fe_field.templates.h
include/deal.II/fe/mapping_internal.h
source/fe/mapping_fe.cc
source/fe/mapping_manifold.cc
source/fe/mapping_q.cc

index c132cced577bf490eef6edbf6c29cad20cac0931..202aba571c6498cf9f3626120d273b30949bad33 100644 (file)
@@ -32,6 +32,7 @@
 #include <deal.II/fe/fe_values.h>
 #include <deal.II/fe/mapping.h>
 #include <deal.II/fe/mapping_fe_field.h>
+#include <deal.II/fe/mapping_internal.h>
 
 #include <deal.II/grid/tria_iterator.h>
 
@@ -2123,27 +2124,14 @@ MappingFEField<dim, spacedim, VectorType>::transform(
     {
       case mapping_covariant_gradient:
         {
-          Assert(data.update_each & update_contravariant_transformation,
+          Assert(data.update_each & update_covariant_transformation,
                  typename FEValuesBase<dim>::ExcAccessToUninitializedField(
                    "update_covariant_transformation"));
 
           for (unsigned int q = 0; q < output.size(); ++q)
-            for (unsigned int i = 0; i < spacedim; ++i)
-              for (unsigned int j = 0; j < spacedim; ++j)
-                for (unsigned int k = 0; k < spacedim; ++k)
-                  {
-                    output[q][i][j][k] = data.covariant[q][j][0] *
-                                         data.covariant[q][k][0] *
-                                         input[q][i][0][0];
-                    for (unsigned int J = 0; J < dim; ++J)
-                      {
-                        const unsigned int K0 = (0 == J) ? 1 : 0;
-                        for (unsigned int K = K0; K < dim; ++K)
-                          output[q][i][j][k] += data.covariant[q][j][J] *
-                                                data.covariant[q][k][K] *
-                                                input[q][i][J][K];
-                      }
-                  }
+            output[q] =
+              internal::apply_covariant_gradient(data.covariant[q], input[q]);
+
           return;
         }
 
index 8a88e7f8c8dd3c86f6bfc5c93d92685e5c0fdd91..f234bd02dec71d6f054cea4b27f3f9383f102ec2 100644 (file)
@@ -27,6 +27,17 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace internal
 {
+  /**
+   * Map the gradient of a covariant vector field. For more information see the
+   * overload of Mapping::transform() which maps 2-differential forms from the
+   * reference cell to the physical cell.
+   */
+  template <int dim, int spacedim, typename Number>
+  Tensor<3, spacedim, Number>
+  apply_covariant_gradient(
+    const DerivativeForm<1, dim, spacedim, Number> &covariant,
+    const DerivativeForm<2, dim, spacedim, Number> &input);
+
   /**
    * Map the Hessian of a contravariant vector field. For more information see
    * the overload of Mapping::transform() which maps 3-differential forms from
@@ -66,6 +77,36 @@ namespace internal
 
 namespace internal
 {
+  template <int dim, int spacedim, typename Number>
+  Tensor<3, spacedim, Number>
+  apply_covariant_gradient(
+    const DerivativeForm<1, dim, spacedim, Number> &covariant,
+    const DerivativeForm<2, dim, spacedim, Number> &input)
+  {
+    Tensor<3, spacedim, Number> output;
+    for (unsigned int i = 0; i < spacedim; ++i)
+      for (unsigned int j = 0; j < spacedim; ++j)
+        {
+          double tmp[dim];
+          for (unsigned int K = 0; K < dim; ++K)
+            {
+              tmp[K] = covariant[j][0] * input[i][0][K];
+              for (unsigned int J = 1; J < dim; ++J)
+                tmp[K] += covariant[j][J] * input[i][J][K];
+            }
+          for (unsigned int k = 0; k < spacedim; ++k)
+            {
+              output[i][j][k] = covariant[k][0] * tmp[0];
+              for (unsigned int K = 1; K < dim; ++K)
+                output[i][j][k] += covariant[k][K] * tmp[K];
+            }
+        }
+
+    return output;
+  }
+
+
+
   template <int dim, int spacedim, typename Number>
   inline Tensor<3, spacedim, Number>
   apply_contravariant_hessian(
index 367c9b287a316c088b5697eb39b22d5f7c04ef2f..85108cb2d80a86156033fff254e623cdbe7ab178 100644 (file)
@@ -2080,28 +2080,14 @@ MappingFE<dim, spacedim>::transform(
     {
       case mapping_covariant_gradient:
         {
-          Assert(data.update_each & update_contravariant_transformation,
+          Assert(data.update_each & update_covariant_transformation,
                  typename FEValuesBase<dim>::ExcAccessToUninitializedField(
                    "update_covariant_transformation"));
 
           for (unsigned int q = 0; q < output.size(); ++q)
-            for (unsigned int i = 0; i < spacedim; ++i)
-              for (unsigned int j = 0; j < spacedim; ++j)
-                {
-                  double tmp[dim];
-                  for (unsigned int K = 0; K < dim; ++K)
-                    {
-                      tmp[K] = data.covariant[q][j][0] * input[q][i][0][K];
-                      for (unsigned int J = 1; J < dim; ++J)
-                        tmp[K] += data.covariant[q][j][J] * input[q][i][J][K];
-                    }
-                  for (unsigned int k = 0; k < spacedim; ++k)
-                    {
-                      output[q][i][j][k] = data.covariant[q][k][0] * tmp[0];
-                      for (unsigned int K = 1; K < dim; ++K)
-                        output[q][i][j][k] += data.covariant[q][k][K] * tmp[K];
-                    }
-                }
+            output[q] =
+              internal::apply_covariant_gradient(data.covariant[q], input[q]);
+
           return;
         }
 
index 26298d42d5272f7244831b39cd6657ebb4122b6f..54b6672b6ebc86a02f3a7d8654f9731c2c8b7072 100644 (file)
@@ -1303,28 +1303,14 @@ MappingManifold<dim, spacedim>::transform(
     {
       case mapping_covariant_gradient:
         {
-          Assert(data.update_each & update_contravariant_transformation,
+          Assert(data.update_each & update_covariant_transformation,
                  typename FEValuesBase<dim>::ExcAccessToUninitializedField(
                    "update_covariant_transformation"));
 
           for (unsigned int q = 0; q < output.size(); ++q)
-            for (unsigned int i = 0; i < spacedim; ++i)
-              for (unsigned int j = 0; j < spacedim; ++j)
-                {
-                  double tmp[dim];
-                  for (unsigned int K = 0; K < dim; ++K)
-                    {
-                      tmp[K] = data.covariant[q][j][0] * input[q][i][0][K];
-                      for (unsigned int J = 1; J < dim; ++J)
-                        tmp[K] += data.covariant[q][j][J] * input[q][i][J][K];
-                    }
-                  for (unsigned int k = 0; k < spacedim; ++k)
-                    {
-                      output[q][i][j][k] = data.covariant[q][k][0] * tmp[0];
-                      for (unsigned int K = 1; K < dim; ++K)
-                        output[q][i][j][k] += data.covariant[q][k][K] * tmp[K];
-                    }
-                }
+            output[q] =
+              internal::apply_covariant_gradient(data.covariant[q], input[q]);
+
           return;
         }
 
index 5b4d7faa1359e80f55ae3045c5a6281001b736be..78be82a0855e77fd9600723930b2efae27eb6695 100644 (file)
@@ -1497,25 +1497,12 @@ MappingQ<dim, spacedim>::transform(
                    "update_covariant_transformation"));
 
           for (unsigned int q = 0; q < output.size(); ++q)
-            for (unsigned int i = 0; i < spacedim; ++i)
-              for (unsigned int j = 0; j < spacedim; ++j)
-                {
-                  double                                 tmp[dim];
-                  const DerivativeForm<1, dim, spacedim> covariant =
-                    data.inverse_jacobians[q].transpose();
-                  for (unsigned int K = 0; K < dim; ++K)
-                    {
-                      tmp[K] = covariant[j][0] * input[q][i][0][K];
-                      for (unsigned int J = 1; J < dim; ++J)
-                        tmp[K] += covariant[j][J] * input[q][i][J][K];
-                    }
-                  for (unsigned int k = 0; k < spacedim; ++k)
-                    {
-                      output[q][i][j][k] = covariant[k][0] * tmp[0];
-                      for (unsigned int K = 1; K < dim; ++K)
-                        output[q][i][j][k] += covariant[k][K] * tmp[K];
-                    }
-                }
+            {
+              const DerivativeForm<1, dim, spacedim> covariant =
+                data.inverse_jacobians[q].transpose();
+              output[q] =
+                internal::apply_covariant_gradient(covariant, input[q]);
+            }
           return;
         }
 

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.