diff --git a/src/CUDA/NVRTCKernels/CudaChiSquareKernel.cu b/src/CUDA/NVRTCKernels/CudaChiSquareKernel.cu index b22cab7..9a7f808 100644 --- a/src/CUDA/NVRTCKernels/CudaChiSquareKernel.cu +++ b/src/CUDA/NVRTCKernels/CudaChiSquareKernel.cu @@ -83,6 +83,56 @@ __device__ double ifld(double t, double alpha, double phi, double nu, double lam return alpha*cos(wt+ph)*exp(-lambdaT*t) + (1.0-alpha)*exp(-lambdaL*t); } +__device__ double ifgk(double t, double alpha, double nu, double sigma, double lambda, double beta) { + double wt = TWO_PI*nu*t; + double rate2 = sigma*sigma*t*t; + double rateL = 0.0; + double result = 0.0; + + // make sure lambda > 0 + if (lambda < 0.0) + return 0.0; + + if (beta < 0.001) { + rateL = 1.0; + } else { + rateL = pow(lambda*t, beta); + } + + if (nu < 0.01) { + result = (1.0-alpha)*exp(-rateL) + alpha*(1.0-rate2)*exp(-0.5*rate2); + } else { + result = (1.0-alpha)*exp(-rateL) + alpha*(cos(wt)-sigma*sigma*t*t/(wt)*sin(wt))*exp(-0.5*rate2); + } + + return result; +} + +__device__ double ifll(double t, double alpha, double nu, double a, double lambda, double beta) { + double wt = TWO_PI*nu*t; + double at = a*t; + double rateL = 0.0; + double result = 0.0; + + // make sure lambda > 0 + if (lambda < 0.0) + return 0.0; + + if (beta < 0.001) { + rateL = 1.0; + } else { + rateL = pow(lambda*t, beta); + } + + if (nu < 0.01) { + result = (1.0-alpha)*exp(-rateL) + alpha*(1.0-at)*exp(-at); + } else { + result = (1.0-alpha)*exp(-rateL) + alpha*(cos(wt)-a/(TWO_PI*nu)*sin(wt))*exp(-at); + } + + return result; +} + __device__ double b(double t, double phi, double nu) { return j0(TWO_PI*nu*t + DEG_TO_RAD*phi); } diff --git a/src/OpenCL/OpenCLKernels/OpenCLChiSquareRuntime.cl b/src/OpenCL/OpenCLKernels/OpenCLChiSquareRuntime.cl index bdc9374..73ddaa4 100644 --- a/src/OpenCL/OpenCLKernels/OpenCLChiSquareRuntime.cl +++ b/src/OpenCL/OpenCLKernels/OpenCLChiSquareRuntime.cl @@ -106,6 +106,56 @@ double ifld(double t, double alpha, double phi, double nu, double lambdaT, doubl return alpha*cos(wt+ph)*exp(-lambdaT*t) + (1.0-alpha)*exp(-lambdaL*t); } +double ifgk(double t, double alpha, double nu, double sigma, double lambda, double beta) { + double wt = TWO_PI*nu*t; + double rate2 = sigma*sigma*t*t; + double rateL = 0.0; + double result = 0.0; + + // make sure lambda > 0 + if (lambda < 0.0) + return 0.0; + + if (beta < 0.001) { + rateL = 1.0; + } else { + rateL = pow(lambda*t, beta); + } + + if (nu < 0.01) { + result = (1.0-alpha)*exp(-rateL) + alpha*(1.0-rate2)*exp(-0.5*rate2); + } else { + result = (1.0-alpha)*exp(-rateL) + alpha*(cos(wt)-sigma*sigma*t*t/(wt)*sin(wt))*exp(-0.5*rate2); + } + + return result; +} + +double ifll(double t, double alpha, double nu, double a, double lambda, double beta) { + double wt = TWO_PI*nu*t; + double at = a*t; + double rateL = 0.0; + double result = 0.0; + + // make sure lambda > 0 + if (lambda < 0.0) + return 0.0; + + if (beta < 0.001) { + rateL = 1.0; + } else { + rateL = pow(lambda*t, beta); + } + + if (nu < 0.01) { + result = (1.0-alpha)*exp(-rateL) + alpha*(1.0-at)*exp(-at); + } else { + result = (1.0-alpha)*exp(-rateL) + alpha*(cos(wt)-a/(TWO_PI*nu)*sin(wt))*exp(-at); + } + + return result; +} + double b(double t, double phi, double nu) { return bessj0(TWO_PI*nu*t + DEG_TO_RAD*phi); }