I'd like to create a square upper triangular matrix that is defined as follows for some float c and some dimension N:
[[1 , c , c^2, ... c^N],
[0, 1 , c, ... c^{N-1}],
[0, 0 , 1, ... c^{N-2}],
.
.
.
[0, 0 , 0, .... 1]]
For concreteness, if N=2, then the matrix should be
[[1, c],
[0, 1]]
And if N=3, then the matrix should be:
[[1, c, c^2],
[0, 1, c],
[0, 0, 1]]
How can I do this?
This is a simple way to do that:
import numpy as np
c = 2
n = 5
r = np.arange(n + 1)
p = r - r[:, np.newaxis]
res = np.triu(c ** p.clip(min=0))
print(res)
# [[ 1 2 4 8 16 32]
# [ 0 1 2 4 8 16]
# [ 0 0 1 2 4 8]
# [ 0 0 0 1 2 4]
# [ 0 0 0 0 1 2]
# [ 0 0 0 0 0 1]]
If you want to make a very big matrix and want to save time and memory you could also do this:
import numpy as np
c = 2
n = 5
b = np.zeros(2 * n + 1, a.dtype)
b[n:] = c ** np.arange(n + 1)
s, = b.strides
res = np.lib.stride_tricks.as_strided(b[n:], shape=(n + 1, n + 1), strides=(-s, s),
writeable=False)
print(res)
# [[ 1 2 4 8 16 32]
# [ 0 1 2 4 8 16]
# [ 0 0 1 2 4 8]
# [ 0 0 0 1 2 4]
# [ 0 0 0 0 1 2]
# [ 0 0 0 0 0 1]]