             
     }
     
+    if (origlhs != lhs) {
+        PERFWARNING_LOG("Scatter the result back");
+        GHOST_CALL_RETURN(origlhs->fromVec(origlhs,lhs,0,0));
+    }
+    
     return GHOST_SUCCESS;
 }
     
     template<typename dt1, typename dt2>
-static ghost_error_t ghost_cu_crsspmmv_cm_tmpl(ghost_sparsemat_t *mat, ghost_densemat_t * lhs, ghost_densemat_t * rhs, ghost_spmv_flags_t options, va_list argp, cusparse_crs_spmmv_cm_kernel_t crskernel)
+static ghost_error_t ghost_cu_crsspmmv_cm_tmpl(ghost_sparsemat_t *mat, ghost_densemat_t * lhs, ghost_densemat_t *origlhs, ghost_densemat_t * rhs, ghost_spmv_flags_t options, va_list argp, cusparse_crs_spmmv_cm_kernel_t crskernel)
 {
     cusparseHandle_t cusparse_handle;
     cusparseMatDescr_t descr;

    return GHOST_SUCCESS;
}

ghost_error_t ghost_cu_crs_spmv_selector(ghost_sparsemat_t *mat, ghost_densemat_t * lhs, ghost_densemat_t * rhs, ghost_spmv_flags_t options, va_list argp)
{
    ghost_error_t ret = GHOST_SUCCESS;
    ghost_densemat_t *lhscompact, *rhscompact;
    
    if (mat->traits->datatype != lhs->traits.datatype) {
        ERROR_LOG("Mixed data types not implemented!");
        ret = GHOST_ERR_NOT_IMPLEMENTED;
        goto err;
    }
    if (lhs->traits.flags & GHOST_DENSEMAT_SCATTERED) {
        INFO_LOG("Cloning (and compressing) lhs before operation");
        GHOST_CALL_GOTO(lhs->clone(lhs,&lhscompact,lhs->traits.nrows,0,lhs->traits.ncols,0),err,ret);
    } else {
        lhscompact = lhs;
    }
    if (rhs->traits.flags & GHOST_DENSEMAT_SCATTERED) {
        INFO_LOG("Cloning (and compressing) v2 before operation");
        GHOST_CALL_GOTO(rhs->clone(rhs,&rhscompact,rhs->traits.nrows,0,rhs->traits.ncols,0),err,ret);
    } else {
        rhscompact = rhs;
    }


    INFO_LOG("Calling cuSparse CRS SpMV %d",options);

    if (lhs->traits.ncols == 1) {
        if (mat->traits->datatype & GHOST_DT_DOUBLE) {
            if (mat->traits->datatype & GHOST_DT_REAL) {
                GHOST_CALL_GOTO((ghost_cu_crsspmv_tmpl<double,double>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmv_kernel_t)cusparseDcsrmv)),err,ret);
            } else {
                GHOST_CALL_GOTO((ghost_cu_crsspmv_tmpl<cuDoubleComplex,complex double>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmv_kernel_t)cusparseZcsrmv)),err,ret);
            }
        } else {
            if (mat->traits->datatype & GHOST_DT_REAL) {
                GHOST_CALL_GOTO((ghost_cu_crsspmv_tmpl<float,float>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmv_kernel_t)cusparseScsrmv)),err,ret);
            } else {
                GHOST_CALL_GOTO((ghost_cu_crsspmv_tmpl<cuFloatComplex,complex float>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmv_kernel_t)cusparseCcsrmv)),err,ret);
            }
        }
    } else if (rhs->traits.storage == GHOST_DENSEMAT_COLMAJOR) {
        INFO_LOG("Calling col-major cuSparse CRS SpMMV");
        if (mat->traits->datatype & GHOST_DT_DOUBLE) {
            if (mat->traits->datatype & GHOST_DT_REAL) {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_cm_tmpl<double,double>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_cm_kernel_t)cusparseDcsrmm)),err,ret);
            } else {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_cm_tmpl<cuDoubleComplex,complex double>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_cm_kernel_t)cusparseZcsrmm)),err,ret);
            }
        } else {
            if (mat->traits->datatype & GHOST_DT_REAL) {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_cm_tmpl<float,float>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_cm_kernel_t)cusparseScsrmm)),err,ret);
            } else {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_cm_tmpl<cuFloatComplex,complex float>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_cm_kernel_t)cusparseCcsrmm)),err,ret);
            }
        }
    } else {
        if (lhs->traits.storage == GHOST_DENSEMAT_ROWMAJOR) {
            PERFWARNING_LOG("Need to memtranspose the result densemat!");
            ghost_densemat_t *lhscompactnew;
            ghost_densemat_traits_t transtraits = lhscompact->traits;
            transtraits.storage = GHOST_DENSEMAT_COLMAJOR;
            GHOST_CALL_GOTO(ghost_densemat_create(&lhscompactnew,lhscompact->context,transtraits),err,ret);
            GHOST_CALL_GOTO(lhscompactnew->fromVec(lhscompactnew,lhscompact,0,0),err,ret);
            if (lhscompact != lhs) {
                lhscompact->destroy(lhscompact);
            }
            lhscompact = lhscompactnew;
        }
        INFO_LOG("Calling row-major cuSparse CRS SpMMV");
        if (mat->traits->datatype & GHOST_DT_DOUBLE) {
            if (mat->traits->datatype & GHOST_DT_REAL) {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_rm_tmpl<double,double>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_rm_kernel_t)cusparseDcsrmm2)),err,ret);
            } else {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_rm_tmpl<cuDoubleComplex,complex double>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_rm_kernel_t)cusparseZcsrmm2)),err,ret);
            }
        } else {
            if (mat->traits->datatype & GHOST_DT_REAL) {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_rm_tmpl<float,float>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_rm_kernel_t)cusparseScsrmm2)),err,ret);
            } else {
                GHOST_CALL_GOTO((ghost_cu_crsspmmv_rm_tmpl<cuFloatComplex,complex float>(mat,lhscompact,rhscompact,options,argp,(cusparse_crs_spmmv_rm_kernel_t)cusparseCcsrmm2)),err,ret);
            }
        }
    }

            
    }
    
    if (origlhs != lhs) {
        PERFWARNING_LOG("Scatter the result back");
        GHOST_CALL_RETURN(origlhs->fromVec(origlhs,lhs,0,0));
    }

    return GHOST_SUCCESS;
}

    template<typename dt1, typename dt2>
static ghost_error_t ghost_cu_crsspmmv_rm_tmpl(ghost_sparsemat_t *mat, ghost_densemat_t * lhs, ghost_densemat_t *origlhs, ghost_densemat_t * rhs, ghost_spmv_flags_t options, va_list argp, cusparse_crs_spmmv_rm_kernel_t crskernel)
{
    cusparseHandle_t cusparse_handle;
    cusparseMatDescr_t descr;
