deal.II version GIT relicensing-3696-g88d7b26363 2025-07-08 15:10:00+00:00
\(\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 - 2025 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
24
26
27namespace internal
28{
29 template <typename VectorType>
31 {
32 public:
33 static void
34 add(const typename VectorType::value_type value,
36 VectorType &V);
37
38 static void
39 set(typename VectorType::value_type value,
41 VectorType &V);
42
43 static typename VectorType::value_type
44 get(const VectorType &V, const types::global_dof_index i);
45 };
46
47
48
49 template <typename VectorType>
50 inline void
51 ElementAccess<VectorType>::add(const typename VectorType::value_type value,
53 VectorType &V)
54 {
55 V(i) += value;
56 }
57
58
59
60 template <typename VectorType>
61 inline void
62 ElementAccess<VectorType>::set(const typename VectorType::value_type value,
64 VectorType &V)
65 {
66 V(i) = value;
67 }
68
69
70
71 template <typename VectorType>
72 inline typename VectorType::value_type
73 ElementAccess<VectorType>::get(const VectorType &V,
75 {
76 return V(i);
77 }
78
79
80
81#ifdef DEAL_II_WITH_TRILINOS
82 template <>
83 inline void
85 const double value,
88 {
89 // Extract local indices in the vector.
90 Epetra_FEVector vector = V.trilinos_vector();
92 vector.Map().LID(static_cast<TrilinosWrappers::types::int_type>(i));
93
94 vector[0][trilinos_i] += value;
95 }
96
97
98
99 template <>
100 inline void
102 const double value,
105 {
106 // Extract local indices in the vector.
107 Epetra_FEVector vector = V.trilinos_vector();
109 vector.Map().LID(static_cast<TrilinosWrappers::types::int_type>(i));
110
111 vector[0][trilinos_i] = value;
112 }
113
114 template <>
115 inline double
119 {
120 // Extract local indices in the vector.
121 Epetra_FEVector vector = V.trilinos_vector();
123 vector.Map().LID(static_cast<TrilinosWrappers::types::int_type>(i));
124
125 return vector[0][trilinos_i];
126 }
127
128
129
130# ifdef DEAL_II_TRILINOS_WITH_TPETRA
131 template <typename NumberType, typename MemorySpace>
133 LinearAlgebra::TpetraWrappers::Vector<NumberType, MemorySpace>>
134 {
135 public:
138 static void
139 add(const typename VectorType::value_type value,
141 VectorType &V);
142
143 static void
146 VectorType &V);
147
148 static typename VectorType::value_type
149 get(const VectorType &V, const types::global_dof_index i);
150 };
151
152
153
154 template <typename NumberType, typename MemorySpace>
155 inline void
158 add(const typename VectorType::value_type value,
161 {
162 // Extract local indices in the vector.
163 auto vector = V.trilinos_vector();
165 vector.getMap()->getLocalElement(
166 static_cast<TrilinosWrappers::types::int_type>(i));
167
168# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
169 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>(
170 Tpetra::Access::ReadWrite);
171# else
172 vector.template sync<Kokkos::HostSpace>();
173 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
174# endif
175 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
176# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
177 // We're going to modify the data on host.
178 vector.template modify<Kokkos::HostSpace>();
179# endif
180 vector_1d(trilinos_i) += value;
181# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
182 vector.template sync<
183 typename Tpetra::Vector<NumberType, int, types::signed_global_dof_index>::
184 device_type::memory_space>();
185# endif
186 }
187
188
189
190 template <typename NumberType, typename MemorySpace>
191 inline void
194 set(const typename VectorType::value_type value,
197 {
198 // Extract local indices in the vector.
199 auto vector = V.trilinos_vector();
201 vector.getMap()->getLocalElement(
202 static_cast<TrilinosWrappers::types::int_type>(i));
203
204# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
205 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>(
206 Tpetra::Access::ReadWrite);
207# else
208 vector.template sync<Kokkos::HostSpace>();
209 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
210# endif
211 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
212 // We're going to modify the data on host.
213# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
214 vector.template modify<Kokkos::HostSpace>();
215# endif
216 vector_1d(trilinos_i) = value;
217# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
218 vector.template sync<
219 typename Tpetra::Vector<NumberType, int, types::signed_global_dof_index>::
220 device_type::memory_space>();
221# endif
222 }
223
224
225
226 template <typename NumberType, typename MemorySpace>
228 MemorySpace>::value_type
233 {
234 // Extract local indices in the vector.
235# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
236 const auto &vector = V.trilinos_vector();
237 auto vector_2d =
238 vector.template getLocalView<Kokkos::HostSpace>(Tpetra::Access::ReadOnly);
239# else
240 auto vector = V.trilinos_vector();
241 vector.template sync<Kokkos::HostSpace>();
242 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
243# endif
244 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
246 vector.getMap()->getLocalElement(
247 static_cast<TrilinosWrappers::types::int_type>(i));
248 return vector_1d(trilinos_i);
249 }
250# endif
251#endif
252} // namespace internal
253
255
256#endif
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
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)