From d6cc28fe60c535c55304ee3c5370f02aba35b088 Mon Sep 17 00:00:00 2001 From: Daniel Arndt Date: Tue, 20 Aug 2024 15:19:39 -0400 Subject: [PATCH] Portable MatrixFree: Allow initializing host vector --- .../matrix_free/portable_matrix_free.h | 4 +-- .../portable_matrix_free.templates.h | 3 +- source/matrix_free/portable_matrix_free.cc | 32 ++++++++++++++++++- .../matrix_free_device_initialize_vector.cc | 7 ++-- ...with_mpi=on.with_p4est=on.mpirun=1.output} | 0 ....with_mpi=on.with_p4est=on.mpirun=2.output | 10 ++++++ ...ith_p4est=on.with_cuda=off.mpirun=1.output | 6 ---- 7 files changed, 49 insertions(+), 13 deletions(-) rename tests/matrix_free_kokkos/{matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.with_cuda=on.mpirun=1.output => matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.mpirun=1.output} (100%) delete mode 100644 tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.with_cuda=off.mpirun=1.output diff --git a/include/deal.II/matrix_free/portable_matrix_free.h b/include/deal.II/matrix_free/portable_matrix_free.h index e336c662b5..bd631563ed 100644 --- a/include/deal.II/matrix_free/portable_matrix_free.h +++ b/include/deal.II/matrix_free/portable_matrix_free.h @@ -380,10 +380,10 @@ namespace Portable * locally owned degrees of freedom and the ghost elements correspond to the * (additional) locally relevant dofs. */ + template void initialize_dof_vector( - LinearAlgebra::distributed::Vector &vec) - const; + LinearAlgebra::distributed::Vector &vec) const; /** * Return the colored graph of locally owned active cells. diff --git a/include/deal.II/matrix_free/portable_matrix_free.templates.h b/include/deal.II/matrix_free/portable_matrix_free.templates.h index f4475662a1..b782081a8c 100644 --- a/include/deal.II/matrix_free/portable_matrix_free.templates.h +++ b/include/deal.II/matrix_free/portable_matrix_free.templates.h @@ -566,9 +566,10 @@ namespace Portable template + template void MatrixFree::initialize_dof_vector( - LinearAlgebra::distributed::Vector &vec) const + LinearAlgebra::distributed::Vector &vec) const { if (partitioner) vec.reinit(partitioner); diff --git a/source/matrix_free/portable_matrix_free.cc b/source/matrix_free/portable_matrix_free.cc index 8422260e60..4290c25876 100644 --- a/source/matrix_free/portable_matrix_free.cc +++ b/source/matrix_free/portable_matrix_free.cc @@ -16,7 +16,7 @@ DEAL_II_NAMESPACE_OPEN - +#ifndef DOXYGEN namespace Portable { @@ -25,6 +25,36 @@ namespace Portable template class MatrixFree<2, double>; template class MatrixFree<3, float>; template class MatrixFree<3, double>; + + template void + MatrixFree<2, float>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const; + template void + MatrixFree<2, float>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const; + template void + MatrixFree<2, double>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const; + template void + MatrixFree<2, double>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) + const; + template void + MatrixFree<3, float>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const; + template void + MatrixFree<3, float>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const; + template void + MatrixFree<3, double>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) const; + template void + MatrixFree<3, double>::initialize_dof_vector( + LinearAlgebra::distributed::Vector &vec) + const; + } // namespace Portable +#endif + DEAL_II_NAMESPACE_CLOSE diff --git a/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.cc b/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.cc index f7e82faadf..d9389e4ce5 100644 --- a/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.cc +++ b/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.cc @@ -34,10 +34,9 @@ #include "../tests.h" -template +template void -check(const LinearAlgebra::distributed::Vector - &vector, +check(const LinearAlgebra::distributed::Vector &vector, const Utilities::MPI::Partitioner &reference_partitioner) { Assert(vector.get_partitioner()->locally_owned_range() == @@ -103,6 +102,8 @@ main(int argc, char **argv) MPILogInitAll mpi_inilog; + test<2, 1, LinearAlgebra::distributed::Vector>(); + test<2, 1, LinearAlgebra::distributed::Vector>(); diff --git a/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.with_cuda=on.mpirun=1.output b/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.mpirun=1.output similarity index 100% rename from tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.with_cuda=on.mpirun=1.output rename to tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.mpirun=1.output diff --git a/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output b/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output index 50a6f81221..5110052051 100644 --- a/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output +++ b/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.mpirun=2.output @@ -1,4 +1,9 @@ +DEAL:0:0::locally owned dofs : +{[0,14]} +DEAL:0:0::locally relevant dofs : +{[0,17], [21,22]} +DEAL:0:0::OK DEAL:0:0::locally owned dofs : {[0,14]} DEAL:0:0::locally relevant dofs : @@ -10,4 +15,9 @@ DEAL:1:1::locally owned dofs : DEAL:1:1::locally relevant dofs : {[2,3], [5,8], 10, [12,24]} DEAL:1:1::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 diff --git a/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.with_cuda=off.mpirun=1.output b/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.with_cuda=off.mpirun=1.output deleted file mode 100644 index 274361dac4..0000000000 --- a/tests/matrix_free_kokkos/matrix_free_device_initialize_vector.with_mpi=on.with_p4est=on.with_cuda=off.mpirun=1.output +++ /dev/null @@ -1,6 +0,0 @@ - -DEAL:0:0::locally owned dofs : -{[0,24]} -DEAL:0:0::locally relevant dofs : -{[0,24]} -DEAL:0:0::OK -- 2.39.5