Cython: SIMD и OpenMP плохо работают вместе - PullRequest
0 голосов
/ 09 апреля 2020

Предполагая, что размер моих 2D-данных кратен 8 в каждом направлении.

Это игрушечный пример, иллюстрирующий проблему для гораздо более сложной проблемы, с которой я на самом деле имею дело. Я знаю, что это можно сделать гораздо проще с помощью Cython или numpy для этого простого примера.

cdef extern from "immintrin.h" nogil:  # in this example, we use AVX    
    ctypedef float  __m256
    __m256 _mm256_loadu_ps  (float *__P) nogil  
    __m256 _mm256_mul_ps    (__m256 __A, __m256 __B) nogil
    void   _mm256_storeu_ps (float *__P, __m256 __A) nogil   


cdef float multiply_sum_SIMD ( float [::1] A, float [::1] B, int m , int simd ) nogil:
    cdef:
        int k, j, l
        float s = 0.0
        float *out_AB  = <float *> malloc( simd * sizeof(float) )        
        __m256  m_A, m_B, m_AB

    for j in range(m):
        k    =  simd * j  
        m_A  = _mm256_loadu_ps( &A[k] )
        m_B  = _mm256_loadu_ps( &B[k] )

        m_AB = _mm256_mul_ps (  m_A , m_B )        
        _mm256_storeu_ps( &out_AB[0], m_AB )
        for l in range(8):
            s = s + out_AB[l]            
    return s

Итак, вот мой cython cpdef для использования вышеупомянутой функции, выполняющей один поток, и результат согласуется с numpy version.

cpdef float Sum_Cython_ST ( float [:,::1] A, float [:,::1] B ):    
    cdef:        
        int   max_row = A.shape[0] , max_col = A.shape[1]
        float s = 0.0
        int   i, j

        # integer variable related to SIMD
        int simd = 8
        int remainder, m = <int>( (max_col - 1) / simd) + 1, chk = max_col % simd

    for i in range(max_row):        
        s = s + multiply_sum_SIMD ( A[i,:], B[i,:] , m, simd)        
    return s

Теперь, если я переключу внешний l oop на использование prange с nogil

for i in prange( max_row, nogil = True):
    s = s + multiply_sum_SIMD( A[i,:], B[i,:], m , simd )  

, тогда у меня начнутся странные проблемы:
1. Он жалуется это 's' является ссылкой до того, как определено (но это было определено выше)
2. Если я переключил s + = умножить ....., то ошибка исчезла и скомпилировалась без ошибки, но результат был совершенно неверным, в мой случай.
3. Если я переключил multiply_sum_SIMD на функцию, которая не использует SIMD и имеет просто l oop внутри. Результаты снова согласуются с версией numpy (однопотоковой или многопоточной).

Я не понимаю, почему хорошо работать с этим однопоточным потоком, но как только я использую prange, появляется эта проблема , Пожалуйста, помогите.

...