static const char* noise_simplex_cl_source =
"#define MAX_RANK 3                                                            \n"
"                                                                              \n"
"float2                                                                        \n"
"philox (uint2 st,                                                             \n"
"        uint k)                                                               \n"
"{                                                                             \n"
"  ulong p;                                                                    \n"
"  int   i;                                                                    \n"
"                                                                              \n"
"  for (i = 0 ; i < 3 ; i += 1)                                                \n"
"    {                                                                         \n"
"      p = st.x * 0xcd9e8d57ul;                                                \n"
"                                                                              \n"
"      st.x = ((uint)(p >> 32)) ^ st.y ^ k;                                    \n"
"      st.y = (uint)p;                                                         \n"
"                                                                              \n"
"      k += 0x9e3779b9u;                                                       \n"
"    }                                                                         \n"
"                                                                              \n"
"  return convert_float2(st) / 2147483648.0f - 1.0f;                           \n"
"}                                                                             \n"
"                                                                              \n"
"__kernel void kernel_noise (__global float *out,                              \n"
"                            const int x_0,                                    \n"
"                            const int y_0,                                    \n"
"                            const uint iterations,                            \n"
"                            const float scale,                                \n"
"                            const uint seed)                                  \n"
"{                                                                             \n"
"  const int gidx = get_global_id(0);                                          \n"
"  const int gidy = get_global_id(1);                                          \n"
"                                                                              \n"
"  float  c, d, m;                                                             \n"
"  float2 p;                                                                   \n"
"  int    j;                                                                   \n"
"                                                                              \n"
"  for (j = 0, m = 0, c = 1, d = scale;                                        \n"
"       j < iterations;                                                        \n"
"       c *= 2, d *= 2, j += 1)                                                \n"
"    {                                                                         \n"
"      float s, t, n;                                                          \n"
"      float2 g[3], u[3], i, di;                                               \n"
"      int k;                                                                  \n"
"                                                                              \n"
"      p = (float2)(gidx + x_0, gidy + y_0) * d;                               \n"
"                                                                              \n"
"      /* Skew the input point and find the lowest corner of the containing    \n"
"         simplex. */                                                          \n"
"                                                                              \n"
"      s = (p.x + p.y) * (sqrt(3.0f) - 1) / 2;                                 \n"
"      i = floor(p + s);                                                       \n"
"                                                                              \n"
"      /* Calculate the (unskewed) distance between the input point and all    \n"
"         simplex corners. */                                                  \n"
"                                                                              \n"
"      s = (i.x + i.y) * (3 - sqrt(3.0f)) / 6;                                 \n"
"      u[0] = p - i + s;                                                       \n"
"                                                                              \n"
"      di = u[0].x >= u[0].y ? (float2)(1, 0) : (float2)(0, 1);                \n"
"                                                                              \n"
"      u[1] = u[0] - di + (3 - sqrt(3.0f)) / 6;                                \n"
"      u[2] = u[0] - 1 + (3 - sqrt(3.0f)) / 3;                                 \n"
"                                                                              \n"
"      /* Calculate gradients for each corner vertex. We convert to            \n"
"       * signed int first to avoid implementation-defined behavior for        \n"
"       * out-of-range values.  See section 6.2.3.3 of the OpenCL              \n"
"       * specification. */                                                    \n"
"                                                                              \n"
"      g[0] = philox(convert_uint2(convert_int2(i)), seed);                    \n"
"      g[1] = philox(convert_uint2(convert_int2(i + di)), seed);               \n"
"      g[2] = philox(convert_uint2(convert_int2(i + 1)), seed);                \n"
"                                                                              \n"
"      for (k = 0, n = 0 ; k < 3 ; k += 1)                                     \n"
"        {                                                                     \n"
"          t = 0.5f - dot(u[k], u[k]);                                         \n"
"                                                                              \n"
"          if (t > 0)                                                          \n"
"            {                                                                 \n"
"              t *= t;                                                         \n"
"              n += t * t * dot(g[k], u[k]);                                   \n"
"            }                                                                 \n"
"        }                                                                     \n"
"                                                                              \n"
"      m += 70 * n / c;                                                        \n"
"    }                                                                         \n"
"                                                                              \n"
"  out[gidy * get_global_size(0) + gidx] = m;                                  \n"
"}                                                                             \n"
;