 }
 
 //--------------------------------------------------------------------------------------------------------------------------------
-__global__ void calcRho_kernel(uint* indexOfIndex,
-                               Real4* sortedPosRad,
+__global__ void calcRho_kernel(Real4* sortedPosRad,
                                Real4* sortedRhoPreMu,
                                Real4* sortedRhoPreMu_old,
                                uint* cellStart,

                        uint* gridMarkerIndex,
                        uint* cellStart,
                        uint* cellEnd,
                        const size_t numFluidMarkers,
                        const size_t numBounMarkers,
                        const size_t numRigidMarkers,
                        volatile bool* isErrorD) {
    uint id = blockIdx.x * blockDim.x + threadIdx.x;
    if (id >= numFluidMarkers + numRigidMarkers)
        return;

    uint index = indexOfIndex[id];
    if (sortedRhoPreMu[index].w > -0.5 && sortedRhoPreMu[index].w < 0.5)
        return;

    // Check the activity of this particle
    uint originalIndex = gridMarkerIndex[index];
    uint activity = activityIdentifierD[originalIndex];
    if(activity == 0)
        return;

    Real hA = sortedPosRad[index].w;

    //------------------------------------------------------------------------
    // thread per particle
    uint numBlocks, numThreads;
    computeGridSize((int)numObjectsH->numAllMarkers, 256, numBlocks, numThreads);

    // Execute the kernel
    thrust::device_vector<Real4> sortedDerivVelRho(numObjectsH->numAllMarkers);
