# Transformation invariant subspace of three-component array

I have a large sparse array A with dimensions $$m \times n \times q$$. I also have two transformation matrices, B ($$n \times n$$) and C ($$q \times q$$) from which I act from the left on A. So, for instance, I might use

TensorTranspose[
B . TensorTranspose[
C . TensorTranspose[A, {2, 3, 1}],
{2, 3, 1}
],
{2, 3, 1}
]


where the final tensor transpose is there to put the indices back in the correct order at the end.

Now, A has the property such that $$m$$ 'rows' (i.e. if I flatten A into an m by n*q matrix) of the transformed matrix are always linear combinations of the original rows. I want to find the subspace of the space spanned by the m 'rows' of the original matrix which are actually invariant under the transformation which I have described above.

My procedure has been as follows:

1. Perform the transformation as described above

2. Flatten the original A and the transformed A into m by n*q matrices, and use LinearSolve to identify the image of each transformed row as a linear combination of the original rows.

3. Form an m by m matrix X using the results of 2, which acts on A from the left and captures the original combined action of B and C on A.

4. Take Z = NullSpace[X - IdentityMatrix[{m,m}]]; this acts as a projector so that Z.A gives the subspace I wanted.

The problem is that the use of LinearSolve in 2) gets extremely slow (these matrices are generally very large); I feel it should be possible to carry out this process directly using just the matrices B and C, since they themselves are sufficient to transform A, but I can't wrap my head around the required linear algebra.

• Are you sure that this is a Mathematica question though? As you said it sounds like the problem is with finding a new linear algebra approach. Perhaps math.SE? Jun 2 at 13:15
• I think the main complexity is to do with handling the various flattenings of the matrices -- the basic linear algebra idea would surely be the same as for the method I describe but it's specifically the implementation of this within Mathematica which is the issue. Jun 2 at 13:18

I do not fully understand what you mean by "the m rows". Here is as far as I can follow you:

You want to calculate:

This can be done using Transpose and Do. As in MMA you should not use capitalized variable names, I am using ma, mb and mc. Note that the Dot product contracts the last index of the first tensor with the first index of the second tensor. Therefore, in a first step we exchange the first and second index of ma to multiply B with A. In a second step we then exchange index 1 and 3 of the result to get C times B times A. Here is an explicit example with m==n==q==2:

The matrices are:

(ma = Array[Subscript[a, #1, #2, #3] &, {2, 2, 2}]) // MatrixForm
(mb = Array[Subscript[b, #1, #2] &, {2, 2}]) // MatrixForm
(mc = Array[Subscript[c, #1, #2] &, {2, 2}]) // MatrixForm


Now we make the different multiplications:

(tmp = mb . Transpose[ma, {3, 2, 1}]) // MatrixForm
(tmp = mc . Transpose[tmp, {2, 1, 3}]) // MatrixForm


Finally we must restore the original order of the indices:

(tmp = Transpose[tmp, {3, 1, 2}]) // MatrixForm


If we now flatten the first level of ma and the result, we get:

Flatten[ma, 1]
Flatten[tmp, 1]


Or if we reshape them as:

ArrayReshape[ma, {2, 2}]
ArrayReshape[tmp, {2, 2}]


Unfortunately, I do not understand what should now be a linear combination of what. Please clarify this.

• Thanks for your help so far. By "the m rows", I mean the m rows of the m by nq matrix one gets if flattening array A appropriately i.e. ArrayReshape[A,{m,nq}]. That is, each such row of the transformed A will (this is a property of my array) be a linear combination of rows of the original A. I want to find the subspace spanned by the rows of the original A which are actually invariant under the transformation. Jun 2 at 14:29
• If your matrix has some has a special form, you must describe this accurately. Did you mean e.g. that {Subscript[c, 1, 1] (Subscript[b, 1, 1] Subscript[a, 1, 1, 1] + Subscript[b, 1, 2] Subscript[a, 1, 1, 2]) + Subscript[c, 1, 2] (Subscript[b, 1, 1] Subscript[a, 1, 2, 1] + Subscript[b, 1, 2] Subscript[a, 1, 2, 2]), Subscript[c, 2, 1] (Subscript[b, 1, 1] Subscript[a, 1, 1, 1] + Subscript[b, 1, 2] Subscript[a, 1, 1, 2]) + Subscript[c, 2, 2] (Subscript[b, 1, 1] Subscript[a, 1, 2, 1] + Subscript[b, 1, 2] Subscript[a, 1, 2, 2])} Jun 2 at 16:02
• is a linear combination of {Subscript[a, 1, 1, 1], Subscript[a, 1, 1, 2]} and {Subscript[a, 1, 2, 1], Subscript[a, 1, 2, 2]} Jun 2 at 16:02
• Yes, this is what I mean by "A has the property such that m 'rows' (i.e. if I flatten A into an m by n*q matrix) of the transformed matrix are always linear combinations of the original rows" in my question; I apologise if this was unclear. I am happy with how to transform A in the way I want/describe-- indeed in my question I give an example of code using TensorTranspose and dot that does this -- but it's specifically the matter of finding the invariant subspace spanned by the 'rows' (in the flattened sense) which puzzles me. Jun 2 at 16:14

I found a solution to this in the end.

Transform A using B and C, and then flatten both A and the transformed A into m by (n*q) matrices. Applying NullSpace to this matrix gives a projector of A onto the invariant subspace.

To transform A, use

ATrans = TensorTranspose[C.TensorTranspose[B.TensorTranspose[A,{3,2,1}],{3,2,1}],{3,2,1}].


Then, flatten both A and ATrans in their second two arguments:

AFlat = ArrayReshape[A,{m,n*q}]
ATransFlat = ArrayReshape[ATrans,{m,n*q}]


Then take the null space of their difference to get the projector onto the invariant subspace.

projector = NullSpace[AFlat-ATransFlat]


Dot this with A to get the invariant subspace.

projector.A

• Including the code you used and the results you obtained would improve this answer. Jun 2 at 19:45