(Numpy C API) Итерации по одному массиву: NpyIter против цикла (с PyArray_DATA) - PullRequest
Новая
       46

(Numpy C API) Итерации по одному массиву: NpyIter против цикла (с PyArray_DATA)

0 голосов
/ 07 января 2019

Я пишу некоторый код расширения C для модуля Python. Функция, которую я хочу написать: (на python) 

output = 1./(1. + input)

, где input - массив пустых фигур любой формы.

Первоначально я использовал NpyIter_MultiNew: 

static PyObject *
helper_calc1(PyObject *self, PyObject *args){

    PyObject * input;
    PyObject * output = NULL;

    if (!PyArg_ParseTuple(args, "O", &input)){
        return NULL;
    }

    // -- input -----------------------------------------------
    PyArrayObject * in_arr;

    in_arr = (PyArrayObject *) PyArray_FROM_OTF(input, NPY_DOUBLE, NPY_ARRAY_IN_ARRAY);
    if (in_arr == NULL){
        goto fail;
    }

    // -- set up iterator -------------------------------------
    PyArrayObject * op[2];
    npy_uint32 op_flags[2];
    npy_uint32 flags;

    op[0] = in_arr;
    op_flags[0] = NPY_ITER_READONLY;

    op[1] = NULL;
    op_flags[1] = NPY_ITER_WRITEONLY | NPY_ITER_ALLOCATE;

    flags = NPY_ITER_EXTERNAL_LOOP | NPY_ITER_BUFFERED | NPY_ITER_GROWINNER;

    NpyIter * iter = NpyIter_MultiNew(2, op,
                            flags, NPY_KEEPORDER, NPY_NO_CASTING,
                            op_flags, NULL);

    if (iter == NULL){
        goto fail;
    };

    NpyIter_IterNextFunc * iternext = NpyIter_GetIterNext(iter, NULL);
    if (iternext == NULL){
        NpyIter_Deallocate(iter);
        goto fail;
    };

    // -- iterate ---------------------------------------------
    npy_intp count;
    char ** dataptr = NpyIter_GetDataPtrArray(iter);
    npy_intp * strideptr = NpyIter_GetInnerStrideArray(iter);
    npy_intp * innersizeptr = NpyIter_GetInnerLoopSizePtr(iter);

    do {
        count = *innersizeptr;

        while (count--){

            *(double *) dataptr[1] = 1. / (1. + *(double *)dataptr[0]);

            dataptr[0] += strideptr[0];
            dataptr[1] += strideptr[1];
        }

    } while (iternext(iter));

    output = NpyIter_GetOperandArray(iter)[1];

    if (NpyIter_Deallocate(iter) != NPY_SUCCEED){
        goto fail;
    }

    Py_DECREF(in_arr);

    return output;

    fail:
        Py_XDECREF(in_arr);
        Py_XDECREF(output);
        return NULL;
}

Однако, так как это всего лишь один массив (т.е. мне не нужно беспокоиться о трансляции нескольких массивов), Есть ли какая-то причина, по которой я не могу перебирать данные самостоятельно, используя PyArray_DATA, цикл for и размер массива? 

static PyObject *
helper_calc2(PyObject *self, PyObject *args){

    PyObject * input;
    PyObject * output = NULL;

    if (!PyArg_ParseTuple(args, "O", & in)){
        return NULL;
    }

    // -- input -----------------------------------------------
    PyArrayObject * in_arr;

    in_arr = (PyArrayObject *) PyArray_FROM_OTF(input, NPY_DOUBLE, NPY_ARRAY_IN_ARRAY);
    if (in_arr == NULL){
        Py_XDECREF(in_arr);
        return NULL;
    }

    int ndim = PyArray_NDIM(in_arr);
    npy_intp * shape = PyArray_DIMS(in_arr);
    int size = (int) PyArray_SIZE(in_arr);

    double * in_data = (double *) PyArray_DATA(in_arr);

    output =  PyArray_SimpleNew(ndim, shape, NPY_DOUBLE);
    double * out_data = (double *) PyArray_DATA((PyArrayObject *) output);

    for (int i = 0; i < size; i++){
        out_data[i] = 1. / (1. + in_data[i]);
    }

    Py_DECREF(in_arr);
    return output;

fail:
    Py_XDECREF(in_arr);
    Py_XDECREF(output);
    return NULL;
}

Эта вторая версия работает быстрее, а код короче.

Есть ли какие-либо опасности, на которые мне нужно обратить внимание при использовании, PyArray_DATA с петлей for вместо NpyIter_MultiNew?

Из PyArray_DATA документации:

Если вы не гарантировали непрерывный и / или выровненный массив, убедитесь, что вы понимаете, как получить доступ к данным в массиве, чтобы избежать проблем с памятью и / или выравниванием.

Но я считаю, что об этом позаботятся PyArray_FROM_OTF с флагом NPY_ARRAY_IN_ARRAY.

...