Repository: scikit-hep/vector

Issue:
`deltaangle` between identical vectors return nan
### Vector Version

1.3.1

### Python Version

3.11.4

### OS / Environment

Kubuntu Linux 22.04 

### Describe the bug

When taking two identical vectors and computing the angle between them, I expect the result to be 0. Instead, I get NaN:
```python
> v = vector.obj(x=1, y=1, z=1)
> v.deltaangle(v)
nan
```
This is also true when using vector.array. I have not tried other backends.

### Any additional but relevant log output

_No response_

Discussion:
codecov[bot]: ## [Codecov](https://app.codecov.io/gh/scikit-hep/vector/pull/465?dropdown=coverage&src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=scikit-hep) Report
All modified and coverable lines are covered by tests :white_check_mark:
> Project coverage is 86.78%. Comparing base [(`fb13df8`)](https://app.codecov.io/gh/scikit-hep/vector/commit/fb13df8a06dc36049211b8a7783adf62ffd8d4c9?dropdown=coverage&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=scikit-hep) to head [(`3ad4e8e`)](https://app.codecov.io/gh/scikit-hep/vector/pull/465?dropdown=coverage&src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=scikit-hep).


<details><summary>Additional details and impacted files</summary>


```diff
@@           Coverage Diff           @@
##             main     #465   +/-   ##
=======================================
  Coverage   86.77%   86.78%           
=======================================
  Files          98       98           
  Lines       11927    11927           
=======================================
+ Hits        10350    10351    +1     
+ Misses       1577     1576    -1     
```



</details>

[:umbrella: View full report in Codecov by Sentry](https://app.codecov.io/gh/scikit-hep/vector/pull/465?dropdown=coverage&src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=scikit-hep).   
:loudspeaker: Have feedback on the report? [Share it here](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=scikit-hep).

Saransh-cpp: A gentle bump on this, @jpivarski.

Patch:
diff --git a/src/vector/_compute/spatial/deltaangle.py b/src/vector/_compute/spatial/deltaangle.py
index 09dea3d..d54e6a0 100644
--- a/src/vector/_compute/spatial/deltaangle.py
+++ b/src/vector/_compute/spatial/deltaangle.py
@@ -34,32 +34,61 @@ from vector._methods import (
 def xy_z_xy_z(lib, x1, y1, z1, x2, y2, z2):
     v1m = mag.xy_z(lib, x1, y1, z1)
     v2m = mag.xy_z(lib, x2, y2, z2)
-    return lib.arccos(dot.xy_z_xy_z(lib, x1, y1, z1, x2, y2, z2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1, lib.minimum(1, dot.xy_z_xy_z(lib, x1, y1, z1, x2, y2, z2) / v1m / v2m)
+        )
+    )
 
 
 def xy_z_xy_theta(lib, x1, y1, z1, x2, y2, theta2):
     v1m = mag.xy_z(lib, x1, y1, z1)
     v2m = mag.xy_theta(lib, x2, y2, theta2)
-    return lib.arccos(dot.xy_z_xy_theta(lib, x1, y1, z1, x2, y2, theta2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_z_xy_theta(lib, x1, y1, z1, x2, y2, theta2) / v1m / v2m
+            ),
+        )
+    )
 
 
 def xy_z_xy_eta(lib, x1, y1, z1, x2, y2, eta2):
     v1m = mag.xy_z(lib, x1, y1, z1)
     v2m = mag.xy_eta(lib, x2, y2, eta2)
-    return lib.arccos(dot.xy_z_xy_eta(lib, x1, y1, z1, x2, y2, eta2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1,
+            lib.minimum(1, dot.xy_z_xy_eta(lib, x1, y1, z1, x2, y2, eta2) / v1m / v2m),
+        )
+    )
 
 
 def xy_z_rhophi_z(lib, x1, y1, z1, rho2, phi2, z2):
     v1m = mag.xy_z(lib, x1, y1, z1)
     v2m = mag.rhophi_z(lib, rho2, phi2, z2)
-    return lib.arccos(dot.xy_z_rhophi_z(lib, x1, y1, z1, rho2, phi2, z2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_z_rhophi_z(lib, x1, y1, z1, rho2, phi2, z2) / v1m / v2m
+            ),
+        )
+    )
 
 
 def xy_z_rhophi_theta(lib, x1, y1, z1, rho2, phi2, theta2):
     v1m = mag.xy_z(lib, x1, y1, z1)
     v2m = mag.rhophi_theta(lib, rho2, phi2, theta2)
     return lib.arccos(
-        dot.xy_z_rhophi_theta(lib, x1, y1, z1, rho2, phi2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.xy_z_rhophi_theta(lib, x1, y1, z1, rho2, phi2, theta2) / v1m / v2m,
+            ),
+        )
     )
 
 
@@ -67,21 +96,39 @@ def xy_z_rhophi_eta(lib, x1, y1, z1, rho2, phi2, eta2):
     v1m = mag.xy_z(lib, x1, y1, z1)
     v2m = mag.rhophi_eta(lib, rho2, phi2, eta2)
     return lib.arccos(
-        dot.xy_z_rhophi_eta(lib, x1, y1, z1, rho2, phi2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_z_rhophi_eta(lib, x1, y1, z1, rho2, phi2, eta2) / v1m / v2m
+            ),
+        )
     )
 
 
 def xy_theta_xy_z(lib, x1, y1, theta1, x2, y2, z2):
     v1m = mag.xy_theta(lib, x1, y1, theta1)
     v2m = mag.xy_z(lib, x2, y2, z2)
-    return lib.arccos(dot.xy_theta_xy_z(lib, x1, y1, theta1, x2, y2, z2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_theta_xy_z(lib, x1, y1, theta1, x2, y2, z2) / v1m / v2m
+            ),
+        )
+    )
 
 
 def xy_theta_xy_theta(lib, x1, y1, theta1, x2, y2, theta2):
     v1m = mag.xy_theta(lib, x1, y1, theta1)
     v2m = mag.xy_theta(lib, x2, y2, theta2)
     return lib.arccos(
-        dot.xy_theta_xy_theta(lib, x1, y1, theta1, x2, y2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.xy_theta_xy_theta(lib, x1, y1, theta1, x2, y2, theta2) / v1m / v2m,
+            ),
+        )
     )
 
 
@@ -89,7 +136,12 @@ def xy_theta_xy_eta(lib, x1, y1, theta1, x2, y2, eta2):
     v1m = mag.xy_theta(lib, x1, y1, theta1)
     v2m = mag.xy_eta(lib, x2, y2, eta2)
     return lib.arccos(
-        dot.xy_theta_xy_eta(lib, x1, y1, theta1, x2, y2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_theta_xy_eta(lib, x1, y1, theta1, x2, y2, eta2) / v1m / v2m
+            ),
+        )
     )
 
 
@@ -97,7 +149,13 @@ def xy_theta_rhophi_z(lib, x1, y1, theta1, rho2, phi2, z2):
     v1m = mag.xy_theta(lib, x1, y1, theta1)
     v2m = mag.rhophi_z(lib, rho2, phi2, z2)
     return lib.arccos(
-        dot.xy_theta_rhophi_z(lib, x1, y1, theta1, rho2, phi2, z2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.xy_theta_rhophi_z(lib, x1, y1, theta1, rho2, phi2, z2) / v1m / v2m,
+            ),
+        )
     )
 
 
@@ -105,7 +163,15 @@ def xy_theta_rhophi_theta(lib, x1, y1, theta1, rho2, phi2, theta2):
     v1m = mag.xy_theta(lib, x1, y1, theta1)
     v2m = mag.rhophi_theta(lib, rho2, phi2, theta2)
     return lib.arccos(
-        dot.xy_theta_rhophi_theta(lib, x1, y1, theta1, rho2, phi2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.xy_theta_rhophi_theta(lib, x1, y1, theta1, rho2, phi2, theta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -113,35 +179,65 @@ def xy_theta_rhophi_eta(lib, x1, y1, theta1, rho2, phi2, eta2):
     v1m = mag.xy_theta(lib, x1, y1, theta1)
     v2m = mag.rhophi_eta(lib, rho2, phi2, eta2)
     return lib.arccos(
-        dot.xy_theta_rhophi_eta(lib, x1, y1, theta1, rho2, phi2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.xy_theta_rhophi_eta(lib, x1, y1, theta1, rho2, phi2, eta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
 def xy_eta_xy_z(lib, x1, y1, eta1, x2, y2, z2):
     v1m = mag.xy_eta(lib, x1, y1, eta1)
     v2m = mag.xy_z(lib, x2, y2, z2)
-    return lib.arccos(dot.xy_eta_xy_z(lib, x1, y1, eta1, x2, y2, z2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1,
+            lib.minimum(1, dot.xy_eta_xy_z(lib, x1, y1, eta1, x2, y2, z2) / v1m / v2m),
+        )
+    )
 
 
 def xy_eta_xy_theta(lib, x1, y1, eta1, x2, y2, theta2):
     v1m = mag.xy_eta(lib, x1, y1, eta1)
     v2m = mag.xy_theta(lib, x2, y2, theta2)
     return lib.arccos(
-        dot.xy_eta_xy_theta(lib, x1, y1, eta1, x2, y2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_eta_xy_theta(lib, x1, y1, eta1, x2, y2, theta2) / v1m / v2m
+            ),
+        )
     )
 
 
 def xy_eta_xy_eta(lib, x1, y1, eta1, x2, y2, eta2):
     v1m = mag.xy_eta(lib, x1, y1, eta1)
     v2m = mag.xy_eta(lib, x2, y2, eta2)
-    return lib.arccos(dot.xy_eta_xy_eta(lib, x1, y1, eta1, x2, y2, eta2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_eta_xy_eta(lib, x1, y1, eta1, x2, y2, eta2) / v1m / v2m
+            ),
+        )
+    )
 
 
 def xy_eta_rhophi_z(lib, x1, y1, eta1, rho2, phi2, z2):
     v1m = mag.xy_eta(lib, x1, y1, eta1)
     v2m = mag.rhophi_z(lib, rho2, phi2, z2)
     return lib.arccos(
-        dot.xy_eta_rhophi_z(lib, x1, y1, eta1, rho2, phi2, z2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.xy_eta_rhophi_z(lib, x1, y1, eta1, rho2, phi2, z2) / v1m / v2m
+            ),
+        )
     )
 
 
@@ -149,7 +245,15 @@ def xy_eta_rhophi_theta(lib, x1, y1, eta1, rho2, phi2, theta2):
     v1m = mag.xy_eta(lib, x1, y1, eta1)
     v2m = mag.rhophi_theta(lib, rho2, phi2, theta2)
     return lib.arccos(
-        dot.xy_eta_rhophi_theta(lib, x1, y1, eta1, rho2, phi2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.xy_eta_rhophi_theta(lib, x1, y1, eta1, rho2, phi2, theta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -157,21 +261,40 @@ def xy_eta_rhophi_eta(lib, x1, y1, eta1, rho2, phi2, eta2):
     v1m = mag.xy_eta(lib, x1, y1, eta1)
     v2m = mag.rhophi_eta(lib, rho2, phi2, eta2)
     return lib.arccos(
-        dot.xy_eta_rhophi_eta(lib, x1, y1, eta1, rho2, phi2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.xy_eta_rhophi_eta(lib, x1, y1, eta1, rho2, phi2, eta2) / v1m / v2m,
+            ),
+        )
     )
 
 
 def rhophi_z_xy_z(lib, rho1, phi1, z1, x2, y2, z2):
     v1m = mag.rhophi_z(lib, rho1, phi1, z1)
     v2m = mag.xy_z(lib, x2, y2, z2)
-    return lib.arccos(dot.rhophi_z_xy_z(lib, rho1, phi1, z1, x2, y2, z2) / v1m / v2m)
+    return lib.arccos(
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.rhophi_z_xy_z(lib, rho1, phi1, z1, x2, y2, z2) / v1m / v2m
+            ),
+        )
+    )
 
 
 def rhophi_z_xy_theta(lib, rho1, phi1, z1, x2, y2, theta2):
     v1m = mag.rhophi_z(lib, rho1, phi1, z1)
     v2m = mag.xy_theta(lib, x2, y2, theta2)
     return lib.arccos(
-        dot.rhophi_z_xy_theta(lib, rho1, phi1, z1, x2, y2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_z_xy_theta(lib, rho1, phi1, z1, x2, y2, theta2) / v1m / v2m,
+            ),
+        )
     )
 
 
@@ -179,7 +302,12 @@ def rhophi_z_xy_eta(lib, rho1, phi1, z1, x2, y2, eta2):
     v1m = mag.rhophi_z(lib, rho1, phi1, z1)
     v2m = mag.xy_eta(lib, x2, y2, eta2)
     return lib.arccos(
-        dot.rhophi_z_xy_eta(lib, rho1, phi1, z1, x2, y2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.rhophi_z_xy_eta(lib, rho1, phi1, z1, x2, y2, eta2) / v1m / v2m
+            ),
+        )
     )
 
 
@@ -187,7 +315,13 @@ def rhophi_z_rhophi_z(lib, rho1, phi1, z1, rho2, phi2, z2):
     v1m = mag.rhophi_z(lib, rho1, phi1, z1)
     v2m = mag.rhophi_z(lib, rho2, phi2, z2)
     return lib.arccos(
-        dot.rhophi_z_rhophi_z(lib, rho1, phi1, z1, rho2, phi2, z2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_z_rhophi_z(lib, rho1, phi1, z1, rho2, phi2, z2) / v1m / v2m,
+            ),
+        )
     )
 
 
@@ -195,7 +329,15 @@ def rhophi_z_rhophi_theta(lib, rho1, phi1, z1, rho2, phi2, theta2):
     v1m = mag.rhophi_z(lib, rho1, phi1, z1)
     v2m = mag.rhophi_theta(lib, rho2, phi2, theta2)
     return lib.arccos(
-        dot.rhophi_z_rhophi_theta(lib, rho1, phi1, z1, rho2, phi2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_z_rhophi_theta(lib, rho1, phi1, z1, rho2, phi2, theta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -203,7 +345,15 @@ def rhophi_z_rhophi_eta(lib, rho1, phi1, z1, rho2, phi2, eta2):
     v1m = mag.rhophi_z(lib, rho1, phi1, z1)
     v2m = mag.rhophi_eta(lib, rho2, phi2, eta2)
     return lib.arccos(
-        dot.rhophi_z_rhophi_eta(lib, rho1, phi1, z1, rho2, phi2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_z_rhophi_eta(lib, rho1, phi1, z1, rho2, phi2, eta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -211,7 +361,13 @@ def rhophi_theta_xy_z(lib, rho1, phi1, theta1, x2, y2, z2):
     v1m = mag.rhophi_theta(lib, rho1, phi1, theta1)
     v2m = mag.xy_z(lib, x2, y2, z2)
     return lib.arccos(
-        dot.rhophi_theta_xy_z(lib, rho1, phi1, theta1, x2, y2, z2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_theta_xy_z(lib, rho1, phi1, theta1, x2, y2, z2) / v1m / v2m,
+            ),
+        )
     )
 
 
@@ -219,7 +375,15 @@ def rhophi_theta_xy_theta(lib, rho1, phi1, theta1, x2, y2, theta2):
     v1m = mag.rhophi_theta(lib, rho1, phi1, theta1)
     v2m = mag.xy_theta(lib, x2, y2, theta2)
     return lib.arccos(
-        dot.rhophi_theta_xy_theta(lib, rho1, phi1, theta1, x2, y2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_theta_xy_theta(lib, rho1, phi1, theta1, x2, y2, theta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -227,7 +391,15 @@ def rhophi_theta_xy_eta(lib, rho1, phi1, theta1, x2, y2, eta2):
     v1m = mag.rhophi_theta(lib, rho1, phi1, theta1)
     v2m = mag.xy_eta(lib, x2, y2, eta2)
     return lib.arccos(
-        dot.rhophi_theta_xy_eta(lib, rho1, phi1, theta1, x2, y2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_theta_xy_eta(lib, rho1, phi1, theta1, x2, y2, eta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -235,7 +407,15 @@ def rhophi_theta_rhophi_z(lib, rho1, phi1, theta1, rho2, phi2, z2):
     v1m = mag.rhophi_theta(lib, rho1, phi1, theta1)
     v2m = mag.rhophi_z(lib, rho2, phi2, z2)
     return lib.arccos(
-        dot.rhophi_theta_rhophi_z(lib, rho1, phi1, theta1, rho2, phi2, z2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_theta_rhophi_z(lib, rho1, phi1, theta1, rho2, phi2, z2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -243,9 +423,17 @@ def rhophi_theta_rhophi_theta(lib, rho1, phi1, theta1, rho2, phi2, theta2):
     v1m = mag.rhophi_theta(lib, rho1, phi1, theta1)
     v2m = mag.rhophi_theta(lib, rho2, phi2, theta2)
     return lib.arccos(
-        dot.rhophi_theta_rhophi_theta(lib, rho1, phi1, theta1, rho2, phi2, theta2)
-        / v1m
-        / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_theta_rhophi_theta(
+                    lib, rho1, phi1, theta1, rho2, phi2, theta2
+                )
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -253,9 +441,15 @@ def rhophi_theta_rhophi_eta(lib, rho1, phi1, theta1, rho2, phi2, eta2):
     v1m = mag.rhophi_theta(lib, rho1, phi1, theta1)
     v2m = mag.rhophi_eta(lib, rho2, phi2, eta2)
     return lib.arccos(
-        dot.rhophi_theta_rhophi_eta(lib, rho1, phi1, theta1, rho2, phi2, eta2)
-        / v1m
-        / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_theta_rhophi_eta(lib, rho1, phi1, theta1, rho2, phi2, eta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -263,7 +457,12 @@ def rhophi_eta_xy_z(lib, rho1, phi1, eta1, x2, y2, z2):
     v1m = mag.rhophi_eta(lib, rho1, phi1, eta1)
     v2m = mag.xy_z(lib, x2, y2, z2)
     return lib.arccos(
-        dot.rhophi_eta_xy_z(lib, rho1, phi1, eta1, x2, y2, z2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1, dot.rhophi_eta_xy_z(lib, rho1, phi1, eta1, x2, y2, z2) / v1m / v2m
+            ),
+        )
     )
 
 
@@ -271,7 +470,15 @@ def rhophi_eta_xy_theta(lib, rho1, phi1, eta1, x2, y2, theta2):
     v1m = mag.rhophi_eta(lib, rho1, phi1, eta1)
     v2m = mag.xy_theta(lib, x2, y2, theta2)
     return lib.arccos(
-        dot.rhophi_eta_xy_theta(lib, rho1, phi1, eta1, x2, y2, theta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_eta_xy_theta(lib, rho1, phi1, eta1, x2, y2, theta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -279,7 +486,13 @@ def rhophi_eta_xy_eta(lib, rho1, phi1, eta1, x2, y2, eta2):
     v1m = mag.rhophi_eta(lib, rho1, phi1, eta1)
     v2m = mag.xy_eta(lib, x2, y2, eta2)
     return lib.arccos(
-        dot.rhophi_eta_xy_eta(lib, rho1, phi1, eta1, x2, y2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_eta_xy_eta(lib, rho1, phi1, eta1, x2, y2, eta2) / v1m / v2m,
+            ),
+        )
     )
 
 
@@ -287,7 +500,15 @@ def rhophi_eta_rhophi_z(lib, rho1, phi1, eta1, rho2, phi2, z2):
     v1m = mag.rhophi_eta(lib, rho1, phi1, eta1)
     v2m = mag.rhophi_z(lib, rho2, phi2, z2)
     return lib.arccos(
-        dot.rhophi_eta_rhophi_z(lib, rho1, phi1, eta1, rho2, phi2, z2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_eta_rhophi_z(lib, rho1, phi1, eta1, rho2, phi2, z2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -295,9 +516,15 @@ def rhophi_eta_rhophi_theta(lib, rho1, phi1, eta1, rho2, phi2, theta2):
     v1m = mag.rhophi_eta(lib, rho1, phi1, eta1)
     v2m = mag.rhophi_theta(lib, rho2, phi2, theta2)
     return lib.arccos(
-        dot.rhophi_eta_rhophi_theta(lib, rho1, phi1, eta1, rho2, phi2, theta2)
-        / v1m
-        / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_eta_rhophi_theta(lib, rho1, phi1, eta1, rho2, phi2, theta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
@@ -305,7 +532,15 @@ def rhophi_eta_rhophi_eta(lib, rho1, phi1, eta1, rho2, phi2, eta2):
     v1m = mag.rhophi_eta(lib, rho1, phi1, eta1)
     v2m = mag.rhophi_eta(lib, rho2, phi2, eta2)
     return lib.arccos(
-        dot.rhophi_eta_rhophi_eta(lib, rho1, phi1, eta1, rho2, phi2, eta2) / v1m / v2m
+        lib.maximum(
+            -1,
+            lib.minimum(
+                1,
+                dot.rhophi_eta_rhophi_eta(lib, rho1, phi1, eta1, rho2, phi2, eta2)
+                / v1m
+                / v2m,
+            ),
+        )
     )
 
 
diff --git a/src/vector/_lib.py b/src/vector/_lib.py
index e0c5702..6299d3c 100644
--- a/src/vector/_lib.py
+++ b/src/vector/_lib.py
@@ -14,11 +14,11 @@ class SympyLib:
     def nan_to_num(self, val: sympy.Expr, **kwargs: typing.Any) -> sympy.Expr:
         return val
 
-    def maximum(self, val1: sympy.Expr, val2: sympy.Expr) -> sympy.Expr:
-        if val2 == 0:
-            return val1
-        else:
-            return sympy.maximum(val1, val2)  # type: ignore[no-untyped-call]
+    def maximum(self, val1: sympy.Expr | int, val2: sympy.Expr | int) -> sympy.Expr:
+        return val1 if isinstance(val1, sympy.Expr) else val2  # type: ignore[return-value]
+
+    def minimum(self, val1: sympy.Expr | int, val2: sympy.Expr | int) -> sympy.Expr:
+        return val1 if isinstance(val1, sympy.Expr) else val2  # type: ignore[return-value]
 
     def arcsin(self, val: sympy.Expr) -> sympy.Expr:
         return sympy.asin(val)


Test patch:
diff --git a/tests/test_compute_features.py b/tests/test_compute_features.py
index 78a660e..e518577 100644
--- a/tests/test_compute_features.py
+++ b/tests/test_compute_features.py
@@ -359,6 +359,7 @@ allowed_lib_functions = [
     "sign",
     "copysign",
     "maximum",
+    "minimum",
     "sqrt",
     "exp",
     "log",
diff --git a/tests/test_issues.py b/tests/test_issues.py
index c1d509e..88c36d7 100644
--- a/tests/test_issues.py
+++ b/tests/test_issues.py
@@ -61,3 +61,10 @@ def test_issue_443():
         {"E": [1], "px": [1], "py": [1], "pz": [1]}, with_name="Momentum4D"
     ) ** 2 == ak.Array([-2])
     assert vector.obj(E=1, px=1, py=1, pz=1) ** 2 == -2
+
+
+def test_issue_463():
+    v = vector.obj(x=1, y=1, z=1)
+    for transform in "xyz", "xytheta", "xyeta", "rhophiz", "rhophitheta", "rhophieta":
+        trv = getattr(v, "to_" + transform)()
+        assert trv.deltaangle(trv) == 0.0