refactor: replace boolean template parameters with enum-based operation dispatch

seiko2plus · seiko2plus · commit e2e1b10a920f · 2025-09-17T13:52:52.000+03:00
- Replace IS_COS boolean template parameters with Operation enum class (kSin/kCos)
- Introduce NPSR_INTRIN macro for consistent intrinsic function attributes
- Update trigonometric function implementations to use new enum-based dispatch
- Improve code readability and maintainability by eliminating boolean parameter confusion
- Update function signatures and documentation to reflect new API design
- Replace HWY_ATTR with HWY_INLINE for better consistency in LUT implementation
diff --git a/npsr/hwy.h b/npsr/hwy.h
@@ -2,7 +2,15 @@
 #define NPSR_HWY_H_
 
 #include <hwy/highway.h>
-
+// This macro is used to define intrinsics that are:
+// NOTE: equals to HWY_API.
+// - always inlined
+// - flattened (no separate stack frame)
+// - marked maybe unused to suppress warnings when they are not used
+// NOTE: we do not need to use HWY_ATTR because we wrap Highway intrinsics in
+// HWY_BEFORE_NAMESPACE()/HWY_AFTER_NAMESPACE()
+// which implies the nessessary target attributes via #pargma.
+#define NPSR_INTRIN static HWY_INLINE HWY_FLATTEN HWY_MAYBE_UNUSED
 #endif  // NPSR_HWY_H_
 
 #if defined(NPSR_HWY_FOREACH_H_) == defined(HWY_TARGET_TOGGLE)  // NOLINT
diff --git a/npsr/lut-inl.h b/npsr/lut-inl.h
@@ -76,7 +76,7 @@ class Lut {
    * @note Index values must be in range [0, kCols)
    */
   template <typename VU, typename... OutV>
-  HWY_ATTR void Load(VU idx, OutV &...out) const {
+  HWY_INLINE void Load(VU idx, OutV &...out) const {
     static_assert(sizeof...(OutV) == kRows,
                   "Number of output vectors must match number of rows in LUT");
     using namespace hn;
@@ -102,7 +102,7 @@ class Lut {
 
   /// Dispatch to optimal row-load implementation based on vector/LUT size
   template <size_t Off = 0, typename VU, typename... OutV>
-  HWY_ATTR void LoadRow_(VU idx, OutV &...out) const {
+  HWY_INLINE void LoadRow_(VU idx, OutV &...out) const {
     using namespace hn;
     using DU = DFromV<VU>;
     const DU du;
@@ -126,7 +126,7 @@ class Lut {
 
   // Load using single table lookup (vector size == table width)
   template <size_t Off = 0, typename VInd, typename OutV0, typename... OutV>
-  HWY_ATTR void LoadX1_(const VInd &ind, OutV0 &out0, OutV &...out) const {
+  HWY_INLINE void LoadX1_(const VInd &ind, OutV0 &out0, OutV &...out) const {
     using namespace hn;
     using D = DFromV<OutV0>;
     const D d;
@@ -141,7 +141,7 @@ class Lut {
 
   // Load using two table lookups (vector size == table width / 2)
   template <size_t Off = 0, typename VInd, typename OutV0, typename... OutV>
-  HWY_ATTR void LoadX2_(const VInd &ind, OutV0 &out0, OutV &...out) const {
+  HWY_INLINE void LoadX2_(const VInd &ind, OutV0 &out0, OutV &...out) const {
     using namespace hn;
     using D = DFromV<OutV0>;
     const D d;
@@ -158,7 +158,7 @@ class Lut {
 
   //  General fallback using gather instructions
   template <size_t Off = 0, typename VU, typename OutV0, typename... OutV>
-  HWY_ATTR void LoadGather_(const VU &idx, OutV0 &out0, OutV &...out) const {
+  HWY_INLINE void LoadGather_(const VU &idx, OutV0 &out0, OutV &...out) const {
     using namespace hn;
     using D = DFromV<OutV0>;
     const D d;
diff --git a/npsr/trig/extended-inl.h b/npsr/trig/extended-inl.h
@@ -7,12 +7,13 @@
 
 #include "npsr/hwy.h"
 #include "npsr/trig/data/data.h"
+#include "npsr/trig/low-inl.h"  // Operation
 
 HWY_BEFORE_NAMESPACE();
 namespace npsr::HWY_NAMESPACE::trig {
 
-template <bool IS_COS, class V>
-HWY_API V Extended(V x) {
+template <Operation OP, class V>
+NPSR_INTRIN V Extended(V x) {
   using namespace hn;
   namespace data = ::npsr::trig::data;
   using hwy::ExponentBits;
@@ -143,9 +144,8 @@ HWY_API V Extended(V x) {
   // PHASE 5: Extract Quotient and Fractional Parts
   // =============================================================================
 
-  // Extract integer quotient
-  constexpr int kQuotientShift =
-      ExponentBits<T>() + 1;  // 9 for F32, 12 for F64
+  // Extract integer quotient. 9 for F32, 12 for F64
+  constexpr int kQuotientShift = ExponentBits<T>() + 1;
   VU u_shifted_n = ShiftRight<kQuotientShift>(u_n_hi);
 
   // fractional shifts derived from magic constants
@@ -255,8 +255,8 @@ HWY_API V Extended(V x) {
   // =============================================================================
 
   // Generated by npsr/trig/data/approx.h.sol
-  const T *table_base =
-      IS_COS ? data::kCosApproxTable<T> : data::kSinApproxTable<T>;
+  const T *table_base = OP == Operation::kCos ? data::kCosApproxTable<T>
+                                              : data::kSinApproxTable<T>;
 
   // Calculate table index
   VU u_n_mask = Set(du, kIsSingle ? 0xFF : 0x1FF);
diff --git a/npsr/trig/high-inl.h b/npsr/trig/high-inl.h
@@ -8,13 +8,14 @@
 #include "npsr/hwy.h"
 #include "npsr/lut-inl.h"
 #include "npsr/trig/data/data.h"
+#include "npsr/trig/low-inl.h"  // Operation
 
 HWY_BEFORE_NAMESPACE();
 
 namespace npsr::HWY_NAMESPACE::trig {
 
-template <bool IS_COS, typename V, HWY_IF_F32(TFromV<V>)>
-HWY_INLINE V High(V x) {
+template <Operation OP, typename V, HWY_IF_F32(TFromV<V>)>
+NPSR_INTRIN V High(V x) {
   using namespace hn;
   namespace data = ::npsr::trig::data;
 
@@ -37,7 +38,7 @@ HWY_INLINE V High(V x) {
   // Transform cosine to sine using identity: cos(x) = sin(x + π/2)
   const V half_pi = Set(d, data::kHalfPi<T>);
   V x_trans = x_abs;
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     x_trans = Add(x_abs, half_pi);
   }
   // check zero input/subnormal for cosine (cos(~0) = 1)
@@ -49,7 +50,7 @@ HWY_INLINE V High(V x) {
   V n = Sub(n_biased, magic_round);
 
   // Adjust quotient for cosine (accounts for π/2 phase shift)
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     // For cosine, we computed N = round((x + π/2)/π) but need N' for x:
     //   N = round((x + π/2)/π) = round(x/π + 0.5)
     // This is often 1 more than round(x/π), so we subtract 0.5:
@@ -83,7 +84,7 @@ HWY_INLINE V High(V x) {
   // Extract octant sign information from quotient and flip the sign bit
   poly = Xor(poly,
              BitCast(d, ShiftLeft<sizeof(T) * 8 - 1>(BitCast(du, n_biased))));
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     poly = IfThenElse(is_cos_near_zero, Set(d, 1.0f), poly);
   } else {
     // Restore original sign for sine (odd function)
@@ -113,8 +114,8 @@ HWY_INLINE V High(V x) {
  *    - cos(x) = cos(n*π/16 + r) = cos(n*π/16)*cos(r) - sin(n*π/16)*sin(r)
  *
  */
-template <bool IS_COS, typename V, HWY_IF_F64(TFromV<V>)>
-HWY_INLINE V High(V x) {
+template <Operation OP, typename V, HWY_IF_F64(TFromV<V>)>
+NPSR_INTRIN V High(V x) {
   using namespace hn;
   namespace data = ::npsr::trig::data;
   using T = TFromV<V>;
@@ -206,7 +207,7 @@ HWY_INLINE V High(V x) {
   // sin(n*π/16 + r) = sin_table + cos_table*remainder (+ corrections)
   // cos(n*π/16 + r) = cos_table - sin_table*remainder (+ corrections)
   V result;
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     // Cosine reconstruction: cos_table - sin_table*remainder
     // Equivalent to: cos(a)*cos(r) - sin(a)*sin(r) but more efficient
     V res_hi = NegMulAdd(r, sin_hi, cos_hi);  // cos_hi - r*sin_hi
diff --git a/npsr/trig/inl.h b/npsr/trig/inl.h
@@ -23,7 +23,6 @@
 HWY_BEFORE_NAMESPACE();
 
 namespace npsr::HWY_NAMESPACE::trig {
-
 /**
  * @brief Unified sine/cosine implementation with configurable precision
  *
@@ -32,13 +31,13 @@ namespace npsr::HWY_NAMESPACE::trig {
  * - Input magnitude (standard vs extended precision for large arguments)
  * - Special case handling (NaN, Inf)
  *
- * @tparam IS_COS   true for cosine, false for sine
+ * @tparam OP       Operation type: kSin or kCos
  * @tparam Prec     Precise configuration class with accuracy/feature flags
  * @tparam V        Highway vector type
  *
  * @param prec      Precise object controlling FP environment and exceptions
  * @param x         Input vector
- * @return          sin(x) or cos(x) depending on IS_COS
+ * @return          sin(x) or cos(x) depending on OP
  *
  * Algorithm selection:
  * 1. If kLowAccuracy: Use Low<> (Cody-Waite with minimal polynomial)
@@ -49,8 +48,8 @@ namespace npsr::HWY_NAMESPACE::trig {
  * - Float: |x| > 10,000 (empirically chosen for accuracy)
  * - Double: |x| > 2^24 (16,777,216 - where 53-bit mantissa loses precision)
  */
-template <bool IS_COS, typename Prec, typename V>
-HWY_API V SinCos(Prec &prec, V x) {
+template <Operation OP, typename Prec, typename V>
+NPSR_INTRIN V Trig(Prec &prec, V x) {
   using namespace hwy::HWY_NAMESPACE;
   constexpr bool kIsSingle = std::is_same_v<TFromV<V>, float>;
   const DFromV<V> d;
@@ -59,19 +58,19 @@ HWY_API V SinCos(Prec &prec, V x) {
   if constexpr (Prec::kLowAccuracy) {
     // Low precision: Cody-Waite reduction with degree-9 polynomial
     // Error: ~2 ULP and 3~ for non-fma
-    ret = Low<IS_COS>(x);
+    ret = Low<OP>(x);
   } else {
     // High precision: π/16 reduction with table lookup + polynomial
     // Error: ~1 ULP
-    ret = High<IS_COS>(x);
+    ret = High<OP>(x);
   }
   // Step 2: Handle special cases (NaN, Inf) if enabled
   auto is_finite = IsFinite(x);
   if constexpr (Prec::kSpecialCases) {
     // IEEE 754 requires: sin(±∞) = NaN, cos(±∞) = NaN
     ret = IfThenElse(is_finite, ret, NaN(d));
     // -0.0 should return -0.0 for sine
-    if constexpr (!IS_COS) {
+    if constexpr (OP == Operation::kSin) {
       ret = IfThenElse(Eq(x, Set(d, 0.0)), x, ret);
     }
   }
@@ -89,7 +88,7 @@ HWY_API V SinCos(Prec &prec, V x) {
     if (HWY_UNLIKELY(!AllFalse(d, has_large_arg))) {
       // Payne-Hanek reduction: Uses ~96-bit (float) or ~192-bit (double)
       // precision for 4/π to maintain accuracy for huge arguments
-      ret = IfThenElse(has_large_arg, Extended<IS_COS>(x), ret);
+      ret = IfThenElse(has_large_arg, Extended<OP>(x), ret);
     }
   }
   // Step 4: Raise invalid operation exception for infinity inputs
@@ -107,41 +106,45 @@ namespace npsr::HWY_NAMESPACE {
 /**
  * @brief Compute sine of vector elements with configurable precision
  *
- * @tparam Prec  Precise configuration (e.g., Precise<kLowAccuracy>)
+ * @tparam Prec  Precise configuration (e.g., Precise{kLowAccuracy})
  * @tparam V     Highway vector type
  * @param prec   Precise object managing FP environment
  * @param x      Input vector
  * @return       sin(x) for each element
  *
  * @example
  * ```cpp
- * Precise<kHighAccuracy> prec;
+ * Precise prec{
+ *  kLowAccuracy, kNoLargeArgument, kNoExceptions, kNoSpecialCases
+ * };
  * auto result = Sin(prec, input_vector);
  * ```
  */
 template <typename Prec, typename V>
-HWY_API V Sin(Prec &prec, V x) {
-  return trig::SinCos<false>(prec, x);
+NPSR_INTRIN V Sin(Prec &prec, V x) {
+  return trig::Trig<trig::Operation::kSin>(prec, x);
 }
 
 /**
  * @brief Compute cosine of vector elements with configurable precision
  *
- * @tparam Prec  Precise configuration (e.g., Precise<kLowAccuracy>)
+ * @tparam Prec  Precise configuration (e.g., Precise{kLowAccuracy})
  * @tparam V     Highway vector type
  * @param prec   Precise object managing FP environment
  * @param x      Input vector
  * @return       cos(x) for each element
  *
  * @example
  * ```cpp
- * Precise<kNoLargeArgument, kNoSpecialCases> prec;
+ * Precise prec{
+ *  kLowAccuracy, kNoLargeArgument, kNoExceptions, kNoSpecialCases
+ * };
  * auto result = Cos(prec, input_vector);
  * ```
  */
 template <typename Prec, typename V>
-HWY_API V Cos(Prec &prec, V x) {
-  return trig::SinCos<true>(prec, x);
+NPSR_INTRIN V Cos(Prec &prec, V x) {
+  return trig::Trig<trig::Operation::kCos>(prec, x);
 }
 
 }  // namespace npsr::HWY_NAMESPACE
diff --git a/npsr/trig/low-inl.h b/npsr/trig/low-inl.h
@@ -12,19 +12,22 @@ HWY_BEFORE_NAMESPACE();
 
 namespace npsr::HWY_NAMESPACE::trig {
 
-template <bool IS_COS, typename V, HWY_IF_F32(TFromV<V>)>
-HWY_API V PolyLow(V r, V r2) {
+enum class Operation { kSin = 0, kCos = 1 };
+
+template <Operation OP, typename V, HWY_IF_F32(TFromV<V>)>
+NPSR_INTRIN V PolyLow(V r, V r2) {
   using namespace hn;
 
   const DFromV<V> d;
-  const V c9 = Set(d, IS_COS ? 0x1.5d866ap-19f : 0x1.5dbdfp-19f);
-  const V c7 = Set(d, IS_COS ? -0x1.9f6d9ep-13 : -0x1.9f6ffep-13f);
-  const V c5 = Set(d, IS_COS ? 0x1.110ec8p-7 : 0x1.110eccp-7f);
+  constexpr bool kCos = OP == Operation::kCos;
+  const V c9 = Set(d, kCos ? 0x1.5d866ap-19f : 0x1.5dbdfp-19f);
+  const V c7 = Set(d, kCos ? -0x1.9f6d9ep-13 : -0x1.9f6ffep-13f);
+  const V c5 = Set(d, kCos ? 0x1.110ec8p-7 : 0x1.110eccp-7f);
   const V c3 = Set(d, -0x1.55554cp-3f);
   V poly = MulAdd(c9, r2, c7);
   poly = MulAdd(r2, poly, c5);
   poly = MulAdd(r2, poly, c3);
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     // Although this path handles cosine, we have already transformed the
     // input using the identity: cos(x) = sin(x + π/2) This means we're no
     // longer directly evaluating a cosine Taylor series; instead, we evaluate
@@ -48,8 +51,8 @@ HWY_API V PolyLow(V r, V r2) {
   return poly;
 }
 
-template <bool IS_COS, typename V, HWY_IF_F64(TFromV<V>)>
-HWY_API V PolyLow(V r, V r2) {
+template <Operation OP, typename V, HWY_IF_F64(TFromV<V>)>
+NPSR_INTRIN V PolyLow(V r, V r2) {
   using namespace hn;
 
   const DFromV<V> d;
@@ -69,8 +72,8 @@ HWY_API V PolyLow(V r, V r2) {
   return poly;
 }
 
-template <bool IS_COS, typename V>
-HWY_API V Low(V x) {
+template <Operation OP, typename V>
+NPSR_INTRIN V Low(V x) {
   using namespace hn;
   using hwy::SignMask;
   namespace data = ::npsr::trig::data;
@@ -87,7 +90,7 @@ HWY_API V Low(V x) {
   // Transform cosine to sine using identity: cos(x) = sin(x + π/2)
   const V half_pi = Set(d, data::kHalfPi<T>);
   V x_trans = x_abs;
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     x_trans = Add(x_abs, half_pi);
   }
   // check zero input/subnormal for cosine (cos(~0) = 1)
@@ -100,7 +103,7 @@ HWY_API V Low(V x) {
   V n = Sub(n_biased, magic_round);
 
   // Adjust quotient for cosine (accounts for π/2 phase shift)
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     // For cosine, we computed N = round((x + π/2)/π) but need N' for x:
     //   N = round((x + π/2)/π) = round(x/π + 0.5)
     // This is often 1 more than round(x/π), so we subtract 0.5:
@@ -124,7 +127,7 @@ HWY_API V Low(V x) {
     r = r_lo;
   }
   V r2 = Mul(r, r);
-  V poly = PolyLow<IS_COS>(r, r2);
+  V poly = PolyLow<OP>(r, r2);
 
   if constexpr (!kIsSingle) {
     V r2_corr = Mul(r2, r_lo);
@@ -134,7 +137,7 @@ HWY_API V Low(V x) {
   // Extract octant sign information from quotient and flip the sign bit
   poly = Xor(poly,
              BitCast(d, ShiftLeft<sizeof(T) * 8 - 1>(BitCast(du, n_biased))));
-  if constexpr (IS_COS) {
+  if constexpr (OP == Operation::kCos) {
     poly = IfThenElse(is_cos_near_zero, Set(d, static_cast<T>(1.0)), poly);
   } else {
     // Restore original sign for sine (odd function)
