Hi,
 * What is the best usage pattern to do a complete for_each over every element of a ublas matrix?  Unlike other foreachs, it is kind of necessary to have the index or iterator itself since position will usually be important.
 * Has anyone written an equivalent to BOOST_FOREACH that would work with the 2 iterators?  Seems to me could be generic given a concept of two iterators named: ::iterator1 and ::iterator2 within the object and not really ublas specific.  Could get fancy later on changing the order for rowmajor vs. columnmajor, but not initially needed.  Or could have a different macro to switch the looping order...
 * I was thinking the notation might be something like:

matrix<double> table(m, n);
...
BOOST_FOREACH2(double d, int i, int j, table) //Would love to have the current index of the iteration as well?
{
            cout << "Table at  " << i << ", " << j << " is " << d << endl;
}


 * Could get fancy later on changing the order for rowmajor vs. columnmajor, but not initially needed.  Or could have a different macro to switch the order...
 * Also, I am a little confused by the column and row iterators of ublas... do they return a ublas::vector?  I also would love to use with a foreach sometimes.
 * Any good documents or examples on these iterators?

And a related quetsion:
* One main use case is for an iterative algorithm.  I have a matrix, I want to do a calculation for each i, j pair within the matrix (execution order doesn't matter) including using the value.  And then I want to put the value back into a next iteration matrix.  What is a clean way to code this?  Using random access, can obviously just do two loops but I imagine it is much more efficient to use iterators.
* Or is this best done with a different algorithm construction like a transform?  But I haven't seen a macro for an STL transform that allows me to keep my logic where I want it, like foreach, instead of in a separate function.
 * The pseudocode would be something like:

matrix A0
matrix A1

foreach(double d, int i, int j)
{
  A1(i, j) = myfunc(A0, d, ...);
}
if norm(A1-A0) < eps ....
else A0 = A1


//or
foreach_transform(double d, A0, A1)
{
 A1 = myfunc(A0, d)
}
...
-Jesse