New distance queries

[radevgit]

Hi,
Do you have any plans to implement distance queries for other geometric objects like: SEG-ARC, ARC-ARC ?

[davideberly]

The distance between 2 segments in any dimension is implemented in DistSegmentSegment.h.

I am currently implementing distance algorithms between line/ray/segment and circle/arc. The line/ray/segment-circle code is finished, including unit testing. I will post this shortly. I need to implement line/ray/segment-arc. This should take me another 1 or 2 days.

[owai1980]

Waw! This is a very good news! Thank you! 👍 Le lun. 31 juil. 2023 à 15:34, David Eberly ***@***.***> a écrit :

…

The distance between 2 segments in any dimension is implemented in DistSegmentSegment.h. I am currently implementing distance algorithms between line/ray/segment and circle/arc. The line/ray/segment-circle code is finished, including unit testing. I will post this shortly. I need to implement line/ray/segment-arc. This should take me another 1 or 2 days. — Reply to this email directly, view it on GitHub <#66 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AG3BAQN3FCBO7KDKTTJLVA3XS6X5LANCNFSM6AAAAAA247QVHU> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[radevgit]

Thanks David,
I checked before the DistSegmentSegment.h, and I was surprised how complex is the algorithm.
My naive thinking was, that is should be simple than intersection algorithm.
Thanks again for the fast response.

[davideberly]

The document Robust Computation of Distance Between Line Segmants has a couple of examples that show the straightforward algorithm, when computed using floating-point arithmetic, can have gross errors because of subtractive cancellation. In particular, see Example 2. If you use the straightforward algorithm computing with rational arithmetic, all is well.

The file I mentioned has a function ComputeRobust(...), where I attempt to avoid the problem with subtractive cancellation. The idea is reasonable in the floating-point era when you do not have a fused-multiply-add instruction. That instruction computes fma(x,y,z) = x*y+z with 1 rounding step. Without such an instruction, w = x*y has 1 rounding step and w+z has 1 rounding step for a total of 2 rounding steps, which generally leads to more rounding error.

With the fma instruction, you can "fix" the subtractive cancellation in expressions x*y-z*w. The implementation is in my Math.h file, functions RobustDOP. The idea is due to Professor William Kahan On the Cost of Floating-Point Computation
Without Extra-Precise Arithmetic and succinctly mentioned by Matt Pharr Accurate Differences of Products with Kahan’s Algorithm. At some time I will revisit the segment-segment distance code and modify it to use this paradigm.

[radevgit]

Thanks for the references. I learned a lot. Regarding the fused-multiply-add instruction, it should be available std::fma. Also, regarding such numerical issues, I noted somewhere in the code the use of dot(v0, v1) and sqrt for the length of a vector. There are claims that hypot uses robust algorithm for the same. It looks from some comments, that the hypot also use fma operations.

[davideberly]

The IEEE Standard for Floating-Point Arithmetic 754-2019 includes the requirement that implementations provide fusedMultiplyAdd. For hypotenuse in 2D (squared length of vectors), you can compute xx+yy using the same idea I mentioned. Implementation is RobustSOP (Robust Sum Of Products) in Math.h. The replacement of DOP and SOP will take some time in my code.

[mgamito]

If I may add something to this discussion, I found at work that although std:fma is always supported by the standard, it may not necessarily map to a machine instruction on older architectures and needs to be emulated by software. This turns out to be woefully slow and definitely slower than just writing y = a*x + b. I found this while working with NanoVdb that uses fused-multiply-adds and we had to create a local patched version that removes those functions, so that the older machines on our render farm could stay efficient.

…

On Tue, Aug 1, 2023 at 4:34 PM David Eberly ***@***.***> wrote: The IEEE Standard for Floating-Point Arithmetic 754-2019 includes the requirement that implementations provide fusedMultiplyAdd. For hypotenuse in 2D (squared length of vectors), you can compute x*x+y*y using the same idea I mentioned. Implementation is RobustSOP (Robust Sum Of Products) in Math.h. The replacement of DOP and SOP will take some time in my code. — Reply to this email directly, view it on GitHub <#66 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AGCZNWNWH5INAS4W7HZ2O4LXTEOXVANCNFSM6AAAAAA247QVHU> . You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>

[davideberly]

I appreciate your input, thank you.

Software implementations of floating-point functions have that curse. To follow the standard, the logic can be quite complicated and expensive to compute without hardware support. My focus lately has been on robust computations, much of this using rational arithmetic (which effectively is floating-point implemented in software), so I suppose much of my code is "woefully slow" when using rational support.

I added to my short-term TODO list the ability to "discard" FMA operations using a preprocessor define GTE_DISCARD_FMA. If not defined, std::fma(x,y,z) is used. If defined, x*y+z is computed instead with 2 floating-point operations. Longer term I could add a FMA software implementation using the ideas of BSNumber<UIntegerFP32> for N sufficiently large and then add a preprocessor define GTE_USE_SOFTWARE_FMA. But of course this will be slow.

Given you are using modern NVIDIA technology and a render farm, I would have thought a primary goal is to replace old hardware with new hardware. It seems strange that the GPU support of NanoVDB accesses GPU hardware (that has fma), but the underlying CPU is too old for fma. However, if you have a large number of old computers, I can see why hardware upgrade is not desireable.

[radevgit]

For hypotenuse in 2D (squared length of vectors), you can compute x_x+y_y using the same idea I mentioned. Implementation is RobustSOP (Robust Sum Of Products) in Math.h. The replacement of DOP and SOP will take some time in my code.

Here, I found the actual implementation in glibc, that seems to use some additional tricks:
hypot

[davideberly]

I have posted 2D distance queries for line-circle, ray-circle, segment-circle, line-arc, ray-arc, segment-arc. I have not yet posted arc-arc queries (need to think about this one a bit).

[davideberly]

I pressed the wrong button, did not mean to close the issue. Reopening it. I am currently working on the 2D distance queries for circle-circle, circle-arc, and arc-arc.