I'll add a more detailed explanation later, but here are all the intermediate forms of 0223202330203313021301011023230323 (skipping stages that don't perform any changes):

0223202330203313021301011023230323 32221321312323323 322213213i3323323 322213213i33i13i1 322-i33-i33i33i13i1 -i12-i33-i33i33i13i1 -i12-i33-i33i33i-i2i1 -i1-2i3-3i33i33-ii2i1 -i1-i23-i33i333-iii21 -i1-i23-i33i333--i21 -i1-i23-ii3i21 -i1-ii1-ii3i21 -i1-ii1-ii3i-i3 -i-1ii-1ii3-ii3 -i-i1i-i1i3-ii3 -i-i1i-i1i3--3 -i-i1i-i1i33 --ii1-ii1i33 --ii1-iii133 ---1--i133 -1i1 -i11 -i -i0 

Explanation

Here is what the individual stages do:

ii - 

Combines all pairs of i into - to reduce the phase characters.

+`(.)\1|0 <empty> 

Now if there are two consecutive identical characters left, it's either -- or two identical matrices. In either case, multiplying them gives the identity. But we don't need identities, so we just remove all of them, and the explicit identities (the 0s) as well. This stage is repeated in itself with + until the result stops changing. This ensures that things like 123321 get resolved completely, such that the next step can assume that all pairs of digits are distinct.

(.)-|(\d)(\d) -$1$3$2 

This is actually two separate transformations in one (for golfitude). Note that if the first alternative matches, $2 and $3 are empty, and if the second one matches $1 is empty. So this can be decomposed into these two steps:

(\d)(\d) -$2$1 

This just swaps all pairs of digits and adds a minus sign. Since we removed all 0s and all identical pairs, this will only match 12, 23, 31, 21, 32, 13. This step may seem weird, but it allows me to only check for half of these cases later on, because the ones I can't process then will be swapped here in the next iteration.

The other part of the above stage was:

(.)- -$1 

This gradually moves - signs all the way to the left (one position per iteration). I do this such that ultimately they are all next to each other and get resolved in the earlier step.

12 i3 23 i1 31 i2 

These three stages now simply resolve the three pairs of products. As I said above, this will only catch half of the relevant cases, but the other half will be taken care of in the next iteration, after the previous step swapped all pairs.

)`(\d)i i$1 

This is the last stage of the loop. It's similar to the one that shifts - to the left, except for i. The main difference is that this one swaps i only with digits. If I used (.)i then in cases where I get a -i or i- the two would be swapped indefinitely and the program wouldn't terminate. So this only swaps them to the right of the - signs. This is sufficient - as long as all - and i appear together at some point, they can be resolved correctly.

^\D*$ $&0 

The final step (outside of the loop). Remember that we always deleted all identities, so if the result is actually the identity (times a phase), then we won't have the required digit in the output any more, so we add it back.

As an example, here are all the intermediate forms of 0223202330203313021301011023230323 (skipping stages that don't perform any changes):

0223202330203313021301011023230323 321321312 # Remove identities -23-31-12-132 # Swap all pairs -23-31-i3-132 # Resolve 12 -i1-31-i3-132 # Resolve 23 -i1-i2-i3-132 # Resolve 31 -i-1i-2i-3-312 # Move - to the left and swap pairs -i-1i-2i-3-3i3 # Resolve 12 -i-i1-i2-3-i33 # Move i to the left -i-i1-i2-3-i # Remove identities --ii-1i-2-3i # Move - to the left --ii-i1-2-i3 # Move i to the left ----i1-2-i3 # Resolve ii i1-2-i3 # Remove identities i-1-2i3 # Move - to the left i-1-i23 # Move i to the left -i-1i-32 # Move - to the left and swap pairs -i-i1-32 # Move i to the left --ii-1-23 # Move - to the left and swap pairs --ii-1-i1 # Resolve 23 ----1-i1 # Resolve ii 1-i1 # Remove identities -1i1 # Move - to the left -i11 # Move i to the left -i # Remove identities. Now the loop can't change this any longer. -i0 # Fix the result by adding in the 0.