Reference documentation for deal.II version 9.6.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
vector_element_access.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2017 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_vector_element_access_h
16#define dealii_vector_element_access_h
17
18
19#include <deal.II/base/config.h>
20
23
25
26namespace internal
27{
28 template <typename VectorType>
30 {
31 public:
32 static void
33 add(const typename VectorType::value_type value,
35 VectorType &V);
36
37 static void
38 set(typename VectorType::value_type value,
40 VectorType &V);
41
42 static typename VectorType::value_type
43 get(const VectorType &V, const types::global_dof_index i);
44 };
45
46
47
48 template <typename VectorType>
49 inline void
50 ElementAccess<VectorType>::add(const typename VectorType::value_type value,
52 VectorType &V)
53 {
54 V(i) += value;
55 }
56
57
58
59 template <typename VectorType>
60 inline void
61 ElementAccess<VectorType>::set(const typename VectorType::value_type value,
63 VectorType &V)
64 {
65 V(i) = value;
66 }
67
68
69
70 template <typename VectorType>
71 inline typename VectorType::value_type
72 ElementAccess<VectorType>::get(const VectorType &V,
74 {
75 return V(i);
76 }
77
78
79
80#ifdef DEAL_II_WITH_TRILINOS
81 template <>
82 inline void
84 const double value,
87 {
88 // Extract local indices in the vector.
89 Epetra_FEVector vector = V.trilinos_vector();
91 vector.Map().LID(static_cast<TrilinosWrappers::types::int_type>(i));
92
93 vector[0][trilinos_i] += value;
94 }
95
96
97
98 template <>
99 inline void
101 const double value,
104 {
105 // Extract local indices in the vector.
106 Epetra_FEVector vector = V.trilinos_vector();
108 vector.Map().LID(static_cast<TrilinosWrappers::types::int_type>(i));
109
110 vector[0][trilinos_i] = value;
111 }
112
113 template <>
114 inline double
118 {
119 // Extract local indices in the vector.
120 Epetra_FEVector vector = V.trilinos_vector();
122 vector.Map().LID(static_cast<TrilinosWrappers::types::int_type>(i));
123
124 return vector[0][trilinos_i];
125 }
126
127
128
129# ifdef DEAL_II_TRILINOS_WITH_TPETRA
130 template <typename NumberType, typename MemorySpace>
132 LinearAlgebra::TpetraWrappers::Vector<NumberType, MemorySpace>>
133 {
134 public:
137 static void
138 add(const typename VectorType::value_type value,
140 VectorType &V);
141
142 static void
145 VectorType &V);
146
147 static typename VectorType::value_type
148 get(const VectorType &V, const types::global_dof_index i);
149 };
150
151
152
153 template <typename NumberType, typename MemorySpace>
154 inline void
157 add(const typename VectorType::value_type value,
160 {
161 // Extract local indices in the vector.
162 auto vector = V.trilinos_vector();
164 vector.getMap()->getLocalElement(
165 static_cast<TrilinosWrappers::types::int_type>(i));
166
167# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
168 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>(
169 Tpetra::Access::ReadWrite);
170# else
171 vector.template sync<Kokkos::HostSpace>();
172 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
173# endif
174 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
175# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
176 // We're going to modify the data on host.
177 vector.template modify<Kokkos::HostSpace>();
178# endif
179 vector_1d(trilinos_i) += value;
180# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
181 vector.template sync<
182 typename Tpetra::Vector<NumberType, int, types::signed_global_dof_index>::
183 device_type::memory_space>();
184# endif
185 }
186
187
188
189 template <typename NumberType, typename MemorySpace>
190 inline void
193 set(const typename VectorType::value_type value,
196 {
197 // Extract local indices in the vector.
198 auto vector = V.trilinos_vector();
200 vector.getMap()->getLocalElement(
201 static_cast<TrilinosWrappers::types::int_type>(i));
202
203# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
204 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>(
205 Tpetra::Access::ReadWrite);
206# else
207 vector.template sync<Kokkos::HostSpace>();
208 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
209# endif
210 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
211 // We're going to modify the data on host.
212# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
213 vector.template modify<Kokkos::HostSpace>();
214# endif
215 vector_1d(trilinos_i) = value;
216# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
217 vector.template sync<
218 typename Tpetra::Vector<NumberType, int, types::signed_global_dof_index>::
219 device_type::memory_space>();
220# endif
221 }
222
223
224
225 template <typename NumberType, typename MemorySpace>
227 MemorySpace>::value_type
232 {
233 // Extract local indices in the vector.
234# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
235 const auto &vector = V.trilinos_vector();
236 auto vector_2d =
237 vector.template getLocalView<Kokkos::HostSpace>(Tpetra::Access::ReadOnly);
238# else
239 auto vector = V.trilinos_vector();
240 vector.template sync<Kokkos::HostSpace>();
241 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
242# endif
243 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
245 vector.getMap()->getLocalElement(
246 static_cast<TrilinosWrappers::types::int_type>(i));
247 return vector_1d(trilinos_i);
248 }
249# endif
250#endif
251} // namespace internal
252
254
255#endif
const Epetra_FEVector & trilinos_vector() const
const TpetraTypes::VectorType< Number, MemorySpace > & trilinos_vector() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:503
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:504
static VectorType::value_type get(const VectorType &V, const types::global_dof_index i)
static void set(typename VectorType::value_type value, const types::global_dof_index i, VectorType &V)
static void add(const typename VectorType::value_type value, const types::global_dof_index i, VectorType &V)