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 template <typename T>
0018 constexpr T abs(const T n) {
0019   if constexpr (std::is_signed_v<T>) {
0020     if (n < 0) {
0021       return -n;
0022     }
0023   }
0024   return n;
0025 }
0026
0027 /// @brief Calculates the ordinary power of the number x.
0028 /// @param x: Number to take the power from
0029 /// @param p: Power to take
0030 template <typename T, std::integral P>
0031 constexpr T pow(T x, P p) {
0032   constexpr T one = 1;
0033   if constexpr (std::is_signed_v<P>) {
0034     if (p < 0 && abs(x) > std::numeric_limits<T>::epsilon()) {
0035       x = one / x;
0036       p = -p;
0037     }
0038   }
0039   using unsigned_p = std::make_unsigned_t<P>;
0040   return p == 0 ? one : x * pow(x, static_cast<unsigned_p>(p) - 1);
0041 }
0042
0043 /// @brief Returns the square of the passed number
0044 template <typename T>
0045 constexpr auto square(T x) {
0046   return x * x;
0047 }
0048
0049 template <typename... T>
0050 constexpr auto hypotSquare(T... args) {
0051   return (square(args) + ...);
0052 }
0053
0054 template <typename... T>
0055 constexpr auto fastHypot(T... args) {
0056   return std::sqrt(hypotSquare(args...));
0057 }
0058
0059 /// @brief Calculates the sum of 1 + 2 + 3+ ... + N using the
0060 ///        Gaussian sum formula
0061 /// @param N: Number until which the sum runs
0062 template <std::integral T>
0063 constexpr T sumUpToN(const T N) {
0064   return N * (N + 1) / 2;
0065 }
0066 /// @brief Calculates the factorial of a number
0067 ///        N!= N*(N-1)....*3*2*1
0068 /// @param upperN: Upper factor until which the factorial is calculated
0069 /// @param lowerN: Optional argument to remove the first factors from the calculation
0070 template <std::integral T>
0071 constexpr T factorial(const T upperN, const T lowerN = 1) {
0072   constexpr T one = 1;
0073   const T& limit = std::max(one, lowerN);
0074   return upperN >= limit ? upperN * factorial(upperN - 1, limit) : one;
0075 }
0076 /// @brief Calculate the binomial coefficient
0077 ///              n        n!
0078 ///                 =  --------
0079 ///              k     k!(n-k)!
0080 template <std::integral T>
0081 constexpr T binomial(const T n, const T k) {
0082   return factorial<T>(n, n - k + 1) / factorial<T>(k);
0083 }
0084
0085 }  // namespace Acts