From 9f8c3bbeb05057bdda96b22a3f208fe0faa9b94a Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Fri, 19 Jul 2019 18:04:16 -0400 Subject: [PATCH] Implement CUDAWrappers::MatrixFree::initialize_dof_vector --- .../deal.II/matrix_free/cuda_matrix_free.h | 22 ++++ .../matrix_free/cuda_matrix_free.templates.h | 28 +++- tests/cuda/matrix_free_initialize_vector.cu | 122 ++++++++++++++++++ ....with_mpi=on.with_p4est=on.mpirun=1.output | 11 ++ ....with_mpi=on.with_p4est=on.mpirun=2.output | 13 ++ 5 files changed, 195 insertions(+), 1 deletion(-) create mode 100644 tests/cuda/matrix_free_initialize_vector.cu create mode 100644 tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=1.output create mode 100644 tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output diff --git a/include/deal.II/matrix_free/cuda_matrix_free.h b/include/deal.II/matrix_free/cuda_matrix_free.h index 304278e6e4..2015252c16 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.h @@ -240,6 +240,23 @@ namespace CUDAWrappers void set_constrained_values(const Number val, VectorType &dst) const; + /** + * Initialize a serial vector. The size corresponds to the number of degrees + * of freedom in the DoFHandler object. + */ + void + initialize_dof_vector( + LinearAlgebra::CUDAWrappers::Vector &vec) const; + + /** + * Initialize a distributed vector. The local elements correspond to the + * locally owned degrees of freedom and the ghost elements correspond to the + * (additional) locally relevant dofs. + */ + void + initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const; + /** * Free all the memory allocated. */ @@ -360,6 +377,11 @@ namespace CUDAWrappers */ ParallelizationScheme parallelization_scheme; + /** + * Total number of degrees of freedom. + */ + types::global_dof_index n_dofs; + /** * Degree of the finite element used. */ 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 19e2f5a3aa..2be980af56 100644 --- a/include/deal.II/matrix_free/cuda_matrix_free.templates.h +++ b/include/deal.II/matrix_free/cuda_matrix_free.templates.h @@ -527,7 +527,8 @@ namespace CUDAWrappers template MatrixFree::MatrixFree() - : constrained_dofs(nullptr) + : n_dofs(0) + , constrained_dofs(nullptr) , padding_length(0) {} @@ -698,6 +699,29 @@ namespace CUDAWrappers + template + void + MatrixFree::initialize_dof_vector( + LinearAlgebra::CUDAWrappers::Vector &vec) const + { + vec.reinit(n_dofs); + } + + + + template + void + MatrixFree::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const + { + if (partitioner) + vec.reinit(partitioner); + else + vec.reinit(n_dofs); + } + + + template unsigned int MatrixFree::get_padding_length() const @@ -783,6 +807,8 @@ namespace CUDAWrappers // TODO: only free if we actually need arrays of different length free(); + n_dofs = dof_handler.n_dofs(); + const FiniteElement &fe = dof_handler.get_fe(); fe_degree = fe.degree; diff --git a/tests/cuda/matrix_free_initialize_vector.cu b/tests/cuda/matrix_free_initialize_vector.cu new file mode 100644 index 0000000000..f919c63b62 --- /dev/null +++ b/tests/cuda/matrix_free_initialize_vector.cu @@ -0,0 +1,122 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// Test CUDAWrappers::MatrixFree::initialize_dof_vector. + +#include + +#include +#include + +#include + +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + + +template +void +check( + const LinearAlgebra::distributed::Vector &vector, + const Utilities::MPI::Partitioner &reference_partitioner) +{ + Assert(vector.get_partitioner()->locally_owned_range() == + reference_partitioner.locally_owned_range(), + ExcInternalError()); + Assert(vector.get_partitioner()->ghost_indices() == + reference_partitioner.ghost_indices(), + ExcInternalError()); +} + +template +void +check(const LinearAlgebra::CUDAWrappers::Vector &vector, + const Utilities::MPI::Partitioner & reference_partitioner) +{ + AssertDimension(vector.size(), reference_partitioner.size()); +} + +template +void +test() +{ + typedef double Number; + + parallel::distributed::Triangulation tria(MPI_COMM_WORLD); + GridGenerator::hyper_cube(tria); + tria.refine_global(2); + + FE_Q fe(fe_degree); + DoFHandler dof(tria); + dof.distribute_dofs(fe); + + IndexSet owned_set = dof.locally_owned_dofs(); + IndexSet relevant_set; + DoFTools::extract_locally_relevant_dofs(dof, relevant_set); + + deallog << "locally owned dofs :" << std::endl; + owned_set.print(deallog.get_file_stream()); + + deallog << "locally relevant dofs :" << std::endl; + relevant_set.print(deallog.get_file_stream()); + + AffineConstraints constraints(relevant_set); + constraints.close(); + + MappingQGeneric mapping(fe_degree); + CUDAWrappers::MatrixFree mf_data; + const QGauss<1> quad(fe_degree + 1); + typename CUDAWrappers::MatrixFree::AdditionalData + additional_data; + mf_data.reinit(mapping, dof, constraints, quad, additional_data); + + VectorType vector; + mf_data.initialize_dof_vector(vector); + + Utilities::MPI::Partitioner reference_partitioner(owned_set, + relevant_set, + MPI_COMM_WORLD); + check(vector, reference_partitioner); + + deallog << "OK" << std::endl; +} + + +int +main(int argc, char **argv) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization( + argc, argv, testing_max_num_threads()); + + unsigned int myid = Utilities::MPI::this_mpi_process(MPI_COMM_WORLD); + deallog.push(Utilities::int_to_string(myid)); + + init_cuda(true); + MPILogInitAll mpi_inilog; + + test<2, 1, LinearAlgebra::distributed::Vector>(); + if (Utilities::MPI::n_mpi_processes(MPI_COMM_WORLD) == 1) + test<2, 1, LinearAlgebra::CUDAWrappers::Vector>(); +} diff --git a/tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=1.output b/tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=1.output new file mode 100644 index 0000000000..a0d4178229 --- /dev/null +++ b/tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=1.output @@ -0,0 +1,11 @@ + +DEAL:0:0::locally owned dofs : +{[0,24]} +DEAL:0:0::locally relevant dofs : +{[0,24]} +DEAL:0:0::OK +DEAL:0:0::locally owned dofs : +{[0,24]} +DEAL:0:0::locally relevant dofs : +{[0,24]} +DEAL:0:0::OK diff --git a/tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output b/tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output new file mode 100644 index 0000000000..50a6f81221 --- /dev/null +++ b/tests/cuda/matrix_free_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output @@ -0,0 +1,13 @@ + +DEAL:0:0::locally owned dofs : +{[0,14]} +DEAL:0:0::locally relevant dofs : +{[0,17], [21,22]} +DEAL:0:0::OK + +DEAL:1:1::locally owned dofs : +{[15,24]} +DEAL:1:1::locally relevant dofs : +{[2,3], [5,8], 10, [12,24]} +DEAL:1:1::OK + -- 2.39.5