 }
 
 
+template<typename T>
+__device__ __host__
+void swap(T& a, T& b)
+{
+    T tmp = b;
+    b = a;
+    a = tmp;
+}
+
+
 // sorts an array of size_t and mirrors the same permutation to two arrays of double
 __device__ __host__
 void sort_sparse_expr(expr_t* const expr)

            first = true;
        }

        std::swap(a, b);
        std::swap(mirror_1_a, mirror_1_b);
        std::swap(mirror_2_a, mirror_2_b);

        i = 0;
        startl = 0;



__global__
void device_layer_create_sparse_exprs(expr_t** expr_array, const double* weights, const double* bias, const size_t* dims, const size_t num_out_neurons, const size_t num_in_neurons)
{
    size_t i = blockIdx.x;

    if(i < num_out_neurons)
    {
        const double* weight_i = weights + i*num_in_neurons;
        const size_t* dims_i = dims + i*num_in_neurons;
        const double bias_i = bias[i];
        expr_array[i] = create_sparse_expr(weight_i, bias_i, dims_i, num_in_neurons);
        sort_sparse_expr(expr_array[i]);
    }
}


void layer_create_sparse_exprs(expr_t** expr_array, const double* weights, const double* bias, const size_t* dims, const size_t num_out_neurons, const size_t num_in_neurons)
{
    double* tmp_weights;
    cudaMalloc((void**) &tmp_weights, num_out_neurons*num_in_neurons*sizeof(double));
    size_t* tmp_dims;
    cudaMalloc((void**) &tmp_dims   , num_out_neurons*num_in_neurons*sizeof(size_t));

    cudaMemcpy(tmp_weights, weights, num_out_neurons*num_in_neurons*sizeof(double), cudaMemcpyHostToDevice);
    cudaMemcpy(tmp_dims   , dims   , num_out_neurons*num_in_neurons*sizeof(size_t), cudaMemcpyHostToDevice);

    double* tmp_bias;
    cudaMalloc((void**) &tmp_bias, num_out_neurons*sizeof(double));
    cudaMemcpy(tmp_bias, bias, num_out_neurons*sizeof(double), cudaMemcpyHostToDevice);

    device_layer_create_sparse_exprs<<<num_out_neurons, 1>>>(expr_array, tmp_weights, tmp_bias, tmp_dims, num_out_neurons, num_in_neurons);

    cudaFree(tmp_weights);
    cudaFree(tmp_dims);
    cudaFree(tmp_bias);
}


