From 87c97f70fb5e7ccd2753990aa043842d08b3d4fe Mon Sep 17 00:00:00 2001
From: Jordan Peck <jordan.me2@gmail.com>
Date: Tue, 13 Oct 2026 21:54:16 +0100
Subject: [PATCH] Fix Log precision and Pow returning 1 for a zero base

- Log multiplied y by the e * -2.12194440e-4 term instead of adding it, which discarded the polynomial correction. Max error vs std::log drops from 615826 ULP to 1 ULP
- Log masks inputs <= 0 to 0, so Pow( 0, p ) returned Exp( 0 ) = 1. Pow now returns 0 for value <= 0
- Add "f32 pow zero base" test
---
 include/FastSIMD/ToolSet/Generic/Functions.h | 5 +++--
 tests/test.inl                               | 8 ++++++++
 2 files changed, 11 insertions(+), 2 deletions(-)

diff --git a/include/FastSIMD/ToolSet/Generic/Functions.h b/include/FastSIMD/ToolSet/Generic/Functions.h
index 217aeb8..ed5ca34 100644
--- a/include/FastSIMD/ToolSet/Generic/Functions.h
+++ b/include/FastSIMD/ToolSet/Generic/Functions.h
@@ -634,7 +634,7 @@ namespace FS
 
             RegisterF xx = value * value;
             y *= xx;
-            y *= e * RegisterF( -2.12194440e-4f );
+            y += e * RegisterF( -2.12194440e-4f );
             y -= xx * RegisterF( 0.5f );
 
             value += y;
@@ -648,12 +648,13 @@ namespace FS
         }
     }
     
+    // value <= 0 returns 0, matching 0^pow for pow > 0
     template<typename T, std::size_t N, FastSIMD::FeatureSet SIMD>
     FS_FORCEINLINE Register<T, N, SIMD> Pow( const Register<T, N, SIMD>& value, const Register<T, N, SIMD>& pow )
     {
         if constexpr( IsNativeV<Register<T, N, SIMD>> )
         {
-            return Exp( pow * Log( value ) );
+            return Masked( value > Register<T, N, SIMD>( 0 ), Exp( pow * Log( value ) ) );
         }
         else
         {
diff --git a/tests/test.inl b/tests/test.inl
index a24187a..716f40f 100644
--- a/tests/test.inl
+++ b/tests/test.inl
@@ -454,6 +454,14 @@ class FastSIMD::DispatchClass<TestFastSIMD<RegisterBytes, Relaxed>, SIMD> : publ
             b = FS::Min( FS::Max( b, TestRegf32( -8.0f ) ), TestRegf32( 8.0f ) );
             return FS::Pow( a, b ); 
         } ).relaxedAccuracy = 16384;
+        RegisterTest( tests, "f32 pow zero base", []( TestRegf32 a, TestRegf32 b ) {
+            // Negative inputs become a zero base, ranges avoid Exp overflow
+            auto zero = a < TestRegf32( 0.0f );
+            a = FS::Min( FS::Max( a, TestRegf32( 1.e-4f ) ), TestRegf32( 1.e+4f ) );
+            a = FS::Select( zero, TestRegf32( 0.0f ), a );
+            b = FS::Min( FS::Max( FS::Abs( b ), TestRegf32( 0.1f ) ), TestRegf32( 8.0f ) );
+            return FS::Pow( a, b );
+        } ).relaxedAccuracy = 16384;
 
         RegisterTest( tests, "i32 convert to f32", []( TestRegi32 a ) { return FS::Convert<float>( a ); } );
         RegisterTest( tests, "f32 convert to i32", []( TestRegf32 a ) { return FS::Convert<int32_t>( FS::Min( FS::Max( a, TestRegf32( -2147483647 - 1 ) ), TestRegf32( 2147483520 ) ) ); } );