Decompose Commutators

Question

A theorem in this paper¹ states that every integral n-by-n matrix M over the integers with trace M = 0 is a commutator, that means there are two integral matrices A,B of the same size as M such that M = AB - BA.

Challenge

Given an integral matrix M with trace M = 0 find some integral matrices A,B such that M = AB - BA.

Details

• Let A,B be two matrices of compatible size, then AB denotes the matrix product of A and B which in general does not commute, so in general AB = BA does not hold. We can "measure" how close two matrices are to commuting with eachother by measuring how close the commutator - which is defined as AB - BA - is to being zero. (The commutator is sometimes also written as [A, B].)
• The trace of a matrix is the sum of the diagonal entries.
• The decomposition is not necessarily unique. You have to provide one decomposition.

Examples

Since the decomposition is not unique, you might get different results that are still valid. In the following we just consider some possible example.

M:       A:      B:
[ 1  0]  [1  1]  [1  0]
[-2 -1]  [1 -1]  [1  1]

In the trivial case of 1 x 1 matrices only M = 0 is possible, and any two integers A,B will solve the problem

M:   A:   B:
[0]  [x]  [y] for any integers x,y

Note that for M = 0 (a zero matrix of arbitrary size) implies AB = BA, so in this case any two matrices (A, B) that commute (for example if they are inverse to each other) solve the problem.

M:     A:     B:
[0 0]  [2 3]  [3   12]
[0 0]  [1 4]  [4   11]

M:           A:         B:
[11 12  12]  [1  1  1]  [1  2  3]
[ 7  4   0]  [0  1  1]  [4  5  6]
[ 0 -7 -15]  [0  0  1]  [7  8  9]

M:                    A:          B:
[-11811 -9700 -2937]  [ 3 14 15]  [ 2 71 82]
[ -7749 -1098  8051]  [92 65 35]  [81 82 84]
[  3292  7731 12909]  [89 79 32]  [59  4 52]

^{1: "Integral similarity and commutators of integral matrices" by Thomas J.Laffey, Robert Reams, 1994}

Answer 1

05AB1E, 36 29 bytes

ÄàD(Ÿsgãsgãã.ΔD2FsN._`øδ*O}-Q

-7 bytes thanks to @ovs in another answer of mine.

Extremely slow since I'm using a brute-force approach with three cartesian products. But performance is irrelevant for code-golf challenges I guess (would love to see this same challenge as fastest-code or fastest-algorithm).

Try it online (only works for matrices of size 2x2, and even then usually times out if the absolute maximum of the integers in the matrix \$M\$ is too large..)
Verify whether the output \$[A,B]\$ results in the input \$M\$.

Explanation:

Assumes all the integers used in the resulting matrices \$[A,B]\$ will be within the range \$[-\max(\lvert M\rvert), \max(\lvert M\rvert)]\$. If you know any way to prove or disprove this assumption, let me know.
Note that there are of course matrices \$[A, B]\$ with integers larger than the integers of the resulting matrix \$M\$. But this \$M\$ will have other possible \$[A, B]\$ which are within this min-max range. I.e. \$[A,B]=\$ [[[1,2],[3,4]], [[-1,-2],[-3,-3]]] will have \$M=\$ [[0,2],[-3,0]], but using this \$M\$ to find a possible \$[A,B]\$ will result in [[[3,2],[3,3]], [[2,2],[3,3]]], which are all within the range \$[-3,3]\$.

Let me start with an explanation of 2FsN._`øδ*O}-, which will calculate \$M = AB-BA\$ of a given matrix-pair \$[A,B]\$ (which is the same code as in the verifier TIO-link above).

2F       # Loop 2 times:
  s      #  Swap the top two values on the stack
         #  (so the matrix-pair is at the top again after the first iteration)
   N._   #  Rotate the pair of matrices the (0-based) loop-index amount of times
         #  (to reverse the order from [A,B] to [B,A] in the second iteration)
      `  #  Push both matrices separated to the stack
  ø      #  Zip/transpose the top matrix; swapping rows/columns
   δ     #  Apply double-vectorized:
    *    #   Multiply
     O   #  And then sum each inner row together
 }       # After the loop:
  -      # Subtract the values at the same positions in the two matrices from each other

As for the rest of the program to create all possible matrices and brute-force over them:

Ä        # Convert each integer in the (implicit) input-matrix `M` to its absolute value
 à       # Pop and push the flattened maximum of these absolute values of `M`
  D(     # Duplicate it, and negate the maximum
    Ÿ    # Pop both and push a list in the range [max, -max]
s        # Swap to get the input matrix again
 g       # Pop and push its length to get the matrix dimension
  ã      # Create the cartesian product of the [max, -max]-list that many times
         # (i.e. if the input-matrix has a dimension of 3, we'll now have all possible 
         #  triplets of integers within the [max, -max]-range)
sgã      # Do the same, to create all possible combinations of these triplets
         # (so we now have all possible matrices of the same dimensions as the input-matrix,
         #  using integers within the [max, -max]-range)
ã        # Cartesian product yet again (defaults to 2) to create all possible pairs of matrices
.Δ       # Then find the first pair of matrices which is truthy for:
  D      #  Duplicate the current pair
   2FsN._`øδ*O}-
         #  Calculate M as mentioned above
    Q    #  And check if its equal to the (implicit) input-matrix
         # (after which this matrix-pair is output implicitly as result)

Just to give an idea of how inefficient this brute-force approach is:

The given example input-matrix \$M=\$ [[11,12,12],[1,1,1],[1,2,3]] of the challenge description, will create all possible matrix-pairs of the same 3x3 dimensions with integers in the range \$[-12,12]\$. This means the list of matrix-pairs will contain \$(((12^2+1)^3)^3)^2=802830827198685151406498569488525390625\$ matrix-pairs.

And the given example input-matrix \$M=\$ [[-11811,-9700,-2937],[3,14,15],[2,71,82]] would result in \$(((11811^2+1)^3)^3)^2=400239694511528196808960294025932709663435907226315389994060779277222573307615163703296535923733538440194289670070277830289553752923910908747710464\$ matrix-pairs..