Добавление сложных 3D-массивов в pyopencl

Я пытаюсь добавить 2 трехмерных сложных массива с помощью библиотеки Python pyopencl. Результат, полученный на графическом процессоре, отличается от эталонного результата, полученного при использовании процессора. Почему код GPU не обеспечивает корректного сложения комплексных чисел? Я ожидаю, что соответствующие переменные res и ref будут равны друг другу. Код, который я использую:

import pyopencl as cl
import numpy as np
from numpy.fft import fftn, ifftn
from numpy.random import random, normal
import os


os.environ['PYOPENCL_COMPILER_OUTPUT'] = '1'

ctx = cl.create_some_context()
queue = cl.CommandQueue(ctx)

Lx = 100
Ly = 100
Lz = 1
L = Lx * Ly * Lz


const = np.zeros(4).astype(np.float32)
const[0] = Lx
const[1] = Ly
const[2] = Lz
const[3] = L

const_buf = cl.Buffer(ctx, cl.mem_flags.READ_ONLY | cl.mem_flags.COPY_HOST_PTR,    hostbuf=const)


arr1 = np.random.rand(Lz, Ly, Lx).astype(np.float32)
arr2 = np.random.rand(Lz, Ly, Lx).astype(np.float32)

F1 = fftn(arr1)
F2 = fftn(arr2)

out = np.zeros((Lz,Ly,Lx)).astype(np.complex64)
out_buf = cl.Buffer(ctx, cl.mem_flags.WRITE_ONLY, out.nbytes)

do_it_code = """
#include <pyopencl-complex.h>
__kernel void doit(__global const float *constants,
__global const cfloat_t *IN1,__global const cfloat_t *IN2,
__global cfloat_t *out)

{   int z = get_global_id(0);
    int y = get_global_id(1);
    int x = get_global_id(2);
    const int len = constants[3];
  const int lx = constants[0];
const int ly = constants[1];
const int lz = constants[2];
const int pl = lx * ly;
int i = x + y * lx + z * pl;

if (x <= lx-1 && y <= ly-1 && z <= lz-1) {
out[i] = cfloat_add(IN1[i],IN2[i]);

};
};
"""
bld_do_it = cl.Program(ctx, do_it_code).build()

def do_it(a,b):
    a_buf = cl.Buffer(ctx, cl.mem_flags.READ_ONLY | cl.mem_flags.COPY_HOST_PTR,     hostbuf=a.astype(np.complex64))
    b_buf = cl.Buffer(ctx, cl.mem_flags.READ_ONLY | cl.mem_flags.COPY_HOST_PTR, hostbuf=b.astype(np.complex64))
    launch = bld_do_it.doit(queue, a.shape, None, const_buf, a_buf, b_buf, out_buf)
    launch.wait()
    cl.enqueue_copy(queue, out, out_buf)

    return out

ref=F1+F2
print(ref)
print("_________")
res=do_it(F1,F2)
print(res)

python complex-numbers pyopencl
person GeorgePhysics    schedule 03.07.2020    source источник
comment
Разница в результатах, кажется, составляет несколько цифр после запятой, поэтому можно предположить, что это может быть связано с разными блоками FPU, которые находятся на ЦП и ускорителе. Подробнее см. здесь   -  person doqtor    schedule 03.07.2020
comment
нет, разница намного больше, чем несколько цифр после запятой.   -  person GeorgePhysics    schedule 03.07.2020
comment
Да, вы правы, есть нечто большее - смотрите мой ответ.   -  person doqtor    schedule 04.07.2020

Ответы (1)


arrow_upward
1
arrow_downward

Чтобы выяснить, что не так, обычно полезно работать с гораздо меньшим набором данных. Ставим Lx=Ly=2 и проверяем, что у нас получилось:

[[[ 3.40425836+0.j -2.09837677+0.j]
  [-0.27708016+0.j  0.60454108+0.j]]]
_________
[[[ 3.4042583 +0.j -0.27708015+0.j]
  [-2.0983768 +0.j  0.60454106+0.j]]]

Цифры почти идентичны, просто не все в тех местах, где они должны были быть. Простое транспонирование может исправить это:

res=np.transpose(res, (0,2,1))

а затем результат соответствует эталонному результату:

[[[ 3.4042583 +0.j -2.0983768 +0.j]
  [-0.27708015+0.j  0.60454106+0.j]]]

Между результатами все еще есть небольшое расхождение - после нескольких цифр после запятой, - но это можно объяснить тем, что между процессором и графическим процессором используются разные блоки FPU (в моем случае). Подробнее об этом можно прочитать, например, здесь.

Для этого не требуется транспонирование - размещение Lx, Ly, Lz должно быть где-то неправильным - я позволю OP разобраться в этом.

person doqtor    schedule 04.07.2020
comment
Спасибо, это имеет смысл! - person GeorgePhysics; 05.07.2020