Repository navigation
Fix _small_angle_of_2_vectors to evaluate Kahan's formula accurately - #1863
Merged
Merged
Conversation
hongyuchen1030
added
improvement
_small_angle_of_2_vectors implements Kahan's angle formula for vectors of arbitrary length, 2*atan2(|v|u| - u|v||, |v|u| + u|v||), but formed each component of v|u| - u|v| from two separately rounded products. For nearly parallel vectors this is a naive difference of products, and the relative error grew as eps/angle (a mean of about 2e-8 at 1e-9 rad). Each component of v|u| - u|v| is now evaluated with diff_of_products. The norms and the sum are still computed directly: compensating the norms only removes a second-order (eps/angle)^2 term below ~1e-7 rad, and the sum's cancellation near antiparallel vectors does not affect the returned angle. The mean relative error is now below 1e-16 down to 1e-7 rad and 3.3e-15 at 1e-9 rad, at a mean of 48.5 ns per call (41.9 ns before). Add tests against an exact rational reference; they fail on the previous implementation.
hongyuchen1030
force-pushed
the
fix/small-angle-kahan-eft
branch
from
October 30, 2026 04:24
0687705 to
470472b
Compare
hongyuchen1030
requested review from
Sevans711,
cmdupuis3 and
rajeeja
and removed request for
Sevans711 and
cmdupuis3
October 30, 2026 04:25
hongyuchen1030
marked this pull request as ready for review
October 30, 2026 04:25
ASV BenchmarkingBenchmark Comparison ResultsBenchmarks that have stayed the same:
Benchmarks that have got worse:
|
Sevans711
requested changes
Oct 30, 2026
Sevans711
left a comment
Sevans711
left a comment
Collaborator
There was a problem hiding this comment.
This looks good overall, and it will be great to have more accurate angle computations! I confirmed locally that tests fail on main but pass here.
Would it make more sense to put these new tests into test/grid/geometry/test_angles.py instead? I do see that there are actually already a couple tests above the ones you added in test_helpers.py which test _angle_of_2_vectors, though. I would slightly prefer to move all of them into test_angles.py, although this would be a slight expansion of scope for the PR and it is just a non-blocking stylistic preference.
github-project-automation
Bot
moved this from 🏗 In progress
to 👀 In review
in UXarray Development
- Add _numba_diff_of_products3(u, a, v, b) to utils/numba_math.py: the component-wise u*a - v*b of 3-vectors via diff_of_products, returning (hi, lo) tuples. - Write _small_angle_of_2_vectors with vector helpers instead of spelling out components; this also drops the float() conversions, which were only there to unpack components. - Move the _angle_of_2_vectors and _small_angle_of_2_vectors tests from test_helpers.py to test/grid/geometry/test_angles.py, and test the new helper against exact rational arithmetic. Accuracy and speed are unchanged.
Sevans711
self-requested a review
October 30, 2026 17:42
Sevans711
approved these changes
Oct 30, 2026
Sevans711
left a comment
Sevans711
left a comment
Collaborator
There was a problem hiding this comment.
All my previous comments have been addressed; everything looks good to me now, happy to approve!
rajeeja
approved these changes
Oct 30, 2026
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Closes #1864
Overview
_small_angle_of_2_vectorsimplements Kahan's angle formula for vectors of arbitrary length,but it formed each component of
v|u| - u|v|from two separately rounded products. For nearly parallel vectors this is a naive difference of two nearly equal products, so the relative error grew as eps/angle (a mean of about 2e-8 at 1e-9 rad). This was found while reviewing #1806.This PR evaluates each component of
v|u| - u|v|with the FMA-baseddiff_of_products(Kahan's 2x2 determinant) fromuxarray/utils/computing.py, through a new vector helper_numba_diff_of_products3inuxarray/utils/numba_math.py. The norms and the sumv|u| + u|v|are still evaluated directly:The function keeps supporting inputs of arbitrary length, which its callers rely on (tangent-plane vectors in
grid/angles.py, chord vectors ingrid/point_in_face.py).Accuracy benchmark
Mean relative error over 200 random pairs per angle:
Speed (mean time per call over 10 runs of 2,000,000 calls inside an
@njitloop, after a warm-up run): main 41.9 ns, this PR 48.5 ns.Setup
u ~ N(0, I_3);v = cos(angle) * u/|u| + sin(angle) * w/|w|withw = u x r,r ~ N(0, I_3); thenuandvare scaled by independent factors fromU(0.1, 10), so the two vectors have different, non-power-of-two lengths. RNG:numpy.random.default_rng(0), the same pairs for both implementations.fractions.Fraction, thenatan(|u x v| / (u . v))is evaluated with a 50-digitdecimalTaylor series (_exact_small_angleintest/grid/geometry/test_angles.py).Benchmark script
Tests
test/grid/geometry/test_angles.pynow holds the angle-function tests: the two existing_angle_of_2_vectorstests (moved fromtest/test_helpers.py) and new tests for_small_angle_of_2_vectors, which previously had none: special values (orthogonal, parallel, antiparallel, zero vector) and nearly parallel vectors of different lengths at 1e-3 to 1e-9 rad against the exact reference above. The accuracy tests fail on main and pass with this change.test/utils/test_numba_math.pyadds a test for_numba_diff_of_products3against exact rational arithmetic.PR Checklist
General
Testing & Benchmarking
Documentation and Examples
docs/api.rst; internal (private) function names start with an underscore (_)AI Disclosure
AI Usage: Claude Opus 5.5,