5

いくつかの変数に型を追加することで、python 関数を cython に相当するものに変換しました。ただし、cython 関数は、元の python 関数とはわずかに異なる出力を生成します。

この投稿でこの違いの理由のいくつかを学びました Cython: unsigned int indexs for numpy array gets different result しかし、この投稿で学んだことを使用しても、まだ cython 関数で同じ結果を生成することはできませんpythonのものとして。

そこで、私が試したことを示す 4 つの関数をまとめました。関数ごとにわずかに異なる結果が得られる理由を誰かが明らかにするのを手伝ってもらえますか? function1とまったく同じ値を返すcython関数を取得する方法は? 以下にコメントします。

%%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は、元の Python 関数で行った操作を表します。tmp1 が結果です。

function2は私の最初の cython バージョンで、わずかに異なる結果を生成します。どうやら、応答配列が型付き変数、unsigned int または int でインデックス付けされている場合、配列の型が np.float32_t であっても、結果は (PyFloat_FromDouble を使用して) double に強制されます。しかし、配列が python int でインデックス付けされている場合、関数 PyObject_GetItem が代わりに使用され、関数 1 で発生する np.float32_t を取得します。したがって、function1 の式は np.float32_t オペランドを使用して計算されますが、function2 の式は double を使用して計算されます。function1 とはわずかに異なる出力が得られます。

function3は、function1 と同じ出力を得ようとする 2 回目の cython の試みです。ここでは unsigned int インデックスを使用して配列応答にアクセスしますが、結果は np.float32_t 中間変数に残され、計算で使用されます。少し異なる結果が得られます。どうやら、print ステートメントは PyFloat_FromDouble を使用するため、np.float32_t を出力できません。

次に、python関数をcython関数に合わせて変更しようとしました。function4は、各式の少なくとも 1 つのオペランドを float に変換することでこれを達成しようとします。そのため、残りのオペランドも Python float に強制されます。これは cython では double であり、式は function2 のように double で計算されます。関数内の出力は function2 とまったく同じですが、返される値はわずかに異なります?!

4

2 に答える 2

2

比較してみましょう:

  • function1float32_tずっと残っています。
  • function2インデックス作成時に に変換しfloat、 で中間ステップを実行してから、最終結果のために にfloat戻します。float32_t
  • function3に変換されfloatますが、すぐに に戻りfloat32_t、中間ステップを実行します。
  • function4に変換しfloat、中間ステップを実行してから、最終結果を として返しますfloat

function4と同じものを出力するのに、異なるものを返す理由についてfunction2は、型を見れば簡単です。値は、同じように発生するほど十分に近いように見えますがprint、同じように十分に近いわけではありませんrepr。彼らが同じタイプではないことを考えると、これは驚くべきことではありません。

于 2013-03-13T03:53:52.840 に答える
2

10 進数で 7.225 桁の精度しかない単精度浮動小数点数を使用している場合、2 倍にすることによる小さな差異が大きな問題になるとは思いません。

の説明を明確にするためにfunction2、オブジェクトでインデックスを作成する場合、Cython はスカラー オブジェクトPyObject_GetItemを取得するために使用します ( C の単なる typedef である ではありません)。代わりにバッファに直接インデックスを付け、Cython がオブジェクトを必要とする場合、Cython は を呼び出します。、、およびを割り当てるオブジェクトが必要です。それらは型付けされていないためです。np.float32np.float32_tfloatPyFloat_FromDoubletmp1tmp2tmp3

一方function3、 では、変数を入力しましたが、出力して結果を返すオブジェクトをtmp作成する必要があります。float代わりに NumPy ndarray(以下を参照) を使用すると、その問題は発生しません。

ところでfunction1、 では、2 で割ったときに結果を に昇格させnp.float64ます。たとえば、次のようになります。

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

対。

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

すべての操作が関数と関数float32の両方にあることを確認したとしても、コンパイルされた拡張モジュールの 2 つの間で最終結果が異なる場合があります。次の例では、中間結果がすべてオブジェクトであることを確認しました。生成された C で、(または同等の typedef)へのキャストがないことを確認しました。それでも、これら 2 つの関数は、わずかに異なる結果を生成します。理由を理解するには、おそらくコンパイルされたアセンブリに飛び込む必要がありますが、単純なものを見落としている可能性があります。defcpdeffunction1np.float32function2double

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)
于 2013-03-13T08:18:30.333 に答える