Showing block diagonal structure of matrix by reordering

There are at least two approaches to this kind of problem, which is sometimes called the "bandwidth minimization problem." You can read about them here:

E. Cuthill and J. McKee, "Reducing the Bandwidth of Sparse Symmetric Matrices," in Proceedings of the 24 National Conference of ACM, New York: ACM Publications, 1969 pp. 157-172.

A. Lim, B. Rodrigues, and F. Xiao, "A Centroid-Based Approach to Solve the Bandwidth Minimization Problem," in Proceedings of the 37 Annual Hawaii International Conference on System Sciences (HICSS'04), Track 3, Vol. 3, Washington D.C.: IEEE Computer Society, 2004 p. 30075.1.

Both are implemented in Mathematica using the command MinimumBandwidthOrdering in the Graph Utilities Package, and the Matlab equivalent is the function symrcm.

First we set $tilde M$ for the matrix obtained by replacing all non-zero elements by 1. Start with a vector $tilde v_1$ having exactly one coordinate equal to $1$ and all other coordinates zero. Consider $v_2=v_1+tilde Mv_1$ and define $tilde v_2$ by replacing all non-zero elements (which are in fact strictly positive integers) of $v_2$ with $1$. Iterate until you get $tilde v_{k+1}=tilde v_k$. The coordinates of $tilde v_k$ form then a block. Removal of this block and iteration of this algorithm gives the block decomposition. (This is more or less one of the methods suggested by Gowers and I guess that it is also related to the Cuthill-McKee algorithm mentioned by Suvrit.)

Correction: The above algorithm does not give decomposition into blocks but a decomposition into blocks on the diagonal plus upper triangular stuff.

One has thus also to check if later appearing blocks interact with blocks considered previously. This can be done by considering an upper triangular matrix with rows and columns indexed by the constructed diagonal blocks. Its coefficients are $1$ on the diagonal and $1$ for indices $i,j$ such that the image of of the $j-$th block involves the subspace spanned by the $i-$th block (where we consider the decomposition $V=V_1oplus V_2oplusdots$ with $V_i$ corresponding to the image spanned by the $i-$th diagonal block).

You might want to have a look at the Cuthill-KcKee algorithm and from there perhaps also chase the other sparse matrix reordering algorithms. The abovecited method seeks to minimize the bandwidth of a sparse matrix---so, although it does not produce a block-diagonal, it produces a banded matrix (with narrow bandwidth), which might be useful to you.

If the matrix corresponds to a $d$-regular graph, then one can make some easy observations. For instance, if it is a block matrix with $k$ components, then there will be a $k$-dimensional eigenspace of vectors with eigenvalue $d$. As you say, this gives significant information about the blocks. You can't say that an arbitrary eigenvector will pick out a component, but you can say that the eigenspace will be generated by characteristic functions of components, and the eigenvectors will be constant on the components. One would have to make the notion of approximation precise, but presumably there are reasonable notions of approximation such that if a matrix is approximately a block matrix with $k$ components and is still $d$-regular, then there are several vectors $v$ that approximately satisfy $Av=dv$ and that are approximately constant on the approximate components. And finally, if all you're interested in is the components, then probably it's not too important what the values of the matrix are (as long as they are clearly zero or clearly non-zero), so perhaps you can reweight them so as to reduce to something like the regular case -- but I'm less sure about that last idea. All this would of course rely on points within components being significantly more connected than points in different components.

Another, related, way of detecting this would be to look at powers of the matrix, which would fill up the approximate components much more quickly than the links between components. In fact, I'd be tempted to pick a point and apply the matrix a few times to that point (considered as a unit vector) in order to identify its "approximate component". In the exact case, the component would consist of all points where the value is non-zero, whereas in the approximate case one would have to go for some kind of cutoff. Then one could repeat for other points. The resulting "components" might overlap a bit, but one could then just arbitrarily assign points in the overlap to some component or other.

I'm not sure what to make of your suggestion that coming up with a precise definition is part of the problem: I think there are probably several reasonable precise definitions, each of which leads to a different version of the problem. So either you should be more specific about why you want to do this or you should just choose your method of proof and define the problem in such a way that the method works. (So, for instance, you might decide on the above cutoff method and then choose a notion of "approximate block matrix" that guarantees that that process approximately identifies -- in some other sense that you make precise -- the approximate blocks.)

This kind of question is of interest in additive combinatorics. If $E$ is a set of integers, one can define a matrix $A$, indexed by $Etimes E$, where $A(x,y)$ is the number of pairs $(z,w)in E^2$ such that $x-y=z-w$. It would be nice to have a good way of splitting $E$ up into "approximate components" of this matrix/weighted graph, which would be "approximate additive components" of $A$. Even nicer would be to be able to say something nice about the additive structure of these components, on the assumption that there are at least $c|E|^3$ quadruples $(x,y,z,w)in E^4$ with $x-y=z-w$.

Given a symmetric matrix $Minmathbb{R}^{ntimes n}$ you are able to visualize the matrix as an adjacency matrix of a graph. Where node $i$ and $j$ have an edge between them if the $(i,j)$ entry of $M$ is non-zero. When the matrix is not symmetric then you end up looking at a directed graph.

A $k$ clustering algorithm on the associated graph should do a good job. Where $k$ is the number of blocks you want. This should also work when the matrix is not completely block diagonal.

If you make believe $M$ is the adjacency matrix of a graph (vertices $i$ and $j$ are adjacent if and only if at least one of the two entries $m_{ij},m_{ji}$ is nonzero), then there are efficient algorithms for finding the connected components of the graph, which should correspond to the blocks of the matrix.