title: Pycuda Blocks and Grids to work with big datas
body: <p>I need help to know the size of my blocks and grids.
  I'm building a python app to perform metric calculations based on scipy as: Euclidean distance, Manhattan, Pearson, Cosine, joined other.</p>
  <p>The project is <a href="https://github.com/vinigracindo/pycudaDistances" rel="nofollow noreferrer">PycudaDistances</a>.</p>
  <p>It seems to work very well with small arrays. When I perform a more exhaustive test, unfortunately it did not work. I downloaded movielens set (<a href="http://www.grouplens.org/node/73" rel="nofollow noreferrer">http://www.grouplens.org/node/73</a>).</p>
  <p>Using <code>Movielens</code> 100k, I declared an array with shape (943, 1682). That is, users are 943 and 1682 films evaluated. The films not by a classifier user I configured the value to 0.</p>
  <p>With a much larger array algorithm no longer works. I face the following error:</p>
  <blockquote>
  <p>pycuda._driver.LogicError: cuFuncSetBlockShape failed: invalid value.</p>
  </blockquote>
  <p>Researching this error, I found an explanation of telling Andrew that supports 512 threads to join and to work with larger blocks it is necessary to work with blocks and grids.</p>
  <p>I wanted a help to adapt the algorithm Euclidean distance arrays to work from small to giant arrays.</p>
  <pre><code>def euclidean_distances(X, Y=None, inverse=True):
      X, Y = check_pairwise_arrays(X,Y)
      rows = X.shape[0]
      cols = Y.shape[0]
      solution = numpy.zeros((rows, cols))
      solution = solution.astype(numpy.float32)
  
      kernel_code_template = &quot;&quot;&quot;
      #include &lt;math.h&gt;
      
      __global__ void euclidean(float *x, float *y, float *solution) {
  
          int idx = threadIdx.x + blockDim.x * blockIdx.x;
          int idy = threadIdx.y + blockDim.y * blockIdx.y;
          
          float result = 0.0;
          
          for(int iter = 0; iter &lt; %(NDIM)s; iter++) {
              
              float x_e = x[%(NDIM)s * idy + iter];
              float y_e = y[%(NDIM)s * idx + iter];
              result += pow((x_e - y_e), 2);
          }
          int pos = idx + %(NCOLS)s * idy;
          solution[pos] = sqrt(result);
      }
      &quot;&quot;&quot;
      kernel_code = kernel_code_template % {
          'NCOLS': cols,
          'NDIM': X.shape[1]
      }
  
      mod = SourceModule(kernel_code)
  
      func = mod.get_function(&quot;euclidean&quot;)
      func(drv.In(X), drv.In(Y), drv.Out(solution), block=(cols, rows, 1))
  
      return numpy.divide(1.0, (1.0 + solution)) if inverse else solution
  </code></pre>
  <p>For more details see: <a href="https://github.com/vinigracindo/pycudaDistances/blob/master/distances.py" rel="nofollow noreferrer">https://github.com/vinigracindo/pycudaDistances/blob/master/distances.py</a></p>