Try http://cap.anu.edu.au/cap/projects/linalg + links. Just a thought: Generally you lose efficiency when opting for a direct solver (scales at best O(n to the 4/3)). Casting your discretization in a form suitable (read: high degree of diagonal dominance) for iterative solvers you can achieve O(nlogn) or even O(n) efficiency, which is a hell of a difference for large-scale problems.