   {
     sum = 0;
     for(j = 0; j <= i; ++j) {
-	  sum += *sigma++ * (data[j] - mean[j]);
+      sum += *sigma++ * (data[j] - mean[j]);
     }
     discrim += sum * sum;
   }

  int data_per = max_threads / params_per;
  // at least 16 data points per block
  while (data_per < 16) {
	params_per /= 2;
	data_per *= 2;
  }

  int half_warp = 16;

  // hide your kids, hide your wife (auto-tuning the GPU)
  while (1) {
	while (compute_shmem(data, params, params_per, data_per) > max_smem) {
	  if (data_per == 0)
		break;

	  if (params_per < half_warp) {
		// no half-warp to protect, maximize data per block
		params_per /= 2;
	  }
	  else {
		// keep the half-warp
		if (data_per > 1) {
		  // TODO: should ask what is the half warp size instead of 16
		  if (params_per > 16) {
			params_per /= 2;
			data_per *= 2;
		  }
		  else {
			// have to do less data
			data_per /= 2;
		  }
		}
		else {
		  params_per /= 2;
		  data_per *= 2;
		}
	  }
	}
	// can't fit max_block_params sets of parameters into the shared memory,
	// uh oh
	if (data_per == 0) {
	  params_per /= 2;

	  // start over the tuning
	  continue;
	}
	else break;
  }

  // possible to squeeze more data?
  while (compute_shmem(data, params, params_per, 2 * data_per) <= max_smem)
	if (2 * data_per * params_per <= max_threads)
	  data_per *= 2;
	else
	  break;

  info->data_per_block = data_per;
  info->params_per_block = params_per;


/*
__device__ void _write_results(PMatrix* data, PMatrix* params,
                               float* output, float* sh_result,
                               int thidx, int thidy,
                               int tid)
{
  // write out in other order to coalesce
  // transpose! to get it to coalesce
