I'm curious if experienced users can comment on CGAL's numerics approach (https://www.cgal.org/exact.html) where CGAL states that they track error bounds explicitly and fall back to high-precision exact arithmetic.
In numerical linear algebra, problems have an intrinsic condition number, and most algorithms admit backward or forward stability in terms of that condition number, implying a kind of Lipschitz continuity in computation: you can't be too far off from the true exact solution, and if your input was rounded, then it's still close enough to the exact solution.
The documentation seems to imply that CG problems don't have this kind of regularity ("no such arsenal"), hence the need for fallback to exact.
But if this regularity doesn't exist in CG (i.e., you can create natural geometric constructs for which low-bit perturbations result in large-scale changes in the output), then why would exact numerics even save us here? Wouldn't it be possible to construct examples where we could make the amount of precision required super-polynomial in the input bitwidth? Or are there different regularity theorems which limit the complexity of the output (but still require, say, polynomially-increased precision).
The numerics in computational geometry are used for making combinatorial decisions (do three points make a left or a right turn in the plane; is point p inside or outside of the circumsphere defined by these three points; etc). These are called predicates in the CG literature. What makes this difficult is that multiple different predicates interact (e.g., take the three points in different order), and the answers they give need to be consistent. If you want some nice examples of how very simple things can go wrong, see L. Kettner, K. Mehlhorn, S. Pion, S. Schirra, and C. Yap, “Classroom examples of robustness problems in geometric computations,” Comput. Geom., vol. 40, no. 1, pp. 61–78, 2008.
Most predicates boil down to a decision of whether some quantity is less than, equal, or greater than zero. CGAL implements filtered predicates — a work of art in my opinion — where they use ordinary computation and if the result is far enough away from zero, return its sign. If not, they switch to higher precision or interval arithmetic. A good explanation of how this works (what "far enough from zero" actually means) is in O. Devillers and S. Pion, “Efficient Exact Geometric Predicates for Delaunay Triangulations,” in Proceedings of the 5th Workshop on Algorithm Engineering and Experiments, pp. 37–44, 2003.
The larger problem of symbolically perturbing input, so that it's in general position (in a sense, what to do if your predicate returns 0) was a very active research topic in the 90s, and both the problem and much of the work is explained really well in R. Seidel, “The Nature and Meaning of Perturbations in Geometric Computing,” Discrete Comput. Geom., vol. 19, no. 1, pp. 1–17, 1998.
Thank you for the detailed reply. "Classroom examples of robustness problems in geometric computations" is great; I look forward to doing a deep read there.
Your reference, "The Nature and Meaning of Perturbations in Geometric Computing" itself refers to "On degeneracy in geometric computations" (https://dl.acm.org/doi/10.5555/314464.314474) which I think points to a less tame reality where the answer to my question (continuity claims about CG outputs) is "it depends on the CG problem".
I think this is an interesting question, and I'd be curious to know the answer, too.
That said, I might push back on the idea that CG problems don't have this kind of regularity. They give the example of a point being either on a line determined by two other points, or to one or the other side of it. They say that you need to know this information exactly, which is combinatorial.
It's not so clear to me that you really do need to know this information. For this example, you could represent the line as a zero level set, and then the query above is answered by the sign of the value of the zero level set function evaluated at this point. There is some latent numerical error here, but this problem is well-conditioned in the sense you mention above.
I think it would be interesting to see whether it's possible to build a library of CG routines based on a "fuzzier" level set conception of things.
I'm not the best person to give counterexamples here, but a classical demo of "lack of regularity" would be the point-in-polygon problem (https://en.wikipedia.org/wiki/Point_in_polygon).
Say you have a spiral-shaped but simple polygon that's wrapped around a query point. Even if your query point is far from the innermost "loop" of the spiral, a raycast could be arbitrarily confused if the sprial has many "thin" loops on the outside. (However, notably, the bits precision required would be a reasonable function of the spiral's specification. My point is just that fuzziness in general can be a compounding problem as a CG algorithm "progresses".)
Not sure where OP learned about this, but I learned about it when using OpenSCAD to programmatically generate designs to 3D print. It's an easier way for a programmer to generate 3D designs than learning something like Blender.
There are also "official" ones there: https://github.com/CGAL/cgal-swig-bindings. It does not cover all components, but is performed on a "when there are enough requests for it" basis.
Several years ago I used FEniCS mshr which wraps CGAL for Python, but I am unsure if is still maintained or if they went for a Gmsh wrapper instead, or abandoned the idea of constructive solid geometry from Python code completely. I have not followed FEniCS closely in the last few years.
It's has an abundantly iterator-based API, which is common in C++ but in my experience CGAL goes overboard with it which makes it hard to learn. I also didn't think it's well-documented when I was working with it between 2016-2018 as a computer graphics researcher. Maybe things improved. Great library overall, but definitely very fragmented and requires fair amount of reading to understand the API.
CGAL is a seminal library, but I agree with the parent's criticism of "fragmented" - for instance, the various sub-libraries have different licenses (GPL, LGPL, others), which may lead to accidental mis-use.
When I needed a few geometric functions back ~2005, I ended up implementing them myself to avoid a huge dependency on CGAL (the library back then was already bigger than my application) and the licensing implications.
One of my proudest professional moments was getting someone with an Image Processing background to hire one of my friends with a Computational Geometry background to work on our problem domain. He was extremely reluctant to hire him, and my friend was a resounding success.
I’ve never had to use CGAL directly, but indirectly it played a role in my PhD thesis, as a backend for a clustering problem called jet finding in high energy nuclear physics :) always get a little ping of happiness when I see these tools pop up.
Thanks. However, while SFCGAL itself is LGPL, it uses parts of CGAL, which are GPLv3+, and because of this, it - and its Rust wrapper crate sfcgal (https://mthh.github.io/sfcgal-rs/sfcgal/), are all effectively covered by GPLv3+.
Therefore, I hope someone re-implements the innards under a liberal (LGPL or Apache2 or BSD or MIT) license.
In numerical linear algebra, problems have an intrinsic condition number, and most algorithms admit backward or forward stability in terms of that condition number, implying a kind of Lipschitz continuity in computation: you can't be too far off from the true exact solution, and if your input was rounded, then it's still close enough to the exact solution.
The documentation seems to imply that CG problems don't have this kind of regularity ("no such arsenal"), hence the need for fallback to exact.
But if this regularity doesn't exist in CG (i.e., you can create natural geometric constructs for which low-bit perturbations result in large-scale changes in the output), then why would exact numerics even save us here? Wouldn't it be possible to construct examples where we could make the amount of precision required super-polynomial in the input bitwidth? Or are there different regularity theorems which limit the complexity of the output (but still require, say, polynomially-increased precision).