-
Notifications
You must be signed in to change notification settings - Fork 15
/
set_positions_rand.cu
78 lines (64 loc) · 3.17 KB
/
set_positions_rand.cu
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
#include <cuda.h>
#include <stdio.h>
#include "datadef.h"
#include "warp_device.cuh"
__global__ void set_positions_rand_kernel(unsigned N , unsigned outer_cell_type, spatial_data * positions_ptr , unsigned * rn_bank , float x_min , float y_min , float z_min , float x_max , float y_max , float z_max ){
int tid = threadIdx.x+blockIdx.x*blockDim.x;
if (tid>=N){return;}
unsigned rn = rn_bank[tid];
float rn1 = get_rand(&rn);
float rn2 = get_rand(&rn);
float rn3 = get_rand(&rn);
float mu = 2.0 * get_rand(&rn)-1.0;
float theta = 6.28318530718 * get_rand(&rn);
// set surface distance very far
positions_ptr[tid].surf_dist = 500000;
// set directions isotropic
positions_ptr[tid].xhat = sqrtf(1.0-mu*mu) * cosf( theta );
positions_ptr[tid].yhat = sqrtf(1.0-mu*mu) * sinf( theta );
positions_ptr[tid].zhat = mu;
// set positions
if(outer_cell_type==0){ // cube
positions_ptr[tid].x = 0.9 * ( ( x_max - x_min ) * rn1 + x_min );
positions_ptr[tid].y = 0.9 * ( ( y_max - y_min ) * rn2 + y_min );
positions_ptr[tid].z = 0.9 * ( ( z_max - z_min ) * rn3 + z_min );
}
else if(outer_cell_type==1){ //cylinder
float t = 6.28318530718 * rn1;
float r = x_max * sqrtf(rn2);
positions_ptr[tid].x = 0.9 * r * cosf(t);
positions_ptr[tid].y = 0.9 * r * sinf(t);
positions_ptr[tid].z = 0.9 * ( ( z_max - z_min ) * rn3 + z_min );
}
else if(outer_cell_type==2){ // hex
positions_ptr[tid].x = 1.4 * ( ( x_max - x_min ) * rn1 + x_min );
positions_ptr[tid].y = 1.4 * ( ( y_max - y_min ) * rn2 + y_min );
positions_ptr[tid].z = 1.4 * ( ( z_max - z_min ) * rn3 + z_min );
}
else if(outer_cell_type==3){ // sphere
float t = 6.28318530718 * rn1;
float p = 2.0 * rn2 - 1.0;
float r = x_max * cbrtf(rn3);
positions_ptr[tid].x = 0.9 * r * sqrt(1.0-p*p) * cos(t);
positions_ptr[tid].y = 0.9 * r * sqrt(1.0-p*p) * sin(t);
positions_ptr[tid].z = 0.9 * r * p;
}
// update rn after using it
rn_bank[tid] = rn;
}
/**
* \brief sets starting cycle points uniformly random
* \details sets starting cycle points uniformly random in the specified volume
*
* @param[in] NUM_THREADS - the number of threads to run per thread block
* @param[in] N - the total number of threads to launch on the grid
* @param[in] outer_cell_type - the outer cell type, sets the shape of the sampling
* @param[in] d_space - device pointer to spatial data array
* @param[in] d_rn_bank - device pointer to random number array
* @param[in] outer_cell_dims - host pointer to array of outer cell extrema
*/
void set_positions_rand( unsigned NUM_THREADS, unsigned N, unsigned outer_cell_type, spatial_data * d_space , unsigned * d_rn_bank, float * outer_cell_dims){
unsigned blks = ( N + NUM_THREADS - 1 ) / NUM_THREADS;
set_positions_rand_kernel <<< blks, NUM_THREADS >>> ( N , outer_cell_type, d_space , d_rn_bank, outer_cell_dims[0], outer_cell_dims[1], outer_cell_dims[2], outer_cell_dims[3], outer_cell_dims[4], outer_cell_dims[5]);
cudaThreadSynchronize();
}