7
$\begingroup$

Is there an iterative solver that can handle potentially zero entries on the central diagonal? I am implementing a polynomial fitting algorithm (up to $10^{th}$-order) and my matrix is a "Vandermonde-like" system that results in zeros for the odd-powered polynomial basis functions when my sample space is symmetric with respect to the vertical axis of the local reference frame in which I'm fitting the polynomial. I have tried Gauss-Seidel and LUSGS, but later realised that they both rely on non-zero diagonal entries.

Edit: I have double-checked that the system is correctly implemented because it solves with Gauss elimination with back-substitution (LUP decomposition). However, this is too slow for my application (the $10^{th}$-order matrix is $121\!\times\! 121$), and would prefer an iterative solver.

$\endgroup$
  • 4
    $\begingroup$ What makes you think that an iterative solver will be faster than LUP for a small dense (and probably poorly conditioned) matrix like this? $\endgroup$ – Brian Borchers Jul 1 at 13:58
  • 2
    $\begingroup$ If the eigenvalues of the matrix are clustered, there is no reason to believe that an iterative solver might not be faster. For example, if there were only 10 distinct eigenvalues, then most iterative solvers will converge in 10 iterations requiring $10n^2$ operations -- still faster than the $O(n^3)$ for an LU decomposition. $\endgroup$ – Wolfgang Bangerth Jul 1 at 17:32
  • 2
    $\begingroup$ @WolfgangBangerth there's no particular reason to believe that eigenvalues will be clustered. $\endgroup$ – Brian Borchers Jul 1 at 18:46
  • 7
    $\begingroup$ @wolfgang bangerth, thinking in terms of complexities won't accurately depict behaviour for such small matrices. A simple LU will probably be orders of magnitude faster here. Especially if using a library like LAPACK. $\endgroup$ – Thijs Steel Jul 1 at 18:47
  • 1
    $\begingroup$ @ThijsSteel -- ah, I had not realized the matrix is small. You're likely right that for a 121x121 matrix, LU decomposition is likely faster than anything else, unless you only care about very relaxed accuracy. $\endgroup$ – Wolfgang Bangerth Jul 2 at 1:31
10
$\begingroup$

Iterative Krylov-subspace solvers generally only require matrix-vector products and don't care whether or where there are zeros in the matrix. In your case, unless you have other information about the matrix (e.g., symmetry), you could for example use GMRES.

What you probably had in mind is the question of preconditioning, and that you can't use things such as Gauss-Seidel or Jacobi preconditioning. That is true, but depending on what else you know of the matrix, you can do something like a "block" version of Jacobi (or Gauss-Seidel, or SSOR, or ...) in which you consider $2\times 2$ blocks of the matrix which might be invertible.

| cite | improve this answer | |
$\endgroup$
  • 1
    $\begingroup$ Okay, thank you for your response! I will implement GMRES and compare speeds. The reason I mentioned that Gauss-Seidel didn't work is because I implemented the algorithm on Wikipedia, which requires division by diagonal entries $A_{ii}$. $\endgroup$ – niran90 Jul 1 at 14:20
  • 1
    $\begingroup$ Iterative solvers generally only require matrix-vector products and don't care whether or where there are zeros in the matrix OP has Gauss-Seidel in mind, which is an iterative solver and does care if there are zeros on the diagonal. $\endgroup$ – Federico Poloni Jul 1 at 19:19
  • $\begingroup$ Ah, fair enough. I had not considered that anyone might still be interested in using fixed point iterations today as actual solvers. I clarified it by adding "Krylov subspace" to "iterative solvers". $\endgroup$ – Wolfgang Bangerth Jul 2 at 1:30
3
$\begingroup$

Conclusion: I have tried CG, BiCG, and GMRES with no pre-conditioner, and these are all at least twice as slow as LU decomposition for the 121x121 matrix. Further, even for convergence tolerances near machine-precision, I do not get the accuracy that I get with LU. So I guess I'm sticking with LU :)

| cite | improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.