Address more places where we call std::pow.

diff --git a/source/base/function_lib.cc b/source/base/function_lib.cc
index 221ceaaa35a92edef7f56366f2bdfd326203d0b5..7e11ff326ab2825a8ec3006ad17dea8a57991618 100644 (file)
--- a/source/base/function_lib.cc
+++ b/source/base/function_lib.cc
@@ -1181,7 +1181,7 @@ namespace Functions
     const double phi       = std::atan2(y, -x) + numbers::PI;
     const double r_squared = x * x + y * y;
 
-    return std::pow(r_squared, 1. / 3.) * std::sin(2. / 3. * phi);
+    return std::cbrt(r_squared) * std::sin(2. / 3. * phi);
   }
 
 
@@ -1206,7 +1206,7 @@ namespace Functions
             const double phi       = std::atan2(y, -x) + numbers::PI;
             const double r_squared = x * x + y * y;
 
-            values[i] = std::pow(r_squared, 1. / 3.) * std::sin(2. / 3. * phi);
+            values[i] = std::cbrt(r_squared) * std::sin(2. / 3. * phi);
           }
       }
   }
@@ -1235,8 +1235,7 @@ namespace Functions
             const double phi       = std::atan2(y, -x) + numbers::PI;
             const double r_squared = x * x + y * y;
 
-            values[i](0) =
-              std::pow(r_squared, 1. / 3.) * std::sin(2. / 3. * phi);
+            values[i](0) = std::cbrt(r_squared) * std::sin(2. / 3. * phi);
           }
       }
   }
@@ -2965,10 +2964,12 @@ namespace Functions
     const double pi_y = numbers::PI * point(1);
     const double pi_t = numbers::PI / T * this->get_time();
 
-    values[0] = -2 * std::cos(pi_t) * std::pow(std::sin(pi_x), 2) *
-                std::sin(pi_y) * std::cos(pi_y);
-    values[1] = +2 * std::cos(pi_t) * std::pow(std::sin(pi_y), 2) *
-                std::sin(pi_x) * std::cos(pi_x);
+    values[0] = -2 * std::cos(pi_t) *
+                Utilities::fixed_power<2>(std::sin(pi_x)) * std::sin(pi_y) *
+                std::cos(pi_y);
+    values[1] = +2 * std::cos(pi_t) *
+                Utilities::fixed_power<2>(std::sin(pi_y)) * std::sin(pi_x) *
+                std::cos(pi_x);
 
     if (dim == 3)
       values[2] = 0;

diff --git a/source/base/quadrature_lib.cc b/source/base/quadrature_lib.cc
index fd1207cff04404a7ff74f63fcfc9310551650eb4..8e640c346373c9b7c7df10fc1b450583fa457535 100644 (file)
--- a/source/base/quadrature_lib.cc
+++ b/source/base/quadrature_lib.cc
@@ -127,8 +127,8 @@ namespace internal
       std::vector<long double> w(q);
 
       const long double factor =
-        std::pow(2., alpha + beta + 1) * gamma(alpha + q) * gamma(beta + q) /
-        ((q - 1) * gamma(q) * gamma(alpha + beta + q + 1));
+        Utilities::pow(2, alpha + beta + 1) * gamma(alpha + q) *
+        gamma(beta + q) / ((q - 1) * gamma(q) * gamma(alpha + beta + q + 1));
       for (unsigned int i = 0; i < q; ++i)
         {
           const long double s =
@@ -1168,22 +1168,19 @@ QTelles<1>::QTelles(const Quadrature<1> &base_quad, const Point<1> &singularity)
   // We need to check if the singularity is at the boundary of the interval.
   if (std::abs(eta_star) <= tol)
     {
-      gamma_bar =
-        std::pow((eta_bar * eta_star + std::abs(eta_star)), 1.0 / 3.0) +
-        std::pow((eta_bar * eta_star - std::abs(eta_star)), 1.0 / 3.0) +
-        eta_bar;
+      gamma_bar = std::cbrt(eta_bar * eta_star + std::abs(eta_star)) +
+                  std::cbrt(eta_bar * eta_star - std::abs(eta_star)) + eta_bar;
     }
   else
     {
-      gamma_bar = (eta_bar * eta_star + std::abs(eta_star)) /
-                    std::abs(eta_bar * eta_star + std::abs(eta_star)) *
-                    std::pow(std::abs(eta_bar * eta_star + std::abs(eta_star)),
-                             1.0 / 3.0) +
-                  (eta_bar * eta_star - std::abs(eta_star)) /
-                    std::abs(eta_bar * eta_star - std::abs(eta_star)) *
-                    std::pow(std::abs(eta_bar * eta_star - std::abs(eta_star)),
-                             1.0 / 3.0) +
-                  eta_bar;
+      gamma_bar =
+        (eta_bar * eta_star + std::abs(eta_star)) /
+          std::abs(eta_bar * eta_star + std::abs(eta_star)) *
+          std::cbrt(std::abs(eta_bar * eta_star + std::abs(eta_star))) +
+        (eta_bar * eta_star - std::abs(eta_star)) /
+          std::abs(eta_bar * eta_star - std::abs(eta_star)) *
+          std::cbrt(std::abs(eta_bar * eta_star - std::abs(eta_star))) +
+        eta_bar;
     }
   for (unsigned int q = 0; q < quadrature_points.size(); ++q)
     {

diff --git a/source/fe/fe_q_hierarchical.cc b/source/fe/fe_q_hierarchical.cc
index 794f0d0e1f998300f32c1a0390fb43c6476b1da2..4c5e51ada7d12e90b5d3f83b3113dc28539084e9 100644 (file)
--- a/source/fe/fe_q_hierarchical.cc
+++ b/source/fe/fe_q_hierarchical.cc
@@ -484,21 +484,17 @@ FE_Q_Hierarchical<dim>::build_dofs_cell(
               // factor == k * (k-1) * ... * (k-j+1) / j! = k! / (k-j)! / j!
               if (c == 0)
                 {
-                  dofs_subcell[c](j, k) =
-                    ((k + j) % 2 == 0) ?
-                      std::pow(.5, static_cast<double>(k)) * factor :
-                      -std::pow(.5, static_cast<double>(k)) * factor;
-                  dofs_cell[c](j, k) =
-                    std::pow(2., static_cast<double>(j)) * factor;
+                  dofs_subcell[c](j, k) = ((k + j) % 2 == 0) ?
+                                            Utilities::pow(.5, k) * factor :
+                                            -Utilities::pow(.5, k) * factor;
+                  dofs_cell[c](j, k)    = Utilities::pow(2, j) * factor;
                 }
               else
                 {
-                  dofs_subcell[c](j, k) =
-                    std::pow(.5, static_cast<double>(k)) * factor;
-                  dofs_cell[c](j, k) =
-                    ((k + j) % 2 == 0) ?
-                      std::pow(2., static_cast<double>(j)) * factor :
-                      -std::pow(2., static_cast<double>(j)) * factor;
+                  dofs_subcell[c](j, k) = Utilities::pow(.5, k) * factor;
+                  dofs_cell[c](j, k)    = ((k + j) % 2 == 0) ?
+                                            Utilities::pow(2, j) * factor :
+                                            -Utilities::pow(2, j) * factor;
                 }
             }
         }
@@ -1082,7 +1078,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                   for (unsigned int i = 2; i < this->n_dofs_per_face(face_no);
                        ++i)
                     {
-                      interpolation_matrix(i, i) = std::pow(0.5, i);
+                      interpolation_matrix(i, i) = Utilities::pow(0.5, i);
                       factorial_i *= i;
                       int factorial_j  = factorial_i;
                       int factorial_ij = 1;
@@ -1096,12 +1092,12 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
 
                           if (((i + j) & 1) != 0u)
                             interpolation_matrix(i, j) =
-                              -1.0 * std::pow(0.5, j) * factorial_j /
+                              -1.0 * Utilities::pow(0.5, j) * factorial_j /
                               (factorial_i * factorial_ij);
 
                           else
                             interpolation_matrix(i, j) =
-                              std::pow(0.5, j) * factorial_j /
+                              Utilities::pow(0.5, j) * factorial_j /
                               (factorial_i * factorial_ij);
                         }
                     }
@@ -1128,7 +1124,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                   for (unsigned int i = 2; i < this->n_dofs_per_face(face_no);
                        ++i)
                     {
-                      interpolation_matrix(i, i) = std::pow(0.5, i);
+                      interpolation_matrix(i, i) = Utilities::pow(0.5, i);
                       factorial_i *= i;
                       int factorial_j  = factorial_i;
                       int factorial_ij = 1;
@@ -1140,7 +1136,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                           factorial_ij *= j - i;
                           factorial_j *= j;
                           interpolation_matrix(i, j) =
-                            std::pow(0.5, j) * factorial_j /
+                            Utilities::pow(0.5, j) * factorial_j /
                             (factorial_i * factorial_ij);
                         }
                     }
@@ -1194,8 +1190,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
 
                   for (unsigned int i = 2; i <= this->degree; ++i)
                     {
-                      double tmp                         = std::pow(0.5, i);
-                      interpolation_matrix(i + 2, i + 2) = tmp;
+                      double tmp = Utilities::pow(0.5, i);
+                      interpolation_matrix(i + 2, i + 2)         = tmp;
                       interpolation_matrix(i + 2 * source_fe.degree,
                                            i + 2 * this->degree) = tmp;
                       tmp *= 0.5;
@@ -1227,7 +1223,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                           interpolation_matrix(i + (j + 2) * source_fe.degree -
                                                  j,
                                                i + (j + 2) * this->degree - j) =
-                            std::pow(0.5, i + j);
+                            Utilities::pow(0.5, i + j);
                           factorial_k *= j;
                           int factorial_kl = 1;
                           int factorial_l  = factorial_k;
@@ -1241,14 +1237,14 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                 interpolation_matrix(
                                   i + (j + 2) * source_fe.degree - j,
                                   i + (k + 2) * this->degree - k) =
-                                  -1.0 * std::pow(0.5, i + k) * factorial_l /
-                                  (factorial_k * factorial_kl);
+                                  -1.0 * Utilities::pow(0.5, i + k) *
+                                  factorial_l / (factorial_k * factorial_kl);
 
                               else
                                 interpolation_matrix(
                                   i + (j + 2) * source_fe.degree - j,
                                   i + (k + 2) * this->degree - k) =
-                                  std::pow(0.5, i + k) * factorial_l /
+                                  Utilities::pow(0.5, i + k) * factorial_l /
                                   (factorial_k * factorial_kl);
                             }
                         }
@@ -1264,8 +1260,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
 
                           if (((i + j) & 1) != 0u)
                             {
-                              tmp = -1.0 * std::pow(0.5, j) * factorial_j /
-                                    (factorial_i * factorial_ij);
+                              tmp = -1.0 * Utilities::pow(0.5, j) *
+                                    factorial_j / (factorial_i * factorial_ij);
                               interpolation_matrix(i + 2, j + 2)         = tmp;
                               interpolation_matrix(i + 2 * source_fe.degree,
                                                    j + 2 * this->degree) = tmp;
@@ -1276,7 +1272,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                   interpolation_matrix(
                                     i + (k + 2) * source_fe.degree - k,
                                     j + (k + 2) * this->degree - k) =
-                                    tmp * std::pow(0.5, k);
+                                    tmp * Utilities::pow(0.5, k);
                                   factorial_k *= k;
                                   int factorial_l  = factorial_k;
                                   int factorial_kl = 1;
@@ -1292,7 +1288,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                         interpolation_matrix(
                                           i + (k + 2) * source_fe.degree - k,
                                           j + (l + 2) * this->degree - l) =
-                                          -1.0 * tmp * std::pow(0.5, l) *
+                                          -1.0 * tmp * Utilities::pow(0.5, l) *
                                           factorial_l /
                                           (factorial_k * factorial_kl);
 
@@ -1300,7 +1296,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                         interpolation_matrix(
                                           i + (k + 2) * source_fe.degree - k,
                                           j + (l + 2) * this->degree - l) =
-                                          tmp * std::pow(0.5, l) * factorial_l /
+                                          tmp * Utilities::pow(0.5, l) *
+                                          factorial_l /
                                           (factorial_k * factorial_kl);
                                     }
                                 }
@@ -1330,7 +1327,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                             }
                           else
                             {
-                              tmp = std::pow(0.5, j) * factorial_j /
+                              tmp = Utilities::pow(0.5, j) * factorial_j /
                                     (factorial_i * factorial_ij);
                               interpolation_matrix(i + 2, j + 2)         = tmp;
                               interpolation_matrix(i + 2 * source_fe.degree,
@@ -1342,7 +1339,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                   interpolation_matrix(
                                     i + (k + 2) * source_fe.degree - k,
                                     j + (k + 2) * this->degree - k) =
-                                    tmp * std::pow(0.5, k);
+                                    tmp * Utilities::pow(0.5, k);
                                   factorial_k *= k;
                                   int factorial_l  = factorial_k;
                                   int factorial_kl = 1;
@@ -1358,7 +1355,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                         interpolation_matrix(
                                           i + (k + 2) * source_fe.degree - k,
                                           j + (l + 2) * this->degree - l) =
-                                          -1.0 * tmp * std::pow(0.5, l) *
+                                          -1.0 * tmp * Utilities::pow(0.5, l) *
                                           factorial_l /
                                           (factorial_k * factorial_kl);
 
@@ -1366,7 +1363,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                         interpolation_matrix(
                                           i + (k + 2) * source_fe.degree - k,
                                           j + (l + 2) * this->degree - l) =
-                                          tmp * std::pow(0.5, l) * factorial_l /
+                                          tmp * Utilities::pow(0.5, l) *
+                                          factorial_l /
                                           (factorial_k * factorial_kl);
                                     }
                                 }
@@ -1440,8 +1438,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
 
                   for (unsigned int i = 2; i <= this->degree; ++i)
                     {
-                      double tmp                         = std::pow(0.5, i + 1);
-                      interpolation_matrix(i + 2, i + 2) = tmp;
+                      double tmp = Utilities::pow(0.5, i + 1);
+                      interpolation_matrix(i + 2, i + 2)                = tmp;
                       interpolation_matrix(i + 2, i + this->degree + 1) = tmp;
                       interpolation_matrix(i + 3 * source_fe.degree - 1,
                                            i + 2 * this->degree)        = tmp;
@@ -1473,7 +1471,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                         {
                           factorial_ij *= j - i;
                           factorial_j *= j;
-                          tmp = std::pow(0.5, j) * factorial_j /
+                          tmp = Utilities::pow(0.5, j) * factorial_j /
                                 (factorial_i * factorial_ij);
                           interpolation_matrix(i + 2 * source_fe.degree,
                                                j + 2 * this->degree) = tmp;
@@ -1484,7 +1482,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                               interpolation_matrix(
                                 i + (k + 2) * source_fe.degree - k,
                                 j + (k + 2) * this->degree - k) =
-                                tmp * std::pow(0.5, k);
+                                tmp * Utilities::pow(0.5, k);
                               factorial_k *= k;
                               int factorial_l  = factorial_k;
                               int factorial_kl = 1;
@@ -1499,7 +1497,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                     interpolation_matrix(
                                       i + (k + 2) * source_fe.degree - k,
                                       j + (l + 2) * this->degree - l) =
-                                      -1.0 * tmp * std::pow(0.5, l) *
+                                      -1.0 * tmp * Utilities::pow(0.5, l) *
                                       factorial_l /
                                       (factorial_k * factorial_kl);
 
@@ -1507,7 +1505,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                     interpolation_matrix(
                                       i + (k + 2) * source_fe.degree - k,
                                       j + (l + 2) * this->degree - l) =
-                                      tmp * std::pow(0.5, l) * factorial_l /
+                                      tmp * Utilities::pow(0.5, l) *
+                                      factorial_l /
                                       (factorial_k * factorial_kl);
                                 }
                             }
@@ -1555,7 +1554,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                           interpolation_matrix(i + (j + 2) * source_fe.degree -
                                                  j,
                                                i + (j + 2) * this->degree - j) =
-                            std::pow(0.5, i + j);
+                            Utilities::pow(0.5, i + j);
                           factorial_k *= j;
                           int factorial_l  = factorial_k;
                           int factorial_kl = 1;
@@ -1569,14 +1568,14 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                 interpolation_matrix(
                                   i + (j + 2) * source_fe.degree - j,
                                   i + (k + 2) * this->degree - k) =
-                                  -1.0 * std::pow(0.5, i + k) * factorial_l /
-                                  (factorial_k * factorial_kl);
+                                  -1.0 * Utilities::pow(0.5, i + k) *
+                                  factorial_l / (factorial_k * factorial_kl);
 
                               else
                                 interpolation_matrix(
                                   i + (j + 2) * source_fe.degree - j,
                                   i + (k + 2) * this->degree - k) =
-                                  std::pow(0.5, i + k) * factorial_l /
+                                  Utilities::pow(0.5, i + k) * factorial_l /
                                   (factorial_k * factorial_kl);
                             }
                         }
@@ -1624,8 +1623,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
 
                   for (unsigned int i = 2; i <= this->degree; ++i)
                     {
-                      double tmp                         = std::pow(0.5, i);
-                      interpolation_matrix(i + 2, i + 2) = tmp;
+                      double tmp = Utilities::pow(0.5, i);
+                      interpolation_matrix(i + 2, i + 2)             = tmp;
                       interpolation_matrix(i + 3 * source_fe.degree - 1,
                                            i + 3 * this->degree - 1) = tmp;
                       tmp *= 0.5;
@@ -1657,7 +1656,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                           interpolation_matrix(i + (j + 2) * source_fe.degree -
                                                  j,
                                                i + (j + 2) * this->degree - j) =
-                            std::pow(0.5, i + j);
+                            Utilities::pow(0.5, i + j);
                           factorial_k *= j;
                           int factorial_kl = 1;
                           int factorial_l  = factorial_k;
@@ -1669,7 +1668,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                               interpolation_matrix(
                                 i + (j + 2) * source_fe.degree - j,
                                 i + (k + 2) * this->degree - k) =
-                                std::pow(0.5, i + k) * factorial_l /
+                                Utilities::pow(0.5, i + k) * factorial_l /
                                 (factorial_k * factorial_kl);
                             }
                         }
@@ -1682,7 +1681,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                         {
                           factorial_ij *= j - i;
                           factorial_j *= j;
-                          tmp = std::pow(0.5, j) * factorial_j /
+                          tmp = Utilities::pow(0.5, j) * factorial_j /
                                 (factorial_i * factorial_ij);
                           interpolation_matrix(i + 2, j + 2) = tmp;
                           tmp *= -1.0;
@@ -1718,7 +1717,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                               interpolation_matrix(
                                 i + (k + 2) * source_fe.degree - k,
                                 j + (k + 2) * this->degree - k) =
-                                tmp * std::pow(0.5, k);
+                                tmp * Utilities::pow(0.5, k);
                               factorial_k *= k;
                               int factorial_l  = factorial_k;
                               int factorial_kl = 1;
@@ -1731,7 +1730,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                   interpolation_matrix(
                                     i + (k + 2) * source_fe.degree - k,
                                     j + (l + 2) * this->degree - l) =
-                                    tmp * std::pow(0.5, l) * factorial_l /
+                                    tmp * Utilities::pow(0.5, l) * factorial_l /
                                     (factorial_k * factorial_kl);
                                 }
                             }
@@ -1790,8 +1789,8 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
 
                   for (unsigned int i = 2; i <= this->degree; ++i)
                     {
-                      double tmp                         = std::pow(0.5, i + 1);
-                      interpolation_matrix(i + 2, i + 2) = tmp;
+                      double tmp = Utilities::pow(0.5, i + 1);
+                      interpolation_matrix(i + 2, i + 2)                = tmp;
                       interpolation_matrix(i + 2, i + this->degree + 1) = tmp;
                       interpolation_matrix(i + 2 * source_fe.degree,
                                            i + 2 * this->degree)        = tmp;
@@ -1822,7 +1821,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                           interpolation_matrix(i + (j + 2) * source_fe.degree -
                                                  j,
                                                i + (j + 2) * this->degree - j) =
-                            std::pow(0.5, i + j);
+                            Utilities::pow(0.5, i + j);
                           factorial_k *= j;
                           int factorial_l  = factorial_k;
                           int factorial_kl = 1;
@@ -1834,7 +1833,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                               interpolation_matrix(
                                 i + (j + 2) * source_fe.degree - j,
                                 i + (k + 2) * this->degree - k) =
-                                std::pow(0.5, i + k) * factorial_l /
+                                Utilities::pow(0.5, i + k) * factorial_l /
                                 (factorial_k * factorial_kl);
                             }
                         }
@@ -1847,7 +1846,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                         {
                           factorial_ij *= j - i;
                           factorial_j *= j;
-                          tmp = std::pow(0.5, j + 1) * factorial_j /
+                          tmp = Utilities::pow(0.5, j + 1) * factorial_j /
                                 (factorial_i * factorial_ij);
                           interpolation_matrix(i + 2, j + 2) = tmp;
                           interpolation_matrix(i + 2, j + this->degree + 1) =
@@ -1868,7 +1867,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                               interpolation_matrix(
                                 i + (k + 2) * source_fe.degree - k,
                                 j + (k + 2) * this->degree - k) =
-                                tmp * std::pow(0.5, k);
+                                tmp * Utilities::pow(0.5, k);
                               factorial_k *= k;
                               int factorial_l  = factorial_k;
                               int factorial_kl = 1;
@@ -1881,7 +1880,7 @@ FE_Q_Hierarchical<dim>::get_subface_interpolation_matrix(
                                   interpolation_matrix(
                                     i + (k + 2) * source_fe.degree - k,
                                     j + (l + 2) * this->degree - l) =
-                                    tmp * std::pow(0.5, l) * factorial_l /
+                                    tmp * Utilities::pow(0.5, l) * factorial_l /
                                     (factorial_k * factorial_kl);
                                 }
                             }

diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc
index f598404ba2cfddb9e8cb54be800fb5e94bb8591a..03a70784104f47a2db7e265b337a86c3924d5aaa 100644 (file)
--- a/source/grid/grid_generator.cc
+++ b/source/grid/grid_generator.cc
@@ -665,7 +665,7 @@ namespace GridGenerator
                 const double x = i * 1 / (1.0 * number_points - 1);
                 const double y_t =
                   5 * t *
-                  (0.2969 * std::pow(x, 0.5) - 0.126 * x -
+                  (0.2969 * std::sqrt(x) - 0.126 * x -
                    0.3516 * Utilities::fixed_power<2>(x) +
                    0.2843 * Utilities::fixed_power<3>(x) -
                    0.1036 * Utilities::fixed_power<4>(
@@ -696,7 +696,7 @@ namespace GridGenerator
 
                 const double y_t =
                   5 * t *
-                  (0.2969 * std::pow(x, 0.5) - 0.126 * x -
+                  (0.2969 * std::sqrt(x) - 0.126 * x -
                    0.3516 * Utilities::fixed_power<2>(x) +
                    0.2843 * Utilities::fixed_power<3>(x) -
                    0.1036 * Utilities::fixed_power<4>(

diff --git a/source/grid/grid_refinement.cc b/source/grid/grid_refinement.cc
index 9d5fd9d40271050811762a426a7504f4e03080c7..eae7ca90e5831907600160deadeb2545982c5188 100644 (file)
--- a/source/grid/grid_refinement.cc
+++ b/source/grid/grid_refinement.cc
@@ -472,14 +472,15 @@ GridRefinement::refine_and_coarsen_optimize(Triangulation<dim, spacedim> &tria,
   double      min_cost = std::numeric_limits<double>::max();
   std::size_t min_arg  = 0;
 
+  const double reduction_factor = (1. - std::pow(2., -1. * order));
   for (std::size_t M = 0; M < criteria.size(); ++M)
     {
-      expected_error_reduction +=
-        (1 - std::pow(2., -1. * order)) * criteria(cell_indices[M]);
+      expected_error_reduction += reduction_factor * criteria(cell_indices[M]);
 
-      const double cost = std::pow(((std::pow(2., dim) - 1) * (1 + M) + N),
-                                   static_cast<double>(order) / dim) *
-                          (original_error - expected_error_reduction);
+      const double cost =
+        std::pow(((Utilities::fixed_power<dim>(2) - 1) * (1 + M) + N),
+                 static_cast<double>(order) / dim) *
+        (original_error - expected_error_reduction);
       if (cost <= min_cost)
         {
           min_cost = cost;

diff --git a/source/grid/grid_tools.cc b/source/grid/grid_tools.cc
index c1d33574e7da582bb4538ed220be326780cf92f6..dad68daa4e2114365264ff6ec1cf051a4e1c2aff 100644 (file)
--- a/source/grid/grid_tools.cc
+++ b/source/grid/grid_tools.cc
@@ -2427,10 +2427,10 @@ namespace GridTools
         double objective = 0;
         for (unsigned int c = 0; c < object->n_children(); ++c)
           for (const unsigned int i : object->child(c)->vertex_indices())
-            objective +=
-              (child_alternating_forms[c][i] -
-               average_parent_alternating_form / std::pow(2., 1. * structdim))
-                .norm_square();
+            objective += (child_alternating_forms[c][i] -
+                          average_parent_alternating_form /
+                            Utilities::fixed_power<structdim>(2))
+                           .norm_square();
 
         return objective;
       }

diff --git a/source/grid/tria.cc b/source/grid/tria.cc
index f1e48f669ca8d4bd84515ff472848eda512988ba..752e9b4ddddedbce5a253b9e1528056a409cf954 100644 (file)
--- a/source/grid/tria.cc
+++ b/source/grid/tria.cc
@@ -1577,7 +1577,8 @@ namespace
         GeometryInfo<dim>::alternating_form_at_vertices(vertices, determinants);
 
         for (const unsigned int i : GeometryInfo<dim>::vertex_indices())
-          if (determinants[i] <= 1e-9 * std::pow(cell->diameter(), 1. * dim))
+          if (determinants[i] <=
+              1e-9 * Utilities::fixed_power<dim>(cell->diameter()))
             {
               distorted_cells.distorted_cells.push_back(cell);
               break;
@@ -1612,7 +1613,7 @@ namespace
 
         for (const unsigned int i : GeometryInfo<dim>::vertex_indices())
           if (determinants[i] <=
-              1e-9 * std::pow(cell->child(c)->diameter(), 1. * dim))
+              1e-9 * Utilities::fixed_power<dim>(cell->child(c)->diameter()))
             return true;
       }
 

diff --git a/source/hp/refinement.cc b/source/hp/refinement.cc
index 36f6851f34dbdede2736702bcad738b7e24839c9..fa3cec1e8d8ac9d9809ecbc576a9e166f2704669 100644 (file)
--- a/source/hp/refinement.cc
+++ b/source/hp/refinement.cc
@@ -633,16 +633,26 @@ namespace hp
           // step 1: exponential decay with p-adaptation
           if (cell->future_fe_index_set())
             {
-              predicted_errors[cell->active_cell_index()] *=
-                std::pow(gamma_p,
-                         int(future_fe_degree) - int(cell->get_fe().degree));
+              if (future_fe_degree > cell->get_fe().degree)
+                predicted_errors[cell->active_cell_index()] *=
+                  Utilities::pow(gamma_p,
+                                 future_fe_degree - cell->get_fe().degree);
+              else if (future_fe_degree < cell->get_fe().degree)
+                predicted_errors[cell->active_cell_index()] /=
+                  Utilities::pow(gamma_p,
+                                 cell->get_fe().degree - future_fe_degree);
+              else
+                {
+                  // The two degrees are the same; we do not need to
+                  // adapt the predicted error
+                }
             }
 
           // step 2: algebraic decay with h-adaptation
           if (cell->refine_flag_set())
             {
               predicted_errors[cell->active_cell_index()] *=
-                (gamma_h * std::pow(.5, future_fe_degree));
+                (gamma_h * Utilities::pow(.5, future_fe_degree));
 
               // predicted error will be split on children cells
               // after adaptation via CellDataTransfer
@@ -650,7 +660,7 @@ namespace hp
           else if (cell->coarsen_flag_set())
             {
               predicted_errors[cell->active_cell_index()] /=
-                (gamma_h * std::pow(.5, future_fe_degree));
+                (gamma_h * Utilities::pow(.5, future_fe_degree));
 
               // predicted error will be summed up on parent cell
               // after adaptation via CellDataTransfer