Комплексное сокращение числа с помощью numba cuda - PullRequest
Новая
       14

Комплексное сокращение числа с помощью numba cuda

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

Я пытаюсь ускорить код на python, используя cuda \ numba. Код работает с большими массивами комплексных, чисел с плавающей точкой и целых чисел. Я включил и версию Python, и версию Numba-Cuda здесь. Версия numba-cuda не компилируется.

Я пытался выполнить вычисление комплексных чисел как отдельных действительных и мнимых чисел, хотя я мог подумать о сложном формате. 

def random_choice_noreplace(m,n, axis=-1):
 # m, n are the number of rows, cols of output
 return np.random.rand(m,n).argsort(axis=axis)

@cuda.jit
def cuda_kernel (d_npart, d_npts, d_data, d_data_index, d_coef, d_datasum, d_tmp):

 row, col = cuda.grid(2)
 if row < d_npart and col < d_npts :
 d_tmp[row, col] = d_data[d_data_index[row, col]]
 d_tmp[row, col] =d_tmp[row, col] * d_coef[row, col]
 # All threads get to this point ===============================
 cuda.syncthreads()
 if row == 0 and col ==0 :
 d_datasum = np.sum(d_tmp, axis=0)

def calculate_cuda (data, data_index, coef):

 npart, npts = data_index.shape
 # arrays to copy to GPU memory =====================================
 d_npart = cuda.to_device(npart)
 d_npts = cuda.to_device(npts)
 d_data = cuda.to_device(data)
 d_data_index = cuda.to_device(data_index)
 d_coef = cuda.to_device(coef)

 d_datasum = cuda.device_array(npts, np.complex64)
 d_tmp = cuda.device_array((npart,npts), np.complex64)

 threadsperblock = (16, 16)
 blockspergrid_x = int(math.ceil(npts / threadsperblock[0]))+1
 blockspergrid_y = int(math.ceil(npart / threadsperblock[1]))+1
 blockspergrid = (blockspergrid_x, blockspergrid_y)
 cuda_kernel[blockspergrid, threadsperblock](d_npart, d_npts, d_data, d_data_index, d_coef, d_datasum, d_tmp)
 # Copy data from GPU to CPU ========================================
 final_data_sum = d_datasum.copy_to_host()
 return final_data_sum



def calculate_python (data, data_index, coef):
 npart, npts = data_index.shape
 data_sum = np.zeros(npts, dtype=np.complex64)
 tmp = np.zeros(npts, dtype=np.complex64)
 print(" Calling python function...")
 start_time = time.time()
 for i in range(npart):
 tmp[:] = data[data_index[i]]
 data_sum += tmp * coef[i]
 return data_sum

if __name__ == '__main__':

 data_size = 1200
 rows = 31
 cols = 1000

 rand_float1 = np.random.randn(data_size)
 rand_float2 = np.random.randn(data_size)

 data = rand_float1 + 1j * rand_float2
 coef = np.random.randn(rows, cols)
 data_index = random_choice_noreplace(rows, cols)

 start_time = time.time()
 gpu_data_sum_python = calculate_python (data, data_index, coef)
 python_time = time.time() - start_time #print("gpu c : ", c_gpu)
 print("---- %s seconds for python ----:" % (python_time))


 start_time = time.time()
 gpu_data_sum = calculate_cuda (data, data_index, coef)
 gpu_time = time.time() - start_time
 print("---- %s seconds for gpu ----:" % (gpu_time))

Когда я запускаю код, я получаю эту ошибку: 

    Calling python function...
---- 0.000344038009644 seconds for python ----:
Traceback (most recent call last):
  File "GPU_Fake_PA_partial.py", line 82, in <module>
    gpu_data_sum = calculate_cuda (data, data_index, coef)
  File "GPU_Fake_PA_partial.py", line 44, in calculate_cuda
    cuda_kernel[blockspergrid, threadsperblock](d_npart, d_npts, d_data, d_data_index, d_coef, d_datasum, d_tmp)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/cuda/compiler.py", line 765, in __call__
    kernel = self.specialize(*args)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/cuda/compiler.py", line 776, in specialize
    kernel = self.compile(argtypes)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/cuda/compiler.py", line 792, in compile
    **self.targetoptions)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler_lock.py", line 32, in _acquire_compile_lock
    return func(*args, **kwargs)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/cuda/compiler.py", line 62, in compile_kernel
    cres = compile_cuda(pyfunc, types.void, args, debug=debug, inline=inline)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler_lock.py", line 32, in _acquire_compile_lock
    return func(*args, **kwargs)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/cuda/compiler.py", line 51, in compile_cuda
    locals={})
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler.py", line 926, in compile_extra
    return pipeline.compile_extra(func)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler.py", line 374, in compile_extra
    return self._compile_bytecode()
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler.py", line 857, in _compile_bytecode
    return self._compile_core()
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler.py", line 844, in _compile_core
    res = pm.run(self.status)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler_lock.py", line 32, in _acquire_compile_lock
    return func(*args, **kwargs)
  File "/disk/home/ajooya/software/venv/lib/python2.7/site-packages/numba/compiler.py", line 255, in run
    raise patched_exception
numba.errors.TypingError: Failed in nopython mode pipeline (step: nopython frontend)
Invalid use of Function(<built-in function lt>) with argument(s) of type(s): (int32, array(int64, 1d, A))
Known signatures:
 * (bool, bool) -> bool
 * (int8, int8) -> bool
 * (int16, int16) -> bool
 * (int32, int32) -> bool
 * (int64, int64) -> bool
 * (uint8, uint8) -> bool
 * (uint16, uint16) -> bool
 * (uint32, uint32) -> bool
 * (uint64, uint64) -> bool
 * (float32, float32) -> bool
 * (float64, float64) -> bool
 * parameterized
In definition 0:
    All templates rejected with literals.
In definition 1:
    All templates rejected without literals.
In definition 2:
    All templates rejected with literals.
In definition 3:
    All templates rejected without literals.
In definition 4:
    All templates rejected with literals.
In definition 5:
    All templates rejected without literals.
In definition 6:
    All templates rejected with literals.
In definition 7:
    All templates rejected without literals.
This error is usually caused by passing an argument of a type that is unsupported by the named function.
[1] During: typing of intrinsic-call at GPU_Fake_PA_partial.py (15)


File "GPU_Fake_PA_partial.py", line 15:
def cuda_kernel (d_npart, d_npts,  d_data, d_data_index, d_coef, d_datasum, d_tmp):
    <source elided>
    row, col = cuda.grid(2)
    if row < d_npart and col < d_npts :

Любая помощь высоко ценится.

1 Ответ

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

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

  1. методы массива numpy (например, np.sum()) нельзя использовать в ядрах numba CUDA .

  2. скалярные величины, передаваемые ядру numba cuda (например, npart, npts), не нуждаются и не должны обрабатываться массивом, как .to_device(). Просто используйте их как есть. Это причина ошибки Python, которую вы показываете: 

    Invalid use of Function(<built-in function lt>) with argument(s) of type(s): (int32, array(int64, 1d, A))

  3. Ваше ядро ​​было слишком сложным. Вы в основном выполняете суммы столбцов матрицы, которая была переставлена ​​в соответствии с шаблоном индекса, умноженного на массив коэффициентов. Мы можем выполнить это с помощью одного цикла на поток.

  4. Для реализации вышеизложенного нам не нужна двумерная сетка потоков.

  5. В коде, который вы разместили, возникли проблемы с отступами.

Для демонстрации я уменьшил размер вашего набора данных с 1000 столбцов до 15 столбцов. Вот пример, который обращается к вышеуказанным пунктам: 

$ cat t31.py
import numba as nb
import numpy as np
from numba import cuda
import time
import math

def random_choice_noreplace(m,n, axis=-1):
 # m, n are the number of rows, cols of output
 return np.random.rand(m,n).argsort(axis=axis)

@cuda.jit
def cuda_kernel (npart, npts, d_data, d_data_index, d_coef, d_datasum):
 col = cuda.grid(1)
 if col < npts:
     temp = 0
     for i in range(npart):
         temp += d_data[d_data_index[i, col]] * d_coef[i, col]
     d_datasum[col] = temp


def calculate_cuda (data, data_index, coef):

 npart, npts = data_index.shape
 # arrays to copy to GPU memory =====================================
 d_data = cuda.to_device(data)
 d_data_imag = cuda.to_device(data_imag)
 d_data_index = cuda.to_device(data_index)
 d_coef = cuda.to_device(coef)

 d_datasum = cuda.device_array(npts, np.complex64)

 threadsperblock = 64
 blockspergrid = int(math.ceil(npts / threadsperblock))+1
 cuda_kernel[blockspergrid, threadsperblock](npart, npts, d_data, d_data_index, d_coef, d_datasum)
 # Copy data from GPU to CPU ========================================
 final_data_sum = d_datasum.copy_to_host()
 return final_data_sum



def calculate_python (data, data_index, coef):
 npart, npts = data_index.shape
 data_sum = np.zeros(npts, dtype=np.complex64)
 tmp = np.zeros(npts, dtype=np.complex64)
 print(" Calling python function...")
 for i in range(npart):
  tmp[:] = data[data_index[i]]
  data_sum += tmp * coef[i]
 return data_sum

if __name__ == '__main__':

 rows = 31
 cols = 15
 data_size = rows * cols

 data_real = np.random.randn(data_size).astype(np.float32)
 data_imag = np.random.randn(data_size).astype(np.float32)

 data = data_real + 1j * data_imag
 coef = np.random.randn(rows, cols)
 data_index = random_choice_noreplace(rows, cols)

 start_time = time.time()
 gpu_data_sum_python = calculate_python (data, data_index, coef)
 python_time = time.time() - start_time #print("gpu c : ", c_gpu)
 print("---- %s seconds for python ----:" % (python_time))
 print(gpu_data_sum_python)

 start_time = time.time()
 gpu_data_sum = calculate_cuda (data, data_index, coef)
 gpu_time = time.time() - start_time
 print("---- %s seconds for gpu ----:" % (gpu_time))
 print(gpu_data_sum)
$ python t31.py
 Calling python function...
---- 0.000281095504761 seconds for python ----:
[-1.10292518+0.90700358j  2.67771578+2.47935939j -5.22553015-2.22675705j
 -3.55810285+2.39755774j  4.11441088-3.89396238j -2.70894790-0.75690132j
  3.24859619+0.65993834j  1.05531025+2.3959775j  -4.27368307+1.6297332j
  0.17896785-7.0437355j  -6.31506491+6.22674656j -1.85534143-6.08459902j
  0.40037563+6.33309507j -1.71916604-0.55055946j  0.49263301+1.08690035j]
---- 0.593510866165 seconds for gpu ----:
[-1.10292506+0.9070037j   2.67771506+2.47935939j -5.22553062-2.22675681j
 -3.55810285+2.39755821j  4.11441135-3.89396238j -2.70894790-0.75690138j
  3.24859619+0.65993822j  1.05531013+2.39597774j -4.27368307+1.62973344j
  0.17896791-7.0437355j  -6.31506491+6.22674656j -1.85534155-6.08459902j
  0.40037528+6.33309603j -1.71916604-0.55055946j  0.49263287+1.08690035j]
$

Обратите внимание, что между результатами вычислений хоста и устройства имеются небольшие различия, в некоторых случаях начиная с шестого знака после запятой. Я связываю это с возможными различиями порядка вычислений между кодом хоста и устройства в сочетании с пределами float32 (или complex64) numpy datatype.

Поскольку у вас есть встроенный тайминг в ваш код, вас может заинтересовать производительность. Для numba python я рекомендую типичную практику бенчмаркинга, которая заключается не в измерении первого прогона, а в измерении второго прогона. Это позволяет избежать однократных накладных расходов при измерении. Кроме того, мы хотели бы выбрать гораздо больший размер набора данных, чем 15 столбцов, чтобы дать графическому процессору достаточно большой объем работы для амортизации различных затрат. С этими модификациями, вот пример, показывающий, что версия GPU в этом коде может быть быстрее, чем версия CPU в этом коде: 

$ cat t31.py
import numba as nb
import numpy as np
from numba import cuda
import time
import math

def random_choice_noreplace(m,n, axis=-1):
 # m, n are the number of rows, cols of output
 return np.random.rand(m,n).argsort(axis=axis)

@cuda.jit
def cuda_kernel (npart, npts, d_data, d_data_index, d_coef, d_datasum):
 col = cuda.grid(1)
 if col < npts:
     temp = 0
     for i in range(npart):
         temp += d_data[d_data_index[i, col]] * d_coef[i, col]
     d_datasum[col] = temp


def calculate_cuda (data, data_index, coef):

 npart, npts = data_index.shape
 # arrays to copy to GPU memory =====================================
 d_data = cuda.to_device(data)
 d_data_imag = cuda.to_device(data_imag)
 d_data_index = cuda.to_device(data_index)
 d_coef = cuda.to_device(coef)

 d_datasum = cuda.device_array(npts, np.complex64)

 threadsperblock = 64
 blockspergrid = int(math.ceil(npts / threadsperblock))+1
 cuda_kernel[blockspergrid, threadsperblock](npart, npts, d_data, d_data_index, d_coef, d_datasum)
 # Copy data from GPU to CPU ========================================
 final_data_sum = d_datasum.copy_to_host()
 return final_data_sum



def calculate_python (data, data_index, coef):
 npart, npts = data_index.shape
 data_sum = np.zeros(npts, dtype=np.complex64)
 tmp = np.zeros(npts, dtype=np.complex64)
 print(" Calling python function...")
 for i in range(npart):
  tmp[:] = data[data_index[i]]
  data_sum += tmp * coef[i]
 return data_sum

if __name__ == '__main__':

 rows = 31
 cols = 1000000
 data_size = rows * cols

 data_real = np.random.randn(data_size).astype(np.float32)
 data_imag = np.random.randn(data_size).astype(np.float32)

 data = data_real + 1j * data_imag
 coef = np.random.randn(rows, cols)
 data_index = random_choice_noreplace(rows, cols)

 gpu_data_sum_python = calculate_python (data, data_index, coef)
 start_time = time.time()
 gpu_data_sum_python = calculate_python (data, data_index, coef)
 python_time = time.time() - start_time #print("gpu c : ", c_gpu)
 print("---- %s seconds for python ----:" % (python_time))
 print(gpu_data_sum_python)

 gpu_data_sum = calculate_cuda (data, data_index, coef)
 start_time = time.time()
 gpu_data_sum = calculate_cuda (data, data_index, coef)
 gpu_time = time.time() - start_time
 print("---- %s seconds for gpu ----:" % (gpu_time))
 print(gpu_data_sum)
$ python t31.py
 Calling python function...
 Calling python function...
---- 0.806931018829 seconds for python ----:
[  6.56164026-7.95271683j  -7.42586899+3.68758106j   3.48999476+3.10376763j
 ...,  13.12746525+4.61855698j   0.08796659+0.9710322j
  -6.54224586+4.89485168j]
---- 0.417661905289 seconds for gpu ----:
[  6.56164074-7.95271683j  -7.42586851+3.68758035j   3.48999500+3.10376763j
 ...,  13.12746525+4.61855745j   0.08796643+0.97103256j
  -6.54224634+4.89485121j]
$

С этими модификациями код графического процессора выглядит примерно в 2 раза быстрее, чем код процессора.

Это было измерено на CUDA 9.2, Fedora 27, Quadro K2000 (относительно небольшой, медленный GPU). Я бы не стал читать слишком много в этом сравнении, так как код процессора почти наверняка также неоптимален, и это все еще сравнительно небольшой объем работы на каждую точку выходных данных, чтобы ускорение графического процессора было интересным.

...