Wow, this email ended up a lot longer than I intended. I'm sorry about that. It's fine if you don't have time to read the whole thing. The TLDR is:
I think with adding some minor generalizations of the matrix-multiply operation, Octave can be a very good tool for performing a BFS.
I've written C++ projects before, but never for public consumption. I have some recent experience writing a project in C++ including using templates and C++-11 rvalue references and I've attached one piece of that project.
----------------------------------------------------------------------
Thank you, Jordi
That's a very good point. Certainly Octave is not a good tool for a depth-first search.
Nonetheless, I still think there are enough applications to warrant such a library. Whenever I confront a graph problem, I usually turn to Octave (if I know a poly-time algorithm) or Clingo (if not or if it's known to be NP-hard).
As for a breadth-first search, if we're just concerned with reachability (eg. for an unweighted graph), this can be done with a number of sparse vector-matrix multiplies linear in depth:
lastreachable = false
reachable = g(from,:)
while any(reachable != lastreachable)
reachable = logical(reachable * g)
if reachable(to)
break
endif
endwhile
res = reachable(to)
endfunctionBetter yet, the whole reachability/transitive-closure matrix can be built with only a number of sparse-matrix multiplies logarithmic in the maximum depth
lastreachable = eye(size(g));
reachable = g | lastreachable;
while any(reachable != lastreachable)
reachable = logical(reachable ^ 2);
endwhile
endfunction
I've actually used this quite a lot. I know of no other language where obtaining the transitive closure is this fast and this easy to write.
Finding min-cost paths on (positive) weighted graphs is a little more tricky. It would be nice to do something like a matrix multiply, but replace the dot product with min(sum(nonzero pairs)). It seems to me that very often, breadth-first type algorithms can be re-written in terms of "for each row of x and each column of y, combine them in some way to get r(x,y)" which is like a generalized matrix multiply.
For instance, one that seems to come up a lot (such as in the max-cliques example) is all(x | ~y)
Using this this operation in place of the dot-product when multiplying a logical matrix A by the transpose of another logical matrix B tells which rows in A are subsets of which rows of B.
I don't know how sparse-matrix-multiplies are done internally (it isn't discussed in the PDF and I still can't really navigate the Octave source very well), but I have a pretty good guess that finishes handling each column with a kind of funky half-merge-sort (because it's already partially ordered). If it's anything like I imagine, then these nonstandard operations shouldn't be too hard to substitute.
I have a decent knowledge of C++. I only recently re-familiarized myself with the language while writing a preprocessor (for grounded Answer Set Programs) that does something analogous to forward-checking to infer binary constraints. (It turns out that for competition problems, the number of binary constraints found by the technique is generally too large to hold in memory so I ultimately abandoned the project. I have another potential usage in mind regarding General Game Playing, but I haven't yet had time to work on that. It would be a very big project). While working on this project, I read "Effective C++" and some of "C++ Templates: The Complete Guide" and I've done extensive online browsing to understand C++-11 rvalue references and perfect forwarding (although I realize Octave doesn't use C++-11 so those won't be available). I also posted a lot of questions on Stack Overflow.
I can make the source-code available if you'd like to look at it. I've attaced a scapegoat-tree class I used with self-reference-counting nodes. The idea is that you should be able to iterate over the elements in the tree and the iterator still behaves reasonably if something erases the element it's is sitting on. I use a path-compression trick like in a union-find data-structure to avoid the potential O(n^2) cost of having each iterator step forward to the next not-deleted node. This compiles with g++-4.6 using --std=c++0x . Also, you might have to add -D__GXX_EXPERIMENTAL_CXX0X__
I only tested it inside my particular application, so it may still be buggy or slow for general use. Of course, when working on Octave, I'll try to be much more thorough since I know I have to try and anticipate everyone's use case, not just my own.
p_tree.h