Block sub diagonals matrix

I have known how to build a block main diagonals matrix. But I also need to fill the sub diagonals.

For example,

A=rand(3,3,5);
B=rand(3,3,5);
C=rand(3,3,5);

and I want to build a matrix D like,

D=[A(:,:,1),B(:,:,1),0,...,0
   C(:,:,1),A(:,:,2),B(:,:,2),...,0
   C(:,:,2),A(:,:,3),B(:,:,3),...,0
   ...
    A(:,:,5)]

D is 15×15

And A(:,:,i) is main diagonal, B and C is sub diagonals.

Kindly help me with this.

Answer 1

[DGM]

If my interpretation is correct, this should be one method:

bksize = [2 2]; % smaller for example

A = ones(bksize);
B = 11*A;
C = 111*A;

blocks = {zeros(size(A)) A B C}; % this makes the blocks indexable
map = toeplitz([2 4 1 1 1],[2 3 1 1 1]); % a map of the block locations

D = cell2mat(reshape(blocks(map),size(map))) % the block array

D = 10×10

1 1 11 11 0 0 0 0 0 0 1 1 11 11 0 0 0 0 0 0 111 111 1 1 11 11 0 0 0 0 111 111 1 1 11 11 0 0 0 0 0 0 111 111 1 1 11 11 0 0 0 0 111 111 1 1 11 11 0 0 0 0 0 0 111 111 1 1 11 11 0 0 0 0 111 111 1 1 11 11 0 0 0 0 0 0 111 111 1 1 0 0 0 0 0 0 111 111 1 1

If the blocks are 3D, that can be done without changing anything.

bksize = [2 2 3]; % 3D blocks

A = ones(bksize);
B = 11*A;
C = 111*A;

blocks = {zeros(size(A)) A B C};
map = toeplitz([2 4 1 1 1],[2 3 1 1 1]);

D = cell2mat(reshape(blocks(map),size(map)))

D = 
D(:,:,1) =

     1     1    11    11     0     0     0     0     0     0
     1     1    11    11     0     0     0     0     0     0
   111   111     1     1    11    11     0     0     0     0
   111   111     1     1    11    11     0     0     0     0
     0     0   111   111     1     1    11    11     0     0
     0     0   111   111     1     1    11    11     0     0
     0     0     0     0   111   111     1     1    11    11
     0     0     0     0   111   111     1     1    11    11
     0     0     0     0     0     0   111   111     1     1
     0     0     0     0     0     0   111   111     1     1


D(:,:,2) =

     1     1    11    11     0     0     0     0     0     0
     1     1    11    11     0     0     0     0     0     0
   111   111     1     1    11    11     0     0     0     0
   111   111     1     1    11    11     0     0     0     0
     0     0   111   111     1     1    11    11     0     0
     0     0   111   111     1     1    11    11     0     0
     0     0     0     0   111   111     1     1    11    11
     0     0     0     0   111   111     1     1    11    11
     0     0     0     0     0     0   111   111     1     1
     0     0     0     0     0     0   111   111     1     1


D(:,:,3) =

     1     1    11    11     0     0     0     0     0     0
     1     1    11    11     0     0     0     0     0     0
   111   111     1     1    11    11     0     0     0     0
   111   111     1     1    11    11     0     0     0     0
     0     0   111   111     1     1    11    11     0     0
     0     0   111   111     1     1    11    11     0     0
     0     0     0     0   111   111     1     1    11    11
     0     0     0     0   111   111     1     1    11    11
     0     0     0     0     0     0   111   111     1     1
     0     0     0     0     0     0   111   111     1     1

[剑豪 戴]

Thanks! Maybe my description of problem is not clear and correct. What I want is:

For 3D matrix A(3,3,5), A(:,:,i) is different from A(:,:,j). And the main diagonal of D will be [A(:,:,1), A(:,:,2), ... ,A(:,:,5)].

The dimension of D is (3*5)×(3*5)×1.

But I follow your interpretation, it work right for me! Thanks!

A = rand(2,2,3);
B = rand(2,2,3);
C = rand(2,2,3);
[k,l,m] = size(A);
A = mat2cell(A,k,l,ones(1,m));
B = mat2cell(B,k,l,ones(1,m));
C = mat2cell(C,k,l,ones(1,m));
blocks = [{zeros(k,l)},A(:)',B(:)',C(:)'];
map = diag(1:m,0)+diag(m+1:(2*m-1),1)+diag((2*m+2):3*m,-1)+1;
D = cell2mat(reshape(blocks(map),size(map)))

D = 6×6

0.3994 0.3321 0.2983 0.3994 0 0 0.8713 0.4807 0.3125 0.3380 0 0 0.6651 0.4984 0.8404 0.2329 0.4791 0.0781 0.4750 0.6980 0.2003 0.1860 0.3050 0.0632 0 0 0.1076 0.2159 0.1666 0.3028 0 0 0.6675 0.6244 0.7260 0.2880

But I worried this is inefficient. Cause the acutal dimension of matrix in my code is 1000×1000. I wonder if it could use some techniques of sparse matrix.

[DGM]

Sorry about the misunderstanding.

I don't know about the efficiency of the method given, but I'd imagine it's fair enough. I don't work with sparse tools, so offhand I don't know if there would be significant advantages or at what scale those advantages manifest. The only way to know would be to test various methods with different size inputs.

EDIT:

This is the above method versus a slightly modified version. On my hardware, the modified version is faster for smaller arrays, but reaches equivalence for about a 6000x6000 output size.

% prepare inputs
A = rand(3,3,100);
B = A*11;
C = A*111;

% original method
a = timeit(@() testA(A,B,C))

a = 0.0095

% avoiding mat2cell calls
b = timeit(@() testB(A,B,C))

b = 0.0040

% time ratio
a/b

ans = 2.3716

function testA(A,B,C)
    [k,l,m] = size(A);
    
    A = mat2cell(A,k,l,ones(1,m));
    B = mat2cell(B,k,l,ones(1,m));
    C = mat2cell(C,k,l,ones(1,m));
    
    blocks = [{zeros(k,l)},A(:)',B(:)',C(:)'];
    
    map = diag(1:m,0)+diag(m+1:(2*m-1),1)+diag((2*m+2):3*m,-1)+1;
    D = cell2mat(reshape(blocks(map),size(map)));
end

function testB(A,B,C)
    [k,l,m] = size(A);
    
    A = reshape(A,k,[]);
    B = reshape(B,k,[]);
    C = reshape(C,k,[]);
    
    blocks = mat2cell([zeros(k,l) A B C],k,l*ones(1,m*3+1));
    
    map = diag(1:m,0)+diag(m+1:(2*m-1),1)+diag((2*m+2):3*m,-1)+1;
    D = cell2mat(reshape(blocks(map),size(map)));
end

[剑豪 戴]

Thanks! Your advanced method worked nice for me! But in the end, I need to do x=D\b in linear equation system. Maybe sparse will be better.