Papers
Chan, Tony, and Henk A. Van der Vorst. "Approximate and Incomplete Factorizations." (PDF)
Benzi, M., G. H. Golub, and J. Liesen. "Numerical Solution of Saddle Point Problems." Acta Numerica 14 (2005): 1-137. (To appear.) (PDF - 1.0 MB)
MATLAB® Documentation
MATLAB® Technical Documentation
Documentation for MathWorks Products, R2006a
MATLAB® Mathematics Documentation
Other Thoughts
Is there a simple way to count the fill-in in elimination?
Run chol (or lu) and use "nnz" (number of non-zeros)
% Generate discrete Laplacian operator
A=delsq(numgrid('S',100));
% Cholesky factorize
R=chol(A);
nnz(R)
ans = 941289
% As above, but permute rows & columns by approximate minimum degree
p=symamd(A);
R=chol(A(p,p));
nnz(R)
ans = 184917
Does min degree just take nodes a row at a time (if it breaks ties correctly)? If that is N^3 - and would sparse backslash do it?
Yes it does, it is a "greedy algorithm." The approximate one might be a little different, but in principle it does the same thing. Sparse backslash does this automatically. If symmetric, then it uses "symmmd" (symmetric minimum degree), then "chol." If non-symmetric, "colamd" (column approximate minimum degree), then umfpack.
This is about as fast as possible except for Fast Poisson Solvers?
Well, it is probably as fast as possible for Gaussian elimination. But iterative solvers do better, I think CG with Incomplete Cholesky should be N^(2.5). Then of course fast solvers, FFT-based and multigrid.
Develop a code for Incomplete LU and try it on the 2D matrix for increasing N. You can make it a general code or specialize it to these particular matrices (and decide which entries to keep in the approximate factors L and U). Try on A2D u = random right side. You could use eig to find the spectral radius of the preconditioned matrix B. That matrix is I - inv(P)A2D where P ≈ LU. Of course if P = A2D then B = 0 and convergence in 1 step.