From 76cef2d8cc36d5ab06cfb9142a6ccaeef5908515 Mon Sep 17 00:00:00 2001 From: Bruno Turcksin Date: Fri, 29 May 2020 19:47:42 +0000 Subject: [PATCH] Fix a bug where only 1 CUDAWrappers::MatrixFree object was valid any given time --- include/deal.II/base/cuda_size.h | 13 + .../deal.II/matrix_free/cuda_fe_evaluation.h | 6 +- .../deal.II/matrix_free/cuda_matrix_free.h | 15 + .../matrix_free/cuda_matrix_free.templates.h | 111 ++++--- .../matrix_free/cuda_tensor_product_kernels.h | 277 ++++++++++-------- tests/cuda/cuda_evaluate_1d_shape.cu | 21 +- tests/cuda/cuda_evaluate_2d_shape.cu | 25 +- tests/cuda/matrix_free_no_index_initialize.cu | 2 +- 8 files changed, 279 insertions(+), 191 deletions(-) diff --git a/include/deal.II/base/cuda_size.h b/include/deal.II/base/cuda_size.h index 09e12eacbe..30a11dc924 100644 --- a/include/deal.II/base/cuda_size.h +++ b/include/deal.II/base/cuda_size.h @@ -38,6 +38,19 @@ namespace CUDAWrappers * Define the number of threads in a warp. */ constexpr int warp_size = 32; + + /** + * Define the largest finite element degree that can be solved using + * CUDAWrappers::MatrixFree. Changing this number will affect the amount of + * constant memory being used. + */ + constexpr unsigned int mf_max_elem_degree = 10; + + /** + * Define the maximum number of valid CUDAWrappers::MatrixFree object. + * Changing this number will affect the amount of constant memory being used. + */ + constexpr unsigned int mf_n_concurrent_objects = 5; } // namespace CUDAWrappers DEAL_II_NAMESPACE_CLOSE diff --git a/include/deal.II/matrix_free/cuda_fe_evaluation.h b/include/deal.II/matrix_free/cuda_fe_evaluation.h index f42b87409a..b065e0964a 100644 --- a/include/deal.II/matrix_free/cuda_fe_evaluation.h +++ b/include/deal.II/matrix_free/cuda_fe_evaluation.h @@ -316,6 +316,7 @@ namespace CUDAWrappers types::global_dof_index *local_to_global; unsigned int n_cells; unsigned int padding_length; + const unsigned int mf_object_id; const unsigned int constraint_mask; @@ -343,6 +344,7 @@ namespace CUDAWrappers SharedData *shdata) : n_cells(data->n_cells) , padding_length(data->padding_length) + , mf_object_id(data->id) , constraint_mask(data->constraint_mask[cell_id]) , use_coloring(data->use_coloring) , values(shdata->values) @@ -426,7 +428,7 @@ namespace CUDAWrappers fe_degree, n_q_points_1d, Number> - evaluator_tensor_product; + evaluator_tensor_product(mf_object_id); if (evaluate_val == true && evaluate_grad == true) { evaluator_tensor_product.value_and_gradient_at_quad_pts(values, @@ -463,7 +465,7 @@ namespace CUDAWrappers fe_degree, n_q_points_1d, Number> - evaluator_tensor_product; + evaluator_tensor_product(mf_object_id); if (integrate_val == true && integrate_grad == true) { evaluator_tensor_product.integrate_value_and_gradient(values, diff --git a/include/deal.II/matrix_free/cuda_matrix_free.h b/include/deal.II/matrix_free/cuda_matrix_free.h index 5b78bb91e6..e621bc1033 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.h @@ -183,6 +183,11 @@ namespace CUDAWrappers */ Number *JxW; + /** + * ID of the associated MatrixFree object. + */ + unsigned int id; + /** * Number of cells. */ @@ -215,6 +220,11 @@ namespace CUDAWrappers */ MatrixFree(); + /** + * Destructor. + */ + ~MatrixFree(); + /** * Return the length of the padding. */ @@ -443,6 +453,11 @@ namespace CUDAWrappers const Number val, LinearAlgebra::CUDAWrappers::Vector &dst) const; + /** + * Unique ID associated with the object. + */ + int my_id; + /** * Parallelization scheme used, parallelization over degrees of freedom or * over cells. diff --git a/include/deal.II/matrix_free/cuda_matrix_free.templates.h b/include/deal.II/matrix_free/cuda_matrix_free.templates.h index 5c2c62620c..c7e81aa35f 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.templates.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.templates.h @@ -48,80 +48,85 @@ namespace CUDAWrappers { namespace internal { - // These variables are stored in the device constant memory. - constexpr unsigned int max_elem_degree = 10; + constexpr unsigned int data_array_size = + (mf_max_elem_degree + 1) * (mf_max_elem_degree + 1); + + // Default initialized to false + std::array used_objects; template - using DataArray = NumberType[(max_elem_degree + 1) * (max_elem_degree + 1)]; + using DataArray = NumberType[data_array_size]; - __constant__ double - global_shape_values_d[(max_elem_degree + 1) * (max_elem_degree + 1)]; - __constant__ float - global_shape_values_f[(max_elem_degree + 1) * (max_elem_degree + 1)]; + // These variables are stored in the device constant memory. + // Shape values + __constant__ double global_shape_values_d[mf_n_concurrent_objects] + [data_array_size]; + __constant__ float global_shape_values_f[mf_n_concurrent_objects] + [data_array_size]; + // Shape gradients + __constant__ double global_shape_gradients_d[mf_n_concurrent_objects] + [data_array_size]; + __constant__ float global_shape_gradients_f[mf_n_concurrent_objects] + [data_array_size]; + // for collocation methods + __constant__ double global_co_shape_gradients_d[mf_n_concurrent_objects] + [data_array_size]; + __constant__ float global_co_shape_gradients_f[mf_n_concurrent_objects] + [data_array_size]; template __host__ __device__ inline DataArray & - get_global_shape_values(); + get_global_shape_values(unsigned int i); template <> __host__ __device__ inline DataArray & - get_global_shape_values() + get_global_shape_values(unsigned int i) { - return global_shape_values_d; + return global_shape_values_d[i]; } template <> __host__ __device__ inline DataArray & - get_global_shape_values() + get_global_shape_values(unsigned int i) { - return global_shape_values_f; + return global_shape_values_f[i]; } - __constant__ double - global_shape_gradients_d[(max_elem_degree + 1) * (max_elem_degree + 1)]; - __constant__ float - global_shape_gradients_f[(max_elem_degree + 1) * (max_elem_degree + 1)]; - template __host__ __device__ inline DataArray & - get_global_shape_gradients(); + get_global_shape_gradients(unsigned int i); template <> __host__ __device__ inline DataArray & - get_global_shape_gradients() + get_global_shape_gradients(unsigned int i) { - return global_shape_gradients_d; + return global_shape_gradients_d[i]; } template <> __host__ __device__ inline DataArray & - get_global_shape_gradients() + get_global_shape_gradients(unsigned int i) { - return global_shape_gradients_f; + return global_shape_gradients_f[i]; } // for collocation methods - __constant__ double global_co_shape_gradients_d[(max_elem_degree + 1) * - (max_elem_degree + 1)]; - __constant__ float global_co_shape_gradients_f[(max_elem_degree + 1) * - (max_elem_degree + 1)]; - template __host__ __device__ inline DataArray & - get_global_co_shape_gradients(); + get_global_co_shape_gradients(unsigned int i); template <> __host__ __device__ inline DataArray & - get_global_co_shape_gradients() + get_global_co_shape_gradients(unsigned int i) { - return global_co_shape_gradients_d; + return global_co_shape_gradients_d[i]; } template <> __host__ __device__ inline DataArray & - get_global_co_shape_gradients() + get_global_co_shape_gradients(unsigned int i) { - return global_co_shape_gradients_f; + return global_co_shape_gradients_f[i]; } template @@ -606,10 +611,19 @@ namespace CUDAWrappers : n_dofs(0) , constrained_dofs(nullptr) , padding_length(0) + , my_id(-1) {} + template + MatrixFree::~MatrixFree() + { + free(); + } + + + template void MatrixFree::reinit(const Mapping & mapping, @@ -661,11 +675,12 @@ namespace CUDAWrappers data_copy.local_to_global = local_to_global[color]; data_copy.inv_jacobian = inv_jacobian[color]; data_copy.JxW = JxW[color]; - data_copy.constraint_mask = constraint_mask[color]; + data_copy.id = my_id; data_copy.n_cells = n_cells[color]; data_copy.padding_length = padding_length; data_copy.row_start = row_start[color]; data_copy.use_coloring = use_coloring; + data_copy.constraint_mask = constraint_mask[color]; return data_copy; } @@ -697,6 +712,9 @@ namespace CUDAWrappers constraint_mask.clear(); Utilities::CUDA::free(constrained_dofs); + + internal::used_objects[my_id].store(false); + my_id = -1; } @@ -882,29 +900,44 @@ namespace CUDAWrappers unsigned int size_co_shape_values = n_q_points_1d * n_q_points_1d * sizeof(Number); + // Check if we already a part of the constant memory allocated to us. If + // not, we try to get a block of memory. + bool found_id = false; + while (!found_id) + { + ++my_id; + Assert( + my_id < mf_n_concurrent_objects, + ExcMessage( + "Maximum number of concurrents MatrixFree objects reached. Increase mf_n_concurrent_objects")); + bool f = false; + found_id = + internal::used_objects[my_id].compare_exchange_strong(f, true); + } + cudaError_t cuda_error = - cudaMemcpyToSymbol(internal::get_global_shape_values(), + cudaMemcpyToSymbol(internal::get_global_shape_values(0), shape_info.data.front().shape_values.data(), size_shape_values, - 0, + my_id * internal::data_array_size * sizeof(Number), cudaMemcpyHostToDevice); AssertCuda(cuda_error); if (update_flags & update_gradients) { cuda_error = - cudaMemcpyToSymbol(internal::get_global_shape_gradients(), + cudaMemcpyToSymbol(internal::get_global_shape_gradients(0), shape_info.data.front().shape_gradients.data(), size_shape_values, - 0, + my_id * internal::data_array_size * sizeof(Number), cudaMemcpyHostToDevice); AssertCuda(cuda_error); cuda_error = - cudaMemcpyToSymbol(internal::get_global_co_shape_gradients(), + cudaMemcpyToSymbol(internal::get_global_co_shape_gradients(0), shape_info_co.data.front().shape_gradients.data(), size_co_shape_values, - 0, + my_id * internal::data_array_size * sizeof(Number), cudaMemcpyHostToDevice); AssertCuda(cuda_error); } diff --git a/include/deal.II/matrix_free/cuda_tensor_product_kernels.h b/include/deal.II/matrix_free/cuda_tensor_product_kernels.h index a5871e1377..fa3d28c2c6 100644 --- a/include/deal.II/matrix_free/cuda_tensor_product_kernels.h +++ b/include/deal.II/matrix_free/cuda_tensor_product_kernels.h @@ -59,7 +59,9 @@ namespace CUDAWrappers int n_q_points_1d, typename Number> struct EvaluatorTensorProduct - {}; + { + const int mf_object_id; + }; @@ -82,7 +84,7 @@ namespace CUDAWrappers Utilities::pow(n_q_points_1d, dim); __device__ - EvaluatorTensorProduct(); + EvaluatorTensorProduct(int mf_object_id); /** * Evaluate the values of a finite element function at the quadrature @@ -147,6 +149,8 @@ namespace CUDAWrappers */ __device__ void integrate_value_and_gradient(Number *u, Number *grad_u[dim]); + + const int mf_object_id; }; @@ -157,7 +161,8 @@ namespace CUDAWrappers dim, fe_degree, n_q_points_1d, - Number>::EvaluatorTensorProduct() + Number>::EvaluatorTensorProduct(int object_id) + : mf_object_id(object_id) {} @@ -256,37 +261,31 @@ namespace CUDAWrappers { case 1: { - values<0, true, false, true>(get_global_shape_values(), - u, - u); + values<0, true, false, true>( + get_global_shape_values(mf_object_id), u, u); break; } case 2: { - values<0, true, false, true>(get_global_shape_values(), - u, - u); + values<0, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<1, true, false, true>(get_global_shape_values(), - u, - u); + values<1, true, false, true>( + get_global_shape_values(mf_object_id), u, u); break; } case 3: { - values<0, true, false, true>(get_global_shape_values(), - u, - u); + values<0, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<1, true, false, true>(get_global_shape_values(), - u, - u); + values<1, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<2, true, false, true>(get_global_shape_values(), - u, - u); + values<2, true, false, true>( + get_global_shape_values(mf_object_id), u, u); break; } @@ -312,37 +311,31 @@ namespace CUDAWrappers { case 1: { - values<0, false, false, true>(get_global_shape_values(), - u, - u); + values<0, false, false, true>( + get_global_shape_values(mf_object_id), u, u); break; } case 2: { - values<0, false, false, true>(get_global_shape_values(), - u, - u); + values<0, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<1, false, false, true>(get_global_shape_values(), - u, - u); + values<1, false, false, true>( + get_global_shape_values(mf_object_id), u, u); break; } case 3: { - values<0, false, false, true>(get_global_shape_values(), - u, - u); + values<0, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<1, false, false, true>(get_global_shape_values(), - u, - u); + values<1, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<2, false, false, true>(get_global_shape_values(), - u, - u); + values<2, false, false, true>( + get_global_shape_values(mf_object_id), u, u); break; } @@ -370,60 +363,68 @@ namespace CUDAWrappers case 1: { gradients<0, true, false, false>( - get_global_shape_gradients(), u, grad_u[0]); + get_global_shape_gradients(mf_object_id), u, grad_u[0]); break; } case 2: { gradients<0, true, false, false>( - get_global_shape_gradients(), u, grad_u[0]); - values<0, true, false, false>(get_global_shape_values(), - u, - grad_u[1]); + get_global_shape_gradients(mf_object_id), u, grad_u[0]); + values<0, true, false, false>( + get_global_shape_values(mf_object_id), u, grad_u[1]); __syncthreads(); - values<1, true, false, true>(get_global_shape_values(), + values<1, true, false, true>(get_global_shape_values( + mf_object_id), grad_u[0], grad_u[0]); gradients<1, true, false, true>( - get_global_shape_gradients(), grad_u[1], grad_u[1]); + get_global_shape_gradients(mf_object_id), + grad_u[1], + grad_u[1]); break; } case 3: { gradients<0, true, false, false>( - get_global_shape_gradients(), u, grad_u[0]); - values<0, true, false, false>(get_global_shape_values(), - u, - grad_u[1]); - values<0, true, false, false>(get_global_shape_values(), - u, - grad_u[2]); + get_global_shape_gradients(mf_object_id), u, grad_u[0]); + values<0, true, false, false>( + get_global_shape_values(mf_object_id), u, grad_u[1]); + values<0, true, false, false>( + get_global_shape_values(mf_object_id), u, grad_u[2]); __syncthreads(); - values<1, true, false, true>(get_global_shape_values(), + values<1, true, false, true>(get_global_shape_values( + mf_object_id), grad_u[0], grad_u[0]); gradients<1, true, false, true>( - get_global_shape_gradients(), grad_u[1], grad_u[1]); - values<1, true, false, true>(get_global_shape_values(), + get_global_shape_gradients(mf_object_id), + grad_u[1], + grad_u[1]); + values<1, true, false, true>(get_global_shape_values( + mf_object_id), grad_u[2], grad_u[2]); __syncthreads(); - values<2, true, false, true>(get_global_shape_values(), + values<2, true, false, true>(get_global_shape_values( + mf_object_id), grad_u[0], grad_u[0]); - values<2, true, false, true>(get_global_shape_values(), + values<2, true, false, true>(get_global_shape_values( + mf_object_id), grad_u[1], grad_u[1]); gradients<2, true, false, true>( - get_global_shape_gradients(), grad_u[2], grad_u[2]); + get_global_shape_gradients(mf_object_id), + grad_u[2], + grad_u[2]); break; } @@ -451,55 +452,61 @@ namespace CUDAWrappers { case 1: { - values<0, true, false, true>(get_global_shape_values(), - u, - u); + values<0, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); gradients<0, true, false, false>( - get_global_co_shape_gradients(), u, grad_u[0]); + get_global_co_shape_gradients(mf_object_id), + u, + grad_u[0]); break; } case 2: { - values<0, true, false, true>(get_global_shape_values(), - u, - u); + values<0, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<1, true, false, true>(get_global_shape_values(), - u, - u); + values<1, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); gradients<0, true, false, false>( - get_global_co_shape_gradients(), u, grad_u[0]); + get_global_co_shape_gradients(mf_object_id), + u, + grad_u[0]); gradients<1, true, false, false>( - get_global_co_shape_gradients(), u, grad_u[1]); + get_global_co_shape_gradients(mf_object_id), + u, + grad_u[1]); break; } case 3: { - values<0, true, false, true>(get_global_shape_values(), - u, - u); + values<0, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<1, true, false, true>(get_global_shape_values(), - u, - u); + values<1, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<2, true, false, true>(get_global_shape_values(), - u, - u); + values<2, true, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); gradients<0, true, false, false>( - get_global_co_shape_gradients(), u, grad_u[0]); + get_global_co_shape_gradients(mf_object_id), + u, + grad_u[0]); gradients<1, true, false, false>( - get_global_co_shape_gradients(), u, grad_u[1]); + get_global_co_shape_gradients(mf_object_id), + u, + grad_u[1]); gradients<2, true, false, false>( - get_global_co_shape_gradients(), u, grad_u[2]); + get_global_co_shape_gradients(mf_object_id), + u, + grad_u[2]); break; } @@ -528,63 +535,73 @@ namespace CUDAWrappers case 1: { gradients<0, false, add, false>( - get_global_shape_gradients(), grad_u[dim], u); + get_global_shape_gradients(mf_object_id), + grad_u[dim], + u); break; } case 2: { gradients<0, false, false, true>( - get_global_shape_gradients(), grad_u[0], grad_u[0]); - values<0, false, false, true>(get_global_shape_values(), + get_global_shape_gradients(mf_object_id), + grad_u[0], + grad_u[0]); + values<0, false, false, true>(get_global_shape_values( + mf_object_id), grad_u[1], grad_u[1]); __syncthreads(); - values<1, false, add, false>(get_global_shape_values(), - grad_u[0], - u); + values<1, false, add, false>( + get_global_shape_values(mf_object_id), grad_u[0], u); __syncthreads(); gradients<1, false, true, false>( - get_global_shape_gradients(), grad_u[1], u); + get_global_shape_gradients(mf_object_id), grad_u[1], u); break; } case 3: { gradients<0, false, false, true>( - get_global_shape_gradients(), grad_u[0], grad_u[0]); - values<0, false, false, true>(get_global_shape_values(), + get_global_shape_gradients(mf_object_id), + grad_u[0], + grad_u[0]); + values<0, false, false, true>(get_global_shape_values( + mf_object_id), grad_u[1], grad_u[1]); - values<0, false, false, true>(get_global_shape_values(), + values<0, false, false, true>(get_global_shape_values( + mf_object_id), grad_u[2], grad_u[2]); __syncthreads(); - values<1, false, false, true>(get_global_shape_values(), + values<1, false, false, true>(get_global_shape_values( + mf_object_id), grad_u[0], grad_u[0]); gradients<1, false, false, true>( - get_global_shape_gradients(), grad_u[1], grad_u[1]); - values<1, false, false, true>(get_global_shape_values(), + get_global_shape_gradients(mf_object_id), + grad_u[1], + grad_u[1]); + values<1, false, false, true>(get_global_shape_values( + mf_object_id), grad_u[2], grad_u[2]); __syncthreads(); - values<2, false, add, false>(get_global_shape_values(), - grad_u[0], - u); + values<2, false, add, false>( + get_global_shape_values(mf_object_id), grad_u[0], u); __syncthreads(); - values<2, false, true, false>(get_global_shape_values(), - grad_u[1], - u); + values<2, false, true, false>( + get_global_shape_values(mf_object_id), grad_u[1], u); __syncthreads(); gradients<2, false, true, false>( - get_global_shape_gradients(), grad_u[2], u); + get_global_shape_gradients(mf_object_id), grad_u[2], u); break; } @@ -613,31 +630,34 @@ namespace CUDAWrappers case 1: { gradients<0, false, true, false>( - get_global_co_shape_gradients(), grad_u[0], u); + get_global_co_shape_gradients(mf_object_id), + grad_u[0], + u); __syncthreads(); - values<0, false, false, true>(get_global_shape_values(), - u, - u); + values<0, false, false, true>( + get_global_shape_values(mf_object_id), u, u); break; } case 2: { gradients<1, false, true, false>( - get_global_co_shape_gradients(), grad_u[1], u); + get_global_co_shape_gradients(mf_object_id), + grad_u[1], + u); __syncthreads(); gradients<0, false, true, false>( - get_global_co_shape_gradients(), grad_u[0], u); + get_global_co_shape_gradients(mf_object_id), + grad_u[0], + u); __syncthreads(); - values<1, false, false, true>(get_global_shape_values(), - u, - u); + values<1, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<0, false, false, true>(get_global_shape_values(), - u, - u); + values<0, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); break; @@ -645,26 +665,29 @@ namespace CUDAWrappers case 3: { gradients<2, false, true, false>( - get_global_co_shape_gradients(), grad_u[2], u); + get_global_co_shape_gradients(mf_object_id), + grad_u[2], + u); __syncthreads(); gradients<1, false, true, false>( - get_global_co_shape_gradients(), grad_u[1], u); + get_global_co_shape_gradients(mf_object_id), + grad_u[1], + u); __syncthreads(); gradients<0, false, true, false>( - get_global_co_shape_gradients(), grad_u[0], u); + get_global_co_shape_gradients(mf_object_id), + grad_u[0], + u); __syncthreads(); - values<2, false, false, true>(get_global_shape_values(), - u, - u); + values<2, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<1, false, false, true>(get_global_shape_values(), - u, - u); + values<1, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); - values<0, false, false, true>(get_global_shape_values(), - u, - u); + values<0, false, false, true>( + get_global_shape_values(mf_object_id), u, u); __syncthreads(); break; diff --git a/tests/cuda/cuda_evaluate_1d_shape.cu b/tests/cuda/cuda_evaluate_1d_shape.cu index e6e0a2fc8f..21f0f63268 100644 --- a/tests/cuda/cuda_evaluate_1d_shape.cu +++ b/tests/cuda/cuda_evaluate_1d_shape.cu @@ -42,14 +42,14 @@ evaluate_tensor_product(double *dst, double *src) M - 1, N, double> - evaluator; + evaluator(0); if (type == 0) evaluator.template values<0, dof_to_quad, add, false>( - CUDAWrappers::internal::get_global_shape_values(), src, dst); + CUDAWrappers::internal::get_global_shape_values(0), src, dst); if (type == 1) evaluator.template gradients<0, dof_to_quad, add, false>( - CUDAWrappers::internal::get_global_shape_values(), src, dst); + CUDAWrappers::internal::get_global_shape_values(0), src, dst); } template @@ -98,16 +98,17 @@ test() unsigned int size_shape_values = M * N * sizeof(double); - cudaError_t cuda_error = cudaMemcpyToSymbol( - CUDAWrappers::internal::get_global_shape_values(), - shape_host.begin(), - size_shape_values, - 0, - cudaMemcpyHostToDevice); + cudaError_t cuda_error = + cudaMemcpyToSymbol(CUDAWrappers::internal::get_global_shape_values( + 0), + shape_host.begin(), + size_shape_values, + 0, + cudaMemcpyHostToDevice); AssertCuda(cuda_error); cuda_error = cudaMemcpyToSymbol( - CUDAWrappers::internal::get_global_shape_gradients(), + CUDAWrappers::internal::get_global_shape_gradients(0), shape_host.begin(), size_shape_values, 0, diff --git a/tests/cuda/cuda_evaluate_2d_shape.cu b/tests/cuda/cuda_evaluate_2d_shape.cu index b424baceb8..48d9172dc8 100644 --- a/tests/cuda/cuda_evaluate_2d_shape.cu +++ b/tests/cuda/cuda_evaluate_2d_shape.cu @@ -42,23 +42,23 @@ evaluate_tensor_product(double *dst, double *src) M - 1, N, double> - evaluator; + evaluator(0); if (type == 0) { evaluator.template values<0, dof_to_quad, false, false>( - CUDAWrappers::internal::get_global_shape_values(), src, src); + CUDAWrappers::internal::get_global_shape_values(0), src, src); __syncthreads(); evaluator.template values<1, dof_to_quad, add, false>( - CUDAWrappers::internal::get_global_shape_values(), src, dst); + CUDAWrappers::internal::get_global_shape_values(0), src, dst); } if (type == 1) { evaluator.template gradients<0, dof_to_quad, false, false>( - CUDAWrappers::internal::get_global_shape_values(), src, src); + CUDAWrappers::internal::get_global_shape_values(0), src, src); __syncthreads(); evaluator.template gradients<1, dof_to_quad, add, false>( - CUDAWrappers::internal::get_global_shape_values(), src, dst); + CUDAWrappers::internal::get_global_shape_values(0), src, dst); } } @@ -123,16 +123,17 @@ test() unsigned int size_shape_values = M * N * sizeof(double); - cudaError_t cuda_error = cudaMemcpyToSymbol( - CUDAWrappers::internal::get_global_shape_values(), - shape_host.begin(), - size_shape_values, - 0, - cudaMemcpyHostToDevice); + cudaError_t cuda_error = + cudaMemcpyToSymbol(CUDAWrappers::internal::get_global_shape_values( + 0), + shape_host.begin(), + size_shape_values, + 0, + cudaMemcpyHostToDevice); AssertCuda(cuda_error); cuda_error = cudaMemcpyToSymbol( - CUDAWrappers::internal::get_global_shape_gradients(), + CUDAWrappers::internal::get_global_shape_gradients(0), shape_host.begin(), size_shape_values, 0, diff --git a/tests/cuda/matrix_free_no_index_initialize.cu b/tests/cuda/matrix_free_no_index_initialize.cu index 0adf975f14..058a65014e 100644 --- a/tests/cuda/matrix_free_no_index_initialize.cu +++ b/tests/cuda/matrix_free_no_index_initialize.cu @@ -47,7 +47,7 @@ public: Number * dst) const { CUDAWrappers::FEEvaluation - fe_eval(cell, gpu_data, shared_data); + fe_eval(0, cell, gpu_data, shared_data); // set to unit vector fe_eval.submit_dof_value(1.); -- 2.39.5