Acts/Utilities/MathHelpers.hpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #pragma once
0010
0011 #include <cmath>
0012
0013 namespace Acts {
0014
0015 /// @brief Returns the absolute of a number
0016 ///        (Can be removed for c++ 23)
0017 /// @param n The number to take absolute value of
0018 /// @return The absolute value of the input
0019 template <typename T>
0020 constexpr T abs(const T n) {
0021   if (std::is_constant_evaluated()) {
0022     if constexpr (std::is_signed_v<T>) {
0023       if (n < 0) {
0024         return -n;
0025       }
0026     }
0027     return n;
0028   } else {
0029     return std::abs(n);
0030   }
0031 }
0032
0033 /// @brief Calculates the ordinary power of the number x.
0034 /// @param x: Number to take the power from
0035 /// @param p: Power to take
0036 /// @return x raised to the power p
0037 template <typename T, std::integral P>
0038 constexpr T pow(T x, P p) {
0039   if (std::is_constant_evaluated()) {
0040     constexpr T one = 1;
0041     if constexpr (std::is_signed_v<P>) {
0042       if (p < 0 && abs(x) > std::numeric_limits<T>::epsilon()) {
0043         x = one / x;
0044         p = -p;
0045       }
0046     }
0047     using unsigned_p = std::make_unsigned_t<P>;
0048     return p == 0 ? one : x * pow(x, static_cast<unsigned_p>(p) - 1);
0049   } else {
0050     return static_cast<T>(std::pow(x, static_cast<T>(p)));
0051   }
0052 }
0053
0054 /// @brief Returns the square of the passed number
0055 /// @param x The number to square
0056 /// @return The square of the input
0057 template <typename T>
0058 constexpr auto square(T x) {
0059   return x * x;
0060 }
0061
0062 /// @brief Calculates the sum of squares of arguments
0063 /// @param args Variable number of arguments to square and sum
0064 /// @return Sum of squares of all arguments
0065 template <typename... T>
0066 constexpr auto hypotSquare(T... args) {
0067   return (square(args) + ...);
0068 }
0069
0070 /// @brief Fast hypotenuse calculation for multiple arguments
0071 /// @param args Variable number of arguments
0072 /// @return Square root of sum of squares of arguments
0073 template <typename... T>
0074 constexpr auto fastHypot(T... args) {
0075   return std::sqrt(hypotSquare(args...));
0076 }
0077
0078 /// @brief Calculates the sum of 1 + 2 + 3+ ... + N using the
0079 ///        Gaussian sum formula
0080 /// @param N: Number until which the sum runs
0081 /// @return Sum of integers from 1 to N
0082 template <std::integral T>
0083 constexpr T sumUpToN(const T N) {
0084   return N * (N + 1) / 2;
0085 }
0086 /// @brief Calculates the factorial of a number
0087 ///        N!= N*(N-1)....*3*2*1
0088 /// @param upperN: Upper factor until which the factorial is calculated
0089 /// @param lowerN: Optional argument to remove the first factors from the calculation
0090 /// @return Factorial result
0091 template <std::integral T>
0092 constexpr T factorial(const T upperN, const T lowerN = 1) {
0093   constexpr T one = 1;
0094   const T& limit = std::max(one, lowerN);
0095   return upperN >= limit ? upperN * factorial(upperN - 1, limit) : one;
0096 }
0097 /// @brief Calculate the binomial coefficient
0098 ///              n        n!
0099 ///                 =  --------
0100 ///              k     k!(n-k)!
0101 /// @param n Upper value in binomial coefficient
0102 /// @param k Lower value in binomial coefficient
0103 /// @return Binomial coefficient n choose k
0104 template <std::integral T>
0105 constexpr T binomial(const T n, const T k) {
0106   return factorial<T>(n, n - k + 1) / factorial<T>(k);
0107 }
0108
0109 }  // namespace Acts