Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Replace std::pow with integer exponent #13321

Closed
bangerth opened this issue Feb 1, 2022 · 7 comments · Fixed by #16449
Closed

Replace std::pow with integer exponent #13321

bangerth opened this issue Feb 1, 2022 · 7 comments · Fixed by #16449

Comments

@bangerth
Copy link
Member

bangerth commented Feb 1, 2022

A possibly surprising realization from geodynamics/aspect#4456 is that std::pow(x,p) with fixed integer exponent p is vastly more expensive than just writing out the product, at least for relatively small integers p. We should take a look at our code base and see whether there is something that's worth replacing. These 108 places might be a good starting point to look at:

> egrep -r 'std::pow.*[2345](\.0)?\)' include/ source/ examples/
include/deal.II/optimization/line_minimization.h:   *     const double f = 100. * std::pow(x, 4) + std::pow(1. - x, 2); // Value
include/deal.II/optimization/line_minimization.h:   *     const double g = 400. * std::pow(x, 3) - 2. * (1. - x); // Gradient
include/deal.II/optimization/line_minimization.h:      std::pow(x2_shift * x3_shift, 2) * (x2_shift - x3_shift);
include/deal.II/optimization/line_minimization.h:      (r1 * std::pow(x3_shift, 2) - r2 * std::pow(x2_shift, 2)) / denom;
include/deal.II/optimization/line_minimization.h:      (r2 * std::pow(x2_shift, 3) - r1 * std::pow(x3_shift, 3)) / denom;
include/deal.II/grid/reference_cell.h:         {{Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
include/deal.II/grid/reference_cell.h:                      +std::pow(1.0 / 3.0, 1.0 / 4.0),
include/deal.II/grid/reference_cell.h:           Point<dim>(-std::pow(1.0 / 3.0, 1.0 / 4.0),
include/deal.II/grid/reference_cell.h:                      +std::pow(1.0 / 3.0, 1.0 / 4.0))}}}};
source/fe/fe_nedelec.cc:        return 1.e-15 * std::exp(std::pow(p, 1.075));
source/base/function_lib_cutoff.cc:                          (-2.0 * r * r / std::pow(-r * r + d * d, 2.0) * d *
source/base/flow_function.cc:            const double rl1 = std::pow(r2, lambda / 2. - .5);
source/base/flow_function.cc:            const double rl1  = std::pow(r2, lambda / 2. - .5);
source/base/data_out_base.cc:        std::sqrt(std::pow(gradient[0], 2.0) + std::pow(gradient[1], 2.0));
source/base/data_out_base.cc:                                    std::pow(nrml[i * d1 + j * d2](2), 2.));
source/base/quadrature_lib.cc:        std::pow((eta_bar * eta_star + std::abs(eta_star)), 1.0 / 3.0) +
source/base/quadrature_lib.cc:        std::pow((eta_bar * eta_star - std::abs(eta_star)), 1.0 / 3.0) +
source/base/quadrature_lib.cc:      double eta   = (std::pow(gamma - gamma_bar, 3.0) +
source/base/function_lib.cc:    return std::pow(r_squared, .25) * std::sin(.5 * phi);
source/base/function_lib.cc:        values[i] = std::pow(r_squared, .25) * std::sin(.5 * phi);
source/base/function_lib.cc:        values[i](0) = std::pow(r_squared, .25) * std::sin(.5 * phi);
source/base/function_lib.cc:    return std::pow(r_squared, .125) * std::sin(.25 * phi);
source/base/function_lib.cc:        values[i] = std::pow(r_squared, .125) * std::sin(.25 * phi);
source/base/function_lib.cc:        values[i](0) = std::pow(r_squared, .125) * std::sin(.25 * phi);
source/base/function_signed_distance.cc:          std::pow(distance, 3);
source/base/function_signed_distance.cc:        val += std::pow((point[d] - center[d]) / radii[d], 2);
source/base/function_signed_distance.cc:          const double ex = (a * a - b * b) * std::pow(std::cos(t), 3) / a;
source/base/function_signed_distance.cc:          const double ey = (b * b - a * a) * std::pow(std::sin(t), 3) / b;
source/numerics/derivative_approximation.cc:        (std::pow(s[0][0], 3) + std::pow(s[1][1], 3) + std::pow(s[2][2], 3) +
source/grid/grid_generator.cc:                  (0.2969 * std::pow(x, 0.5) - 0.126 * x -
source/grid/grid_generator.cc:                   0.3516 * std::pow(x, 2) + 0.2843 * std::pow(x, 3) -
source/grid/grid_generator.cc:                   0.1036 * std::pow(x, 4)); // half thickness at a position x
source/grid/grid_generator.cc:                  (x <= p) ? m / std::pow(p, 2) * (2 * p * x - std::pow(x, 2)) :
source/grid/grid_generator.cc:                             m / std::pow(1 - p, 2) *
source/grid/grid_generator.cc:                               ((1 - 2 * p) + 2 * p * x - std::pow(x, 2));
source/grid/grid_generator.cc:                                      2 * m / std::pow(p, 2) * (p - x) :
source/grid/grid_generator.cc:                                      2 * m / std::pow(1 - p, 2) * (p - x);
source/grid/grid_generator.cc:                  (0.2969 * std::pow(x, 0.5) - 0.126 * x -
source/grid/grid_generator.cc:                   0.3516 * std::pow(x, 2) + 0.2843 * std::pow(x, 3) -
source/grid/grid_generator.cc:                   0.1036 * std::pow(x, 4)); // half thickness at a position x
source/grid/manifold_lib.cc:  double w     = std::sqrt(std::pow(y - std::sin(phi) * R, 2.0) +
source/grid/manifold_lib.cc:                       std::pow(x - std::cos(phi) * R, 2.0) + z * z) /
examples/step-85/step-85.cc:                  std::pow(error_at_point, 2) * fe_values->JxW(q);
examples/step-47/step-47.cc:        return 4 * std::pow(PI, 4.0) * std::sin(PI * p[0]) *
examples/step-12/step-12.cc:        std::pow(cell->diameter(), 1 + 1.0 * dim / 2);
examples/step-9/step-9.cc:      (std::pow(cell->diameter(), 1 + 1.0 * dim / 2) * gradient.norm());
examples/step-25/step-25.cc:    system_matrix.add(std::pow(time_step * theta, 2), laplace_matrix);
examples/step-25/step-25.cc:    system_matrix.add(std::pow(time_step * theta, 2), tmp_matrix);
examples/step-25/step-25.cc:    system_rhs.add(std::pow(time_step * theta, 2), tmp_vector);
examples/step-25/step-25.cc:    system_rhs.add(std::pow(time_step, 2) * theta * (1 - theta), tmp_vector);
examples/step-25/step-25.cc:    system_rhs.add(std::pow(time_step, 2) * theta, tmp_vector);
examples/step-82/step-82.cc:        return_value = 24.0 * std::pow(p(1) * (1.0 - p(1)), 2) +
examples/step-82/step-82.cc:                       +24.0 * std::pow(p(0) * (1.0 - p(0)), 2) +
examples/step-82/step-82.cc:          24.0 * std::pow(p(1) * (1.0 - p(1)) * p(2) * (1.0 - p(2)), 2) +
examples/step-82/step-82.cc:          24.0 * std::pow(p(0) * (1.0 - p(0)) * p(2) * (1.0 - p(2)), 2) +
examples/step-82/step-82.cc:          24.0 * std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2) +
examples/step-82/step-82.cc:            std::pow(p(2) * (1.0 - p(2)), 2) +
examples/step-82/step-82.cc:            std::pow(p(1) * (1.0 - p(1)), 2) +
examples/step-82/step-82.cc:            std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:        return_value = std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(2) - 6.0 * std::pow(p(2), 2) + 4.0 * std::pow(p(2), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:                               std::pow(p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3));
examples/step-82/step-82.cc:                               std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          std::pow(p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(2) - 6.0 * std::pow(p(2), 2) + 4.0 * std::pow(p(2), 3)) *
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(2) - 6.0 * std::pow(p(2), 2) + 4.0 * std::pow(p(2), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:              1.0 / std::pow(face->diameter(), 3); // ĥ_e^{-3}
examples/step-82/step-82.cc:              1.0 / std::pow(face->diameter(), 3); // h^{-3}
examples/step-82/step-82.cc:                                std::pow(u_exact_q - solution_values[q], 2) *
examples/step-82/step-82.cc:                                std::pow(u_exact_q - solution_values[q], 2) *
examples/step-62/step-62.cc:                   std::exp(-(std::pow(p[0] - force_center[0], 2) /
examples/step-62/step-62.cc:                                (2 * std::pow(force_sigma_x, 2)) +
examples/step-62/step-62.cc:                              std::pow(p[1] - force_center[1], 2) /
examples/step-62/step-62.cc:                                (2 * std::pow(force_sigma_y, 2))));
examples/step-62/step-62.cc:                      matrix_sum += -std::pow(omega, 2) *
examples/step-55/step-55.cc:                0.1 * std::pow(-std::sqrt(25.0 + 4 * pi2) + 5.0, 2) *
examples/step-55/step-55.cc:                0.05 * std::pow(-std::sqrt(25.0 + 4 * pi2) + 5.0, 3) *
examples/step-71/step-71.cc:        std::pow(1.0 / std::cosh(two_h_dot_h_div_h_sat_squ), 2.0);
examples/step-71/step-71.cc:                                std::pow(det_F, -2.0 / dim - 2.0) * ddet_F_dC) +
examples/step-4/step-4.cc:    return_value += 4.0 * std::pow(p(i), 4.0);
examples/step-30/step-30.cc:        std::pow(cell->diameter(), 1 + 1.0 * dim / 2);
examples/step-76/step-76.cc:  const double courant_number = 0.15 / std::pow(fe_degree, 1.5);
examples/step-12b/step-12b.cc:        std::pow(cell->diameter(), 1 + 1.0 * dim / 2);
examples/step-43/step-43.cc:    const double denominator = std::pow(temp, 2.0);
examples/step-44/step-44.cc:              std::pow((det_F_qp - J_tilde_qp), 2);
examples/step-67/step-67.cc:  const double courant_number = 0.15 / std::pow(fe_degree, 1.5);
examples/step-53/step-53.cc:      std::atan2(x(2) + ep * ep * b * std::pow(std::sin(th), 3),
examples/step-53/step-53.cc:                  (ellipticity * ellipticity * R * std::pow(std::cos(th), 3))));
@Rombur
Copy link
Member

Rombur commented Feb 1, 2022

Which optimization flags are you passing to the compiler? With -O3, gcc generates the same code but the code is different if you use -O0 https://godbolt.org/z/bz3bsn33K

@bangerth
Copy link
Member Author

bangerth commented Feb 1, 2022

I used whatever release mode gave me. But you already get different answers if p=3.

@bangerth
Copy link
Member Author

bangerth commented Feb 1, 2022

Here is a better regex:

> egrep -r 'std::pow.*, [2345](\.0)?\)' include/ source/ examples/
include/deal.II/optimization/line_minimization.h:   *     const double f = 100. * std::pow(x, 4) + std::pow(1. - x, 2); // Value
include/deal.II/optimization/line_minimization.h:   *     const double g = 400. * std::pow(x, 3) - 2. * (1. - x); // Gradient
include/deal.II/optimization/line_minimization.h:      std::pow(x2_shift * x3_shift, 2) * (x2_shift - x3_shift);
include/deal.II/optimization/line_minimization.h:      (r1 * std::pow(x3_shift, 2) - r2 * std::pow(x2_shift, 2)) / denom;
include/deal.II/optimization/line_minimization.h:      (r2 * std::pow(x2_shift, 3) - r1 * std::pow(x3_shift, 3)) / denom;
source/base/function_lib_cutoff.cc:                          (-2.0 * r * r / std::pow(-r * r + d * d, 2.0) * d *
source/base/data_out_base.cc:        std::sqrt(std::pow(gradient[0], 2.0) + std::pow(gradient[1], 2.0));
source/base/quadrature_lib.cc:      double eta   = (std::pow(gamma - gamma_bar, 3.0) +
source/base/function_signed_distance.cc:          std::pow(distance, 3);
source/base/function_signed_distance.cc:        val += std::pow((point[d] - center[d]) / radii[d], 2);
source/base/function_signed_distance.cc:          const double ex = (a * a - b * b) * std::pow(std::cos(t), 3) / a;
source/base/function_signed_distance.cc:          const double ey = (b * b - a * a) * std::pow(std::sin(t), 3) / b;
source/numerics/derivative_approximation.cc:        (std::pow(s[0][0], 3) + std::pow(s[1][1], 3) + std::pow(s[2][2], 3) +
source/grid/grid_generator.cc:                   0.3516 * std::pow(x, 2) + 0.2843 * std::pow(x, 3) -
source/grid/grid_generator.cc:                   0.1036 * std::pow(x, 4)); // half thickness at a position x
source/grid/grid_generator.cc:                  (x <= p) ? m / std::pow(p, 2) * (2 * p * x - std::pow(x, 2)) :
source/grid/grid_generator.cc:                             m / std::pow(1 - p, 2) *
source/grid/grid_generator.cc:                               ((1 - 2 * p) + 2 * p * x - std::pow(x, 2));
source/grid/grid_generator.cc:                                      2 * m / std::pow(p, 2) * (p - x) :
source/grid/grid_generator.cc:                                      2 * m / std::pow(1 - p, 2) * (p - x);
source/grid/grid_generator.cc:                   0.3516 * std::pow(x, 2) + 0.2843 * std::pow(x, 3) -
source/grid/grid_generator.cc:                   0.1036 * std::pow(x, 4)); // half thickness at a position x
source/grid/manifold_lib.cc:  double w     = std::sqrt(std::pow(y - std::sin(phi) * R, 2.0) +
source/grid/manifold_lib.cc:                       std::pow(x - std::cos(phi) * R, 2.0) + z * z) /
examples/step-85/step-85.cc:                  std::pow(error_at_point, 2) * fe_values->JxW(q);
examples/step-47/step-47.cc:        return 4 * std::pow(PI, 4.0) * std::sin(PI * p[0]) *
examples/step-25/step-25.cc:    system_matrix.add(std::pow(time_step * theta, 2), laplace_matrix);
examples/step-25/step-25.cc:    system_matrix.add(std::pow(time_step * theta, 2), tmp_matrix);
examples/step-25/step-25.cc:    system_rhs.add(std::pow(time_step * theta, 2), tmp_vector);
examples/step-25/step-25.cc:    system_rhs.add(std::pow(time_step, 2) * theta * (1 - theta), tmp_vector);
examples/step-25/step-25.cc:    system_rhs.add(std::pow(time_step, 2) * theta, tmp_vector);
examples/step-82/step-82.cc:        return_value = 24.0 * std::pow(p(1) * (1.0 - p(1)), 2) +
examples/step-82/step-82.cc:                       +24.0 * std::pow(p(0) * (1.0 - p(0)), 2) +
examples/step-82/step-82.cc:          24.0 * std::pow(p(1) * (1.0 - p(1)) * p(2) * (1.0 - p(2)), 2) +
examples/step-82/step-82.cc:          24.0 * std::pow(p(0) * (1.0 - p(0)) * p(2) * (1.0 - p(2)), 2) +
examples/step-82/step-82.cc:          24.0 * std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2) +
examples/step-82/step-82.cc:            std::pow(p(2) * (1.0 - p(2)), 2) +
examples/step-82/step-82.cc:            std::pow(p(1) * (1.0 - p(1)), 2) +
examples/step-82/step-82.cc:            std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:        return_value = std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(2) - 6.0 * std::pow(p(2), 2) + 4.0 * std::pow(p(2), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:                               std::pow(p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3));
examples/step-82/step-82.cc:                               std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          std::pow(p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(0) - 6.0 * std::pow(p(0), 2) + 4.0 * std::pow(p(0), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(2) - 6.0 * std::pow(p(2), 2) + 4.0 * std::pow(p(2), 3)) *
examples/step-82/step-82.cc:          std::pow(p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(2) * (1.0 - p(2)), 2);
examples/step-82/step-82.cc:          (2.0 * p(1) - 6.0 * std::pow(p(1), 2) + 4.0 * std::pow(p(1), 3)) *
examples/step-82/step-82.cc:          (2.0 * p(2) - 6.0 * std::pow(p(2), 2) + 4.0 * std::pow(p(2), 3)) *
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)), 2);
examples/step-82/step-82.cc:          std::pow(p(0) * (1.0 - p(0)) * p(1) * (1.0 - p(1)), 2);
examples/step-82/step-82.cc:              1.0 / std::pow(face->diameter(), 3); // ĥ_e^{-3}
examples/step-82/step-82.cc:              1.0 / std::pow(face->diameter(), 3); // h^{-3}
examples/step-82/step-82.cc:                                std::pow(u_exact_q - solution_values[q], 2) *
examples/step-82/step-82.cc:                                std::pow(u_exact_q - solution_values[q], 2) *
examples/step-62/step-62.cc:                   std::exp(-(std::pow(p[0] - force_center[0], 2) /
examples/step-62/step-62.cc:                                (2 * std::pow(force_sigma_x, 2)) +
examples/step-62/step-62.cc:                              std::pow(p[1] - force_center[1], 2) /
examples/step-62/step-62.cc:                                (2 * std::pow(force_sigma_y, 2))));
examples/step-62/step-62.cc:                      matrix_sum += -std::pow(omega, 2) *
examples/step-55/step-55.cc:                0.1 * std::pow(-std::sqrt(25.0 + 4 * pi2) + 5.0, 2) *
examples/step-55/step-55.cc:                0.05 * std::pow(-std::sqrt(25.0 + 4 * pi2) + 5.0, 3) *
examples/step-71/step-71.cc:        std::pow(1.0 / std::cosh(two_h_dot_h_div_h_sat_squ), 2.0);
examples/step-4/step-4.cc:    return_value += 4.0 * std::pow(p(i), 4.0);
examples/step-43/step-43.cc:    const double denominator = std::pow(temp, 2.0);
examples/step-44/step-44.cc:              std::pow((det_F_qp - J_tilde_qp), 2);
examples/step-53/step-53.cc:      std::atan2(x(2) + ep * ep * b * std::pow(std::sin(th), 3),
examples/step-53/step-53.cc:                  (ellipticity * ellipticity * R * std::pow(std::cos(th), 3))));

@Rombur
Copy link
Member

Rombur commented Feb 1, 2022

If you use p >= 3, you need to use ffast-math. It works at least up to pow(x, 8). It might have unwanted effect though.

@blaisb
Copy link
Member

blaisb commented Feb 1, 2022

That's funny, I thought this was automatic (at least for n =2,3,4). I thought this was only an issue if you made the mistake of writing std::pow(x,2.) instead of std::pow(x,2)...

@marcfehling
Copy link
Member

IIRC, we have Utilities::fixed_power (7c9f785) and Utilities::pow (#9764) that should be faster than std::pow.

@bangerth
Copy link
Member Author

bangerth commented Jan 9, 2024

Current state: The egrep expression only triggers for two more places. One is in step-4, where I would like to avoid the complication of having to explain why we use the function; the other is in step-82, and that change is currently stuck behind #16440.

So the only thing left to finish this is the change to step-82, which will require the following regex replace:

perl -pi -e 's/std::pow\(([^,]*), ([2345])(\.0?)?\)/Utilities::fixed_power<\2>(\1)/g;' examples/step-82/step-82.cc 

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging a pull request may close this issue.

4 participants