2
$\begingroup$

I need to find inverse of a large sparse matrix with small diagonal element. I follow this question - Is there any way to obtain an approximate inverse for very large sparse matrices? and try to use Daniel Lichtblau's answer. But the problem is it is not converging for small diagonal element. For example

n = 40; (*Dimension of matrix*)
X = 10.; (*Diagonal element*)
s = SparseArray[{
    Band[{1, 1}] -> X,
    Band[{1, 2}, {20, 21}] -> {1, 1, 1, 1, 1, 1, 0},
    Band[{2, 1}, {21, 20}] -> {1, 1, 1, 1, 1, 1, 0},
    Band[{1, 8}, {14, 21}] -> 1,
    Band[{8, 1}, {21, 14}] -> 1
    }, {n, n}]

f = LinearSolve[s];
aInv = f[SparseArray[{Band[{1, 1}] -> 1.}, {n, n}, 0.]];

inv1 = SparseArray[Band[{1, 1}] -> 1/Normal[Diagonal[s]], {n, n},0.];
sparseIden = SparseArray[Band[{1, 1}] -> 1., {n, n}, 0.];

niter = 5;(*number of iteration*)
Do[
  residual1 = inv1.s - sparseIden;
  delta1 = -residual1.inv1;
  inv1 = inv1 + delta1;
  Max[Abs[residual1]]//Print;
, {i, niter}]

gives an output

$0.1 \\ 0.04 \\ 0.0034 \\ 0.00003332 \\ 4.11443 \times 10^{-9}$

which shows a fair convergence. But for smaller X it is not converging at all. I try to check it over a range of of X and it looks like it is converging only for X>5.

enter image description here

Since it is a small matrix I use Inverse which gives an exact inverse (Max[Abs[residual1]] ~ $10^{-16}$ for all X). But it is not an good option if n->10000.

What would be a suitable method to find inverse in such cases?

$\endgroup$
  • $\begingroup$ Of course, your "approximate inverse" may well be much denser than your original one. $\endgroup$ – J. M. is away Jul 30 '17 at 3:42
  • $\begingroup$ @J.M., my main concern is the applicability of the method for small diagonal elements. I didn't find any method which can handle such large sparse matrices. Later I tried to get some selective elements (How to find selected elements of inverse of a banded matrix without inverting it?) near diagonal, but that didn't go well either. $\endgroup$ – Sumit Jul 30 '17 at 13:30

Your Answer

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

Browse other questions tagged or ask your own question.