corecel/math/Turn.hh

0001 //------------------------------- -*- C++ -*- -------------------------------//
0002 // Copyright Celeritas contributors: see top-level COPYRIGHT file for details
0003 // SPDX-License-Identifier: (Apache-2.0 OR MIT)
0004 //---------------------------------------------------------------------------//
0005 //! \file corecel/math/Turn.hh
0006 //---------------------------------------------------------------------------//
0007 #pragma once
0008
0009 #include <cmath>
0010 #include <cstdlib>
0011 #include <type_traits>
0012
0013 #include "corecel/Constants.hh"
0014 #include "corecel/Types.hh"
0015
0016 #include "Algorithms.hh"
0017 #include "Quantity.hh"
0018
0019 namespace celeritas
0020 {
0021 //---------------------------------------------------------------------------//
0022 //! Unit for 2*pi radians
0023 struct TwoPi
0024 {
0025     static CELER_CONSTEXPR_FUNCTION Constant value()
0026     {
0027         return 2 * constants::pi;
0028     }
0029     //! Text label for output
0030     static char const* label() { return "tr"; }
0031 };
0032
0033 //---------------------------------------------------------------------------//
0034 //! Unit for pi/2 radians
0035 struct HalfPi
0036 {
0037     static CELER_CONSTEXPR_FUNCTION Constant value()
0038     {
0039         return constants::pi / 2;
0040     }
0041     //! Text label for output
0042     static char const* label() { return "qtr"; }
0043 };
0044
0045 //---------------------------------------------------------------------------//
0046 /*!
0047  * Quantity denoting a full turn.
0048  *
0049  * Turns are a useful way of representing angles without the historical
0050  * arbitrariness of degrees or the roundoff errors of radians. See, for
0051  * example, https://www.computerenhance.com/p/turns-are-better-than-radians .
0052  */
0053 template<class T>
0054 using Turn_t = Quantity<TwoPi, T>;
0055
0056 //! Turn with default precision (DEPRECATEDish)
0057 using Turn = Turn_t<real_type>;
0058 //! Turn with default precision
0059 using RealTurn = Turn_t<real_type>;
0060
0061 //! Create a turn using template deduction
0062 template<class T>
0063 CELER_CONSTEXPR_FUNCTION Turn_t<T> make_turn(T value)
0064 {
0065     static_assert(std::is_floating_point_v<T>,
0066                   "turn type must be floating point");
0067     return Turn_t<T>{value};
0068 }
0069
0070 //---------------------------------------------------------------------------//
0071 //! Quantity for an integer number of turns for axis swapping (DEPRECATEDish)
0072 using QuarterTurn = Quantity<HalfPi, int>;
0073 //! Quantity for an integer number of turns for axis swapping
0074 using IntQuarterTurn = Quantity<HalfPi, int>;
0075
0076 //---------------------------------------------------------------------------//
0077 //!@{
0078 //! Special overrides for math functions for more precise arithmetic
0079 template<class T>
0080 CELER_FORCEINLINE_FUNCTION T sin(Turn_t<T> r)
0081 {
0082     return sinpi(r.value() * 2);
0083 }
0084
0085 template<class T>
0086 CELER_FORCEINLINE_FUNCTION T cos(Turn_t<T> r)
0087 {
0088     return cospi(r.value() * 2);
0089 }
0090
0091 template<class T>
0092 CELER_FORCEINLINE_FUNCTION T tan(Turn_t<T> r)
0093 {
0094     return std::tan(native_value_from(r));
0095 }
0096
0097 template<class T>
0098 CELER_FORCEINLINE_FUNCTION void sincos(Turn_t<T> r, T* sinv, T* cosv)
0099 {
0100     return sincospi(r.value() * 2, sinv, cosv);
0101 }
0102
0103 CELER_CONSTEXPR_FUNCTION int cos(IntQuarterTurn r)
0104 {
0105     // Cosine is symmetric and periodic: modulo with 4
0106     // (note: avoid std::abs which isn't constexpr till C++23)
0107     auto i = (r.value() > 0 ? r.value() : -r.value()) % 4;
0108     // Map to {1, 0, -1, 0}[i] by encoding the 2-bit values into an int
0109     // offset by one: {2, 1, 0, 1}. Since ints are written big endian,
0110     // the bits below are:   { 1,0,1,2}
0111     constexpr int valbits = 0b01000110;
0112
0113     // Select the two bits we care about
0114     auto result_plus_one = (valbits >> (i << 1)) & 0b11;
0115     return result_plus_one - 1;
0116 }
0117
0118 CELER_CONSTEXPR_FUNCTION int sin(IntQuarterTurn r)
0119 {
0120     // Define in terms of the symmetric "cos": sin(x) = cos(x - pi/2)
0121     return cos(r - IntQuarterTurn{1});
0122 }
0123
0124 CELER_CONSTEXPR_FUNCTION void sincos(IntQuarterTurn r, int* sinv, int* cosv)
0125 {
0126     *sinv = sin(r);
0127     *cosv = cos(r);
0128 }
0129 //!@}
0130
0131 //! Math functions that return turns from types
0132 template<class T>
0133 CELER_FORCEINLINE_FUNCTION Turn_t<T> atan2turn(T y, T x)
0134 {
0135     // TODO: some hardware/libraries have this natively; use that instead
0136     return native_value_to<Turn_t<T>>(std::atan2(y, x));
0137 }
0138
0139 //---------------------------------------------------------------------------//
0140 }  // namespace celeritas