trig improved

Author: SkibidiProduction

numpower.c

@@ -474,6 +474,11 @@ static NDArray *ndarray_promote_and_op(zend_uchar opcode, NDArray *nda,
                                         NDArray *ndb,
                                         const char **result_type_out);
 
+/* Forward declaration: the strict numeric-string dtype inferrer is defined
+   further down; the shared string-scalar validator below (used by both the
+   binary and unary intakes) needs it before its definition. */
+static const char *ndarray_infer_dtype_from_string(const char *str, size_t len);
+
 /**
  * @brief Widen an integer dtype to floating point for binary ops that always
  *        return a float (true division and atan2).
@@ -602,6 +607,50 @@ static NDArray *ndarray_make_typed_scalar(zval *z, const char *target_dt)
     return r;
 }
 
+/**
+ * @brief Infer the dtype of a numeric-string scalar, throwing a precise
+ *        diagnostic on a malformed / empty / whitespace-only literal.
+ *
+ * Shared by the binary (`ndarray_arith_resolve_operand`) and unary
+ * (`ndarray_resolve_unary_input`) string-scalar intakes so both reject
+ * non-numeric input identically — the bare `strto*` parsers behind
+ * `NDArray_EncodeZvalToDtype` ignore trailing junk and read "abc" as 0,
+ * which would otherwise let `arctan2($t, "abc")` silently compute
+ * `atan2(x, 0)`. The three failure modes are split out so a typo is easy
+ * to spot.
+ *
+ * @param[in] value PHP zval of type IS_STRING.
+ * @return Canonical dtype string on success; NULL after throwing a PHP
+ *         exception when @p value is not a valid numeric literal.
+ */
+static const char *ndarray_string_scalar_dtype_or_throw(zval *value)
+{
+    const char *p  = Z_STRVAL_P(value);
+    size_t      n  = Z_STRLEN_P(value);
+    const char *dt = ndarray_infer_dtype_from_string(p, n);
+    if (dt != NULL) {
+        return dt;
+    }
+    size_t ws = 0;
+    while (ws < n && (p[ws] == ' ' || p[ws] == '\t' ||
+                      p[ws] == '\n' || p[ws] == '\r')) {
+        ws++;
+    }
+    if (n == 0) {
+        zend_throw_error(NULL, "Numeric string expected, got an empty value.");
+    } else if (ws == n) {
+        zend_throw_error(NULL,
+            "Numeric string expected, got a whitespace-only value.");
+    } else {
+        /* Length-bounded (`%.*s`) so an embedded NUL doesn't truncate the
+           offending literal in the diagnostic. */
+        zend_throw_error(NULL,
+            "Numeric string expected, got malformed literal: \"%.*s\".",
+            (int)n, p);
+    }
+    return NULL;
+}
+
 /**
  * @brief Resolve a PHP operand zval to an NDArray, honouring weak-scalar
  *        promotion against the peer operand's dtype.
@@ -630,6 +679,16 @@ static NDArray *ndarray_arith_resolve_operand(zval *value, NDArray *other,
 {
     *is_owned = 0;
     if (other != NULL && ndarray_is_promotable_scalar(value)) {
+        /* A string operand must be a syntactically valid numeric literal:
+           validate with the same strict inferrer the unary intake uses so
+           malformed / empty / whitespace input throws instead of being
+           silently coerced to 0. The value still adopts the peer's dtype
+           (below) for loss-free float128 / uint64 intake — we only borrow the
+           inferrer's syntax check, not its inferred dtype. */
+        if (Z_TYPE_P(value) == IS_STRING &&
+            ndarray_string_scalar_dtype_or_throw(value) == NULL) {
+            return NULL;
+        }
         const char *target = ndarray_pick_scalar_dtype(NDArray_TYPE(other), value);
         NDArray *r = ndarray_make_typed_scalar(value, target);
         if (r == NULL) {
@@ -812,12 +871,12 @@ static NDArray *ndarray_promote_and_op(zend_uchar opcode, NDArray *nda, NDArray
 
     if (both_gpu) {
         /* GPU stays on GPU for every supported dtype. We promote types, cast
-           on GPU via NDArray_AsType (now GPU-aware), call the typed GPU binop,
+           on GPU via NDArray_AsType (GPU-aware), call the typed GPU binop,
            then cast back. No CPU round-trip for float32, float64, float16,
-           int8..uint64, and (via dd kernels) float128.
-           float4/float8 fall back to CPU because there are no native CUDA
-           intrinsics and they're 1-byte values (we go through NDArray_AsType
-           which already routes them through CPU for those source/target types). */
+           int8..uint64, float4/float8 (cast on GPU via cuda_cast_fp4/fp8_*),
+           and (via dd kernels) float128. The CPU fallthrough below the
+           NDArray_TypedBinOp_GPU call is a defensive guard for any AsType cast
+           a given build cannot keep on the device. */
         const char *gpu_result_type =
             ndarray_binop_result_type(opcode, NDArray_TYPE(nda), NDArray_TYPE(ndb));
         const char *gpu_comp_type = compute_dtype_for_arithmetic(gpu_result_type);
@@ -3723,32 +3782,8 @@ static const char *ndarray_infer_dtype_from_string(const char *str, size_t len)
 static NDArray *ndarray_resolve_unary_input(zval *array, int *owned)
 {
     if (Z_TYPE_P(array) == IS_STRING) {
-        const char *dt = ndarray_infer_dtype_from_string(
-            Z_STRVAL_P(array), Z_STRLEN_P(array));
+        const char *dt = ndarray_string_scalar_dtype_or_throw(array);
         if (dt == NULL) {
-            /* Differentiate the three failure modes so callers can spot
-               typos quickly: empty literal, whitespace-only, or
-               syntactically malformed (non-empty, non-whitespace). */
-            const char *p = Z_STRVAL_P(array);
-            size_t      n = Z_STRLEN_P(array);
-            size_t      ws = 0;
-            while (ws < n && (p[ws] == ' ' || p[ws] == '\t' ||
-                              p[ws] == '\n' || p[ws] == '\r')) {
-                ws++;
-            }
-            if (n == 0) {
-                zend_throw_error(NULL,
-                    "Numeric string expected, got an empty value.");
-            } else if (ws == n) {
-                zend_throw_error(NULL,
-                    "Numeric string expected, got a whitespace-only value.");
-            } else {
-                /* Length-bounded (`%.*s`) so an embedded NUL doesn't
-                   truncate the offending literal in the diagnostic. */
-                zend_throw_error(NULL,
-                    "Numeric string expected, got malformed literal: \"%.*s\".",
-                    (int)n, p);
-            }
             return NULL;
         }
         NDArray *nda = ndarray_make_typed_scalar(array, dt);

src/ndmath/arithmetics.c

@@ -2614,14 +2614,35 @@ NDArray* NDArray_Mod_Float128(NDArray* a, NDArray* b) {
    numerator / y-coordinate, b the denominator / x-coordinate). */
 #define DEFINE_ATAN2_FLOAT_CPU(NAME, T, DT_CONST, FN)                              \
 NDArray* NAME(NDArray* a, NDArray* b) {                                            \
+    /* Expand a 0-D scalar operand to the peer's shape first, matching            \
+       NDArray_Add_Double et al. AND the GPU path: NumPy/PyTorch broadcast        \
+       arctan2(0-D, shape-(n,)) to (n,), never to a 0-D scalar. Without this,     \
+       a 0-D vs numel-1-array pair would tie on element count, fall through to    \
+       the `else` below, and take the 0-D operand's rank — yielding a 0-D CPU     \
+       result while the GPU (which broadcasts the scalar) returns (n,).           \
+       NDArray_Broadcast replicates the single element for any dtype. */          \
+    NDArray *sa = NULL, *sb = NULL;                                               \
+    if (NDArray_NDIM(a) == 0 && NDArray_NDIM(b) > 0) {                            \
+        sa = NDArray_Broadcast(a, b); if (sa == NULL) return NULL; a = sa;        \
+    } else if (NDArray_NDIM(b) == 0 && NDArray_NDIM(a) > 0) {                     \
+        sb = NDArray_Broadcast(b, a); if (sb == NULL) return NULL; b = sb;        \
+    }                                                                             \
     NDArray *broadcasted = NULL, *a_broad, *b_broad;                              \
     if (NDArray_NUMELEMENTS(a) < NDArray_NUMELEMENTS(b)) {                        \
         broadcasted = NDArray_Broadcast(a, b);                                    \
-        if (broadcasted == NULL) return NULL;                                     \
+        if (broadcasted == NULL) {                                                \
+            if (sa) NDArray_FREE(sa);                                             \
+            if (sb) NDArray_FREE(sb);                                             \
+            return NULL;                                                          \
+        }                                                                         \
         a_broad = broadcasted; b_broad = b;                                       \
     } else if (NDArray_NUMELEMENTS(b) < NDArray_NUMELEMENTS(a)) {                 \
         broadcasted = NDArray_Broadcast(b, a);                                    \
-        if (broadcasted == NULL) return NULL;                                     \
+        if (broadcasted == NULL) {                                                \
+            if (sa) NDArray_FREE(sa);                                             \
+            if (sb) NDArray_FREE(sb);                                             \
+            return NULL;                                                          \
+        }                                                                         \
         b_broad = broadcasted; a_broad = a;                                       \
     } else { a_broad = a; b_broad = b; }                                          \
     int ndim   = NDArray_NDIM(a_broad);                                           \
@@ -2630,15 +2651,19 @@ NDArray* NAME(NDArray* a, NDArray* b) {
     else          shape[0] = 1;                                                   \
     NDArray *result = NDArray_Empty(shape, ndim, DT_CONST, NDARRAY_DEVICE_CPU);   \
     if (result == NULL) {                                                         \
-        if (broadcasted) NDArray_FREE(broadcasted);                              \
+        if (broadcasted) NDArray_FREE(broadcasted);                               \
+        if (sa) NDArray_FREE(sa);                                                 \
+        if (sb) NDArray_FREE(sb);                                                 \
         return NULL;                                                              \
-    }                                                                            \
+    }                                                                             \
     T       *rd = (T *)NDArray_DATA(result);                                      \
     const T *ad = (const T *)NDArray_DATA(a_broad);                               \
     const T *bd = (const T *)NDArray_DATA(b_broad);                               \
     long n = NDArray_NUMELEMENTS(result);                                         \
     for (long i = 0; i < n; i++) rd[i] = FN(ad[i], bd[i]);                        \
     if (broadcasted) NDArray_FREE(broadcasted);                                   \
+    if (sa) NDArray_FREE(sa);                                                     \
+    if (sb) NDArray_FREE(sb);                                                     \
     return result;                                                                \
 }
 DEFINE_ATAN2_FLOAT_CPU(NDArray_Arctan2_Float,  float,  NDARRAY_TYPE_FLOAT32, atan2f)
@@ -2928,9 +2953,9 @@ ndarray_int_binop_cpu(NDArray *a, NDArray *b, int opcode) {
  * Both operands must already share the same dtype — the calling
  * dispatcher (`ndarray_promote_and_op`) handles the cast through
  * `NDArray_AsType` before reaching this entry point. Falls back to
- * NULL with a PHP error for opcodes outside the supported set; `/` is
- * already promoted to a float dtype by `ndarray_div_promote` and never
- * reaches here.
+ * NULL with a PHP error for opcodes outside the supported set; `/` (and
+ * `arctan2`) are already promoted to a float dtype by
+ * `ndarray_widen_int_to_float` and never reach here.
  *
  * @param[in] opcode ZEND_ADD / SUB / MUL / MOD / POW.
  * @param[in] a, b   Same-dtype operands on CPU; one may be 0-D.

src/ndmath/arithmetics.h

@@ -86,7 +86,7 @@ NDArray* NDArray_TypedBinOp_GPU(int opcode, NDArray* a, NDArray* b);
       and GPU once intermediates spilled past 2^53).
    Both operands must be on CPU and share the dtype. One operand may be
    a 0-D scalar. opcode is ZEND_ADD/SUB/MUL/MOD/POW; ZEND_DIV is promoted
-   to float upstream by ndarray_div_promote and never reaches here. */
+   to float upstream by ndarray_widen_int_to_float and never reaches here. */
 NDArray* NDArray_TypedBinOp_CPU_Int(int opcode, NDArray* a, NDArray* b);
 
 NDArray* NDArray_Add(NDArray* a, NDArray* b);

tests/math/120-arctan2-all-dtypes.phpt

@@ -93,17 +93,33 @@ $qy = NumPower::array(['1.0', '0.0',  '1.0'], 'float128');
 $rq = NumPower::arctan2($qx, $qy);
 ok(dt($rq) === 'float128', 'fp128 result dtype');
 $rqa = $rq->toArray();
-/* π/4 and π/2 to > 30 digits (decimal reference). */
-$PI_4 = '0.78539816339744830961566084581988';
-$PI_2 = '1.5707963267948966192313216916398';
-ok(strncmp((string)$rqa[0], $PI_4, 30) === 0, "fp128 atan2(1,1)=π/4 got={$rqa[0]}");
-ok(strncmp((string)$rqa[1], $PI_2, 30) === 0, "fp128 atan2(1,0)=π/2 got={$rqa[1]}");
-ok(near((float)$rqa[2], -M_PI / 4, 1e-13), 'fp128 atan2(-1,1)=-π/4');
-
-/* String scalar adopts the fp128 peer dtype → full precision. */
+$PI_4 = '0.78539816339744830961566084581988';   /* decimal reference, 32 digits */
+/* Value checks use fp64 tolerance so they are PORTABLE: with libquadmath the
+   fp128 atan2 is full 113-bit, but on the double-double fallback build
+   (macOS / non-x86) fp128 *transcendentals* compute at fp64 precision (the
+   NDARRAY_FP128_ATAN2 macro routes through atan2(double)). A >fp64 digit-prefix
+   assertion would pass on Linux and FAIL on a DD build — so only assert it
+   behind a runtime libquadmath probe below. */
+ok(dt($rq) === 'float128', 'fp128 result dtype');
+ok(near((float)$rqa[0],  M_PI / 4, 1e-13), 'fp128 atan2(1,1)=π/4 (fp64 tol)');
+ok(near((float)$rqa[1],  M_PI / 2, 1e-13), 'fp128 atan2(1,0)=π/2 (fp64 tol)');
+ok(near((float)$rqa[2], -M_PI / 4, 1e-13), 'fp128 atan2(-1,1)=-π/4 (fp64 tol)');
+
+/* Probe for full-precision fp128 transcendentals (libquadmath): if the result
+   already agrees with π/4 to 20 sig digits it cannot be the fp64 (~16-digit)
+   DD result, so the build has libquadmath and we can assert the 30-digit
+   prefix. On a DD build this probe is false and the strict check is skipped
+   (the fp64-tolerance checks above already cover correctness there). */
+$has_quadmath = (strncmp((string)$rqa[0], $PI_4, 20) === 0);
+if ($has_quadmath) {
+    ok(strncmp((string)$rqa[0], $PI_4, 30) === 0, 'fp128 atan2(1,1)=π/4 full precision (libquadmath)');
+}
+
+/* String scalar adopts the fp128 peer dtype (intake is loss-free on every
+   build); the atan2 result value is checked at fp64 tolerance for portability. */
 $rqs = NumPower::arctan2(NumPower::array(['1.0', '1.0'], 'float128'), '1.0');
 ok(dt($rqs) === 'float128', 'fp128 + string scalar dtype');
-ok(strncmp((string)$rqs->toArray()[0], $PI_4, 30) === 0, 'fp128 + string scalar value');
+ok(near((float)$rqs->toArray()[0], M_PI / 4, 1e-13), 'fp128 + string scalar value');
 
 /* uint64 string intake keeps a > 2^53 magnitude loss-free in the denominator
    (result ~0 because the numerator 1 is tiny next to it, but the point is the
@@ -156,8 +172,11 @@ ok(is_float($s) && near($s, M_PI / 4, 1e-14), '0-D float64 → PHP float');
    arithmetic operators — still a PHP float, at float32 precision. */
 $sbare = NumPower::arctan2(1.0, 1.0);
 ok(is_float($sbare) && near($sbare, M_PI / 4, 1e-6), '0-D bare scalar → PHP float');
+/* The important contract here is that a 0-D float128 result returns as a PHP
+   *string* (not a lossy float); the value is checked at fp64 tolerance so the
+   assertion is portable across the libquadmath and double-double builds. */
 $sq = NumPower::arctan2(NumPower::array('1.0', 'float128'), NumPower::array('1.0', 'float128'));
-ok(is_string($sq) && strncmp($sq, $PI_4, 30) === 0, '0-D float128 → PHP string');
+ok(is_string($sq) && near((float)$sq, M_PI / 4, 1e-13), '0-D float128 → PHP string');
 
 /* ── Edge values: ±inf, NaN, signed zero ────────────────────────────────── */
 $ex = NumPower::array([INF,  INF, 1.0, NAN, 1.0], 'float64');
@@ -188,6 +207,40 @@ ok(abs((float)$pr->toArray()[0] - $want) < 1e-15, 'precision guard: full float64
 $e0 = NumPower::arctan2(NumPower::zeros([0, 4], 'float64'), NumPower::zeros([0, 4], 'float64'));
 ok($e0->shape() === [0, 4] && dt($e0) === 'float64', 'empty (0,4) shape+dtype');
 
+/* ── 0-D scalar broadcasts to a numel-1 array's shape (regression: CPU must
+      not collapse to a 0-D scalar — it has to match the GPU / NumPy result
+      shape (1,), see the DEFINE_ATAN2_FLOAT_CPU 0-D expansion) ────────────── */
+$sc1 = NumPower::arctan2(NumPower::array(1.0, 'float64'), NumPower::array([2.0], 'float64'));
+ok(is_object($sc1) && $sc1->shape() === [1], '0-D numerator + (1,) denominator -> shape (1,)');
+ok(near($sc1->toArray()[0], atan2(1.0, 2.0), 1e-14), '0-D + (1,) value');
+$sc2 = NumPower::arctan2(NumPower::array([2.0], 'float64'), NumPower::array(1.0, 'float64'));
+ok(is_object($sc2) && $sc2->shape() === [1], '(1,) numerator + 0-D denominator -> shape (1,)');
+/* float128 takes the same 0-D-expansion path */
+$sc3 = NumPower::arctan2(NumPower::array('1.0', 'float128'), NumPower::array(['2.0'], 'float128'));
+ok(is_object($sc3) && $sc3->shape() === [1], 'fp128 0-D + (1,) -> shape (1,)');
+/* a genuine 0-D pair still collapses to a PHP scalar */
+ok(is_float(NumPower::arctan2(NumPower::array(1.0, 'float64'), NumPower::array(1.0, 'float64'))),
+   '0-D + 0-D -> PHP scalar');
+
+/* ── String-operand validation: malformed literals throw, they are NOT
+      silently coerced to 0 (regression: the binary dispatch must reject
+      garbage like the unary path does) ─────────────────────────────────────── */
+$peer = NumPower::array(['1.0'], 'float128');
+foreach (['abc', '', '   ', '1.5.5', '0xff', '1,5'] as $bad) {
+    $threw = false;
+    try { NumPower::arctan2($peer, $bad); } catch (\Throwable $e) { $threw = true; }
+    ok($threw, "arctan2(fp128, malformed '" . addslashes($bad) . "') throws");
+}
+/* valid numeric strings (incl. inf / nan / exponent / sign) are still accepted.
+   Use a float64 peer so the 0-D result returns as a PHP float (PHP's
+   (float)"nan" cast yields 0.0, not NAN — testing via the float64 toArray
+   element is the reliable check). */
+$peerf = NumPower::array([1.0], 'float64');
+ok(near(NumPower::arctan2($peerf, 'inf')->toArray()[0], 0.0, 1e-30), "string 'inf' -> atan2(1,inf)=0");
+ok(is_nan(NumPower::arctan2($peerf, 'nan')->toArray()[0]), "string 'nan' -> NaN");
+ok(near(NumPower::arctan2($peerf, '-2')->toArray()[0], atan2(1.0, -2.0), 1e-12), "string '-2' accepted");
+ok(near(NumPower::arctan2($peerf, '1e3')->toArray()[0], atan2(1.0, 1000.0), 1e-12), "string '1e3' accepted");
+
 echo $FAILS === 0 ? "ALL CHECKS PASSED\n" : "TOTAL FAILURES: $FAILS\n";
 ?>
 --EXPECT--

tests/math/121-arctan2-cpu-gpu-parity.phpt

@@ -95,6 +95,18 @@ for ($i = 0; $i < 2; $i++) for ($j = 0; $j < 3; $j++)
     if (abs($cpu_m[$i][$j] - $gpu_m[$i][$j]) > 1e-12) $m_ok = false;
 ok($m_ok, 'broadcast row-vector → matrix CPU/GPU parity');
 
+/* ── 0-D scalar + numel-1 array: CPU and GPU must agree on BOTH shape and
+      value (CPU must broadcast the scalar, not collapse to 0-D) ───────────── */
+$z0 = NumPower::array(1.0, 'float64');
+$z1 = NumPower::array([2.0], 'float64');
+$cpu_z = NumPower::arctan2($z0, $z1);
+$gpu_z = NumPower::arctan2($z0->gpu(), $z1->gpu());
+ok($gpu_z->isGPU(), '0-D+(1,) stays on GPU');
+ok(is_object($cpu_z) && $cpu_z->shape() === $gpu_z->cpu()->shape(),
+   '0-D+(1,) CPU/GPU shape match');
+ok(abs($cpu_z->toArray()[0] - $gpu_z->cpu()->toArray()[0]) < 1e-12,
+   '0-D+(1,) CPU/GPU value match');
+
 /* ── Mixed device: a CPU scalar operand migrates to the GPU array's device ─ */
 $ga = NumPower::array([1.0, -1.0, 2.0, -2.0], 'float64')->gpu();
 $mixed = NumPower::arctan2(2.0, $ga);          /* CPU scalar numerator + GPU array */