Most programming languages have a built-in implementation of exponentiation for integers and reals only.



Task

Demonstrate how to implement matrix exponentiation as an operator.




F matrix_mul(m1, m2)
   assert(m1[0].len == m2.len)
   V r = [[0] * m2[0].len] * m1.len
   L(j) 0 .< m1.len
      L(i) 0 .< m2[0].len
         V s = 0
         L(k) 0 .< m2.len
            s += m1[j][k] * m2[k][i]
         r[j][i] = s
   R r
 
F identity(size)
   V rsize = 0 .< size
   R rsize.map(j -> @rsize.map(i -> Int(i == @j)))
 
F matrixExp(m, pow)
   assert(pow >= 0 & Int(pow) == pow, ‘Only non-negative, integer powers allowed’)
   V accumulator = identity(m.len)
   L(i) 0 .< pow
      accumulator = matrix_mul(accumulator, m)
   R accumulator
 
F printtable(data)
   L(row) data
      print(row.map(cell -> ‘#<5’.format(cell)).join(‘ ’))
 
V m = [[3, 2], [2, 1]]
L(i) 5
   print("\n#.:".format(i))
   printtable(matrixExp(m, i))
 
print("\n10:")
printtable(matrixExp(m, 10))