     // free memory on device
     cudaFree(d_src_data);
     cudaFree(d_dst_data);
-    cudaFree(d_pixel_equ);
+    // cudaFree(d_pixel_equ);
     cudaFree(d_block_equ);
     cudaFree(d_grid_equ);
 

    // allocate memory on device
    uint8_t* d_src_data = NULL;
    uint8_t* d_dst_data = NULL;
    uint8_t* d_pixel_equ = NULL;
    uint8_t* d_block_equ = NULL;
    uint8_t* d_grid_equ = NULL;

    unsigned int data_size = (*src_bitmap)->width * (*src_bitmap)->height * sizeof(uint8_t);
    unsigned int pixel_equ_size = ((*src_bitmap)->width - padding.left - padding.right) * ((*src_bitmap)->height - padding.top - padding.bottom) * sizeof(uint8_t);
    unsigned int block_equ_size = grid_dim.x * grid_dim.y * sizeof(uint8_t);
    unsigned int grid_equ_size = 1 * sizeof(uint8_t);

    cudaError d_src_malloc_success = cudaMalloc((void**) &d_src_data, data_size);
    cudaError d_dst_malloc_success = cudaMalloc((void**) &d_dst_data, data_size);
    cudaError d_pixel_equ_malloc_success = cudaMalloc((void**) &d_pixel_equ, pixel_equ_size);
    cudaError d_block_equ_malloc_success = cudaMalloc((void**) &d_block_equ, block_equ_size);
    cudaError d_grid_equ_malloc_success = cudaMalloc((void**) &d_grid_equ, grid_equ_size);

    assert((d_src_malloc_success == cudaSuccess) && "Error: could not allocate memory for d_src_data");
    assert((d_dst_malloc_success == cudaSuccess) && "Error: could not allocate memory for d_dst_data");
    assert((d_pixel_equ_malloc_success == cudaSuccess) && "Error: could not allocate memory for d_pixel_equ");
    assert((d_block_equ_malloc_success == cudaSuccess) && "Error: could not allocate memory for d_block_equ");
    assert((d_grid_equ_malloc_success == cudaSuccess) && "Error: could not allocate memory for d_grid_equ");


    do {
        // 2D grid & 2D block
        skeletonize_pass<<<grid_dim, block_dim>>>(d_src_data, d_dst_data, (*src_bitmap)->width, padding);
        // pixel_equality<<<grid_dim, block_dim>>>(d_src_data, d_dst_data, d_pixel_equ, (*src_bitmap)->width, padding);

        // 1D grid & 1D block reduction
        // and_reduction(grid_dim, block_dim, d_block_equ, d_grid_equ, block_equ_size);

        // ===================================

        // First reduction from d_src_data and d_dst_data to d_block_equ
        and_reduction_2D<<<grid_dim, block_dim, block_dim.x * block_dim.y * sizeof(uint8_t)>>>(d_src_data, d_dst_data, d_block_equ, (*src_bitmap)->width, data_size, padding);

        // iterative reductions of block_equ: if the number of blocks in the grid
        // exceeds the number of threads in a block, then we cannot go to the "leaf"
        // reduction where 1 block is only running in the grid, and must perform
        // another multi-block reduction.

        // used for reduction operation, since we have to modify the grid sizes
        dim3 reduction_grid_dim(grid_dim.x, grid_dim.y);

        while ((reduction_grid_dim.x * reduction_grid_dim.y) > (block_dim.x * block_dim.y)) {
            reduction_grid_dim.x = (unsigned int) ceil(reduction_grid_dim.x / ((double) block_dim.x));
            reduction_grid_dim.y = (unsigned int) ceil(reduction_grid_dim.y / ((double) block_dim.y));

            and_reduction_1D<<<reduction_grid_dim.x * reduction_grid_dim.y, block_dim.x * block_dim.y, block_dim.x * block_dim.y * sizeof(uint8_t)>>>(d_block_equ, d_block_equ, block_equ_size);
        }

        and_reduction_1D<<<1, block_dim.x * block_dim.y, block_dim.x * block_dim.y * sizeof(uint8_t)>>>(d_block_equ, d_grid_equ, block_equ_size);

        // ================

        // bring d_grid_equ back from device
        cudaMemcpy(&grid_equ, d_grid_equ, grid_equ_size, cudaMemcpyDeviceToHost);
