title: How to structure this OpenCL brute-force code
body: <p>I'm just starting to play about with OpenCL, and I'm stuck with how to structure the program in a reasonably efficient manner (mainly avoiding lots of transferring of data to/from the GPU or wherever the work is being done)</p>
  <p>What I'm trying to do is, given:</p>
  <pre><code>v = r*i + b*j + g*k
  </code></pre>
  <p>..I know <code>v</code> for various values of <code>r</code>, <code>g</code> and <code>b</code>, but <code>i</code>, <code>j</code> and <code>k</code> are unknown. I want to calculate reasonable values for <code>i</code>/<code>j</code>/<code>k</code> via brute force</p>
  <p>In other words, I have a bunch of &quot;raw&quot; RGB pixel values, and I have a desaturated version of these colours. I do not know the weightings (i/j/k) used calculate the desaturated values.</p>
  <p>My initial plan was to:</p>
  <ol>
  <li><p>load the data into a CL buffer (so the input r/g/b values, and the output)</p>
  </li>
  <li><p>have a kernel which takes the three possible matrix values, and the various pixel-data buffers.</p>
  <p>It then performs <code>v = r*i + b*j + g*k</code>, and subtracts the value of <code>v</code> to the known value, and stores this in a &quot;score&quot; buffer</p>
  </li>
  <li><p>Another kernel calculates the RMS error for that value (if the difference is zero for all input values, the values for i/j/k are &quot;correct&quot;)</p>
  </li>
  </ol>
  <p>I have this working (written using Python and PyCL, <a href="https://gist.github.com/1176640/ad53c61e1b016dd4b138bb8b1c701e53085f93c8" rel="nofollow noreferrer">the code is here</a>), but I'm wondering how I can parallelise this chunk of work more (by try multiple i/j/k values at once)</p>
  <p>I issue is, I have the 4 read-only buffers (3 for the input values, 1 for the expected values), but I need a separate &quot;score&quot; buffer for every combination of i/j/k</p>
  <p>Another issue is the RMS calculation is the slowest part, since it's effectively single-threaded (total up all the values in &quot;score&quot; and sqrt() the total)</p>
  <p>Basically, I'm wondering if there's a sensible way to structure such a program.</p>
  <p>It seems like a task well-suited to OpenCL - hopefully the description of my goal wasn't too convoluted! As mentioned, <a href="https://gist.github.com/1176640/ad53c61e1b016dd4b138bb8b1c701e53085f93c8" rel="nofollow noreferrer">my current code is here</a>, and in case it is clearer, this is the Python version of what I'm trying to do:</p>
  <pre><code>import sys
  import math
  import random
  
  
  def make_test_data(w = 128, h = 128):
      in_r, in_g, in_b = [], [], []
  
      print &quot;Make raw data&quot;
      for x in range(w):
          for y in range(h):
              in_r.append(random.random())
              in_g.append(random.random())
              in_b.append(random.random())
  
      # the unknown values
      mtx = [random.random(), random.random(), random.random()]
  
      print &quot;Secret numbers were: %s&quot; % mtx
  
      out_r = [(r*mtx[0] + g*mtx[1] + b*mtx[2]) for (r, g, b) in zip(in_r, in_g, in_b)]
  
      return {'in_r': in_r, 'in_g': in_g, 'in_b': in_b,
              'expected_r': out_r}
  
  
  def score_matrix(ir, ig, ib, expected_r, mtx):
      ms = 0
      for i in range(len(ir)):
          val = ir[i] * mtx[0] + ig[i] * mtx[1] + ib[i] * mtx[2]
          ms += abs(val - expected_r[i]) ** 2
      rms = math.sqrt(ms / float(len(ir)))
      return rms
  
  
  # Make random test data
  test_data = make_test_data(16, 16)
  
  
  lowest_rms = sys.maxint
  closest = []
  
  divisions = 10
  for possible_r in range(divisions):
      for possible_g in range(divisions):
          for possible_b in range(divisions):
  
              pr, pg, pb = [x / float(divisions-1) for x in (possible_r, possible_g, possible_b)]
  
              rms = score_matrix(
                  test_data['in_r'], test_data['in_g'], test_data['in_b'], 
                  test_data['expected_r'],
                  mtx = [pr, pg, pb])
  
              if rms &lt; lowest_rms:
                  closest = [pr, pg, pb]
                  lowest_rms = rms
  
  print closest
  </code></pre>