     const size_t last_y = blockIdx.y;
     const size_t last_z = blockIdx.z;
 
+    const size_t inp_z = threadIdx.x;
+    const size_t y_shift = threadIdx.y;
+    const size_t x_shift = threadIdx.z;
+
     const size_t n = last_x*gridDim.y*gridDim.z + last_y*gridDim.z + last_z;
 
     const long int min_out_x = offset_x + last_x*shift_x;


elina_manager_t* fppoly_manager_alloc()
{
    std::cout << "This is the GPU version of fppoly!" << std::endl;
    results_calculated = false;
    output_counter = 1;


{
    const size_t n = blockIdx.x;

    const size_t inp_z = threadIdx.x;
    const size_t y_shift = threadIdx.y;
    const size_t x_shift = threadIdx.z;

    float_type tmp1, tmp2;
    float_type maxRes, maxMul;

    if(x_shift*filter_size_y*input_size_z + y_shift*input_size_z + inp_z < filter_size_x*filter_size_y*input_size_z)
    {
        for(size_t out_x = 0; out_x < output_size_x; out_x++)
        {
            for(size_t out_y = 0; out_y < output_size_y; out_y++)
            {
                const long int x_val = out_x*stride_x + x_shift - pad_x;
                const long int y_val = out_y*stride_y + y_shift - pad_y;

                if(!((y_val < 0) || (y_val >= (long int)input_size_y)))
                {
                    if(!((x_val < 0) || (x_val >= (long int)input_size_x)))
                    {
                        const size_t mat_in = x_val*input_size_y*input_size_z + y_val*input_size_z + inp_z;
                        const size_t b = n*input_size_x*input_size_y*input_size_z + mat_in;

                        float_type inf_coeff = res_inf_coeff[b];
                        float_type sup_coeff = res_sup_coeff[b];

                        for(size_t out_z = 0; out_z < output_size_z; out_z++)
                        {
                            const size_t mat_out = out_x*output_size_y*output_size_z + out_y*output_size_z + out_z;

                            const size_t a = n*output_size_x*output_size_y*output_size_z + mat_out;

                            const float_type prev_inf_coeff = expr_inf_coeff[a];
                            const float_type prev_sup_coeff = expr_sup_coeff[a];

                            if((prev_inf_coeff != 0) || (prev_sup_coeff != 0))
                            {
                                const size_t filter_index = out_z*filter_size_x*filter_size_y*input_size_z + x_shift*filter_size_y*input_size_z + y_shift*input_size_z + inp_z;

                                const float_type aux_coeff = aux_coeffs[filter_index];
                                elina_double_interval_mul_expr_coeff(&tmp1, &tmp2, prev_inf_coeff, prev_sup_coeff, -aux_coeff, aux_coeff);

                                maxRes = fmax(fabs(inf_coeff), fabs(sup_coeff));
                                maxMul = fmax(fabs(tmp1), fabs(tmp2));

                                inf_coeff = inf_coeff + tmp1 + (maxRes + maxMul)*ulp;
                                sup_coeff = sup_coeff + tmp2 + (maxRes + maxMul)*ulp;
                            }
                        }

                        res_inf_coeff[b] = inf_coeff;
                        res_sup_coeff[b] = sup_coeff;
                    }
                }

                __syncthreads();
            }
        }
    }
