Using python/numpy to create a complex matrix

Question

Using python/numpy, I would like to create a 2D matrix M whose components are:

I know I can do this with a bunch of for loops but is there a better way to do this by using numpy (not using for loops)?

---

This is how I tried, which end up giving me a value error.

I tried to first define a function that takes the sum over k:

define sum_function(i,j):
    initial_array = np.arange(g(i,j),h(i,j)+1)
    applied_array = f(i,j,initial_array)
    return applied_array.sum()

then I tried to create the M matrix with np.mgrid as follows:

ii, jj = np.mgrid(start:fin, start:fin)
M_matrix = sum_function(ii,jj)

--

(Edited) Let me write down the concrete form of a matrix as an example: M_{i,j} = \sum_{k=min(i,j)}^{i+j}\sin{\left( (i+j)^k \right)}

if i,j = 0,1, then this matrix is 2 by 2 and it's form will be \bigl(\begin{smallmatrix} \sin(0) & \sin(1) \ \sin(1)& \sin(2)+\sin(4) \end{smallmatrix}\bigr)

Now if the matrix gets really big, how would I create this matrix without using for loops?

Answer 1

To simplify thinking, lets ravel the i,j dimensions to one, ij dimension. Can we evaluate 3 arrays:

G = g(ij)   # for all ij values
H = h(ij)
F = f(ij, kk)  # for all ij, and all kk

In other words, can g,h,f be evaluated at multiple values, to produce whole-arrays?

If the G and H values were the same for all ij, or subsets (preferably slices), then

F[:, G:H].sum(axis=1)

would be the value for all ij.

If the H-G difference, the size of each slice, was the same, then we can construct a 2d indexing array, GH such that

F[:, GH].sum(axis=1)

In other words we are summing constant size windows of the F rows.

But if the H-G differences vary across ij, I think we are stuck with doing the sum for each ij element separately - with Python level loops, or ones complied with numba or cython.