Applying BGL to Computational Geometry

Introduction

The Boost libraries [1] extend the C++ Standard library by providing new functionality (e.g., threads [2] and smart pointers [3]) and establishing an extensive set of language idioms. Not surprisingly, Boost will have an impact on the emerging C++ Library Extensions Technical Report (see Herb Sutter’s discussion in [4]). Among Boost’s several dozen libraries, the BGL (Boost Graph Library) stands as a prime example of modern C++ design applied to advanced algorithmic concepts. You do not need to be doing scientific computing for a living to appreciate graph-theoretic concepts. For instance, your build system analyzes the directed (and, hopefully, acyclic!) graph of file dependencies to control compilation. Many network routing problems are solved in the graph framework.

In this article, I apply BGL to the domain of computational geometry (see sidebar). First, I formulate a concrete problem in graph terms. Second, I develop a way to transform the output of an existing algorithm into an appropriate BGL data structure. Finally, I implement two new algorithms for my BGL graph. The first algorithm gets the job done, but could have been written in any programming language. The second algorithm, however, shows the power of BGL’s generic programming approach.

A Brief Overview of BGL

BGL originated at the Laboratory for Scientific Computing at the University of Notre Dame. The initiative started around 1998 and was initially termed the GGCL (Generic Graph Component Library). GGCL joined forces with Boost two years later and became BGL at the end of 2000. As of this writing, the latest version of Boost is version 1.27.0. You can master BGL through online documentation [5] or a recently published book [6].

BGL is not about data structures and algorithms, although it contains both. It is about writing interfaces in an efficient and extensible way. In doing so, BGL tries to correct several key problems that plague scientific software. First, the majority of such software does not separate algorithms from underlying data types, guaranteeing that the next person to come along will have to modify both to solve their particular problem. Second, and this applies particularly to problems in graph domains, an application often already has its own data representation and converting that data structure into some graph library format may severely degrade performance. For example, in the domain of geographic information systems, one might want to compute object mobility on a particular kind of terrain. This can easily be done using Dijkstra’s single-source shortest-paths algorithm. The problem is that DEMs (digital elevation maps) can take up megabytes (and even gigabytes) of storage, and you would not want to copy such a behemoth into a new graph data structure (represented as, say, a linked list) just to run one algorithm. Instead, you should superimpose the graph structure on the existing data in an implicit fashion by providing the right interface and then run your graph algorithm against that interface.

BGL achieves these and several other goals through consistent application of generic programming principles influenced by the C++ Standard library. BGL clearly separates containers (graphs) from algorithms (graph-traversal, etc.) by writing the latter against the graph_traits<> interface [7] instead of hard-coded member functions. Although BGL provides several reference graph representations such as adjacency_list and adjacency_matrix, it supplies examples of interfacing to external graphs (such as LEDA graphs) and encourages you to write your own graph adapters to suit your problem domain.

The unifying implementation theme of BGL is template metaprogramming [8]. For example, in Listing 1, I pass a list of visitors (functors) to a graph-traversal algorithm. Do not search for vector<Visitor> in Listing 1 — it is not there. Instead the list of visitors consists of nested std::pairs and is treated as a type known at compile time. Similarly, should you look at the graph-traversal algorithm implementation, do not try to find a switch statement that selects the right visitor at the right time. Again, the decision of which visitor to apply is made at compile time via the compile-time if idiom. I cannot even begin to list the benefits of compile-time checks for static program validation, but I will briefly mention the performance gain. On a lark, I wrote a benchmark to compare element access using a compile-time list data structure coupled with the compile-time if idiom and a dynamic container (std::vector) coupled with an inheritance-based idiom. The compile-time approach outperformed the other by a factor of 20. When I introduced inlining, the Microsoft Visual C++ compiler increased the performance gap to a factor of 100. After that, I stopped doubting BGL creators’ claim that generic programming matches hand-optimized math libraries.

BGL Meets Legacy Code

My favorite book on computational geometry by O’Rourke [9] refers to a very efficient plane sweep algorithm for constructing Voronoi diagrams published by Steve Fortune in 1987. I located the implementation [10], which was in plain C and produced only ASCII output (i.e., not an in-memory graph data structure). Here is how I “boosted” it to C++.

As Fortune’s algorithm sweeps across the plane, it emits events corresponding to newly constructed portions of the Voronoi diagram. Each event is output in ASCII in one of the four formats:

1. s x y indicates an input point at coordinates (x,y).
2. l a b c indicates a line 1 with equation ax + by = c.
3. v x y indicates a vertex at (x,y).
4. e l v1 v2 indicates a Voronoi segment, which is a subsegment of line l above, with endpoints numbered v1 and v2. If v1 or v2 is -1, the line extends to infinity.

I don’t care about output records of type 1 since they represent the input data. I can also deduce that an output record of type 4 has to refer to valid values of v1, v2, and 1, and therefore must be emitted only after the corresponding records of type 2 and 3 have been emitted.

I am going to make one simplification. Since the final result I want is the medial axis, I know I will have to prune unbounded edges. Therefore, I ignore records of type 2 and type 4 when the line extends to infinity.

To define the BGL graph, I chose the adjacency_list representation over adjacency_matrix because the graph will be sparse and will change over time.

typedef
adjacency_list<vecS, vecS,
   undirectedS,VertexProp,
   no_property> MedialAxisGraph;

This defines an undirected graph that has no properties attached to edges, but has data associated with vertices. Furthermore, I use the vecS selector tags to pick std::vector as the underlying container to hold vertices and edges. This choice will simplify some of the code in Listing 2.

The vertex property (VertexProp) ties in the abstract graph with the concrete geometric problem I am trying to solve. You need to follow a three-step process to attach properties to either vertices or edges. First, you need to decide what you want stored at each vertex. This can be an object of either a fundamental (e.g., bool) or a user-defined type (e.g., Point). In this example, I use both. Second, you need to define a property type tag. BGL’s internal machinery needs this tag type to distinguish properties. Finally, you define the property type by constructing a compile-time property list. I show these steps in Listing 3 by attaching a boolean and a Point vertex property. Listing 2 shows how the graph is populated during the run of Fortune’s Voronoi algorithm.

Here is what’s going on. When the record type indicates a Voronoi vertex, the code adds a new graph vertex and then sets its two properties. Notice how type tags allow overload resolution (once again, a compile-time decision). The coord_t property gets the x,y coordinates. The insidecurve_t property gets true as a default value. The case is trickier for handling Voronoi segments. All I want to do is add an edge, but what are the corresponding graph vertices? Because my Voronoi vertices are zero-based, and because I chose std::vector as the underlying storage for vertices, I can type cast the former into the latter.

BGL Graph Processing (A Non-BGL Way)

I skip the details of iterating through the graph and setting the insidecurve_t property to the correct boolean value. Testing the location of a point with respect to a closed polygon is one of the fundamental topics in computational geometry, and you have a choice of many algorithms. The next step, however, is more interesting because it illustrates Boost graph editing at run time. I want to prune the graph, removing vertices whose insidecurve_t property is false. Since the BGL interface is iterator based, I have to make sure that I do not accidentally invalidate my iterators while removing vertices. A quick look through online documentation reveals that I am out of luck: adjacency_list<vecS,vecS,...> invalidates all iterators on vertex removal. Hence the double for-loop in Listing 4. Notice the call to clear_vertex prior to remove_vertex. The code transforms the graph in Figure 3 into the graph in Figure 4.

BGL All the Way

My goal is to compute the shape structure graph shown in Figure 5. Intuitively this entails iterating through the medial axis graph and marking each vertex whose degree is not equal to two. Then for each marked vertex you can search for the nearest marked vertex and add the appropriate edge to the shape structure graph. Such a naive implementation is opaque and has O(N2) complexity.

A better approach takes advantage of two BGL features: generic graph-traversal algorithms and visitors. Generic graph traversals, such as depth_first_search and breadth_first_search parallel std::transform from the C++ Standard library. Since graph traversal involves backtracking and other “interesting” events, BGL lets you pass a list of visitors rather than a single functor to the traversal function.

I chose the depth_first_search for graph traversal and a single visitor, DiscoverVertex, that processes on_discover_vertex events (i.e., the first time the traversal algorithm encounters the vertex). Listing 1 shows key elements of DiscoverVertex. The typedef sets up the event filter so that DiscoverVertex::operator() is called only once per vertex. The depth_first_search algorithm expects a visitor list, which I create using std::make_pair. Since DiscoverVertex builds the shape structure graph along the way, it stores a pointer to the (initially empty) graph in sGraph_. You can see how it all fits together in the source code available from CUJ’s website (<www.cuj.com/code>).

Conclusions

BGL provides an elegant and extensible framework for building graph data structures and algorithms. Whether you are retrofitting an old but reliable algorithm or are developing a brand-new one, the extra effort to “boostify” your code will pay off. However, the beauty of BGL comes at the price of complexity, both to the programmer and the compiler. Unless you dream of partially instantiated templates at night, expect to spend some time getting used to the Boost design philosophy. Likewise, compilers have a rough time handling some of the code’s idioms. I primarily use Microsoft Visual C++ v6.0 and have stared at the “Internal Compiler Error” message many times while working with BGL. (Things are better, but still not perfect, with GCC 3.0.x.) Most of the time (but not always), you can push forward by breaking up complex expressions into multiple statements or rearranging the code in some other way. The end result is certainly worth the investment.

Acknowledgements

This work was supported by the Air Force Research Labs, Eglin AFB, under Contract F08630-99-C-0022. In addition, I want to thank A. J. Fox for helpful suggestions.

References

[1] Boost. <www.boost.org>.

[2] Bill Kempf. “The Boost.Threads Library,” C/C++ Users Journal, May 2002.

[3] Bjorn Karlsson. “Smart Pointers in Boost,” C/C++ Users Journal, April 2002.

[4] Herb Sutter. “The New C++: The Group of Seven — Extensions under Consideration for the C++ Standard Library,” C/C++ Users Journal C++ Experts Forum, April 2002, <www.cuj.com/experts/2004/sutter.htm>.

[5] BGL. <www.boost.org/libs/graph/doc/table_of_contents.html>.

[6] Jeremy G. Siek, Lie-Quan Lee, and Andrew Lumsdaine. The Boost Graph Library: User Guide and Reference Manual (Addison-Wesley, 2002).

[7] Thomas Becker. “STL & Generic Programming: Traits Classes,” C/C++ Users Journal, December 2001.

[8] For more information on template metaprogramming, see Thomas Becker’s column in this issue. (Thomas Becker. “STL & Generic Programming: C++ Template Metaprogramming,” C/C++ Users Journal, August 2002.)

[9] Joseph O’Rourke. Computational Geometry in C (Cambridge University Press, 1998).

[10] Steven Fortune’s Voronoi diagram code. <http://cm.bell-labs.com/who/sjf>.

Vitaly Ablavsky holds a B.A. in Mathematics from Brandeis University and an M.S. in Computer Science from UMass Amherst. He is a senior software engineer at Charles River Analytics in Cambridge, MA. His interests are computer vision and pattern recognition. He can be contacted at [email protected].