How can I direct sum matrices into the middle of one another another?

Question

I would like to execute the mathematical operation of the direct sum of matrices in the case where the matrices are not appended one after the other along the diagonal, but instead mixed among one another. This must be generalizable to any permutation of mixing indices along the diagonal.

For example, say we have matrix A, B and C $$A=\begin{equation} \begin{bmatrix} \cos (t) & \sin(t) \\ -\sin(t) & \cos(t) \end{bmatrix} \end{equation}$$

$$B=\begin{equation} \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} \end{equation}$$

$$C=\begin{equation} \begin{bmatrix} \cos (u) & \sin(u) \\ -\sin(u) & \cos(u) \end{bmatrix} \end{equation}$$

$A ? B ? C=\begin{equation} \begin{bmatrix} \cos (t) & \sin(t) & 0 & 0 & 0 & 0\\ -\sin(t) & \cos(t) & 0 & 0 & 0 & 0\\ 0 & 0 & 1 & 0 & 0 & 0\\ 0 & 0 & 0 & \cos(u) & 0 & \sin(u)\\ 0 & 0 & 0 & 0 & 1 & 0\\ 0 & 0 & 0 & -\sin(u) & 0 & \cos(u) \end{bmatrix} \end{equation}$

How would this be done in a generalizable way in Mathematica? One could use a permutation matrix, $P=\begin{equation} \begin{bmatrix} 1 & 0 & 0 & 0 & 0 & 0\\ 0 & 1 & 0 & 0 & 0 & 0\\ 0 & 0 & 1 & 0 & 0 & 0\\ 0 & 0 & 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 & 0 & 0\\ 0 & 0 & 0 & 0 & 0 & 1 \end{bmatrix} \end{equation}$ where, $A?B?C=P[A⊕B⊕C]P^T$, but how could I create such a matrix P in Mathematica knowing only what the mapping is? In this case, the index mapping being: (1,2,3,4,5,6)->(1,2,3,5,4,6)

Answer 1

ClearAll[r1, t, u, a, b, c, abc]

r1 = {Cos[t], Sin[t]};

a = {r1, Cross[r1]};
b = IdentityMatrix[2];
c = a /. t -> u;

abc = SparseArray[Band[{1, 1}] -> {a, b, c}];

abc // MatrixForm // TeXForm

$\left( \begin{array}{cccccc} \cos (t) & \sin (t) & 0 & 0 & 0 & 0 \\ -\sin (t) & \cos (t) & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & \cos (u) & \sin (u) \\ 0 & 0 & 0 & 0 & -\sin (u) & \cos (u) \\ \end{array} \right)$

perm = {1, 2, 3, 5, 4, 6};

abc[[perm, perm]] // MatrixForm // TeXForm

$\left( \begin{array}{cccccc} \cos (t) & \sin (t) & 0 & 0 & 0 & 0 \\ -\sin (t) & \cos (t) & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & \cos (u) & 0 & \sin (u) \\ 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & -\sin (u) & 0 & \cos (u) \\ \end{array} \right)$

Alternatively,

p = IdentityMatrix[6][[perm]];

p // MatrixForm // TeXForm

$\left( \begin{array}{cccccc} 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 1 \\ \end{array} \right)$

p.abc.Transpose[p] // MatrixForm // TeXForm

$\left( \begin{array}{cccccc} \cos (t) & \sin (t) & 0 & 0 & 0 & 0 \\ -\sin (t) & \cos (t) & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & \cos (u) & 0 & \sin (u) \\ 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & -\sin (u) & 0 & \cos (u) \\ \end{array} \right)$

Note: An alternative way to get abc is to use SparseArray`SparseBlockMatrix:

abc2 = SparseArray`SparseBlockMatrix[MapIndexed[#2[[{1, 1}]] -> # &, {a, b, c}]];

abc2 == abc

True