Are iterative methods always better than direct methods for solving large linear systems?

I am trying to solve a very sparse linear system Ax = b. A is very sparse - if A is of size N^2 x N^2, all nontrivial elements for any row k are located in between (k-2N, k+2N). The overall the number of nontrivial elements in A is bounded by 13*N^2 (for a matrix with N^4 elements)

Now back to the original question: Currently I am using mldivide to solve the system. On a 1e6 x 1e6 matrix, the program takes around 60-75s 64 bit machine. Now it is possible to beat this performance using iterative methods? The ones I have tried (built into MATLAB) do not seem to offer any advantages; however, I think that may be due to the fact that I am not using a preconditioner. The problem is, how do I effectively get a preconditioner? I tried using ilu, but that is also pretty slow (and there is no guarantee that it will do the trick).

Thanks, Peter

Answer 1

[Adrien Leygue]

From my personal experience, there is only one case where I have been able to outperform mldivide through the use of one of the Matlab-provided iterative solvers. In this particular application I had to solve many linear systems (several hundredth, size > 1e5 unknowns), all the systems were different (left and right hand side) but each system to solve was very close to the previous one. I could therefore provide a very good initial guess for the solution. There PCG with incomplete factorization provided some speed improvement as only a few iterations were needed.

The problem is that mldivide is so efficient that unless you have memory issues (or if you have a very good estimate of your solution), the provided iterative solvers are no match.

Hope this helps

A.

[Peter]

Hi Adrien,

Sorry for the late response, but with my recent work I seem to be coming to the same response. I recently tried to use iterative solvers again as I am looking at some very big matrices (4e6 x 4e6) for a finite difference scheme. The direct solver with mldivide can invert the system in ~ 40 minutes. The iterative solvers though are not that promising. I think the issue is that I cannot find a good preconditioner (I tried Jacobi and ILU). They don't seem to affect any of the methods except GMRES. The issue with GMRES is that while 10 iterations gives a rel residual of ~ 1e-2 - 1e-3, it seems to stangate as one approaches the 1e-6 tolerance

Answer 2

[John D'Errico]

The use case which some may miss here, is the very large array, where a dense solver fails due to memory, but even the sparse one fails due to fill-in, BUT, the problem is solvable using an iterative solver. This may require a sufficiently well-posed array that the iterative solver has no problem in converging. In this example, I'll make up a case where the matrix is extremely well-posed, just so I can solve it trivially online.

A = sprandn(1e5,1e5,0.001)*0.001 + speye(1e5);
max(sum(abs(A),2))

ans = sparse double
   (1,1)       1.1255

Clearly A is moderately sparse, but it is far too non-sparse that a sparse solver would fail miserably due to fill-in. And a dense solver would be a complete waste of time. As well, I've constructed A to be quite well-posed. Any iterative solver worth its salt will converge easily enough.

X0 = ones(1e5,1);
b = A*X0;
tic,X = lsqr(A,b); toc

lsqr converged at iteration 4 to a solution with relative residual 4.1e-08.
Elapsed time is 0.155120 seconds.

norm(X - X0)

ans = 1.2897e-05

On a somewhat worse system, now we could employ an incomplete preconditioner,

Have I ever had a real world case where the sparse solvers were the only recourse, and I needed to employ a preconditioner? YES.

Answer 3

[Wolfgang Schwanghart]

Why is performance a problem? Do you have to solve the system repeatedly or do you want to solve larger systems? If former is the case, then take a look at Tim Davis' factorize

http://www.mathworks.com/matlabcentral/fileexchange/24119-dont-let-that-inv-go-past-your-eyes-to-solve-that-system-factorize

If latter is the case, then I hope someone else will be able to provide an answer.

HTH, W.