def heisenberg_mult(v1, v2, p): i1, j1, k1 = v1 i2, j2, k2 = v2 i3 = (i1 + i2)%p j3 = (j1 + j2)%p k3 = (-i2*j1 + k1 + k2)%p return(i3, j3, k3) def heisenberg_inv(v1, p): i1, j1, k1 = v1 i2 = p-i1 j2 = p-j1 k2 = p-k1 k2 = (k2 - i2*j2)%p return (i2, j2, k2) def heisenberg_power(v1, n, p): i1, j1, k1 = v1 if n < 0: return heisenberg_inv(heisenberg_power((i1, j1, k1), -n, p), p) if n == 0: return (0, 0, 0) return heisenberg_mult((i1, j1, k1), heisenberg_power((i1, j1, k1), n-1, p), p)