Reference documentation for deal.II version Git 191d06ed00 2021-05-11 21:15:49 -0400
\(\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\}}\)
copy.cc
Go to the documentation of this file.
1 //-----------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2019 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 #include <deal.II/sundials/copy.h>
17 
18 #ifdef DEAL_II_WITH_SUNDIALS
19 
21 namespace SUNDIALS
22 {
23  namespace internal
24  {
32  inline std::size_t
33  N_Vector_length(const N_Vector &vec)
34  {
35  const N_Vector_ID id = N_VGetVectorID(vec);
36  long int length = -1;
37  switch (id)
38  {
39  case SUNDIALS_NVEC_SERIAL:
40  {
41  length = NV_LENGTH_S(vec);
42  break;
43  }
44 # ifdef DEAL_II_WITH_MPI
45  case SUNDIALS_NVEC_PARALLEL:
46  {
47  length = NV_LOCLENGTH_P(vec);
48  break;
49  }
50 # endif
51  default:
52  Assert(false, ExcNotImplemented());
53  }
54 
55  Assert(length >= 0, ExcInternalError());
56  return static_cast<std::size_t>(length);
57  }
58 
59 
60  void
62  {
63  const IndexSet is = dst.locally_owned_elements();
64  const std::size_t N = is.n_elements();
67  for (std::size_t i = 0; i < N; ++i)
68  {
69 # ifdef DEAL_II_WITH_MPI
70  dst.local_element(i) = NV_Ith_P(src, i);
71 # else
72  dst.local_element(i) = NV_Ith_S(src, i);
73 # endif
74  }
76  }
77 
78  void
80  {
81  const IndexSet is = src.locally_owned_elements();
82  const std::size_t N = is.n_elements();
84  for (std::size_t i = 0; i < N; ++i)
85  {
86 # ifdef DEAL_II_WITH_MPI
87  NV_Ith_P(dst, i) = src.local_element(i);
88 # else
89  NV_Ith_S(dst, i) = src.local_element(i);
90 # endif
91  }
92  }
93 
94  void
96  const N_Vector & src)
97  {
98  const IndexSet is = dst.locally_owned_elements();
99  const std::size_t N = is.n_elements();
101  dst.zero_out_ghost_values();
102  for (std::size_t i = 0; i < N; ++i)
103  {
104 # ifdef DEAL_II_WITH_MPI
105  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
106 # else
107  dst[is.nth_index_in_set(i)] = NV_Ith_S(src, i);
108 # endif
109  }
111  }
112 
113  void
114  copy(N_Vector & dst,
116  {
117  IndexSet is = src.locally_owned_elements();
118  const std::size_t N = is.n_elements();
120  for (std::size_t i = 0; i < N; ++i)
121  {
122 # ifdef DEAL_II_WITH_MPI
123  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
124 # else
125  NV_Ith_S(dst, i) = src[is.nth_index_in_set(i)];
126 # endif
127  }
128  }
129 
130 # ifdef DEAL_II_WITH_MPI
131 
132 # ifdef DEAL_II_WITH_TRILINOS
133 
134 
135  void
136  copy(TrilinosWrappers::MPI::Vector &dst, const N_Vector &src)
137  {
138  const IndexSet is = dst.locally_owned_elements();
139  const std::size_t N = is.n_elements();
141  for (std::size_t i = 0; i < N; ++i)
142  {
143  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
144  }
146  }
147 
148  void
149  copy(N_Vector &dst, const TrilinosWrappers::MPI::Vector &src)
150  {
151  const IndexSet is = src.locally_owned_elements();
152  const std::size_t N = is.n_elements();
154  for (std::size_t i = 0; i < N; ++i)
155  {
156  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
157  }
158  }
159 
160  void
161  copy(TrilinosWrappers::MPI::BlockVector &dst, const N_Vector &src)
162  {
163  const IndexSet is = dst.locally_owned_elements();
164  const std::size_t N = is.n_elements();
166  for (std::size_t i = 0; i < N; ++i)
167  {
168  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
169  }
171  }
172 
173  void
174  copy(N_Vector &dst, const TrilinosWrappers::MPI::BlockVector &src)
175  {
176  IndexSet is = src.locally_owned_elements();
177  const std::size_t N = is.n_elements();
179  for (std::size_t i = 0; i < N; ++i)
180  {
181  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
182  }
183  }
184 
185 # endif // DEAL_II_WITH_TRILINOS
186 
187 # ifdef DEAL_II_WITH_PETSC
188 # ifndef PETSC_USE_COMPLEX
189 
190  void
191  copy(PETScWrappers::MPI::Vector &dst, const N_Vector &src)
192  {
193  const IndexSet is = dst.locally_owned_elements();
194  const std::size_t N = is.n_elements();
196  for (std::size_t i = 0; i < N; ++i)
197  {
198  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
199  }
201  }
202 
203  void
204  copy(N_Vector &dst, const PETScWrappers::MPI::Vector &src)
205  {
206  const IndexSet is = src.locally_owned_elements();
207  const std::size_t N = is.n_elements();
209  for (std::size_t i = 0; i < N; ++i)
210  {
211  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
212  }
213  }
214 
215  void
216  copy(PETScWrappers::MPI::BlockVector &dst, const N_Vector &src)
217  {
218  const IndexSet is = dst.locally_owned_elements();
219  const std::size_t N = is.n_elements();
221  for (std::size_t i = 0; i < N; ++i)
222  {
223  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
224  }
226  }
227 
228  void
229  copy(N_Vector &dst, const PETScWrappers::MPI::BlockVector &src)
230  {
231  const IndexSet is = src.locally_owned_elements();
232  const std::size_t N = is.n_elements();
234  for (std::size_t i = 0; i < N; ++i)
235  {
236  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
237  }
238  }
239 
240 # endif // PETSC_USE_COMPLEX
241 # endif // DEAL_II_WITH_PETSC
242 
243 # endif // mpi
244 
245  void
246  copy(BlockVector<double> &dst, const N_Vector &src)
247  {
248  const std::size_t N = dst.size();
250  for (std::size_t i = 0; i < N; ++i)
251  {
252  dst[i] = NV_Ith_S(src, i);
253  }
254  }
255 
256  void
257  copy(N_Vector &dst, const BlockVector<double> &src)
258  {
259  const std::size_t N = src.size();
261  for (std::size_t i = 0; i < N; ++i)
262  {
263  NV_Ith_S(dst, i) = src[i];
264  }
265  }
266 
267  void
268  copy(Vector<double> &dst, const N_Vector &src)
269  {
270  const std::size_t N = dst.size();
272  for (std::size_t i = 0; i < N; ++i)
273  {
274  dst[i] = NV_Ith_S(src, i);
275  }
276  }
277 
278  void
279  copy(N_Vector &dst, const Vector<double> &src)
280  {
281  const std::size_t N = src.size();
283  for (std::size_t i = 0; i < N; ++i)
284  {
285  NV_Ith_S(dst, i) = src[i];
286  }
287  }
288  } // namespace internal
289 } // namespace SUNDIALS
291 
292 #endif
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1622
virtual void compress(::VectorOperation::values operation) override
void compress(::VectorOperation::values operation)
void copy(TrilinosWrappers::MPI::Vector &dst, const N_Vector &src)
Definition: copy.cc:136
Number local_element(const size_type local_index) const
IndexSet locally_owned_elements() const
void compress(const VectorOperation::values operation)
virtual void compress(::VectorOperation::values operation) override
#define Assert(cond, exc)
Definition: exceptions.h:1465
void compress(::VectorOperation::values operation)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
std::size_t N_Vector_length(const N_Vector &vec)
Definition: copy.cc:33
IndexSet locally_owned_elements() const
virtual ::IndexSet locally_owned_elements() const override
virtual ::IndexSet locally_owned_elements() const override
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:395
static const char N
IndexSet locally_owned_elements() const
static ::ExceptionBase & ExcNotImplemented()
size_type n_elements() const
Definition: index_set.h:1832
size_type nth_index_in_set(const size_type local_index) const
Definition: index_set.h:1880
static ::ExceptionBase & ExcInternalError()