Reference documentation for deal.II version 9.0.0
copy.cc
1 //-----------------------------------------------------------
2 //
3 // Copyright (C) 2017 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 at
12 // the top level of the deal.II distribution.
13 //
14 //-----------------------------------------------------------
15 
16 #include <deal.II/sundials/copy.h>
17 
18 #ifdef DEAL_II_WITH_SUNDIALS
19 
20 DEAL_II_NAMESPACE_OPEN
21 namespace SUNDIALS
22 {
23  namespace internal
24  {
32  inline
33  std::size_t
34  N_Vector_length(const N_Vector &vec)
35  {
36  const N_Vector_ID id = N_VGetVectorID(vec);
37  long int length = -1;
38  switch (id)
39  {
40  case SUNDIALS_NVEC_SERIAL:
41  {
42  length = NV_LENGTH_S(vec);
43  break;
44  }
45 #ifdef DEAL_II_WITH_MPI
46  case SUNDIALS_NVEC_PARALLEL:
47  {
48  length = NV_LOCLENGTH_P(vec);
49  break;
50  }
51 #endif
52  default:
53  Assert(false, ExcNotImplemented());
54  }
55 
56  Assert(length >= 0, ExcInternalError());
57  return static_cast<std::size_t>(length);
58  }
59 
60 #ifdef DEAL_II_WITH_MPI
61 
62 #ifdef DEAL_II_WITH_TRILINOS
63 
64 
65  void copy(TrilinosWrappers::MPI::Vector &dst, const N_Vector &src)
66  {
67  const IndexSet is = dst.locally_owned_elements();
68  const size_t N = is.n_elements();
69  AssertDimension(N, N_Vector_length(src));
70  for (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 copy(N_Vector &dst, const TrilinosWrappers::MPI::Vector &src)
78  {
79  const IndexSet is = src.locally_owned_elements();
80  const size_t N = is.n_elements();
81  AssertDimension(N, N_Vector_length(dst));
82  for (size_t i=0; i<N; ++i)
83  {
84  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
85  }
86  }
87 
88  void copy(TrilinosWrappers::MPI::BlockVector &dst, const N_Vector &src)
89  {
90  const IndexSet is = dst.locally_owned_elements();
91  const size_t N = is.n_elements();
92  AssertDimension(N, N_Vector_length(src));
93  for (size_t i=0; i<N; ++i)
94  {
95  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
96  }
98  }
99 
100  void copy(N_Vector &dst, const TrilinosWrappers::MPI::BlockVector &src)
101  {
102  IndexSet is = src.locally_owned_elements();
103  const size_t N = is.n_elements();
104  AssertDimension(N, N_Vector_length(dst));
105  for (size_t i=0; i<N; ++i)
106  {
107  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
108  }
109  }
110 
111 #endif //DEAL_II_WITH_TRILINOS
112 
113 #ifdef DEAL_II_WITH_PETSC
114 #ifndef PETSC_USE_COMPLEX
115 
116  void copy(PETScWrappers::MPI::Vector &dst, const N_Vector &src)
117  {
118  const IndexSet is = dst.locally_owned_elements();
119  const size_t N = is.n_elements();
120  AssertDimension(N, N_Vector_length(src));
121  for (size_t i=0; i<N; ++i)
122  {
123  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
124  }
126  }
127 
128  void copy(N_Vector &dst, const PETScWrappers::MPI::Vector &src)
129  {
130  const IndexSet is = src.locally_owned_elements();
131  const size_t N = is.n_elements();
132  AssertDimension(N, N_Vector_length(dst));
133  for (size_t i=0; i<N; ++i)
134  {
135  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
136  }
137  }
138 
139  void copy(PETScWrappers::MPI::BlockVector &dst, const N_Vector &src)
140  {
141  const IndexSet is = dst.locally_owned_elements();
142  const size_t N = is.n_elements();
143  AssertDimension(N, N_Vector_length(src));
144  for (size_t i=0; i<N; ++i)
145  {
146  dst[is.nth_index_in_set(i)] = NV_Ith_P(src, i);
147  }
149  }
150 
151  void copy(N_Vector &dst, const PETScWrappers::MPI::BlockVector &src)
152  {
153  const IndexSet is = src.locally_owned_elements();
154  const size_t N = is.n_elements();
155  AssertDimension(N, N_Vector_length(dst));
156  for (size_t i=0; i<N; ++i)
157  {
158  NV_Ith_P(dst, i) = src[is.nth_index_in_set(i)];
159  }
160  }
161 
162 #endif //PETSC_USE_COMPLEX
163 #endif //DEAL_II_WITH_PETSC
164 
165 #endif //mpi
166 
167  void copy(BlockVector<double> &dst, const N_Vector &src)
168  {
169  const size_t N = dst.size();
170  AssertDimension(N_Vector_length(src), N);
171  for (size_t i=0; i<N; ++i)
172  {
173  dst[i] = NV_Ith_S(src, i);
174  }
175  }
176 
177  void copy(N_Vector &dst, const BlockVector<double> &src)
178  {
179  const size_t N = src.size();
180  AssertDimension(N_Vector_length(dst), N);
181  for (size_t i=0; i<N; ++i)
182  {
183  NV_Ith_S(dst, i) = src[i];
184  }
185  }
186 
187  void copy(Vector<double> &dst, const N_Vector &src)
188  {
189  const size_t N = dst.size();
190  AssertDimension(N_Vector_length(src), N);
191  for (size_t i=0; i<N; ++i)
192  {
193  dst[i] = NV_Ith_S(src, i);
194  }
195  }
196 
197  void copy(N_Vector &dst, const Vector<double> &src)
198  {
199  const size_t N = src.size();
200  AssertDimension(N_Vector_length(dst), N);
201  for (size_t i=0; i<N; ++i)
202  {
203  NV_Ith_S(dst, i) = src[i];
204  }
205  }
206  }
207 }
208 DEAL_II_NAMESPACE_CLOSE
209 
210 #endif
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1248
void compress(::VectorOperation::values operation)
size_type nth_index_in_set(const unsigned int local_index) const
Definition: index_set.h:1769
IndexSet locally_owned_elements() const
void compress(const VectorOperation::values operation)
#define Assert(cond, exc)
Definition: exceptions.h:1142
void compress(::VectorOperation::values operation)
IndexSet locally_owned_elements() const
size_type size() const
IndexSet locally_owned_elements() const
static ::ExceptionBase & ExcNotImplemented()
size_type n_elements() const
Definition: index_set.h:1717
static ::ExceptionBase & ExcInternalError()