I have a really big symmetric 7.000.000 X 7.000.000 matrix that i would like to invert. The matrix is extremely sparse and it can be rearranged as to become a block diagonal matrix. The biggest blocks are around 1500 X 1500 but most blocks are much smaller like 4 X 4. What is an efficient way to find the inverse? scipy.sparse.linalg.inv
in Python gives me a memory error.
-
3$\begingroup$ What’s the application for this inverse? Do you want to solve a linear system? There are very few cases that you actually needs an inverse. $\endgroup$– Mithridates the GreatJan 17, 2020 at 1:43
-
$\begingroup$ Do you mean to say that it is block diagonal / upper triangular? Or just that it could be in principle divided in blocks, like every matrix can? $\endgroup$– Federico PoloniJan 17, 2020 at 7:59
-
1$\begingroup$ I mean block diagonal. I could cut the matrix up in pieses, invert them, and then paste them back together. $\endgroup$– CasparJan 17, 2020 at 12:31
-
2$\begingroup$ I think it’s a bad design of your software that only accepts inverse matrix. In my experience, it’s impossible to find the inverse of 7Mx7M matrix in a reasonable time. $\endgroup$– Mithridates the GreatJan 18, 2020 at 21:33
-
1$\begingroup$ You can likely speed things up a lot by looking at using conjugate gradient methods. These are perfect for very sparse matrices where the matrix-vector products can be computed efficiently, which I suspect is the case here $\endgroup$– whpowell96Jan 18, 2020 at 23:17
2 Answers
I think you already know the answer- find the inverse of each of the small diagonal blocks. The details of how you might break the matrix up into the block in Python are something out of scope for this stackexchange group.
As others have suggested you most likely don't really need the inverse of this matrix and could probably do everything you need with a Cholesky factorization of the individual blocks.
I'd encourage you to step back and ask a broader question that includes relevant background and explains what problem you are actually trying to solve.
-
1$\begingroup$ He never said his matrix was SPD, so LU of each block may be a better choice instead of Cholesky $\endgroup$– vibeJan 20, 2020 at 23:30
I believe breaking your monstrous matrix into smaller blocks won't help you that much cause as I mentioned in my comment, using inverse of a huge sparse matrix is a really bad idea cause generally inverse of a huge sparse matrix could be a dense matrix, which probably you can't afford to even store it in your memory. But, in order to show you how you can find the inverse of a general matrix in terms of inverse of blocks (keep in mind that you need to break up your matrix into blocks that actually each block should be invertible.), I would suggest this mechanism, which I doubt would be useful at all here. Let's say you can break up your $\mathbf{X}$ matrix into this block form:
$$\mathbf{X} = \begin{bmatrix} \mathbf{A} & \mathbf{B} \\ \mathbf{C} & \mathbf{D} \end{bmatrix}$$
Remember $\mathbf{A}$ and $\mathbf{D}$ must be square matrices and be invertible and also $\mathbf{D} - \mathbf{C} \mathbf{A}^{-1} \mathbf{B}$ must be invertible too. So, you see it's not really easy to satisfy these conditions, specifically the last one, but if you satisfy them, you find the inverse of $\mathbf{X}$ as:
$$\mathbf{X}^{-1} = \begin{bmatrix} \mathbf{A}^{-1} + \mathbf{A}^{-1}\mathbf{B}(\mathbf{D} - \mathbf{C} \mathbf{A}^{-1} \mathbf{B})^{-1} \mathbf{C} \mathbf{A}^{-1} & -\mathbf{A}^{-1} \mathbf{B} (\mathbf{D} - \mathbf{C} \mathbf{A}^{-1} \mathbf{B})^{-1} \\ -(\mathbf{D} - \mathbf{C} \mathbf{A}^{-1} \mathbf{B})^{-1} \mathbf{C} \mathbf{A}^{-1} & (\mathbf{D} - \mathbf{C} \mathbf{A}^{-1} \mathbf{B})^{-1} \end{bmatrix}$$
Also, somewhere in your comments, you mentioned that your matrix is block diagonal. So, if it is the case, instead of using the above mechanism, you have this for $\mathbf{X}$:
$$\mathbf{X} = \begin{bmatrix} \mathbf{A}_1 & 0 & \cdots & 0 \\ 0 & \mathbf{A}_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & \mathbf{A}_n \end{bmatrix}$$
Where $\mathbf{A}_{i}$ are block diagonals. So, you find the inverse much easier but still you need the inverse of each block diagonal, as:
$$\mathbf{X}^{-1} = \begin{bmatrix} \mathbf{A}_{1}^{-1} & 0 & \cdots & 0 \\ 0 & \mathbf{A}_{2}^{-1} & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & \mathbf{A}_{n}^{-1} \end{bmatrix}$$
It looks easier, but the challenge is that can you guarantee that each block diagonal $\mathbf{A}_{i}$ is invertible? I doubt you can, but give it a try!
Finally, I second Brian Borchers suggestion and recommend you to think if it is possible to avoid using inverse, which makes your life much much easier.