title: OpenCL nested loop produces unexpected results
body: <p>I'm a bit new to OpenCL/C99, but can't understand why the the two kernels below give different results.  X has been initialized to zeros, but requires "re-zeroing" at each outer loop step else incorrect results are obtained (see plot).  Notice that I have not invoked any parallelism here; to the best of my knowledge, this is entirely serial:</p>
  
  <p>content of <code>kernels.cl</code>:</p>
  
  <pre><code>#pragma OPENCL EXTENSION cl_khr_fp64: enable
  #define __CL_ENABLE_EXCEPTIONS
  
  __kernel void nested_sum_succeeds(
      int L,
      __global __read_only double* X,
      __global double* Y
  )
  {
      for (int k=0; k&lt;=L; k++) {
          Y[k] = 0;
          for (int j=0; j&lt;=k; j++) {
              Y[k] += X[j];
          }
      }
  }
  
  __kernel void nested_sum_fails(
      int L,
      __global __read_only double* X,
      __global double* Y
  )
  {
      for (int k=0; k&lt;=L; k++) {
          // Y[k] = 0;
          for (int j=0; j&lt;=k; j++) {
              Y[k] += X[j];
          }
      }
  }
  </code></pre>
  
  <p>content of <code>script.py</code>:</p>
  
  <pre><code>import numpy as np
  import pyopencl as cl
  import pyopencl.array as cl_array
  import matplotlib.pyplot as plt
  
  with open("./kernels.cl") as fp:
      prog_str = fp.read()
  ctx = cl.create_some_context()
  queue = cl.CommandQueue(ctx)
  prog = cl.Program(ctx, prog_str).build()
  
  L = 1000
  X = np.linspace(0, 10, L)
  X_dev = cl_array.to_device(queue, X)
  Y_succeeds_dev = cl_array.to_device(queue, np.zeros(shape=X.shape, dtype=np.float64))
  Y_fails_dev = cl_array.to_device(queue, np.zeros(shape=X.shape, dtype=np.float64))
  
  
  nested_sum_succeeds = prog.nested_sum_succeeds
  nested_sum_succeeds.set_scalar_arg_dtypes([
      np.int64,
      None,
      None,
  ])
  
  nested_sum_succeeds(
      queue,
      (len(X),),
      None,
      L,
      X_dev.data,
      Y_succeeds_dev.data,
  )
  
  
  nested_sum_fails = prog.nested_sum_fails
  nested_sum_fails.set_scalar_arg_dtypes([
      np.int64,
      None,
      None,
  ])
  
  nested_sum_fails(
      queue,
      (len(X),),
      None,
      L,
      X_dev.data,
      Y_fails_dev.data,
  )
  
  np.allclose(Y_succeeds_dev.get(), Y_fails_dev.get()) #False
  
  plt.ion()
  plt.plot(Y_succeeds_dev.get())
  plt.plot(Y_fails_dev.get())
  </code></pre>
  
  <p>Results:</p>
  
  <p><a href="https://i.stack.imgur.com/zgU0m.png" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/zgU0m.png" alt="enter image description here"></a></p>