Reference documentation for deal.II version 9.2.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\}}\)
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 # ifdef DEAL_II_WITH_MPI
60 
61 # ifdef DEAL_II_WITH_TRILINOS
62 
63 
64  void
65  copy(TrilinosWrappers::MPI::Vector &dst, const N_Vector &src)
66  {
67  const IndexSet is = dst.locally_owned_elements();
68  const std::size_t N = is.n_elements();
70  for (std::size_t i = 0; i < N; ++i)
71  {
72  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
73  }
75  }
76 
77  void
78  copy(N_Vector &dst, const TrilinosWrappers::MPI::Vector &src)
79  {
80  const IndexSet is = src.locally_owned_elements();
81  const std::size_t N = is.n_elements();
83  for (std::size_t i = 0; i < N; ++i)
84  {
85  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
86  }
87  }
88 
89  void
90  copy(TrilinosWrappers::MPI::BlockVector &dst, const N_Vector &src)
91  {
92  const IndexSet is = dst.locally_owned_elements();
93  const std::size_t N = is.n_elements();
95  for (std::size_t i = 0; i < N; ++i)
96  {
97  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
98  }
100  }
101 
102  void
103  copy(N_Vector &dst, const TrilinosWrappers::MPI::BlockVector &src)
104  {
105  IndexSet is = src.locally_owned_elements();
106  const std::size_t N = is.n_elements();
108  for (std::size_t i = 0; i < N; ++i)
109  {
110  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
111  }
112  }
113 
114 # endif // DEAL_II_WITH_TRILINOS
115 
116 # ifdef DEAL_II_WITH_PETSC
117 # ifndef PETSC_USE_COMPLEX
118 
119  void
120  copy(PETScWrappers::MPI::Vector &dst, const N_Vector &src)
121  {
122  const IndexSet is = dst.locally_owned_elements();
123  const std::size_t N = is.n_elements();
125  for (std::size_t i = 0; i < N; ++i)
126  {
127  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
128  }
130  }
131 
132  void
133  copy(N_Vector &dst, const PETScWrappers::MPI::Vector &src)
134  {
135  const IndexSet is = src.locally_owned_elements();
136  const std::size_t N = is.n_elements();
138  for (std::size_t i = 0; i < N; ++i)
139  {
140  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
141  }
142  }
143 
144  void
145  copy(PETScWrappers::MPI::BlockVector &dst, const N_Vector &src)
146  {
147  const IndexSet is = dst.locally_owned_elements();
148  const std::size_t N = is.n_elements();
150  for (std::size_t i = 0; i < N; ++i)
151  {
152  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
153  }
155  }
156 
157  void
158  copy(N_Vector &dst, const PETScWrappers::MPI::BlockVector &src)
159  {
160  const IndexSet is = src.locally_owned_elements();
161  const std::size_t N = is.n_elements();
163  for (std::size_t i = 0; i < N; ++i)
164  {
165  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
166  }
167  }
168 
169 # endif // PETSC_USE_COMPLEX
170 # endif // DEAL_II_WITH_PETSC
171 
172 # endif // mpi
173 
174  void
175  copy(BlockVector<double> &dst, const N_Vector &src)
176  {
177  const std::size_t N = dst.size();
179  for (std::size_t i = 0; i < N; ++i)
180  {
181  dst[i] = NV_Ith_S(src, i);
182  }
183  }
184 
185  void
186  copy(N_Vector &dst, const BlockVector<double> &src)
187  {
188  const std::size_t N = src.size();
190  for (std::size_t i = 0; i < N; ++i)
191  {
192  NV_Ith_S(dst, i) = src[i];
193  }
194  }
195 
196  void
197  copy(Vector<double> &dst, const N_Vector &src)
198  {
199  const std::size_t N = dst.size();
201  for (std::size_t i = 0; i < N; ++i)
202  {
203  dst[i] = NV_Ith_S(src, i);
204  }
205  }
206 
207  void
208  copy(N_Vector &dst, const Vector<double> &src)
209  {
210  const std::size_t N = src.size();
212  for (std::size_t i = 0; i < N; ++i)
213  {
214  NV_Ith_S(dst, i) = src[i];
215  }
216  }
217  } // namespace internal
218 } // namespace SUNDIALS
220 
221 #endif
IndexSet
Definition: index_set.h:74
TrilinosWrappers::MPI::Vector
Definition: trilinos_vector.h:400
VectorOperation::insert
@ insert
Definition: vector_operation.h:49
BlockVector< double >
StandardExceptions::ExcNotImplemented
static ::ExceptionBase & ExcNotImplemented()
BlockVectorBase::size
std::size_t size() const
PETScWrappers::VectorBase::compress
void compress(const VectorOperation::values operation)
Definition: petsc_vector_base.cc:361
TrilinosWrappers::MPI::Vector::locally_owned_elements
IndexSet locally_owned_elements() const
SUNDIALS::internal::copy
void copy(TrilinosWrappers::MPI::Vector &dst, const N_Vector &src)
Definition: copy.cc:65
TrilinosWrappers::MPI::BlockVector
Definition: trilinos_parallel_block_vector.h:75
Vector::size
size_type size() const
SUNDIALS
Definition: arkode.h:58
BlockVectorBase::compress
void compress(::VectorOperation::values operation)
copy.h
IndexSet::n_elements
size_type n_elements() const
Definition: index_set.h:1833
TrilinosWrappers::MPI::Vector::compress
void compress(::VectorOperation::values operation)
Definition: trilinos_vector.cc:583
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
AssertDimension
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1579
SUNDIALS::internal::N_Vector_length
std::size_t N_Vector_length(const N_Vector &vec)
Definition: copy.cc:33
PETScWrappers::MPI::BlockVector
Definition: petsc_block_vector.h:61
PETScWrappers::MPI::Vector
Definition: petsc_vector.h:158
StandardExceptions::ExcInternalError
static ::ExceptionBase & ExcInternalError()
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
IndexSet::nth_index_in_set
size_type nth_index_in_set(const size_type local_index) const
Definition: index_set.h:1881
PETScWrappers::VectorBase::locally_owned_elements
IndexSet locally_owned_elements() const
BlockVectorBase::locally_owned_elements
IndexSet locally_owned_elements() const
internal
Definition: aligned_vector.h:369
LAPACKSupport::N
static const char N
Definition: lapack_support.h:159
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
Vector< double >