From 24f394c6938598c385b027937fb46405d27ca302 Mon Sep 17 00:00:00 2001 From: Uldis Locans Date: Mon, 24 Apr 2017 10:44:41 +0200 Subject: [PATCH] add enableRutherfordScattering option to OPALs collimatorPhysics GPU version --- src/Algorithms/CollimatorPhysics.h | 3 ++- src/CUDA/CudaCollimatorPhysics.cu | 26 ++++++++++++++++---------- src/CUDA/CudaCollimatorPhysics.cuh | 2 +- src/DKSBase.cpp | 6 ++++-- src/DKSBase.h | 3 ++- src/MIC/MICCollimatorPhysics.cpp | 4 +++- src/MIC/MICCollimatorPhysics.h | 3 ++- src/OpenCL/OpenCLCollimatorPhysics.cpp | 2 +- src/OpenCL/OpenCLCollimatorPhysics.h | 3 ++- 9 files changed, 33 insertions(+), 19 deletions(-) diff --git a/src/Algorithms/CollimatorPhysics.h b/src/Algorithms/CollimatorPhysics.h index b7e8190..212957b 100644 --- a/src/Algorithms/CollimatorPhysics.h +++ b/src/Algorithms/CollimatorPhysics.h @@ -19,7 +19,8 @@ public: virtual ~DKSCollimatorPhysics() { } - virtual int CollimatorPhysics(void *mem_ptr, void *par_ptr, int numpartices) = 0; + virtual int CollimatorPhysics(void *mem_ptr, void *par_ptr, int numpartices, + bool enableRutherfordScattering = true) = 0; virtual int CollimatorPhysicsSoA(void *label_ptr, void *localID_ptr, void *rx_ptr, void *ry_ptr, void *rz_ptr, diff --git a/src/CUDA/CudaCollimatorPhysics.cu b/src/CUDA/CudaCollimatorPhysics.cu index 495ad59..ed9cd09 100644 --- a/src/CUDA/CudaCollimatorPhysics.cu +++ b/src/CUDA/CudaCollimatorPhysics.cu @@ -23,6 +23,7 @@ #define X0_M 9 #define I_M 10 #define DT_M 11 +#define LOWENERGY_THR 12 #define BLOCK_SIZE 128 #define NUMPAR 12 @@ -81,7 +82,7 @@ __device__ inline void energyLoss(double &Eng, bool &pdead, curandState &state, Eng = Eng + delta_E / 1E3; } - pdead = ((Eng<1E-4) || (dEdx>0)); + pdead = ( (Eng < par[LOWENERGY_THR]) || (dEdx > 0) ); } @@ -117,7 +118,9 @@ __device__ inline void Rot(double &px, double &pz, double &x, double &z, double pz = -pxz*sin(Psixz)*sin(thetacou) + pxz*cos(Psixz)*cos(thetacou); } -__device__ inline void coulombScat(double3 &R, double3 &P, curandState &state, double* par) { +__device__ inline void coulombScat(double3 &R, double3 &P, curandState &state, + double* par, bool enableRutherfordScattering) +{ double Eng = sqrt(dot(P, P) + 1.0) * M_P - M_P; double gamma = (Eng + M_P) / M_P; @@ -145,7 +148,7 @@ __device__ inline void coulombScat(double3 &R, double3 &P, curandState &state, d Rot(P.x, P.z, R.x, R.z, xplane, normP, thetacou, deltas, 1, par); double P2 = curand_uniform_double(&state);//gsl_rng_uniform(rGen_m); - if(P2 < 0.0047) { + if( (P2 < 0.0047) && enableRutherfordScattering) { double P3 = curand_uniform_double(&state);//gsl_rng_uniform(rGen_m); double thetaru = 2.5 * sqrt(1 / P3) * sqrt(2.0) * theta0; double P4 = curand_uniform_double(&state);//gsl_rng_uniform(rGen_m); @@ -171,7 +174,7 @@ __device__ inline void coulombScat(double3 &R, double3 &P, curandState &state, d Rot(P.y,P.z,R.y,R.z, yplane, normP, thetacou, deltas, 2, par); P2 = curand_uniform_double(&state);//gsl_rng_uniform(rGen_m); - if(P2 < 0.0047) { + if( (P2 < 0.0047) && enableRutherfordScattering) { double P3 = curand_uniform_double(&state);//gsl_rng_uniform(rGen_m); double thetaru = 2.5 * sqrt(1 / P3) * sqrt(2.0) * theta0; double P4 = curand_uniform_double(&state);//gsl_rng_uniform(rGen_m); @@ -185,7 +188,7 @@ __device__ inline void coulombScat(double3 &R, double3 &P, curandState &state, d template __global__ void kernelCollimatorPhysics(T *data, double *par, curandState *state, - int numparticles) + int numparticles, bool enableRutherfordScattering) { //get global id and thread id @@ -227,7 +230,7 @@ __global__ void kernelCollimatorPhysics(T *data, double *par, curandState *state P.x = P.x * ptot / sq; P.y = P.y * ptot / sq; P.z = P.z * ptot / sq; - coulombScat(R[tid], P, s, p); + coulombScat(R[tid], P, s, p, enableRutherfordScattering); data[idx].Pincol = P; } else { @@ -250,7 +253,8 @@ __global__ void kernelCollimatorPhysics(T *data, double *par, curandState *state } __global__ void kernelCollimatorPhysics2(CUDA_PART2_SMALL data, double *par, - curandState *state, int numparticles) + curandState *state, int numparticles, + bool enableRutherfordScattering) { //get global id and thread id @@ -288,7 +292,7 @@ __global__ void kernelCollimatorPhysics2(CUDA_PART2_SMALL data, double *par, P.x = P.x * ptot / sq; P.y = P.y * ptot / sq; P.z = P.z * ptot / sq; - coulombScat(R[tid], P, s, p); + coulombScat(R[tid], P, s, p, enableRutherfordScattering); data.Pincol[idx] = P; } else { @@ -611,7 +615,8 @@ struct less_then } }; -int CudaCollimatorPhysics::CollimatorPhysics(void *mem_ptr, void *par_ptr, int numparticles) +int CudaCollimatorPhysics::CollimatorPhysics(void *mem_ptr, void *par_ptr, int numparticles, + bool enableRutherfordScattering) { int threads = BLOCK_SIZE; @@ -624,7 +629,8 @@ int CudaCollimatorPhysics::CollimatorPhysics(void *mem_ptr, void *par_ptr, int n kernelCollimatorPhysics<<>>((CUDA_PART_SMALL*)mem_ptr, (double*)par_ptr, m_base->cuda_getCurandStates(), - numparticles); + numparticles, + enableRutherfordScattering); cudaError_t err = cudaGetLastError(); if (err != cudaSuccess) diff --git a/src/CUDA/CudaCollimatorPhysics.cuh b/src/CUDA/CudaCollimatorPhysics.cuh index 9808f33..9f417d8 100644 --- a/src/CUDA/CudaCollimatorPhysics.cuh +++ b/src/CUDA/CudaCollimatorPhysics.cuh @@ -110,7 +110,7 @@ public: * */ int CollimatorPhysics(void *mem_ptr, void *par_ptr, - int numpartices); + int numpartices, bool enableRutherfordScattering = true); int CollimatorPhysicsSoA(void *label_ptr, void *localID_ptr, void *rx_ptr, void *ry_ptr, void *rz_ptr, diff --git a/src/DKSBase.cpp b/src/DKSBase.cpp index 8ea5f58..b72767e 100644 --- a/src/DKSBase.cpp +++ b/src/DKSBase.cpp @@ -746,11 +746,13 @@ int DKSBase::callCollimatorPhysics(void *mem_ptr, void *par_ptr, } -int DKSBase::callCollimatorPhysics2(void *mem_ptr, void *par_ptr, int numparticles) +int DKSBase::callCollimatorPhysics2(void *mem_ptr, void *par_ptr, int numparticles, + bool enableRutherfordScattering) { if (apiCuda()) - return CUDA_SAFECALL( ccol->CollimatorPhysics(mem_ptr, par_ptr, numparticles) ); + return CUDA_SAFECALL( ccol->CollimatorPhysics(mem_ptr, par_ptr, numparticles, + enableRutherfordScattering) ); else if (apiOpenMP()) return MIC_SAFECALL( miccol->CollimatorPhysics(mem_ptr, par_ptr, numparticles) ); diff --git a/src/DKSBase.h b/src/DKSBase.h index 79e23c8..fe0e70c 100644 --- a/src/DKSBase.h +++ b/src/DKSBase.h @@ -1041,7 +1041,8 @@ public: * For specifics check OPAL docs and CudaCollimatorPhysics class documentation. * TODO: opencl and mic implementations. */ - int callCollimatorPhysics2(void *mem_ptr, void *par_ptr, int numparticles); + int callCollimatorPhysics2(void *mem_ptr, void *par_ptr, int numparticles, + bool enableRutherfordScattering = true); /** * Monte carlo code for the degrader from OPAL classic/5.0/src/Solvers/CollimatorPhysics.cpp on device. diff --git a/src/MIC/MICCollimatorPhysics.cpp b/src/MIC/MICCollimatorPhysics.cpp index 6a1b937..72b7028 100644 --- a/src/MIC/MICCollimatorPhysics.cpp +++ b/src/MIC/MICCollimatorPhysics.cpp @@ -368,7 +368,9 @@ void energyLoss(double &Eng, double &dEdx, double *par, double *randv, int ri) { } -int MICCollimatorPhysics::CollimatorPhysics(void *mem_ptr, void *par_ptr, int numparticles) { +int MICCollimatorPhysics::CollimatorPhysics(void *mem_ptr, void *par_ptr, + int numparticles, boll enableRutherfordScattering) +{ //cast device memory pointers to appropriate types MIC_PART_SMALL *data = (MIC_PART_SMALL*) mem_ptr; diff --git a/src/MIC/MICCollimatorPhysics.h b/src/MIC/MICCollimatorPhysics.h index 0795779..d006136 100644 --- a/src/MIC/MICCollimatorPhysics.h +++ b/src/MIC/MICCollimatorPhysics.h @@ -40,7 +40,8 @@ public: ~MICCollimatorPhysics() { }; - int CollimatorPhysics(void *mem_ptr, void *par_ptr, int numparticles); + int CollimatorPhysics(void *mem_ptr, void *par_ptr, int numparticles, + bool enableRutherfordScattering = true); int CollimatorPhysicsSoA(void *label_ptr, void *localID_ptr, void *rx_ptr, void *ry_ptr, void *rz_ptr, diff --git a/src/OpenCL/OpenCLCollimatorPhysics.cpp b/src/OpenCL/OpenCLCollimatorPhysics.cpp index 46d8b24..da1be1f 100644 --- a/src/OpenCL/OpenCLCollimatorPhysics.cpp +++ b/src/OpenCL/OpenCLCollimatorPhysics.cpp @@ -34,7 +34,7 @@ TODO: 2. boost.compute sort for user defined structure crashes */ int OpenCLCollimatorPhysics::CollimatorPhysics(void *mem_ptr, void *par_ptr, - int numparticles) + int numparticles, bool enableRutherfordScattering) { /* //set number of total threads, and number threads per block diff --git a/src/OpenCL/OpenCLCollimatorPhysics.h b/src/OpenCL/OpenCLCollimatorPhysics.h index 7b532ff..3a299cb 100644 --- a/src/OpenCL/OpenCLCollimatorPhysics.h +++ b/src/OpenCL/OpenCLCollimatorPhysics.h @@ -52,7 +52,8 @@ public: } /* execute degrader code on device */ - int CollimatorPhysics(void *mem_ptr, void *par_ptr, int numparticles); + int CollimatorPhysics(void *mem_ptr, void *par_ptr, int numparticles, + bool enableRutherfordScattering = true); int CollimatorPhysicsSoA(void *label_ptr, void *localID_ptr, void *rx_ptr, void *ry_ptr, void *rz_ptr,