Cython Boundscheck & nonecheck - PullRequest
Новая
       7

Cython Boundscheck & nonecheck

1 голос
/ 10 апреля 2019

Я несколько дней борюсь со сценарием, дающим неожиданные результаты.Сегодня я просто понимаю, что если бы я использовал функцию cython с или без boundscheck и nonecheck декораторов, я не получил бы те же результаты!

Вот пример: 

import numpy as np
cimport numpy as np
cimport cython

cdef double[4] c
c[0] = 0.1
c[1] = 0.2
c[2] = 0.3
c[3] = 0.4


def cp1(double[:,::1] u, double[:,::1] K, int ixmin, int ixmax, int izmin, int izmax):
    cpc1(u, K, ixmin, ixmax, izmin, izmax)


def cp2(double[:,::1] u, double[:,::1] K, int ixmin, int ixmax, int izmin, int izmax):
    cpc2(u, K, ixmin, ixmax, izmin, izmax)


@cython.boundscheck(False)
@cython.nonecheck(False)
cdef void cpc1(double[:,::1] u, double[:,::1] K, int ixmin, int ixmax, int izmin, int izmax) nogil:

    cdef Py_ssize_t ix, iz
    cdef double dpu, dmu

    for ix in range(ixmin+2, ixmax-1):
        for iz in range(izmin, izmax):
            dpu = c[0]*u[ix-1, iz] + c[1]*u[ix, iz] + c[2]*u[ix+1, iz]
            dmu = c[1]*u[ix-1, iz] + c[2]*u[ix, iz] + c[3]*u[ix+1, iz]
            K[ix, iz] = 0.5*dpu - 0.5*dmu

@cython.boundscheck(True)
@cython.nonecheck(True)
cdef void cpc2(double[:,::1] u, double[:,::1] K, int ixmin, int ixmax, int izmin, int izmax) nogil:

    cdef Py_ssize_t ix, iz
    cdef double dpu, dmu

    for ix in range(ixmin+2, ixmax-1):
        for iz in range(izmin, izmax):
            dpu = c[0]*u[ix-1, iz] + c[1]*u[ix, iz] + c[2]*u[ix+1, iz]
            dmu = c[1]*u[ix-1, iz] + c[2]*u[ix, iz] + c[3]*u[ix+1, iz]
            K[ix, iz] = 0.5*dpu - 0.5*dmu

Если я запускаю эти строки: 

u = np.random.rand(256, 256)
K1 = np.zeros_like(u)
K2 = np.zeros_like(u)

cp1(u, K1, 100, 150, 100, 150)
cp2(u, K2, 100, 150, 100, 150)

, инструкция np.all(K1 == K2) возвращает False.Разница между двумя массивами близка к точности станка (около 5e-17), но использования этой функции в тысячу раз достаточно, чтобы получить большие различия в конечных результатах.

Теперь, если удалить nogil инструкции в cpc1 и cpc2, и я заменяю cdef double[4] c на c = np.zeros(4), обе функции cp1 и cp2 возвращают одинаковые результаты.Проблема в том, что я теряю около 50% производительности, используя ndarray вместо c array.

. Думаю, проблема в точности массива c, но почему значения boundscheck и nonecheck имеютвлияние на результаты в этом случае (нет доступа за пределы, нет, нет, ...)

Есть ли способ решить эту проблему?

РЕДАКТИРОВАТЬ

Как подчеркнул ead, если я скомпилирую код без -03 --ffast-math -march=native, оба cp1 и cp2 вернут одинаковые результаты!Но за счет удвоения времени исполнения!Я могу более или менее понять, почему O3 и ffast-math приводят к неожиданным результатам при агрессивной оптимизации, но я не понимаю, почему march=native также нарушает код.

Есть ли способ сохранить производительность и точность?

...