# Creating large matrices from sub blocks

I am trying to generate a matrix from square blocks. Effectively, I have three blocks, call them $A_{diag}$, $A_{up}$, $A_{down}$, and a block of zeroes, which I'll just call $0$. I have some integer $L$ such that if $L=1$ then I am generating the matrix

$$M = \begin{pmatrix} A_{diag} \end{pmatrix}$$

If $L=2$ then I should generate

$$M = \begin{pmatrix} A_{diag} & A_{up} \\ A_{down} & A_{diag} \end{pmatrix}$$

If $L=3$ then

$$M = \begin{pmatrix} A_{diag} & A_{up} & 0 \\ A_{down} & A_{diag} & A_{up} \\ 0 & A_{down} & A_{diag} \end{pmatrix}$$

If $L=4$ then

$$M = \begin{pmatrix} A_{diag} & A_{up} & 0 & 0\\ A_{down} & A_{diag} & A_{up} & 0 \\ 0 & A_{down} & A_{diag} & A_{up} \\ 0 & 0 & A_{down} & A_{diag}\end{pmatrix}$$

and so forth. I suspect this is relatively simple, I found other questions which were similar but more complicated and specific. Any help appreciated, thank you.

While the other posted methods work, this is the one I wound up using, based on corey979's suggestion:

If[L == 1,
mup   = SparseArray[Band[{1, 2}] -> 1, {L, L}];
mdown = SparseArray[Band[{2, 1}] -> 1, {L, L}];
mdiag = SparseArray[Band[{1, 1}] -> 1, {L, L}];
A[r_] :=
KroneckerProduct[mup  ,   Aup[r]] +
];


It has the advantage of using SparseArray's, which is desirable if L is large. One can always use Normal to convert the SparseArray's into regular matrices for small L. For my application for example I found that using normal matrices was faster for L < 3 but using SparseArray's was faster for L > 3, with both being equally fast for L = 3.

• Wow, this one is really fast. Respect. – Henrik Schumacher Mar 11 '18 at 22:46
• For your mdiag, I would instead use IdentityMatrix[L, SparseArray]. – J. M. will be back soon Mar 17 '18 at 11:34
• Thanks, I don't notice much change in speed but it certainly makes more sense. – Kai Mar 18 '18 at 5:38

You can use ArrayFlatten and DiagonalMatrix:

mat[diag_, up_, down_, n_] := ArrayFlatten @ Block[{Identity},
DiagonalMatrix[ConstantArray[Identity[diag], n]]+
DiagonalMatrix[ConstantArray[Identity[up], n-1], 1]+
DiagonalMatrix[ConstantArray[Identity[down], n-1], -1]
]


I used a blocked Identity head so that DiagonalMatrix didn't complain.

Example:

mat[Array[a, {2,2}], Array[b, {2,2}], Array[c, {2,2}], 3] //TeXForm


$\left( \begin{array}{cccccc} a(1,1) & a(1,2) & b(1,1) & b(1,2) & 0 & 0 \\ a(2,1) & a(2,2) & b(2,1) & b(2,2) & 0 & 0 \\ c(1,1) & c(1,2) & a(1,1) & a(1,2) & b(1,1) & b(1,2) \\ c(2,1) & c(2,2) & a(2,1) & a(2,2) & b(2,1) & b(2,2) \\ 0 & 0 & c(1,1) & c(1,2) & a(1,1) & a(1,2) \\ 0 & 0 & c(2,1) & c(2,2) & a(2,1) & a(2,2) \\ \end{array} \right)$

• Thanks, the only downside to this is it doesn't work when n=1 – Kai Feb 1 '18 at 1:16
• I simply added If[n>1, (your code) , Adiag] which works perfectly it seems – Kai Feb 1 '18 at 1:31

The "SparseArray" context contains many hidden gems, for example this one:

SparseArraySparseBlockMatrix[{
Band[{1, 1}] -> SparseArray@Array[a, {2, 2}],
Band[{1, 2}] -> SparseArray@Array[b, {2, 2}],
Band[{2, 1}] -> SparseArray@Array[c, {2, 2}]
}, {3, 3}] // MatrixForm


$$\left( \begin{array}{cccccc} a(1,1) & a(1,2) & b(1,1) & b(1,2) & 0 & 0 \\ a(2,1) & a(2,2) & b(2,1) & b(2,2) & 0 & 0 \\ c(1,1) & c(1,2) & a(1,1) & a(1,2) & b(1,1) & b(1,2) \\ c(2,1) & c(2,2) & a(2,1) & a(2,2) & b(2,1) & b(2,2) \\ 0 & 0 & c(1,1) & c(1,2) & a(1,1) & a(1,2) \\ 0 & 0 & c(2,1) & c(2,2) & a(2,1) & a(2,2) \\ \end{array} \right)$$

It's not very efficient, though...

There it is,

M[n_] := Adiag IdentityMatrix[n] +
Table[Aup KroneckerDelta[i + 1, j], {i, 1, n}, {j, 1, n}] +
Table[Adown KroneckerDelta[i - 1, j], {i, 1, n}, {j, 1, n}];

MatrixForm[M]


to elaborate more, you can have Aup[i, j] and Adown[i, j] and Adiag[i, i].