Cómo explicar la matriz contigua a la columna al ampliar numpy con C

Question

Tengo una función C para normalizar las filas de una matriz en el espacio de registro (esto evita el flujo inferior numérico).

El prototipo de mi C-función es la siguiente:

void normalize_logspace_matrix(size_t nrow, size_t ncol, double* mat); 

Se puede ver que se necesita un puntero a un array y lo modifica en su lugar. Por supuesto, el código C supone que los datos se guardan como una matriz contigua a C, es decir, contigua a la fila.

envuelvo la función como sigue usando Cython (importaciones y cdef extern from omitidas):

def normalize_logspace(np.ndarray[np.double_t, ndim=2] mat): 
    cdef Py_ssize_t n, d 
    n = mat.shape[0] 
    d = mat.shape[1] 
    normalize_logspace_matrix(n, d, <double*> mat.data) 
    return mat 

mayoría de las veces numpy-arrays son fila contigua y la función funciona bien. Sin embargo, si una matriz numpy se ha transpuesto previamente, los datos no se copian, sino que se devuelve una nueva vista de los datos. En este caso, mi función falla porque la matriz ya no es fila contigua.

puedo evitar esto mediante la definición de la matriz para tener orden Fortran contiguos, de tal manera que después de la transposición será C contiguas:

A = np.array([some_func(d) for d in range(D)], order='F').T 
A = normalize_logspace(A) 

Obviamente eso es muy propenso a errores y el usuario tiene que tomar se preocupa de que la matriz esté en el orden correcto, que es algo que el usuario no debería tener que preocuparse en Python.

¿Cuál es la mejor manera de cómo puedo hacer que esto funcione tanto con matrices contiguas como con filas y columnas? Supongo que algún tipo de comprobación de orden de matriz en Cython es el camino a seguir. Por supuesto, preferiría una solución que no requiera copiar los datos en una nueva matriz, pero casi asumo que es necesario.

Answer 1

Si desea admitir matrices en orden C y Fortran sin tener que copiar, su función C debe ser lo suficientemente flexible como para admitir ambas ordenaciones. Esto puede lograrse haciendo pasar los avances de la matriz NumPy a la función C: Cambiar el prototipo de

void normalize_logspace_matrix(size_t nrow, size_t ncol, 
           size_t rowstride, size_t colstride, 
           double* mat); 

y la llamada Cython a

def normalize_logspace(np.ndarray[np.double_t, ndim=2] mat): 
    cdef Py_ssize_t n, d, rowstride, colstride 
    n = mat.shape[0] 
    d = mat.shape[1] 
    rowstride = mat.strides[0] // mat.itemsize 
    colstride = mat.strides[1] // mat.itemsize 
    normalize_logspace_matrix(n, d, rowstride, colstride, <double*> mat.data) 
    return mat 

A continuación, reemplazar cada aparición de mat[row*ncol + col] en su C código por mat[row*rowstride + col*colstride].

Answer 2

En este caso realmente desea crear una copia de la matriz de entrada (que puede ser un punto de vista sobre un verdadero arsenal) con el orden de las filas contiguas garantizada. Usted puede lograr esto con algo como esto:

a = numpy.array(A, copy=True, order='C') 

Además, considere echar un vistazo a la exacta array interface de Numpy (no es parte C también).

Answer 3

1 a Sven, cuya respuesta resuelve el Gotcha (bueno, me consiguió) que dstack devuelve una matriz F_contiguous?!

# don't use dstack to stack a,a,a -> rgb for a C func 

import sys 
import numpy as np 

h = 2 
w = 4 
dim = 3 
exec("\n".join(sys.argv[1:])) # run this.py h= ... 

a = np.arange(h*w, dtype=np.uint8) .reshape((h,w)) 
rgb = np.empty((h,w,dim), dtype=np.uint8) 
rgb[:,:,0] = rgb[:,:,1] = rgb[:,:,2] = a 
print "rgb:", rgb 
print "rgb.flags:", rgb.flags # C_contiguous 
print "rgb.strides:", rgb.strides # (12, 3, 1) 

dstack = np.dstack((a, a, a)) 
print "dstack:", dstack 
print "dstack.flags:", dstack.flags # F_contiguous 
print "dstack.strides:", dstack.strides # (1, 2, 8)