     }
 }
 
+#define TASKCFG static __inline__ __device__
+
+struct doubleSum {
+  TASKCFG double init() {return 0.0;}
+  TASKCFG double reduce(double x, double acc){return x + acc;}
+  TASKCFG double f(cuDoubleComplex c){return c.x;}
+};
+
+template <unsigned int blockSize, bool nIsPow2>
+__device__ __inline__
+void reduceSinglePass_dev(const cuDoubleComplex *g_idata, double *g_odata, unsigned int n){
+  reduceSinglePass_devGen<blockSize, nIsPow2, double, cuDoubleComplex, doubleSum>(g_idata, g_odata, n);
+}
+
 #endif // #ifndef _REDUCE_KERNEL_H_

    See the CUDA SDK "reduction" sample for more information.
*/

template <unsigned int blockSize>
__device__ __inline__ void
reduceBlock(volatile double *sdata, double mySum, const unsigned int tid)
{
    sdata[tid] = mySum;
    __syncthreads();

    // do reduction in shared mem
    if (blockSize >= 512)
    {
        if (tid < 256) sdata[tid] = mySum = mySum + sdata[tid + 256];
        __syncthreads();
    }

    if (blockSize >= 256)
    {
        if (tid < 128) sdata[tid] = mySum = mySum + sdata[tid + 128];
        __syncthreads();
    }

    if (blockSize >= 128)
    {
        if (tid <  64) sdata[tid] = mySum = mySum + sdata[tid +  64];
        __syncthreads();
    }

    if (tid < 32)
    {
        sdata[tid] = mySum = mySum + sdata[tid + 32];
        sdata[tid] = mySum = mySum + sdata[tid + 16];
        sdata[tid] = mySum = mySum + sdata[tid +  8];
        sdata[tid] = mySum = mySum + sdata[tid +  4];
        sdata[tid] = mySum = mySum + sdata[tid +  2];
        sdata[tid] = mySum = mySum + sdata[tid +  1];
    }
}

template <unsigned int blockSize, bool nIsPow2>
__device__ __inline__ void
reduceBlocks(const cuDoubleComplex *g_idata, double *g_odata, unsigned int n)
{
    extern __shared__ double sdata[];

    // perform first level of reduction,
    // reading from global memory, writing to shared memory
    unsigned int tid = threadIdx.x;
    unsigned int i = blockIdx.x*(blockSize*2) + threadIdx.x;
    unsigned int gridSize = blockSize*2*gridDim.x;
    double mySum = 0;

    // we reduce multiple elements per thread.  The number is determined by the
    // number of active thread blocks (via gridDim).  More blocks will result
    // in a larger gridSize and therefore fewer elements per thread
    while (i < n)
    {
        mySum += g_idata[i].x;

        // ensure we don't read out of bounds -- this is optimized away for powerOf2 sized arrays
        if (nIsPow2 || i + blockSize < n)
            mySum += g_idata[i+blockSize].x;

        i += gridSize;
    }

    // do reduction in shared mem
    reduceBlock<blockSize>(sdata, mySum, tid);

    // write result for this block to global mem
    if (tid == 0) g_odata[blockIdx.x] = sdata[0];

//
// For more details on the reduction algorithm (notably the multi-pass approach), see
// the "reduction" sample in the CUDA SDK.
template <
    unsigned int blockSize
  , bool nIsPow2
  , typename typ
  , typename styp
  , class task
  >
__device__  __inline__
void reduceSinglePass_devGen(const styp *g_idata, typ *g_odata, unsigned int n)
{

    //
    // PHASE 1: Process all inputs assigned to this block
    //

    reduceBlocksGen<blockSize, nIsPow2, typ, styp, task>(g_idata, g_odata, n);

    //
    // PHASE 2: Last block finished will process all partial sums
