extern "C" {
    __global__ void
    butterfly(long* param, float* x, float* w)
    {
        unsigned long gid = blockIdx.x*blockDim.x+threadIdx.x; // (unsigned long)s->get_param(0);
        unsigned long nid = blockIdx.y*blockDim.y+threadIdx.y; // (unsigned long)s->get_param(1);
        
        long n = param[0];
        long direction_flag = param[1];
        long iter = param[2];
        
        int butterflySize = 1 << (iter-1);
        int butterflyGrpDist = 1 << iter;
        int butterflyGrpNum = n >> iter;
        int butterflyGrpBase = (gid >> (iter-1))*(butterflyGrpDist);
        int butterflyGrpOffset = gid & (butterflySize-1);
        
        int a = nid * n + butterflyGrpBase + butterflyGrpOffset;
        int b = a + butterflySize;
        
        int l = butterflyGrpNum * butterflyGrpOffset;
        
        float xa[2], xb[2], xbxx[2], xbyy[2], wab[2], wayx[2], wbyx[2], resa[2], resb[2];
        
        xa[0] = x[2*a];
        xa[1] = x[2*a+1];
        xb[0] = x[2*b];
        xb[1] = x[2*b+1];
        xbxx[0] = xbxx[1] = xb[0];
        xbyy[0] = xbyy[1] = xb[1];
        
        wab[0] = w[2*l];
        if(direction_flag == 0x80000000) {
            wab[1] = -w[2*l+1];
        } else {
            wab[1] = w[2*l+1];
        }
        
        wayx[0] = -wab[1];
        wayx[1] = wab[0];
        
        wbyx[0] = wab[1];
        wbyx[1] = -wab[0];
        
        resa[0] = xa[0] + xbxx[0]*wab[0] + xbyy[0]*wayx[0];
        resa[1] = xa[1] + xbxx[1]*wab[1] + xbyy[1]*wayx[1];
        
        resb[0] = xa[0] - xbxx[0]*wab[0] + xbyy[0]*wbyx[0];
        resb[1] = xa[1] - xbxx[1]*wab[1] + xbyy[1]*wbyx[1];

        x[2*a] = resa[0];
        x[2*a+1] = resa[1];
        x[2*b] = resb[0];
        x[2*b+1] = resb[1];
    }
}