function f = amdSolveF(la)
% function f = amdSolveF(la)
%
% la should be a laplacian, or other
% diag dom matrix.  applied the amd order,
% and produces a solver for la
%
% exmple:
% a = grid2(50);
% la = lap(a) + speye(2500)/10^9;
% f = amdSolveF(la);
% b = randn(2500,1);
% norm(la*f(b) - b)


p = symamd(la);
plap = la(p,p);
[l,u] = lu(plap);

f = lupSolveF(l,u,p);