From 3822e8f8c6174710244f783e08026e80de6894b1 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Sun, 28 Jul 2019 23:34:30 -0400 Subject: [PATCH] Add AffineConstraints::set_zero for MemorySpace::CUDA --- .../lac/affine_constraints.templates.h | 59 +++++++++++- source/lac/CMakeLists.txt | 1 + source/lac/affine_constraints.cu | 30 +++++++ tests/cuda/affine_constraints_set_zero.cu | 90 +++++++++++++++++++ ...ffine_constraints_set_zero.mpirun=2.output | 31 +++++++ 5 files changed, 208 insertions(+), 3 deletions(-) create mode 100644 source/lac/affine_constraints.cu create mode 100644 tests/cuda/affine_constraints_set_zero.cu create mode 100644 tests/cuda/affine_constraints_set_zero.mpirun=2.output diff --git a/include/deal.II/lac/affine_constraints.templates.h b/include/deal.II/lac/affine_constraints.templates.h index f9b30da3ac..5aad3b34a6 100644 --- a/include/deal.II/lac/affine_constraints.templates.h +++ b/include/deal.II/lac/affine_constraints.templates.h @@ -16,6 +16,7 @@ #ifndef dealii_affine_constraints_templates_h #define dealii_affine_constraints_templates_h +#include #include #include #include @@ -1883,6 +1884,58 @@ namespace internal vec.zero_out_ghosts(); } +#ifdef DEAL_II_COMPILER_CUDA_AWARE + template + __global__ void + set_zero_kernel(const size_type * constrained_dofs, + const unsigned int n_constrained_dofs, + Number * dst) + { + const unsigned int index = threadIdx.x + blockDim.x * blockIdx.x; + if (index < n_constrained_dofs) + dst[constrained_dofs[index]] = 0; + } + + template + void + set_zero_parallel( + const std::vector & cm, + LinearAlgebra::distributed::Vector &vec, + size_type shift = 0) + { + Assert(shift == 0, ExcNotImplemented()); + const int n_constraints = cm.size(); + std::vector constrained_local_dofs_host(n_constraints); + + std::transform(cm.begin(), + cm.end(), + constrained_local_dofs_host.begin(), + [&vec](const size_type global_index) { + return vec.get_partitioner()->global_to_local( + global_index); + }); + + size_type *constrained_local_dofs_device; + Utilities::CUDA::malloc(constrained_local_dofs_device, n_constraints); + Utilities::CUDA::copy_to_dev(constrained_local_dofs_host, + constrained_local_dofs_device); + + const int n_blocks = 1 + (n_constraints - 1) / CUDAWrappers::block_size; + set_zero_kernel<<>>( + constrained_local_dofs_device, cm.size(), vec.get_values()); +# ifdef DEBUG + // Check that the kernel was launched correctly + AssertCuda(cudaGetLastError()); + // Check that there was no problem during the execution of the kernel + AssertCuda(cudaDeviceSynchronize()); +# endif + + Utilities::CUDA::free(constrained_local_dofs_device); + + vec.zero_out_ghosts(); + } +#endif + template void set_zero_in_parallel(const std::vector &cm, @@ -1947,7 +2000,7 @@ template void AffineConstraints::set_zero(VectorType &vec) const { - // since we lines is a private member, we cannot pass it to the functions + // since lines is a private member, we cannot pass it to the functions // above. therefore, copy the content which is cheap std::vector constrained_lines(lines.size()); for (unsigned int i = 0; i < lines.size(); ++i) @@ -2087,7 +2140,7 @@ namespace internal // this is an operation that is different for all vector types and so we // need a few overloads #ifdef DEAL_II_WITH_TRILINOS - void + inline void import_vector_with_ghost_elements( const TrilinosWrappers::MPI::Vector &vec, const IndexSet & /*locally_owned_elements*/, @@ -2110,7 +2163,7 @@ namespace internal #endif #ifdef DEAL_II_WITH_PETSC - void + inline void import_vector_with_ghost_elements( const PETScWrappers::MPI::Vector &vec, const IndexSet & locally_owned_elements, diff --git a/source/lac/CMakeLists.txt b/source/lac/CMakeLists.txt index 79aef23709..4371f69da3 100644 --- a/source/lac/CMakeLists.txt +++ b/source/lac/CMakeLists.txt @@ -65,6 +65,7 @@ IF(DEAL_II_WITH_CUDA) SET(_separate_src ${_separate_src} vector_memory.cu + affine_constraints.cu ) ENDIF() diff --git a/source/lac/affine_constraints.cu b/source/lac/affine_constraints.cu new file mode 100644 index 0000000000..c686dc4fb2 --- /dev/null +++ b/source/lac/affine_constraints.cu @@ -0,0 +1,30 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +#include + +DEAL_II_NAMESPACE_OPEN + +template void +AffineConstraints::set_zero< + LinearAlgebra::distributed::Vector>( + LinearAlgebra::distributed::Vector &) const; + +template void +AffineConstraints::set_zero< + LinearAlgebra::distributed::Vector>( + LinearAlgebra::distributed::Vector &) const; + +DEAL_II_NAMESPACE_CLOSE diff --git a/tests/cuda/affine_constraints_set_zero.cu b/tests/cuda/affine_constraints_set_zero.cu new file mode 100644 index 0000000000..53c13c39d9 --- /dev/null +++ b/tests/cuda/affine_constraints_set_zero.cu @@ -0,0 +1,90 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// check AffineConstraints::set_zero(Vector) for parallel distributed +// vectors + +#include + +#include + +#include + +#include + +#include "../tests.h" + + +__global__ void +initialize_vector(double *vector, int local_size, int offset) +{ + const int index = threadIdx.x + blockIdx.x * blockDim.x; + if (index < local_size) + vector[index] = 1.0 + index + offset; +} + + +void +test() +{ + unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + unsigned int numproc = Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD); + + IndexSet local_active; + local_active.set_size(2 * numproc); + local_active.add_range(myid * numproc, (myid + 1) * numproc); + + LinearAlgebra::distributed::Vector v; + v.reinit(local_active, complete_index_set(2 * numproc), MPI_COMM_WORLD); + + const int n_blocks = 1 + v.size() / CUDAWrappers::block_size; + initialize_vector<<>>(v.get_values(), + numproc, + myid * numproc); + v.compress(VectorOperation::insert); + + AffineConstraints cm; + cm.add_line(numproc * myid + 1); + cm.close(); + + deallog << "vector before:" << std::endl; + v.print(deallog.get_file_stream()); + + deallog << std::endl; + deallog << "CM:" << std::endl; + cm.print(deallog.get_file_stream()); + + cm.set_zero(v); + + deallog << "vector after:" << std::endl; + v.print(deallog.get_file_stream()); + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + MPILogInitAll log; + + init_cuda(); + + test(); + return 0; +} diff --git a/tests/cuda/affine_constraints_set_zero.mpirun=2.output b/tests/cuda/affine_constraints_set_zero.mpirun=2.output new file mode 100644 index 0000000000..15ecc47505 --- /dev/null +++ b/tests/cuda/affine_constraints_set_zero.mpirun=2.output @@ -0,0 +1,31 @@ + +DEAL:0::vector before: +Process #0 +Local range: [0, 2), global size: 4 +Vector data: +1.000e+00 2.000e+00 +DEAL:0:: +DEAL:0::CM: + 1 = 0 +DEAL:0::vector after: +Process #0 +Local range: [0, 2), global size: 4 +Vector data: +1.000e+00 0.000e+00 +DEAL:0::OK + +DEAL:1::vector before: +Process #1 +Local range: [2, 4), global size: 4 +Vector data: +3.000e+00 4.000e+00 +DEAL:1:: +DEAL:1::CM: + 3 = 0 +DEAL:1::vector after: +Process #1 +Local range: [2, 4), global size: 4 +Vector data: +3.000e+00 0.000e+00 +DEAL:1::OK + -- 2.39.5