From: Daniel Arndt Date: Wed, 17 May 2023 20:55:21 +0000 (-0400) Subject: Fix Tpetra for newer Trilinos X-Git-Tag: v9.5.0-rc1~212^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=208e853cf7cb39babfe25c25c93b01dbab16fb32;p=dealii.git Fix Tpetra for newer Trilinos --- diff --git a/include/deal.II/lac/trilinos_tpetra_vector.templates.h b/include/deal.II/lac/trilinos_tpetra_vector.templates.h index c937f9292a..ca0c37c0d1 100644 --- a/include/deal.II/lac/trilinos_tpetra_vector.templates.h +++ b/include/deal.II/lac/trilinos_tpetra_vector.templates.h @@ -208,17 +208,26 @@ namespace LinearAlgebra tpetra_export.getSourceMap()); { +# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) + auto x_2d = source_vector.template getLocalView( + Tpetra::Access::ReadWrite); +# else source_vector.template sync(); auto x_2d = source_vector.template getLocalView(); +# endif auto x_1d = Kokkos::subview(x_2d, Kokkos::ALL(), 0); +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) source_vector.template modify(); +# endif const size_t localLength = source_vector.getLocalLength(); auto values_it = V.begin(); for (size_t k = 0; k < localLength; ++k) x_1d(k) = *values_it++; +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) source_vector.template sync< typename Tpetra::Vector:: device_type::memory_space>(); +# endif } if (operation == VectorOperation::insert) vector->doExport(source_vector, tpetra_export, Tpetra::REPLACE); @@ -331,18 +340,27 @@ namespace LinearAlgebra { AssertIsFinite(a); +# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) + auto vector_2d = vector->template getLocalView( + Tpetra::Access::ReadWrite); +# else vector->template sync(); auto vector_2d = vector->template getLocalView(); +# endif auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0); +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) vector->template modify(); +# endif const size_t localLength = vector->getLocalLength(); for (size_t k = 0; k < localLength; ++k) { vector_1d(k) += a; } +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) vector->template sync< typename Tpetra::Vector:: device_type::memory_space>(); +# endif } @@ -643,9 +661,14 @@ namespace LinearAlgebra else out.setf(std::ios::fixed, std::ios::floatfield); +# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) + auto vector_2d = vector->template getLocalView( + Tpetra::Access::ReadOnly); +# else vector->template sync(); auto vector_2d = vector->template getLocalView(); - auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0); +# endif + auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0); const size_t local_length = vector->getLocalLength(); if (across) diff --git a/include/deal.II/lac/vector_element_access.h b/include/deal.II/lac/vector_element_access.h index 62332e8a04..901d16b715 100644 --- a/include/deal.II/lac/vector_element_access.h +++ b/include/deal.II/lac/vector_element_access.h @@ -163,15 +163,24 @@ namespace internal vector.getMap()->getLocalElement( static_cast(i)); +# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) + auto vector_2d = vector.template getLocalView( + Tpetra::Access::ReadWrite); +# else vector.template sync(); auto vector_2d = vector.template getLocalView(); +# endif auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0); +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) // We're going to modify the data on host. vector.template modify(); +# endif vector_1d(trilinos_i) += value; +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) vector.template sync< typename Tpetra::Vector:: device_type::memory_space>(); +# endif } @@ -190,15 +199,24 @@ namespace internal vector.getMap()->getLocalElement( static_cast(i)); +# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) + auto vector_2d = vector.template getLocalView( + Tpetra::Access::ReadWrite); +# else vector.template sync(); auto vector_2d = vector.template getLocalView(); +# endif auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0); // We're going to modify the data on host. +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) vector.template modify(); +# endif vector_1d(trilinos_i) = value; +# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) vector.template sync< typename Tpetra::Vector:: device_type::memory_space>(); +# endif } @@ -216,10 +234,14 @@ namespace internal vector.getMap()->getLocalElement( static_cast(i)); +# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) + auto vector_2d = + vector.template getLocalView(Tpetra::Access::ReadOnly); +# else vector.template sync(); auto vector_2d = vector.template getLocalView(); +# endif auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0); - // We're going to modify the data on host. return vector_1d(trilinos_i); } # endif