Re: [boost] Preview 2 of a template Geometry Library

From: Hervé Brönnimann (hervebronnimann_at_[hidden])
Date: 2008-03-28 09:46:58

Noah, Fernando:  Apart from the simple known expressions (such as  
eliminating a sqrt from the comparison), symbolic manipulation won't  
be able to do that much, although it does have real benefits for  
automatic code generation.  Here is what I mean:

I had that idea early on, and advised a BSc Thesis on the topic  
(Allen Clement, "MGK: A geometric meta kernel", May 1, 2000,  
Princeton Univerisity Dept of CS).  I was trying to see if I could,  
from the simplest computations, compute (complicated) expressions and  
then use a CAS (computer algebra system, in this case Maple, which  
has a function to_C() and functions factor()/simplify() to reduce  
expressions and output them as a C program).

Specifically, I wanted to use bootstrap and derive algebraic  
expressions for complex primitives in terms of a small core of very  
simple ones (e.g., line through two points, line intersection, perp.  
line passing through a point).  I created a CGAL "symbolic" kernel,  
with a number type that created an expression tree.  I would then  
implement and use the simplest constructions, e.g.:
      squared_point_line_distance(x, y, a, b, c) {
         Point p(x,y);
         Line l(a,b,c);  // ax + by + c = 0
         Line m = construct_perp_line_passing_through(p, l);
         Point q = construct_intersection(l, m);
         return point_point_squared_distance(x, y);
     }
Similar geometric constructions for circumcenter(p, q, r) as  
intersection of perp_bisector(p,q) and perp_bisector(p,r).  You can  
generalize this to 3D, polar coordinates, etc. as long as you've  
implemented your small core of simple functions.

line_point_distance as above results in a "naive" expression using 9  
additions, 20 mult., and 2 divisions.  After Maple is done, you get  
back the well known (x * a + y * b + c)^2 / (a^2 + b^2), which is 3  
additions, 5 multiplications, and 1 division.

circumcenter(p,q,r) construction as above gives 18 adds, 16 mults, 1  
div, and after simplification you get the well-known formula with 12  
adds, 12 mults, 1 div.

The approach isn't without benefits and can automate the construction  
of a suite of C++ code for computing common geometric primitives.   
But beyond those who were rather well-known such as above, we didn't  
witness that much improvement. Where we thought we would gain more  
was in lesser studied representations (lat-lon, polar coordinates)  
where getting efficient symbolic expressions demands more computation.

The approach doesn't really pan out for predicates, due to the  
branching (if statements, mostly signs of algebraic expressions)  
which are usually not well handled by a CAS.

It worked up to a point, but after that, not much.  The final code  
(efficient) is also unreadable, but its correctness stems from the  
correctness/simplicity of the original construction and of the CAS  
expression manipulation.  We could compute equivalent expressions and  
reduce the number of primitive operations by some quantity. However,  
there was also no guarantee that the generated expressions would be  
numerically stable.

See Allen's thesis for details, and his conclusion for why this  
system might be interesting and its limitations.
Maybe some specialized CAS engine for simplification tailored to the  
application domain can help.  I played around with Yacas, but with no  
results...
--
Hervé Brönnimann
hervebronnimann_at_[hidden]
[1] Allen Clement, "MGK: A geometric meta kernel", May 1, 2000,  
Princeton University, Dept of CS.  Available at
http://photon.poly.edu/~hbr/publi/allen_bsthesis.pdf
On Mar 28, 2008, at 3:04 AM, Noah Stein wrote:
> Fernando Cacciola wrote:
>
>> 5) There is a reason why CGAL doesn't provide a distance function  
>> at all
> but
>> only a squared_distance, and I agree with that reason: sqrt() kills
>> robustness and distances are usually used for comparisons, so you  
>> seldom
>> need that.
>
> This is one of the things I think proto can help with. I'm playing  
> around
> with a simple vector/point/matrix library just to get a feel for  
> it. proto
> should be able to transform expressions of the sort "length(a) >  
> length(b)"
> into "dot(a,a) > dot(b,b)". I think proper use of proto can lead to  
> high
> levels of expressivity combined with good optimization.
>
>
> -- Noah
>
>
> No virus found in this outgoing message.
> Checked by AVG.
> Version: 7.5.519 / Virus Database: 269.22.1/1346 - Release Date:  
> 3/27/2008
> 10:03 AM
>
>
> _______________________________________________
> Unsubscribe & other changes: http://lists.boost.org/mailman/ 
> listinfo.cgi/boost