Reference documentation for deal.II version 9.5.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// Copyright (C) 2017 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_vector_element_access_h
17#define dealii_vector_element_access_h
18
19
20#include <deal.II/base/config.h>
21
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>
132 struct ElementAccess<LinearAlgebra::TpetraWrappers::Vector<NumberType>>
133 {
134 public:
136 static void
137 add(const typename VectorType::value_type value,
139 VectorType & V);
140
141 static void
144 VectorType & V);
145
146 static typename VectorType::value_type
147 get(const VectorType &V, const types::global_dof_index i);
148 };
149
150
151
152 template <typename NumberType>
153 inline void
155 const typename VectorType::value_type value,
158 {
159 // Extract local indices in the vector.
160 Tpetra::Vector<NumberType, int, types::signed_global_dof_index> vector =
161 V.trilinos_vector();
163 vector.getMap()->getLocalElement(
164 static_cast<TrilinosWrappers::types::int_type>(i));
165
166# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
167 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>(
168 Tpetra::Access::ReadWrite);
169# else
170 vector.template sync<Kokkos::HostSpace>();
171 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
172# endif
173 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
174# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
175 // We're going to modify the data on host.
176 vector.template modify<Kokkos::HostSpace>();
177# endif
178 vector_1d(trilinos_i) += value;
179# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
180 vector.template sync<
181 typename Tpetra::Vector<NumberType, int, types::signed_global_dof_index>::
182 device_type::memory_space>();
183# endif
184 }
185
186
187
188 template <typename NumberType>
189 inline void
191 const typename VectorType::value_type value,
194 {
195 // Extract local indices in the vector.
196 Tpetra::Vector<NumberType, int, types::signed_global_dof_index> vector =
197 V.trilinos_vector();
199 vector.getMap()->getLocalElement(
200 static_cast<TrilinosWrappers::types::int_type>(i));
201
202# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
203 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>(
204 Tpetra::Access::ReadWrite);
205# else
206 vector.template sync<Kokkos::HostSpace>();
207 auto vector_2d = vector.template getLocalView<Kokkos::HostSpace>();
208# endif
209 auto vector_1d = Kokkos::subview(vector_2d, Kokkos::ALL(), 0);
210 // We're going to modify the data on host.
211# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
212 vector.template modify<Kokkos::HostSpace>();
213# endif
214 vector_1d(trilinos_i) = value;
215# if !DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
216 vector.template sync<
217 typename Tpetra::Vector<NumberType, int, types::signed_global_dof_index>::
218 device_type::memory_space>();
219# endif
220 }
221
222
223
224 template <typename NumberType>
229 {
230 // Extract local indices in the vector.
231 Tpetra::Vector<NumberType, int, types::signed_global_dof_index> vector =
232 V.trilinos_vector();
234 vector.getMap()->getLocalElement(
235 static_cast<TrilinosWrappers::types::int_type>(i));
236
237# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
238 auto vector_2d =
239 vector.template getLocalView<Kokkos::HostSpace>(Tpetra::Access::ReadOnly);
240# else
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);
245 return vector_1d(trilinos_i);
246 }
247# endif
248#endif
249} // namespace internal
250
252
253#endif
const Epetra_FEVector & trilinos_vector() const
const Tpetra::Vector< Number, int, types::signed_global_dof_index > & trilinos_vector() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
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)