How to declare lists of lists in cython

Question

I have the following .pyx code:

import cython
@cython.boundscheck(False)
@cython.cdivision(True)
@cython.wraparound(False)
def f(m):
  cdef int n = len(m)/2
  cdef int j, k
  z = [[0]*(n+1) for _ in range(n*(2*n-1))]
  for j in range(1, 2*n):
    for k in range(j):
      z[j*(j-1)/2+k][0] = m[j][k]
  return solve(z, 2*n, 1, [1] + [0]*n, n)


cdef solve(b, int s, int w, g, int n):
  cdef complex h
  cdef int u,v,j,k
  if s == 0:
    return w*g[n]
  c = [b[(j+1)*(j+2)/2+k+2][:] for j in range(1, s-2) for k in range(j)]
  h = solve(c, s-2, -w, g, n)
  e = g[:]
  for u in range(n):
    for v in range(n-u):
      e[u+v+1] += g[u]*b[0][v]
  for j in range(1, s-2):
    for k in range(j):
      for u in range(n):
        for v in range(n-u):
          c[j*(j-1)/2+k][u+v+1] += b[(j+1)*(j+2)/2][u]*b[(k+1)*(k+2)/2+1][v] + b[(k+1)*(k+2)/2][u]*b[(j+1)*(j+2)/2+1][v]
  return h + solve(c, s-2, w, e, n)

I don't know how to declare the list of lists in cython to speed the code up.

For example, the variable m is a matrix represented as a list of lists of floating point numbers. The variable z is a list of lists of floating point numbers too. What should the line def f(m) look like for example?

---

Following the advice in the answer by @DavidW here is my latest version.

import cython
import numpy as np
def f(complex[:,:] m):
  cdef int n = len(m)/2
  cdef int j, k
  cdef complex[:,:] z = np.zeros((n*(2*n-1), n+1), dtype = complex)
  for j in range(1, 2*n):
    for k in range(j):
      z[j*(j-1)/2+k, 0] = m[j, k]
  return solve(z, 2*n, 1, [1] + [0]*n, n)


cdef solve(complex[:,:] b, int s, int w, g, int n):
  cdef complex h
  cdef int u,v,j,k
  cdef complex[:,:] c
  if s == 0:
    return w*g[n]
  c = [b[(j+1)*(j+2)/2+k+2][:] for j in range(1, s-2) for k in range(j)]
  print("c stats:", len(c), [len(c[i]) for i in len(c)]) 
  h = solve(c, s-2, -w, g, n)
  e = g[:]
  for u in range(n):
    for v in range(n-u):
      e[u+v+1] = e[u+v+1] + g[u]*b[0][v]
  for j in range(1, s-2):
    for k in range(j):
      for u in range(n):
        for v in range(n-u):
          c[j*(j-1)/2+k][u+v+1] = c[j*(j-1)/2+k][u+v+1] + b[(j+1)*(j+2)/2][u]*b[(k+1)*(k+2)/2+1][v] + b[(k+1)*(k+2)/2][u]*b[(j+1)*(j+2)/2+1][v]
  return h + solve(c, s-2, w, e, n)

The main problem now is how to declare c as it is currently a list of lists.

Answer 1

A list-of-lists is not a structure that can get much speed-up from Cython. The structure you should use is a 2D typed memoryview:

def f(double[:,:] m):
    # ...

These are indexed as m[j,k] rather than m[j][k]. You can pass them any suitably shaped object that exposes the Python buffer protocol. Most frequency this will be a Numpy array.

---

You should also avoid using decorators like @cython.boundscheck(False) and @cython.wraparound(False) unless you understand what they do and have considered whether they're appropriate for your function. For your current version (where you're using lists) they actually do nothing and suggest you have copied them without understanding. They do speed up the indexing of memoryviews (at the cost of some safety).

---

Edit: In terms of initializing c you have two options.

1. Initialize a numpy array with the list of lists. This probably isn't hugely quick (but if other steps are slower then it may not matter):

   c = np.array([b[(j+1)*(j+2)/2+k+2,:] for j in range(1, s-2) for k in range(j)], dtype=complex)
   # note that I've changed the indexing of b slightly
   

2. Set up c with an appropriately sized np.zeros array. Swap the list comprehension for two loops. It isn't 100% obvious to me exactly what this translates to, but it's something like

   c = np.zeros("some size you'll have to work out",dtype=complex)
   for k in range(j):
       for j in range(1,s-2):
           c["some function of j and k",:] = b["some function of j and k",:] 
   

You'll also want to replace len(c) with c.shape[0], etc.

Answer 2

Like C:

cdef float z[100][100]

For more detail please refer to this link.