4
$\begingroup$

I have a $50000\times 50000$ matrix $A$ sparse matrix containing only 5 non-zero elements in each row. Now the problem is that the diagonal elements and the constants (in $B$ matrix such that $AX=B$) get updated after every iteration. If I go by the usual Matlab function inv or via Gauss Elimination, it takes around 130 seconds for the solution to be computed for a single iteration and the number of iterations required by the problem is somewhat of the order 100-500 to compute the final solution $X$. Need some suggestions on this one if I am not wanting to use parallel computation.

Improve this question
$\endgroup$
9
  • 5
    $\begingroup$ Inverting your sparse matrix will inevitably yield a dense one, so not a good plan. You can use an iterative method instead. Speaking of which, you haven't mentioned if your matrix is symmetric (positive definite) or not... $\endgroup$
    – J. M.
    Jun 19, 2013 at 14:57
  • $\begingroup$ Do you update both A and B? $\endgroup$
    – Jan
    Jun 19, 2013 at 15:36
  • 3
    $\begingroup$ Welcome to SciComp.SE! Just to make sure: You have a matrix $A\in\mathbb{R}^{n\times n}$ and a matrix $B\in\mathbb{R}^{n\times m}$, and are interested in computing the matrix $X\in\mathbb{R}^{n\times n}$ such that $AX=B$ (that is, all you really want is $X$, and not the inverse $A^{-1}$), and you are using Matlab? In this case, have you tried X=A\B (never use inv for solving linear systems, that's not what it is for!), which is the most efficient "black-box" way of solving your problem in Matlab? $\endgroup$
    – Christian Clason
    Jun 19, 2013 at 17:05
  • 1
    $\begingroup$ Also, an important information is the size of B (same size as A or much smaller)? $\endgroup$
    – Christian Clason
    Jun 19, 2013 at 17:07
  • 1
    $\begingroup$ @clipper: Good point. I think without the help of the OP, there is little else we can suggest. I'm startled that none of our questions have gotten any answers -- does he no longer care about the question? $\endgroup$
    – Wolfgang Bangerth
    Jun 20, 2013 at 6:52

2 Answers 2

Reset to default
3
$\begingroup$

The first idea that comes to my mind is to use an iterative method for solving the linear system. This might speed-up the computation of a single iteration (depending on the method (GMRES, conjugate gradient, BiCGStab,... and the size of B) and it will take advantage of the sparsity of the matrix A.

Moreover, iterative methods let you define a preconditioner. If the diagonal coefficient does not change too much between the iterations, it might be worth spending some time at the begining to build a good preconditioner (e.g. ILU) and reuse it for all the iterations.

Hope it helps!

Edit: The answers to this question might be useful.

Improve this answer
$\endgroup$
1
$\begingroup$

You may use a sparse factorization algorithm, it means computing matrices $P$, $L$, $U$, such that $M = PLU$ where $P$ is a permutation matrix, $L$ a sparse lower triangular matrix and $U$ a sparse upper triangular matrix. The permutation matrix is there and computed in such a way that $L$ and $U$ remain reasonably sparse (without it it is not the case in general, as indicated in the other answer). Then when $M$ is factorized in this form, it is trivial to solve a linear system with an arbitrary rhs. If $M$ is symmetric definite, then there is a sparse cholesky factorization ($M = PLL^t$), if it is symmetric only, then it can be decomosed as $M = PLDL^t$, with $D$ a diagonal matrix.

There are several available implementations of sparse factorization: SuperLU, Choldmod, Mumps, TAUCS (depending on whether you need $LU$, $LL^t$ or $LDL^t$). There probably exists MATLAB bindings for most of them.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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