OmniSciDB  16c4e035a1
ExampleTableFunctions.cpp
Go to the documentation of this file.
1 /*
2  * Copyright 2021 OmniSci, Inc.
3  *
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
9  *
10  * Unless required by applicable law or agreed to in writing, software
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16
17 #include "ExampleTableFunctions.h"
18
19 #ifndef __CUDACC__
20
21 #ifdef HAVE_TBB
22 #include <tbb/parallel_for.h>
23 #endif
24 #endif
25
26 namespace Mandelbrot {
27 // Modified from original code by Akanksha Jolly,
28 // originally posted at https://www.geeksforgeeks.org/fractals-in-cc/,
29 // Author profile: https://auth.geeksforgeeks.org/user/akankshajolly/articles
30 // Used under Creative Commons License as allowed by GeeksForGeeks,
32 template <typename T>
34  const T cy,
35  const int32_t max_iterations) {
36  // z_real
37  T zx = 0;
38  // z_imaginary
39  T zy = 0;
40  int32_t num_iterations = 0;
41  // Calculate whether c(c_real + c_imaginary) belongs
42  // to the Mandelbrot set or not and draw a pixel
43  // at coordinates (x, y) accordingly
44  // If you reach the Maximum number of iterations
45  // and If the distance from the origin is
46  // greater terthan 2 exit the loop
47  while ((zx * zx + zy * zy < 4) && (num_iterations < max_iterations)) {
48  // Calculate Mandelbrot function
49  // z = z*z + c where z is a complex number
50  // tempx = z_real*_real - z_imaginary*z_imaginary + c_real
51  const T temp_x = zx * zx - zy * zy + cx;
52  // 2*z_real*z_imaginary + c_imaginary
53  zy = 2 * zx * zy + cy;
54  // Updating z_real = tempx
55  zx = temp_x;
56
57  // Increment counter
58  ++num_iterations;
59  }
60  return num_iterations;
61 }
62
63 DEVICE inline double get_scale(const double domain_min,
64  const double domain_max,
65  const int32_t num_bins) {
66  return (domain_max - domain_min) / num_bins;
67 }
68
69 template <typename T>
71 #ifdef _WIN32
73 #else
74  __attribute__((__used__))
75 #endif
76  void
77  mandelbrot_impl(const int32_t x_pixels,
78  const int32_t y_begin,
79  const int32_t y_end,
80  const T x_min,
81  const T y_min,
82  const T x_scale,
83  const T y_scale,
84  const int32_t max_iterations,
85  Column<T>& output_x,
86  Column<T>& output_y,
87  Column<int32_t>& output_num_iterations) {
88  // scanning every point in that rectangular area.
89  // each point represents a complex number (x + yi).
90  // iterate that complex number
91
92  for (int32_t y = y_begin; y < y_end; ++y) {
93  // c_imaginary
94  const T cy = y * y_scale + y_min;
95  for (int32_t x = 0; x < x_pixels; ++x) {
96  // c_real
97  const T cx = x * x_scale + x_min;
98  const int32_t output_pixel = y * x_pixels + x;
99  output_x[output_pixel] = cx;
100  output_y[output_pixel] = cy;
101  output_num_iterations[output_pixel] = mandelbrot_pixel(cx, cy, max_iterations);
102  }
103  }
104 }
105
106 #ifndef __CUDACC__
107
108 template <typename T>
110 #ifdef _WIN32
112 #else
113  __attribute__((__used__))
114 #endif
115  int32_t
117  const int32_t x_pixels,
118  const int32_t y_pixels,
119  const T x_min,
120  const T x_max,
121  const T y_min,
122  const T y_max,
123  const int32_t max_iterations,
124  Column<T>& output_x,
125  Column<T>& output_y,
126  Column<int32_t>& output_num_iterations) {
127  if (max_iterations < 1) {
128  return mgr.ERROR_MESSAGE("max_iterations must be a positive integer");
129  }
130  if (x_pixels < 1 || y_pixels < 1) {
131  return mgr.ERROR_MESSAGE("x_pixels and y_pixels must be positive integers");
132  }
133  if (x_max <= x_min || y_max <= y_min) {
134  return mgr.ERROR_MESSAGE(
135  "Max x and y bounds must be greater than min x and y bounds");
136  }
137
138  const T x_scale = get_scale(x_min, x_max, x_pixels);
139  const T y_scale = get_scale(y_min, y_max, y_pixels);
140
141  const int32_t num_pixels = x_pixels * y_pixels;
142  mgr.set_output_row_size(num_pixels);
143 #ifdef HAVE_TBB
144  tbb::parallel_for(tbb::blocked_range<int32_t>(0, y_pixels),
145  [&](const tbb::blocked_range<int32_t>& y_itr) {
146  const int32_t y_begin = y_itr.begin();
147  const int32_t y_end = y_itr.end();
148 #else
149  const int32_t y_begin = 0;
150  const int32_t y_end = y_pixels;
151 #endif
152  mandelbrot_impl(x_pixels,
153  y_begin,
154  y_end,
155  x_min,
156  y_min,
157  x_scale,
158  y_scale,
159  max_iterations,
160  output_x,
161  output_y,
162  output_num_iterations);
163 #ifdef HAVE_TBB
164  });
165 #endif
166  return num_pixels;
167 }
168
169 #else // #ifndef __CUDACC__
170
171 template <typename T>
173 #ifdef _WIN32
175 #else
176  __attribute__((__used__))
177 #endif
178  int32_t
179  mandelbrot_cuda_template__gpu_(const int32_t x_pixels,
180  const int32_t y_pixels,
181  const T x_min,
182  const T x_max,
183  const T y_min,
184  const T y_max,
185  const int32_t max_iterations,
186  Column<T>& output_x,
187  Column<T>& output_y,
188  Column<int32_t>& output_num_iterations) {
189  const T x_scale = get_scale(x_min, x_max, x_pixels);
190  const T y_scale = get_scale(y_min, y_max, y_pixels);
191  const int32_t num_pixels = x_pixels * y_pixels;
192
193  int32_t start = threadIdx.x + blockDim.x * blockIdx.x;
194  int32_t step = blockDim.x * gridDim.x;
195
196  for (int32_t output_pixel = start; output_pixel < num_pixels; output_pixel += step) {
197  const int32_t y = output_pixel / x_pixels;
198  const int32_t x = output_pixel % x_pixels;
199  const T cy = y * y_scale + y_min;
200  const T cx = x * x_scale + x_min;
201  T zx = 0;
202  T zy = 0;
203  int32_t num_iterations = 0;
204  for (; num_iterations < max_iterations; ++num_iterations) {
205  const T temp_x = zx * zx - zy * zy + cx;
206  zy = 2 * zx * zy + cy;
207  zx = temp_x;
208  if (zx * zx + zy * zy > 4.0) {
209  break;
210  }
211  }
212  output_x[output_pixel] = cx;
213  output_y[output_pixel] = cy;
214  output_num_iterations[output_pixel] = num_iterations;
215  }
216  return output_x.size();
217 }
218 #endif // #ifndef __CUDACC__
219
220 } // namespace Mandelbrot
221
222 #ifndef __CUDACC__
223
225 #ifdef _WIN32
227 #else
228 __attribute__((__used__))
229 #endif
231  const int32_t x_pixels,
232  const int32_t y_pixels,
233  const double x_min,
234  const double x_max,
235  const double y_min,
236  const double y_max,
237  const int32_t max_iterations,
238  Column<double>& output_x,
239  Column<double>& output_y,
240  Column<int32_t>& output_num_iterations) {
241  return Mandelbrot::mandelbrot_cpu_template<double>(mgr,
242  x_pixels,
243  y_pixels,
244  x_min,
245  x_max,
246  y_min,
247  y_max,
248  max_iterations,
249  output_x,
250  output_y,
251  output_num_iterations);
252 }
253
255 #ifdef _WIN32
257 #else
258 __attribute__((__used__))
259 #endif
261  const int32_t x_pixels,
262  const int32_t y_pixels,
263  const float x_min,
264  const float x_max,
265  const float y_min,
266  const float y_max,
267  const int32_t max_iterations,
268  Column<float>& output_x,
269  Column<float>& output_y,
270  Column<int32_t>& output_num_iterations) {
271  return Mandelbrot::mandelbrot_cpu_template<float>(mgr,
272  x_pixels,
273  y_pixels,
274  x_min,
275  x_max,
276  y_min,
277  y_max,
278  max_iterations,
279  output_x,
280  output_y,
281  output_num_iterations);
282 }
283
284 #else // #ifndef __CUDACC__
285
287 #ifdef _WIN32
289 #else
290 __attribute__((__used__))
291 #endif
292 int32_t tf_mandelbrot_cuda__gpu_(const int32_t x_pixels,
293  const int32_t y_pixels,
294  const double x_min,
295  const double x_max,
296  const double y_min,
297  const double y_max,
298  const int32_t max_iterations,
299  Column<double>& output_x,
300  Column<double>& output_y,
301  Column<int32_t>& output_num_iterations) {
302  return Mandelbrot::mandelbrot_cuda_template__gpu_(x_pixels,
303  y_pixels,
304  x_min,
305  x_max,
306  y_min,
307  y_max,
308  max_iterations,
309  output_x,
310  output_y,
311  output_num_iterations);
312 }
313
315 #ifdef _WIN32
317 #else
318 __attribute__((__used__))
319 #endif
320 int32_t tf_mandelbrot_cuda_float__gpu_(const int32_t x_pixels,
321  const int32_t y_pixels,
322  const float x_min,
323  const float x_max,
324  const float y_min,
325  const float y_max,
326  const int32_t max_iterations,
327  Column<float>& output_x,
328  Column<float>& output_y,
329  Column<int32_t>& output_num_iterations) {
330  return Mandelbrot::mandelbrot_cuda_template__gpu_(x_pixels,
331  y_pixels,
332  x_min,
333  x_max,
334  y_min,
335  y_max,
336  max_iterations,
337  output_x,
338  output_y,
339  output_num_iterations);
340 }
341
342 #endif // #ifndef __CUDACC__
void set_output_row_size(int64_t num_rows)
Definition: OmniSciTypes.h:329
DEVICE int64_t size() const
Definition: OmniSciTypes.h:243
#define TEMPLATE_NOINLINE
Definition: OmniSciTypes.h:36
DEVICE double get_scale(const double domain_min, const double domain_max, const int32_t num_bins)
#define DEVICE
TEMPLATE_NOINLINE void mandelbrot_impl(const int32_t x_pixels, const int32_t y_begin, const int32_t y_end, const T x_min, const T y_min, const T x_scale, const T y_scale, const int32_t max_iterations, Column< T > &output_x, Column< T > &output_y, Column< int32_t > &output_num_iterations)
#define EXTENSION_NOINLINE
Definition: OmniSciTypes.h:34
FORCE_INLINE T __attribute__((__may_alias__))*may_alias_ptr(T *ptr)
Definition: TypePunning.h:32
EXTENSION_NOINLINE int32_t tf_mandelbrot__cpu_(TableFunctionManager &mgr, const int32_t x_pixels, const int32_t y_pixels, const double x_min, const double x_max, const double y_min, const double y_max, const int32_t max_iterations, Column< double > &output_x, Column< double > &output_y, Column< int32_t > &output_num_iterations)
EXTENSION_NOINLINE int32_t tf_mandelbrot_float__cpu_(TableFunctionManager &mgr, const int32_t x_pixels, const int32_t y_pixels, const float x_min, const float x_max, const float y_min, const float y_max, const int32_t max_iterations, Column< float > &output_x, Column< float > &output_y, Column< int32_t > &output_num_iterations)
void parallel_for(const blocked_range< Int > &range, const Body &body, const Partitioner &p=Partitioner())
TEMPLATE_INLINE int32_t mandelbrot_pixel(const T cx, const T cy, const int32_t max_iterations)
TEMPLATE_NOINLINE int32_t mandelbrot_cpu_template(TableFunctionManager &mgr, const int32_t x_pixels, const int32_t y_pixels, const T x_min, const T x_max, const T y_min, const T y_max, const int32_t max_iterations, Column< T > &output_x, Column< T > &output_y, Column< int32_t > &output_num_iterations)
#define TEMPLATE_INLINE
Definition: OmniSciTypes.h:35