Inverse of a large sparse Hermitian block matrix

Question

I am looking for a method (if it exists) for the inverse of a large sparse Hermitian block matrix. The off diagonal sparse matrices, named δ are 4x4, and they have only one non null element on the diagonal part. The matrices M are Hermitian and also 4x4. I need a way to perform the inverse of M3 or the determinant. I am interested in a generalization of this structure for a matrix of higher rank with the same "telescopic" structure. Here is the code

dimension        = 4;
mat              = ( {{m, a, b, b},{Conjugate[a], m, b, b},{Conjugate[b], Conjugate[b], m, a},{Conjugate[b], Conjugate[b], Conjugate[a], m}});
Subscript[O, 4]  = ( {{0, 0, 0, 0},{0, 0, 0, 0},{0, 0, 0, 0},{0, 0, 0, 0}} );
Subscript[δ, 1]  = ( {{I ρ, 0, 0, 0},{0, 0, 0, 0},{0, 0, 0, 0},{0, 0, 0, 0}} );
Subscript[δ, 2]  = ( {{0, 0, 0, 0},{0, I ρ, 0, 0},{0, 0, 0, 0},{0, 0, 0, 0}} );
Subscript[δ, 3]  = ( {{0, 0, 0, 0},{0, 0, 0, 0},{0, 0, I ρ, 0},{0, 0, 0, 0}} );
Subscript[δ, 4]  = ( {{0, 0, 0, 0},{0, 0, 0, 0},{0, 0, 0, 0},{0, 0, 0, I ρ}} );
Subscript[δ, 22] = ArrayFlatten[( {{Subscript[O, 4], Subscript[δ, 2]}, {Subscript[δ, 2], Subscript[O, 4]}})];
Subscript[δ, 33] = ArrayFlatten[( {
    {Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[δ, 3]},
    {Subscript[O, 4], Subscript[O, 4], Subscript[δ, 3], Subscript[O, 4]},
    {Subscript[O, 4], Subscript[δ, 3], Subscript[O, 4], Subscript[O, 4]},
    {Subscript[δ, 3], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4]}
} )];
Subscript[δ, 44] = ArrayFlatten[( {
    {Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[δ, 4]},
    {Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[δ, 4], Subscript[O, 4]},
    {Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[δ, 4], Subscript[O, 4], Subscript[O, 4]}, 
    {Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[δ, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4]}, 
    {Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[δ, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4]},
    {Subscript[O, 4], Subscript[O, 4], Subscript[δ, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4]},
    {Subscript[O, 4], Subscript[δ, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4]},
    {Subscript[δ, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4], Subscript[O, 4]}
} )];

OutM[M_, δ_] := ArrayFlatten[{{M, δ},{-δ, M}}];
M1 = OutM[mat, Subscript[δ, 1]];
M2 = FullSimplify[OutM[M1, Subscript[δ, 22]]];
M3 = FullSimplify[OutM[M2, Subscript[δ, 33]]];
M4 = FullSimplify[OutM[M3, Subscript[δ, 44]]];

And I need the inverse of M4. The structure is telescopic in some sense.

Even better, could you provide code to construct it so we don't have to re-write it? — acl, May 28 '12 at 11:11
You are aware that inverting a large, sparse matrix invariably results in a large, dense matrix, right? Or do you have any guarantee that the inverse has a known/predictable sparsity pattern? — J. M.'s missing motivation, May 28 '12 at 11:20
If the matrix isn't too big, what exactly is the problem with using Inverse[]/Det[] directly? — J. M.'s missing motivation, May 28 '12 at 12:12
i wrote the structure again. When I raise up the rank of this kind of matrices mathematica is not able to perform the inverse. I need some tip in order to use the sparsity.... — Raffaele Carlone, May 28 '12 at 12:22
I think there's an error in your definition of OutM. A bracket and a comma are missing. I presume it should be OutM[M_, \[Delta]_] := ArrayFlatten[( {{M, \[Delta]},{-\[Delta], M}} )]; — Sjoerd C. de Vries, May 28 '12 at 13:34
I don't see the problem. On my laptop MMA finds the Inverse of M1, M2, M3, M4 in 0.01, 0.4, 5, and 66 seconds, respectively. — Sjoerd C. de Vries, May 28 '12 at 13:43
@SjoerdC.deVries, interesting - MMA 8.0.4 on my Core i5 machine takes 0.04s for Inverse[M1], 9 seconds for Inverse[M2] and I aborted Inverse[M3] after 5 minutes. — Simon Woods, May 28 '12 at 17:57
I have a Core i7, and I see MMA use 3 or 4 of the 8 hyperthreaded cores. That might make a difference. Did you use the same correction as I did (see comment above)? — Sjoerd C. de Vries, May 28 '12 at 18:34
OutM[M_, [Delta]_] := ArrayFlatten[( { {M, [Delta]}, {-[Delta], M} } )]; — Raffaele Carlone, May 28 '12 at 21:06
@SjoerdC.deVries if you can please send me the computation of the inverse at my address ozellus@gmail.com. Thank you — Raffaele Carlone, May 29 '12 at 09:07
@SimonWoods and Raffaele: I was trying to repeat what I did before and it didn't work out the same way. I had been working on this question which resulted in a and b already being defined as a numerical value. Now, with a and b symbolic, MMA indeed runs into problems. Sorry for the confusion. — Sjoerd C. de Vries, May 29 '12 at 10:15
in fact i am using your same cpu configuration and the result is not coming out. — Raffaele Carlone, May 29 '12 at 10:17
Here's a related question. The approach to block-inversion discussed in Daniel Lichtblau's answer should work here. — Eli Lansey, May 29 '12 at 21:36

Inverse of a large sparse Hermitian block matrix

0 Answers0