Search code examples
matrixmathematical-optimizationfractals

Directly compute element in recursion using bitwise XOR operator


Let's have a matrix M X N matrix A[i][j] given partially as :(starting with row=0 column =0):

1) for all 1<=i<=N A[0][i]=0

2) for all 0<=j<=M A[j][0]=1

The matrix constructs A[i][j] further as:

for all 1<=i<=N and 1<=j<=M, A[i][j]=A[i-1][j]^A[i][j-1], where ^ denotes the bitwise XOR operation.

If I want to get some value A[i][j], how can I get that directly (without actually calculating A[i][j] for all 1...j-1 and 1...i-1)? Is there a pattern? Given M and N are too large.


Solution

  • Let's plot out the first 16 rows and columns of the matrix A:

    1000000000000000
    1111111111111111
    1010101010101010
    1100110011001100
    1000100010001000
    1111000011110000
    1010000010100000
    1100000011000000
    1000000010000000
    1111111100000000
    1010101000000000
    1100110000000000
    1000100000000000
    1111000000000000
    1010000000000000
    1100000000000000
    

    The answer to your question "is there a pattern" is pretty clearly "Yes!" The 8x8 submatrix in the top left is copied directly below itself and also directly to the right, except the copy to the right has a 0 in its top-left corner instead of a 1. The bottom-right 8x8 submatrix is all 0's, except it has a 1 in the top-left corner. This exact same pattern emerges if we study just the first 8 rows and columns:

    10000000
    11111111
    10101010
    11001100
    10001000
    11110000
    10100000
    11000000
    

    The top-left 4x4 submatrix is copied directly below itself and also directly to the right, except the copy to the right has a 0 in its top-left corner instead of a 1. The bottom-right 4x4 submatrix is all 0's, except it has a 1 in the top-left corner.

    This recursive self-similarity makes the matrix A a fractal, quite similar to the Sierpinski triangle.

    The recursive self-similarity also enables us to easily compute A[i][j] using the binary representations of i and j. Let B be the highest-order bit set in the binary representation of either i or j. Then the following procedure will return the correct value of A[i][j]:

    • for b = B to 0:
      • if j=0 when limited to bits 0 through b, then return 1 (we are on the first column, which is all 1's)
      • else if i=0 when limited to bits 0 through b, then return 0 (we are on the first row)
      • else if i and j both have a 1 stored in bit b, then return 0 unless all remaining bits in each are 0, in which case return 1
      • else continue (we recurse to a smaller sub-matrix)

    This algorithm runs in O(log(max(i,j))) runtime, which is much faster than the O(ij) runtime of the naive approach.

    Let's see this algorithm in action with a few examples:

    • A[9][9] = 0 from the matrix above. In binary representation, i = 1001 and j = 1001. Both have a 1 in the most significant bit and have some bits set after that one, so by the rules above we return 0.

    • A[9][3] = 1 from the matrix above. In binary representation, i = 1001 and j = 0011. In the leftmost (most significant bit), i has a 1 and j has a 0. Therefore we proceed to the next bit (recurse), where both have a 0. We again proceed to the next bit, where i has a 0 and j has a 1. We therefore proceed to the last bit. Both have a 1 in the bit and all following bits are 0's (this is trivially true since there are no following bits), so we return 1.

    • A[9][4] = 1 from the matrix above. In binary representation, i = 1001 and j = 0100. In the leftmost (most significant bit), i has a 1 and j has a 0. Therefore we proceed to the next bit (recurse), where i has a 0 and j has a 1. We again proceed to the next bit, where both have a 0. j has a 0 in this and all following bits, so we return 1.

    One interesting implication of this algorithm is that each row k (for k > 0) has a repeating pattern with period 2^B, where B is the most significant bit of the row number's binary representation. For instance, row 3 (binary representation 11) repeats 1100, while row 9 (binary representation 1001) repeats 1111111100000000. This means the full repeating series of row k > 0 can be represented in O(k) storage and can be computed in O(k log k) time by separately computing A[k,j] for j = 0, 1, ..., 2^ceiling(log_2 k)-1.