Выход функции cython немного отличается от вывода функции python

Question

Я преобразовал функцию python в эквивалент cython, добавив типы к некоторым переменным. Однако функция cython создает несколько отличный результат, чем исходная функция python.

Я узнал некоторые из причин этой разницы в этом посте Cython: unsigned int indices for numpy arrays gives different result Но даже с тем, что я узнал в этом посте я до сих пор не могу получить функцию Cython производить те же результаты, что и питон один ,

Итак, я собрал 4 функции, иллюстрирующие то, что я пробовал. Может ли кто-нибудь помочь раскрывать, почему я получаю несколько разные результаты для каждой функции? и как получить функцию cython, которая возвращает те же точные значения, что и function1? Я делаю некоторые комментарии ниже:

%%cython 
import numpy as np 
cimport numpy as np  

def function1(response, max_loc):  
    x, y = int(max_loc[0]), int(max_loc[1]) 

    tmp1 = (response[y,x+1] - response[y,x-1])/2*(response[y,x] - min(response[y,x-1], response[y,x+1])) 
    tmp2 = (response[y,x+1] - response[y,x-1]) 
    tmp3 = 2*(response[y,x] - min(response[y,x-1], response[y,x+1])) 

    print tmp1, tmp2, tmp3   
    return tmp1, tmp2, tmp3 

cpdef function2(np.ndarray[np.float32_t, ndim=2] response, np.ndarray[np.float64_t, ndim=1] max_loc): 
    cdef unsigned int x, y 
    x, y = int(max_loc[0]), int(max_loc[1]) 

    tmp1 = (response[y,x+1] - response[y,x-1])/2*(response[y,x] - min(response[y,x-1], response[y,x+1]))   
    tmp2 = (response[y,x+1] - response[y,x-1]) 
    tmp3 = 2*(response[y,x] - min(response[y,x-1], response[y,x+1]))  

    print tmp1, tmp2, tmp3   
    return tmp1, tmp2, tmp3 


cpdef function3(np.ndarray[np.float32_t, ndim=2] response, np.ndarray[np.float64_t, ndim=1] max_loc):  
    cdef unsigned int x, y 
    x, y = int(max_loc[0]), int(max_loc[1]) 

    cdef np.float32_t tmp1, tmp2, tmp3 
    cdef np.float32_t r1 =response[y,x+1] 
    cdef np.float32_t r2 =response[y,x-1] 
    cdef np.float32_t r3 =response[y,x] 
    cdef np.float32_t r4 =response[y,x-1] 
    cdef np.float32_t r5 =response[y,x+1]  

    tmp1 = (r1 - r2)/2*(r3 - min(r4, r5)) 
    tmp2 = (r1 - r2) 
    tmp3 = 2*(r3 - min(r4, r5)) 

    print tmp1, tmp2, tmp3   
    return tmp1, tmp2, tmp3 

def function4(response, max_loc):  
    x, y = int(max_loc[0]), int(max_loc[1]) 

    tmp1 = (float(response[y,x+1]) - response[y,x-1])/2*(float(response[y,x]) - min(response[y,x-1], response[y,x+1])) 
    tmp2 = (float(response[y,x+1]) - response[y,x-1]) 
    tmp3 = 2*(float(response[y,x]) - min(response[y,x-1], response[y,x+1])) 

    print tmp1, tmp2, tmp3   
    return tmp1, tmp2, tmp3 

max_loc = np.asarray([ 15., 25.], dtype=np.float64) 
response = np.zeros((49,49), dtype=np.float32)  
x, y = int(max_loc[0]), int(max_loc[1]) 

response[y,x] = 0.959878861904 
response[y,x-1] = 0.438348740339 
response[y,x+1] = 0.753262758255 

result1 = function1(response, max_loc) 
result2 = function2(response, max_loc) 
result3 = function3(response, max_loc) 
result4 = function4(response, max_loc) 
print result1 
print result2 
print result3 
print result4 

и результаты:

0.0821185777156 0.314914 1.04306030273 
0.082118573023 0.314914017916 1.04306024313 
0.0821185708046 0.314914017916 1.04306030273 
0.082118573023 0.314914017916 1.04306024313 
(0.082118577715618812, 0.31491402, 1.043060302734375) 
(0.08211857302303427, 0.3149140179157257, 1.0430602431297302) 
(0.08211857080459595, 0.3149140179157257, 1.043060302734375) 
(0.082118573023034269, 0.31491401791572571, 1.0430602431297302) 

function1 представляет операции, которые я сделал в моей первоначальной функции питона. Результат tmp1.

function2 - моя первая версия cython, которая дает несколько разные результаты. По-видимому, если массив ответов индексируется с типизированной переменной, unsigned int или int, результат принуждается к двойному (используя PyFloat_FromDouble), даже если тип массива - np.float32_t. Но если массив индексируется с помощью python int, вместо него используется функция PyObject_GetItem, и я получаю np.float32_t, что и происходит в функции1. Таким образом, выражения в функции 1 вычисляются с использованием np.float32_t операндов, тогда как выражения в функции2 вычисляются с использованием удвоений. Я получаю немного другой отпечаток, чем в функции1.

function3 - моя вторая попытка cython, пытающаяся получить тот же результат, что и функция1. Здесь я использую unsigned int индексы для доступа к ответу массива, но результаты остаются на промежуточных переменных np.float32_t, которые затем используются в вычислении. Я получаю немного другой результат. Очевидно, оператор print будет использовать PyFloat_FromDouble, поэтому он не сможет напечатать файл np.float32_t.

Затем я попытался изменить функцию python в соответствии с cython. function4 пытается достичь этого путем преобразования в поплавок по меньшей мере одного операнда в каждом выражении, чтобы остальные операнды были принудительно применены к float python, который является двойным в цитоне, а выражение вычисляется с удвоением, как в функции2. Отпечатки внутри функции точно такие же, как функция2, но возвращаемые значения немного отличаются ?!

Answer 1

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

Для уточнения вашего описания function2, если индекс с объектом, Cython использует PyObject_GetItem получить np.float32 скалярный объект (не np.float32_t, который просто ЬурейиЙ для C float). Если вы указали непосредственно в буфер, а Cython нужен объект, он вызывает PyFloat_FromDouble. Ему нужны объекты для назначения tmp1, tmp2 и tmp3, так как они не печатаются.

В function3, с другой стороны, вы ввели tmp переменные, но она по-прежнему необходимо создать float объекты для печати и возвращают результаты. Если вместо того, чтобы использовать NumPy ndarray (смотри ниже), вы не будете иметь эту проблему:

В function1, кстати, вы продвигаете результат np.float64 когда вы делите на 2. Например:

>>> type(np.float32(1)/2) 
<type 'numpy.float64'> 

против

>>> type(np.float32(1)/np.float32(2)) 
<type 'numpy.float32'> 

Даже если вы уверены, что все операции float32 в обеих функциях def и cpdef, конечный результат все еще может варьироваться между двумя в ком сложенный модуль расширения. В следующем примере я проверил, что промежуточные результаты в function1 - все объекты np.float32. В сгенерированном C function2 я проверил, что не было отлито от double (или эквивалентного typedef). Однако эти две функции по-прежнему дают несколько иные результаты. Я, вероятно, должен был бы погрузиться в сборку сборки, чтобы понять, почему, но, возможно, я пропустил что-то простое.

def function1(response, max_loc):  
    tmp = np.zeros(3, dtype=np.float32) 
    x, y = int(max_loc[0]), int(max_loc[1]) 
    tmp[0] = (((response[y,x+1] - response[y,x-1])/np.float32(2)) * 
      (response[y,x] - min(response[y,x-1], response[y,x+1]))) 
    tmp[1] = response[y,x+1] - response[y,x-1] 
    tmp[2] = 2*(response[y,x] - min(response[y,x-1], response[y,x+1])) 

    print tmp[0], tmp[1], tmp[2] 
    return tmp 

cpdef function2(np.ndarray[np.float32_t, ndim=2] response, max_loc): 
    cdef np.ndarray[np.float32_t, ndim=1] tmp = np.zeros(3, dtype=np.float32) 
    cdef unsigned int x, y 
    x, y = int(max_loc[0]), int(max_loc[1]) 
    tmp[0] = (((response[y,x+1] - response[y,x-1])/<np.float32_t>2) * 
      (response[y,x] - min(response[y,x-1], response[y,x+1]))) 
    tmp[1] = response[y,x+1] - response[y,x-1] 
    tmp[2] = 2*(response[y,x] - min(response[y,x-1], response[y,x+1])) 

    print tmp[int(0)], tmp[int(1)], tmp[int(2)] 
    return tmp 

Сравнение:

>>> function1(response, max_loc) 
0.0821186 0.314914 1.04306 
array([ 0.08211858, 0.31491402, 1.0430603 ], dtype=float32) 

>>> function2(response, max_loc) 
0.0821186 0.314914 1.04306 
array([ 0.08211857, 0.31491402, 1.0430603 ], dtype=float32) 

Answer 2

Давайте сравним:

• function1 остается float32_t весь путь до конца.
• function2 преобразует в float при индексировании, выполняет ли промежуточные шаги с float, а затем возвращается к float32_t для окончательных результатов.
• function3 преобразует в float, но затем сразу возвращается к float32_t, выполняя промежуточные шаги с этим.
• function4 преобразует в float, выполняет ли промежуточные шаги, а затем возвращает окончательные результаты как float.

А почему function4 печатает то же самое, как function2, но возвращает что-то другое: Если вы посмотрите на типы, это просто. Значения, по-видимому, достаточно близки, что они происходят с print таким же образом, но не настолько близки к repr таким же образом. Что неудивительно, учитывая, что они не одного типа.