gpu_mandelbrot_between.cu 4.8 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114
  1. #include <stdio.h>
  2. #include <stdint.h>
  3. #include <starpu.h>
  4. #define THREADS_PER_BLOCK 64
  5. static inline long long jlstarpu_max(long long a, long long b)
  6. {
  7. return (a > b) ? a : b;
  8. }
  9. static inline long long jlstarpu_interval_size(long long start, long long step, long long stop)
  10. {
  11. if (stop >= start){
  12. return jlstarpu_max(0, (stop - start + 1) / step);
  13. } else {
  14. return jlstarpu_max(0, (stop - start - 1) / step);
  15. }
  16. }
  17. __device__ static inline long long jlstarpu_max__device(long long a, long long b)
  18. {
  19. return (a > b) ? a : b;
  20. }
  21. __device__ static inline long long jlstarpu_interval_size__device(long long start, long long step, long long stop)
  22. {
  23. if (stop >= start){
  24. return jlstarpu_max__device(0, (stop - start + 1) / step);
  25. } else {
  26. return jlstarpu_max__device(0, (stop - start - 1) / step);
  27. }
  28. }
  29. __global__ void mandelbrot(int64_t kernel_ids__start_1, int64_t kernel_ids__step_1, int64_t kernel_ids__dim_1, int64_t kernel_ids__start_2,
  30. int64_t kernel_ids__step_2, int64_t kernel_ids__dim_2, double* ptr_hF6lCYyJ, int64_t local_width,
  31. int64_t* ptr_qoUGBRtY, int64_t local_height, double conv_limit, int64_t* ptr_A5zD9sJZ,
  32. uint32_t ld_A5zD9sJZ)
  33. {
  34. int64_t THREAD_ID = (int64_t) ((((blockIdx).x) * ((blockDim).x)) + ((threadIdx).x));
  35. if ((THREAD_ID) >= (((1) * (kernel_ids__dim_2)) * (kernel_ids__dim_1)))
  36. {
  37. return ;
  38. };
  39. int64_t kernel_ids__index_1 = (int64_t) (((THREAD_ID) / ((1) * (kernel_ids__dim_2))) % (kernel_ids__dim_1));
  40. int64_t kernel_ids__index_2 = (int64_t) (((THREAD_ID) / (1)) % (kernel_ids__dim_2));
  41. int64_t x = (int64_t) ((kernel_ids__start_1) + ((kernel_ids__index_1) * (kernel_ids__step_1)));
  42. int64_t y = (int64_t) ((kernel_ids__start_2) + ((kernel_ids__index_2) * (kernel_ids__step_2)));
  43. double max_iterations = (double) (ptr_hF6lCYyJ[(5) - (1)]);
  44. double zoom = (double) ((ptr_hF6lCYyJ[(3) - (1)]) * (0.25296875));
  45. int64_t X = (int64_t) ((x) + ((local_width) * ((ptr_qoUGBRtY[(2) - (1)]) - (1))));
  46. int64_t Y = (int64_t) ((y) + ((local_height) * ((ptr_qoUGBRtY[(1) - (1)]) - (1))));
  47. double cr = (double) ((ptr_hF6lCYyJ[(1) - (1)]) + (((X) - ((ptr_hF6lCYyJ[(3) - (1)]) / (2))) / (zoom)));
  48. double zr = (double) (cr);
  49. double ci = (double) ((ptr_hF6lCYyJ[(2) - (1)]) + (((Y) - ((ptr_hF6lCYyJ[(4) - (1)]) / (2))) / (zoom)));
  50. double zi = (double) (ci);
  51. int64_t n = (int64_t) (0);
  52. int64_t b1 = (int64_t) (((n) < (max_iterations)) + ((((zr) * (zr)) + ((zi) * (zi))) < ((conv_limit) * (conv_limit))));
  53. while ((b1) >= (2))
  54. {
  55. double tmp = (double) ((((zr) * (zr)) - ((zi) * (zi))) + (cr));
  56. zi = ((2) * (zr) * (zi)) + (ci);
  57. zr = tmp;
  58. n = (n) + (1);
  59. b1 = ((n) <= (max_iterations)) + ((((zr) * (zr)) + ((zi) * (zi))) <= ((conv_limit) * (conv_limit)));
  60. }
  61. ;
  62. if ((n) < (max_iterations))
  63. {
  64. ptr_A5zD9sJZ[((y) + (((x) - (1)) * (ld_A5zD9sJZ))) - (1)] = ((255) * (n)) / (max_iterations);
  65. } else
  66. {
  67. ptr_A5zD9sJZ[((y) + (((x) - (1)) * (ld_A5zD9sJZ))) - (1)] = 0;
  68. }
  69. ;
  70. }
  71. extern "C" void CUDA_mandelbrot(void** buffers_uwrYFDVe, void* cl_arg_uwrYFDVe)
  72. {
  73. uint32_t ld_A5zD9sJZ = (uint32_t) (STARPU_MATRIX_GET_LD(buffers_uwrYFDVe[(1) - (1)]));
  74. int64_t* ptr_A5zD9sJZ = (int64_t*) (STARPU_MATRIX_GET_PTR(buffers_uwrYFDVe[(1) - (1)]));
  75. double* ptr_hF6lCYyJ = (double*) (STARPU_VECTOR_GET_PTR(buffers_uwrYFDVe[(2) - (1)]));
  76. int64_t* ptr_qoUGBRtY = (int64_t*) (STARPU_VECTOR_GET_PTR(buffers_uwrYFDVe[(3) - (1)]));
  77. int64_t local_width = (int64_t) (STARPU_MATRIX_GET_NY(buffers_uwrYFDVe[(1) - (1)]));
  78. int64_t local_height = (int64_t) (STARPU_MATRIX_GET_NX(buffers_uwrYFDVe[(1) - (1)]));
  79. double conv_limit = (double) (2.0);
  80. int64_t kernel_ids__start_1 = (int64_t) (1);
  81. int64_t kernel_ids__step_1 = (int64_t) (1);
  82. int64_t kernel_ids__dim_1 = (int64_t) (jlstarpu_interval_size(kernel_ids__start_1, kernel_ids__step_1, local_width));
  83. int64_t kernel_ids__start_2 = (int64_t) (1);
  84. int64_t kernel_ids__step_2 = (int64_t) (1);
  85. int64_t kernel_ids__dim_2 = (int64_t) (jlstarpu_interval_size(kernel_ids__start_2, kernel_ids__step_2, local_height));
  86. int64_t nthreads = (int64_t) (((1) * (kernel_ids__dim_1)) * (kernel_ids__dim_2));
  87. int64_t nblocks = (int64_t) ((((nthreads) + (THREADS_PER_BLOCK)) - (1)) / (THREADS_PER_BLOCK));
  88. mandelbrot
  89. <<< nblocks, THREADS_PER_BLOCK, 0, starpu_cuda_get_local_stream()
  90. >>> (kernel_ids__start_1, kernel_ids__step_1, kernel_ids__dim_1, kernel_ids__start_2,
  91. kernel_ids__step_2, kernel_ids__dim_2, ptr_hF6lCYyJ, local_width,
  92. ptr_qoUGBRtY, local_height, conv_limit, ptr_A5zD9sJZ,
  93. ld_A5zD9sJZ);
  94. ;
  95. cudaStreamSynchronize(starpu_cuda_get_local_stream());
  96. }