マンデルブロ集合をOpenCLで描く（Mathematica） http://blog.unfindable.net/archives/5772

opencl-mandelbrot.m

Needs["OpenCLLink`"]

kernelCode = "
# ifdef USING_DOUBLE _PRECISIONQ
# pragma OPENCL EXTENSION cl_khr _fp64 : enable
# pragma OPENCL EXTENSION cl_amd _fp64 : enable
# endif /* USING_DOUBLE _PRECISIONQ */

int check(Real_t x, Real_t y, int maxSteps) {
  Real_t re = 0, im = 0;
  int i = 0;
  while (i++ < maxSteps) {
    Real_t newRe = re * re - im * im + x;
    Real_t newIm = 2 * re * im + y;
    if (4 < newRe * newRe + newIm * newIm) return i;
    re = newRe;
    im = newIm;
  }
  return -1;
}
 
__kernel void mandelbrot(__global int* result, int maxSteps,
    int width, int height,
    Real_t xStart, Real_t xEnd, Real_t yStart, Real_t yEnd) {
  int xIndex = get_global_id(0);
  int yIndex = get_global_id(1);
  if (xIndex < width && yIndex < height) {
    Real_t x = xStart + (xEnd - xStart) * xIndex / width;
    Real_t y = yStart + (yEnd - yStart) * yIndex / height;
    result[xIndex + yIndex * width] = check(x, y, maxSteps);
  }
}";

kernel = OpenCLFunctionLoad[kernelCode, 
   "mandelbrot", {{"Integer32", _, "Output"}, "Integer32", "Integer32", 
    "Integer32", _Real, _Real, _Real, _Real}, {16, 16}, 
   ShellOutputFunction -> Print];

xMin = -2;
xMax = 1;
yMin = -3/2;
yMax = 3/2;
width = 300;
height = 300;

buffer = OpenCLMemoryAllocate["Integer32", {height, width}];
kernel[buffer, 100, width, height, N@xMin, N@xMax, N@yMin, N@yMax];
result = OpenCLMemoryGet[buffer];
OpenCLMemoryUnload[buffer];
ArrayPlot[result]