# Programming for matrix manipulation

I want to write a programme which operates on nd x nd matrices. The programme looks the matrix as a n x n block matrix where each block is of size d x d. It applies a map, say f on each of these d x d matrices and output the global matrix. As a sample, I am giving a code, which is on 64 x 64 matrices and the map is applying on the 4 x 4 block matrices.

Big[A0_] :=
Module[{A = A0},
For[i = 1, i <= 64, i = i + 4,
For[j = 1, j <= 64, j = j + 4,
A[[i ;; i + 3, j ;; j + 3]] = f[A[[i ;; i + 3, j ;; j + 3]]]
]
];
A]


I then used this simple code to define the function f and a list L.

X = Table[Subscript[x, i, j], {i, 64}, {j, 64}];
L = {a, b, c, d}
f[X_] := DiagonalMatrix[(m + 1)*Diagonal[X] +
L*RotateLeft[Diagonal[X], 1]];
Big[X] // Simplify // MatrixForm


I think, I am forcing Mathematica to do some unwanted calculations (in those for loops), which do not reflect in the output matrix and are unnecessary. So question, is it possible to make it simpler.

Further, I am putting the dimension of the larger matrix by hand (here it is 64) and also the break up of n x d block matrices. If I want to apply this map on dimension 4 to a matrix of order 16 x 16, viewed as a 4 x 4 block matrix, we need to access the main programme and change the limits accordingly by hand. Can it be done in a better way, so that I can write only one programme for any arbitrary dimensions.

Note: My initial programme was more complicated. I just now realised that It can be simplified. So I have modified accordingly.

• Could you add some minimal examples of desired input/output? – Yves Klett Oct 24 '14 at 10:31
• I used this matrix for checking. X = Table[Subscript[x, i, j], {i, 64}, {j, 64}]; Big[X,m,p,q,r,s] This is because of the program which we have written earlier. This makes a 64 x 64 matrix, which can be seen as a 16 x 16 block matrix, where each block is a diagonal matrix. Of course, my actual program is to write a more complicated operation instead of the simple one. This is just a test cast – RSG Oct 24 '14 at 10:49

I am going to demonstrate this with a smaller matrix to enable us to more readily view the results.

A0 = Table[Subscript[a, i, j], {i, 8}, {j, 8}];
A0//MatrixForm


Step 1. Partition A0 into 4x4 blocks

A0byBlock = Partition[A0, {4, 4}];
A0byBlock // MatrixForm


Step 2. Map f onto A0 at level 2

A0out = Map[f, A0byBlock, {2}];
A0out // MatrixForm


Step 3. Repack A0out into 8x8 matrix

A0outF = ArrayFlatten[A0out];
A0outF // MatrixForm


Now these steps can be bundled up into a function. I will define the function name to be big2 as a substitute for your Big (I use small letters for user defined function names)

big2[A0_, f_, nd_] := Module[
{
A0byBlock,
A0Out
},

A0byBlock = Partition[A0, {nd, nd}];
A0Out = Map[f, A0byBlock, {2}];

(* The output is given by flattening A0out *)
ArrayFlatten[A0out]
]


I have tested this on X (your 64x64 matrix) and it worked fine.

• ArrayFlatten[] is quite useful for your third step. – J. M. will be back soon Jul 4 '15 at 0:37
• @Guess who it is - Nice tip. I was unfamiliar with ArrayFlatten. I have modified the answer to use your suggestion. – Jack LaVigne Jul 4 '15 at 11:58