f. c.
f. c.

Reputation: 1135

Calculate the inverse matrix of a diagonal blockwise matrix in matlab

I have a large matrix M like this

M=[A1, Z,  Z,  Z,  Z,  Z ; 
   Z,  A2, Z,  Z,  Z,  Z ; 
   Z,  Z,  A3, Z,  Z,  Z ; 
   Z,  Z,  Z,  A4, Z,  Z ; 
   Z,  Z,  Z,  Z,  A5, Z ; 
   Z,  Z,  Z,  Z,  Z,  A6];

A1,A2,A3,A4,A5,A6 are 4×4 real symmetric matrices, and Z=zeros(4,4).

How to calculate the inverse of M when there are millions of A in the matrix A1,A2,A3,..., An?

I know that I can simplify the inverse matrix to

invM=[B1, Z,  Z,  Z,  Z,  Z 
      Z,  B2, Z,  Z,  Z,  Z 
      Z,  Z,  B3, Z,  Z,  Z 
      Z,  Z,  Z,  B4, Z,  Z 
      Z,  Z,  Z,  Z,  B5, Z 
      Z,  Z,  Z,  Z,  Z,  B6];

B1,B2,B3,B4,B5,B6 are inverse matrices of A1,A2,A3,A4,A5,A6. But when there are many many B, how to do batch processing?

Thank you in advance.

Upvotes: 3

Views: 3762

Answers (2)

Rody Oldenhuis
Rody Oldenhuis

Reputation: 38032

Frankly, I wouldn't bother about the inverse. Chances are you don't need the inverses at all, but rather, you need the products

x(n) = inv(A(n))*b(n)

where b is the solution vector in the equation Ax = b.

Here's why that is important:

clc

N = 4;    % Size of each A
m = 300;  % Amount of matrices

% Create sparse, block diagonal matrix, and a solution vector
C = cellfun(@(~)rand(4), cell(m,1), 'UniformOutput', false);
A = sparse(blkdiag(C{:}));
b = randn(m*N,1);


% Block processing: use inv()
tic
for ii = 1:1e3
    for jj = 1:m
        inds = (1:4) + 4*(jj-1);
        inv(A(inds,inds)) * b(inds); %#ok<MINV>
    end    
end
toc

% Block processing: use mldivide()
tic
for ii = 1:1e3
    for jj = 1:m
        inds = (1:4) + 4*(jj-1);
        A(inds,inds) \ b(inds);
    end
end
toc

% All in one go: use inv()
tic
for ii = 1:1e3
    inv(A)*b;
end
toc

% All in one go: use mldivide()
tic
for ii = 1:1e3
    A\b;
end
toc

Results:

Elapsed time is 4.740024 seconds.  % Block processing, with inv()
Elapsed time is 3.587495 seconds.  % Block processing, with mldivide()
Elapsed time is 69.482007 seconds. % All at once, with inv()
Elapsed time is 0.319414 seconds.  % All at once, with mldivide()

Now, my PC is a little different than most PCs out there, so you might want to re-do this test. But the ratios will be roughly the same -- computing explicit inverses simply takes a lot more time than computing the product x = inv(A)*b.

If you run out of memory when using mldivide, you shouldn't loop over all the individual matrices, but rather split the problem up in larger chunks. Something like this:

chunkSize = N*100;
x = zeros(size(A,1),1);
for ii = 1:N*m/chunkSize
    inds = (1:chunkSize) + chunkSize*(ii-1);
    x(inds) = A(inds,inds) \ b(inds);
end

Upvotes: 1

Vuwox
Vuwox

Reputation: 2359

The inverse of a diagonal matrix is only B = 1/A.

Here the proof : http://www.proofwiki.org/wiki/Inverse_of_Diagonal_Matrix

Upvotes: 0

Related Questions