[Trilinos-Users] Best practice to modify elements of a Tpetra::Vector with openmp?

[Wen Yan]

Thanks!

Could you please also confirm the following?

1. x_2d is a 2D Kokkos::View with Kokkos::LayoutLeft, and the entries x_2d(i,c) are continuous in memory for continuous index i. Here i ranges from 0 to local length - 1, and c ranges from 0 to the number of vectors in the Tpetra::MultiVector - 1

2. x_2d is safe to use with openmp, like:
  for (int c = 0; c < x_2d.dimension_1(); c++) {
#pragma omp parallel for
    for (int i = 0; i < x_2d.dimension_0(); i++) {
      const double temp = x_2d(i, c);
      x_2d(i, c) = temp < 0 ? 0 : temp;
    }
  }

3. For read-only access there is no need to call
modify<Kokkos::HostSpace>()

Thanks,
Wen Yan




> On Jun 8, 2017, at 4:48 PM, Hoemmen, Mark <mhoemme at sandia.gov> wrote:
> 
> 
> Please refer to the following, commented example:
> 
> Trilinos/packages/tpetra/core/example/Lesson02-Map-Vector/lesson02_read_modify_vec_kokkos.cpp
> 
> mfh
> 
> 
> 
> On 6/6/17, 3:11 PM, "Trilinos-Users on behalf of trilinos-users-request at trilinos.org" <trilinos-users-bounces at trilinos.org on behalf of trilinos-users-request at trilinos.org> wrote:
> 
>> Date: Fri, 12 May 2017 17:36:42 -0400
>> From: Wen Yan <wenyan4work at gmail.com>
>> To: "trilinos-users at trilinos.org" <trilinos-users at trilinos.org>
>> Subject: [Trilinos-Users] Best practice to modify elements of a
>> 	Tpetra::Vector with openmp?
>> Message-ID: <D1AADFE6-3AB8-4259-8C00-040F2BE02EDF at gmail.com>
>> Content-Type: text/plain; charset=us-ascii
>> 
>> Hi Trilinos,
>> 
>> I was wondering what is the best practice to modify the elements of a Tpetra::Vector or a Tpetra::MultiVector. The document mentioned that ArrayRCP<T> is not thread safe, and suggested using Kokkos::View, but not much details are discussed. I am using Kokkos::OpenMP as the execution backend of Tpetra.
>> 
>> Assume xRcp is Teuchos::Rcp<Tpetra::MultiVector>
>> 
>> To use Kokkos::View, I can do this, but it is awkward and the Kokkos documentation does not include much useful information about these:
>> 
>> //   auto x_2d = xRcp->getLocalView<Kokkos::HostSpace>();
>> //   auto x_1d = Kokkos::subview(x_2d, Kokkos::ALL(), 0);
>> //   xRcp->modify<Kokkos::HostSpace>();
>> //   const size_t localLength = vecRcp->getLocalLength();
>> //   #pragma omp parallel for
>> //   for (size_t k = 0; k < localLength; ++k) {
>> //        x_1d[k]=k;
>> //   }
>> 
>> 
>> The ArrayRCP interface is more straightforward, but it is mentioned in the document that its reference counting is not compatible with multi threading. I was wondering if it is unsafe to write a simple for loop like this to modify the content of a Tpetra::Vector
>> 
>> // auto xArrayRcp = xRcp->getDataNonConst();
>> // const int localLength = xRcp->getLocalLength();
>> // #pragma omg parallel for
>> // for (int i = 0; i < localLength; i++) {
>> //     xArrayRcp[i]=I; 
>> // }
>> 
>> 
>> Thanks a lot!
>> Wen Yan
>> 
>>