fix: float32 unwrap bit-exact match with numpy using floor-based modulo

- Replace round-based ddmod with numpy's exact path: floor((dd+p2)/period)*period
- Add boundary_ambiguous handling: when dd>0 and ddmod==-p2, use +p2
- Restore float32 large-value unwrap test coverage (now bit-exact)
- All 475 tests pass

Author: peng.li24

.gitignore

@@ -1,6 +1,7 @@
 # Build directories
 build/
 cmake-build-*/
+issue/
 
 # Compiled output
 *.o

numpy/core.h

@@ -720,17 +720,26 @@ inline double safe_divide(double a, double b, double default_val = 0.0) {
 
 /// numpy.unwrap(p, discont=None, axis=-1)
 /// 1D only: unwrap phase angles by correcting jumps > discont.
+/// Uses numpy's exact computation path: (dd+period/2)%period - period/2
+/// to guarantee bit-identical float32 results.
 template<typename T>
 inline void unwrap(const T* src, T* dst, size_t n, T discont = T(M_PI)) {
     if (n == 0) return;
     T period = T(2) * discont;
+    T p2 = period / T(2);
     dst[0] = src[0];
     T cum_correct = T(0);
     for (size_t i = 1; i < n; ++i) {
         T dd = src[i] - src[i - 1];
         T ph_correct = T(0);
         if (std::abs(dd) >= discont) {
-            T ddmod = dd - std::round(dd / period) * period;
+            // numpy: ddmod = (dd + period/2) % period - period/2
+            // Python-style mod using floor division (numpy's mod):
+            T val = dd + p2;
+            T val_mod = val - std::floor(val / period) * period;
+            T ddmod = val_mod - p2;
+            // boundary_ambiguous: when dd > 0 and ddmod == -period/2, use +period/2
+            if (dd > T(0) && ddmod == -p2) ddmod = p2;
             ph_correct = ddmod - dd;
         }
         cum_correct += ph_correct;

tests/test_all.py

@@ -736,10 +736,10 @@ def test_unwrap(cpp):
     for dt in [np.float64, np.float32]:
         a = np.array([0.0, 0.5, 0.8, -0.9, -0.5, 0.2], dtype=dt)
         assert_bit_aligned(cpp.unwrap(a), np.unwrap(a), f"unwrap_{dt}")
-    # Large values: numpy uses float64 π internally even for float32 input,
-    # so float32 unwrap is not bit-exact on the correction path. Test float64 only.
-    a2 = np.array([0.0, 2.5, 5.0, -2.5, -5.0], dtype=np.float64) * np.pi
-    assert_bit_aligned(cpp.unwrap(a2), np.unwrap(a2), "unwrap_large")
+    # Large values — bit-exact for both float64 and float32
+    for dt in [np.float64, np.float32]:
+        a2 = np.array([0.0, 2.5, 5.0, -2.5, -5.0], dtype=dt) * np.pi
+        assert_bit_aligned(cpp.unwrap(a2), np.unwrap(a2), f"unwrap_large_{dt.__name__}")
 
 def test_cumsum(cpp):
     for dt in [np.float64, np.float32]: