  When solving linear systems arising from PDE discretizations, iterative
methods (such as Conjugate Gradient, GMRES, or MINRES) are often the only
practical choice. To converge in a small number of iterations, however, they
have to be coupled with an efficient preconditioner. The efficiency of the
preconditioner depends largely on its accuracy on the eigenvectors
corresponding to small eigenvalues, and unfortunately, black-box methods
typically cannot guarantee sufficient accuracy on these eigenvectors. Thus,
constructing the preconditioner becomes a problem-dependent task. However, for
a large class of problems, including many elliptic equations, the eigenvectors
corresponding to small eigenvalues are smooth functions of the PDE grid. In
this paper, we describe a hierarchical approximate factorization approach which
focuses on improving accuracy on the smooth eigenvectors. The improved accuracy
is achieved by preserving the action of the factorized matrix on piecewise
polynomial functions of the grid. Based on the factorization, we propose a
family of sparse preconditioners with $O(n)$ or $O(n \log{n})$ construction
complexities. Our methods exhibit the optimal $O(n)$ solution times in
benchmarks run on large elliptic problems of different types, arising for
example in flow or mechanical simulations. In the case of the linear elasticity
equation the preconditioners are exact on the near-kernel rigid body modes.
