OpenCL中的复杂数字支持
我知道OpenCL不支持复杂的数字,而且我读过的这个功能不会很快出现。 尽管如此,一些例子在OpenCL内核中使用了复数(例如FFT实现)。
有没有人有这方面的经验? 在OpenCL中启用复杂数字支持的“最佳”方法是什么? 我假设使用float2来包含实部和虚部,但是我应该写一组宏还是内联函数更好? 有没有人知道一组函数/宏是否已经存在用于此目的?
所以,因为我需要一组函数来处理OpenCL中的复数,我最终实现了一组函数。 具体而言,我需要总和和减法(平凡,可以用标准向量操作完成),乘法,除法,得到复数模,参数(或角度)和平方根。
相关维基百科文章:
http://en.wikipedia.org/wiki/Complex_number#Absolute_value_and_argument
http://en.wikipedia.org/wiki/Square_root#Principal_square_root_of_a_complex_number
这大部分是微不足道的,但这确实需要一些时间,所以希望这次可以节省一些人,这里是:
//2 component vector to hold the real and imaginary parts of a complex number:
typedef float2 cfloat;
#define I ((cfloat)(0.0, 1.0))
/*
* Return Real (Imaginary) component of complex number:
*/
inline float real(cfloat a){
return a.x;
}
inline float imag(cfloat a){
return a.y;
}
/*
* Get the modulus of a complex number (its length):
*/
inline float cmod(cfloat a){
return (sqrt(a.x*a.x + a.y*a.y));
}
/*
* Get the argument of a complex number (its angle):
* http://en.wikipedia.org/wiki/Complex_number#Absolute_value_and_argument
*/
inline float carg(cfloat a){
if(a.x > 0){
return atan(a.y / a.x);
}else if(a.x < 0 && a.y >= 0){
return atan(a.y / a.x) + M_PI;
}else if(a.x < 0 && a.y < 0){
return atan(a.y / a.x) - M_PI;
}else if(a.x == 0 && a.y > 0){
return M_PI/2;
}else if(a.x == 0 && a.y < 0){
return -M_PI/2;
}else{
return 0;
}
}
/*
* Multiply two complex numbers:
*
* a = (aReal + I*aImag)
* b = (bReal + I*bImag)
* a * b = (aReal + I*aImag) * (bReal + I*bImag)
* = aReal*bReal +I*aReal*bImag +I*aImag*bReal +I^2*aImag*bImag
* = (aReal*bReal - aImag*bImag) + I*(aReal*bImag + aImag*bReal)
*/
inline cfloat cmult(cfloat a, cfloat b){
return (cfloat)( a.x*b.x - a.y*b.y, a.x*b.y + a.y*b.x);
}
/*
* Divide two complex numbers:
*
* aReal + I*aImag (aReal + I*aImag) * (bReal - I*bImag)
* ----------------- = ---------------------------------------
* bReal + I*bImag (bReal + I*bImag) * (bReal - I*bImag)
*
* aReal*bReal - I*aReal*bImag + I*aImag*bReal - I^2*aImag*bImag
* = ---------------------------------------------------------------
* bReal^2 - I*bReal*bImag + I*bImag*bReal -I^2*bImag^2
*
* aReal*bReal + aImag*bImag aImag*bReal - Real*bImag
* = ---------------------------- + I* --------------------------
* bReal^2 + bImag^2 bReal^2 + bImag^2
*
*/
inline cfloat cdiv(cfloat a, cfloat b){
return (cfloat)((a.x*b.x + a.y*b.y)/(b.x*b.x + b.y*b.y), (a.y*b.x - a.x*b.y)/(b.x*b.x + b.y*b.y));
}
/*
* Square root of complex number.
* Although a complex number has two square roots, numerically we will
* only determine one of them -the principal square root, see wikipedia
* for more info:
* http://en.wikipedia.org/wiki/Square_root#Principal_square_root_of_a_complex_number
*/
inline cfloat csqrt(cfloat a){
return (cfloat)( sqrt(cmod(a)) * cos(carg(a)/2), sqrt(cmod(a)) * sin(carg(a)/2));
}
PyOpenCL在OpenCL中有一个更完整和强大的复数实现:
https://github.com/pyopencl/pyopencl/blob/master/pyopencl/cl/pyopencl-complex.h
链接地址: http://www.djcxy.com/p/24373.html