From: Marc Fehling <mafehling.git@gmail.com>
Date: Sat, 25 Dec 2021 07:50:57 +0000 (+0100)
Subject: Calculate cross product once.
X-Git-Tag: v9.4.0-rc1~710^2
X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F13121%2Fhead;p=dealii.git

Calculate cross product once.
---

diff --git a/include/deal.II/grid/manifold_lib.h b/include/deal.II/grid/manifold_lib.h
index 989bc3a235..06d606c835 100644
--- a/include/deal.II/grid/manifold_lib.h
+++ b/include/deal.II/grid/manifold_lib.h
@@ -472,6 +472,11 @@ private:
    * Relative tolerance to measure zero distances.
    */
   const double tolerance;
+
+  /**
+   * The direction vector perpendicular to both direction and normal_direction.
+   */
+  const Tensor<1, spacedim> dxn;
 };
 
 /**
diff --git a/source/grid/manifold_lib.cc b/source/grid/manifold_lib.cc
index bc78bc58b5..ad4ad875a9 100644
--- a/source/grid/manifold_lib.cc
+++ b/source/grid/manifold_lib.cc
@@ -1066,6 +1066,7 @@ CylindricalManifold<dim, spacedim>::CylindricalManifold(
   , direction(direction / direction.norm())
   , point_on_axis(point_on_axis)
   , tolerance(tolerance)
+  , dxn(cross_product_3d(this->direction, normal_direction))
 {
   // do not use static_assert to make dimension-independent programming
   // easier.
@@ -1150,9 +1151,8 @@ CylindricalManifold<dim, spacedim>::push_forward(
          ExcMessage("CylindricalManifold can only be used for spacedim==3!"));
 
   // Rotate the orthogonal direction by the given angle
-  const double sine_r           = std::sin(chart_point(1)) * chart_point(0);
-  const double cosine_r         = std::cos(chart_point(1)) * chart_point(0);
-  const Tensor<1, spacedim> dxn = cross_product_3d(direction, normal_direction);
+  const double sine_r   = std::sin(chart_point(1)) * chart_point(0);
+  const double cosine_r = std::cos(chart_point(1)) * chart_point(0);
   const Tensor<1, spacedim> intermediate =
     normal_direction * cosine_r + dxn * sine_r;
 
@@ -1175,7 +1175,6 @@ CylindricalManifold<dim, spacedim>::push_forward_gradient(
   // Rotate the orthogonal direction by the given angle
   const double              sine   = std::sin(chart_point(1));
   const double              cosine = std::cos(chart_point(1));
-  const Tensor<1, spacedim> dxn = cross_product_3d(direction, normal_direction);
   const Tensor<1, spacedim> intermediate =
     normal_direction * cosine + dxn * sine;