Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
BelosXpetraAdapterMultiVector_MP_Vector.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// MueLu: A package for multigrid based preconditioning
6// Copyright 2012 Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact
39// Jonathan Hu (jhu@sandia.gov)
40// Andrey Prokopenko (aprokop@sandia.gov)
41// Ray Tuminaro (rstumin@sandia.gov)
42//
43// ***********************************************************************
44//
45// @HEADER
46#ifndef BELOS_XPETRA_ADAPTER_MULTIVECTOR_MP_VECTOR_HPP
47#define BELOS_XPETRA_ADAPTER_MULTIVECTOR_MP_VECTOR_HPP
48
49#include "BelosXpetraAdapterMultiVector.hpp"
52
53#ifdef HAVE_XPETRA_TPETRA
54
55namespace Belos { // should be moved to Belos or Xpetra?
56
57 using Teuchos::RCP;
58 using Teuchos::rcp;
59
61 //
62 // Implementation of the Belos::MultiVecTraits for Xpetra::MultiVector.
63 //
65
70 template<class Storage, class LO, class GO, class Node>
71 class MultiVecTraits<typename Storage::value_type,
72 Xpetra::MultiVector<Sacado::MP::Vector<Storage>,
73 LO,GO,Node> > {
74 public:
75 typedef typename Storage::ordinal_type s_ordinal;
76 typedef typename Storage::value_type BaseScalar;
78 typedef typename Tpetra::MultiVector<Scalar,LO,GO,Node>::dot_type dot_type;
79 typedef typename Tpetra::MultiVector<Scalar,LO,GO,Node>::mag_type mag_type;
80
81 private:
82
83 typedef Xpetra::TpetraMultiVector<Scalar,LO,GO,Node> TpetraMultiVector;
84 typedef MultiVecTraits<dot_type,Tpetra::MultiVector<Scalar,LO,GO,Node> > MultiVecTraitsTpetra;
85
86 public:
87
88#ifdef HAVE_BELOS_XPETRA_TIMERS
89 static RCP<Teuchos::Time> mvTimesMatAddMvTimer_, mvTransMvTimer_;
90#endif
91
92 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > Clone( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const int numvecs )
93 {
94 if (mv.getMap()->lib() == Xpetra::UseTpetra)
95 return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::Clone(toTpetra(mv), numvecs)));
96 }
97
98 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneCopy( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
99 {
100 if (mv.getMap()->lib() == Xpetra::UseTpetra)
101 return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv))));
102 }
103
104 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneCopy( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
105 {
106 if (mv.getMap()->lib() == Xpetra::UseTpetra)
107 return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
108 }
109
110 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
111 CloneCopy (const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
112 const Teuchos::Range1D& index)
113 {
114 if (mv.getMap()->lib() == Xpetra::UseTpetra)
115 return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneCopy(toTpetra(mv), index)));
116 }
117
118 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneViewNonConst( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
119 {
120 if (mv.getMap()->lib() == Xpetra::UseTpetra)
121 return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
122 }
123
124 static RCP<Xpetra::MultiVector<Scalar,LO,GO,Node> >
125 CloneViewNonConst(Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
126 const Teuchos::Range1D& index)
127 {
128 if (mv.getMap()->lib() == Xpetra::UseTpetra)
129 return rcp(new TpetraMultiVector(MultiVecTraitsTpetra::CloneViewNonConst(toTpetra(mv), index)));
130 }
131
132 static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> > CloneView(const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<int>& index )
133 {
134 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
135 //TODO: double check if the const_cast is safe here.
136 RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
137 return rcp(new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
138 }
139 }
140
141 static RCP<const Xpetra::MultiVector<Scalar,LO,GO,Node> >
142 CloneView (const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv,
143 const Teuchos::Range1D& index)
144 {
145 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
146 //TODO: double check if the const_cast is safe here.
147 RCP<const Tpetra::MultiVector<Scalar,LO,GO,Node> > r = MultiVecTraitsTpetra::CloneView(toTpetra(mv), index);
148 return rcp(new TpetraMultiVector(Teuchos::rcp_const_cast<Tpetra::MultiVector<Scalar,LO,GO,Node> >(r)));
149 }
150 }
151
152 static ptrdiff_t GetGlobalLength( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
153 {
154 if (mv.getMap()->lib() == Xpetra::UseTpetra)
155 return MultiVecTraitsTpetra::GetGlobalLength(toTpetra(mv));
156 }
157
158 static int GetNumberVecs( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
159 {
160 if (mv.getMap()->lib() == Xpetra::UseTpetra)
161 return MultiVecTraitsTpetra::GetNumberVecs(toTpetra(mv));
162 }
163
164 static bool HasConstantStride( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
165 {
166 if (mv.getMap()->lib() == Xpetra::UseTpetra)
167 return MultiVecTraitsTpetra::HasConstantStride(toTpetra(mv));
168 }
169
170 static void MvTimesMatAddMv( dot_type alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
171 const Teuchos::SerialDenseMatrix<int,dot_type>& B,
172 dot_type beta, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
173 {
174#ifdef HAVE_BELOS_XPETRA_TIMERS
175 Teuchos::TimeMonitor lcltimer(*mvTimesMatAddMvTimer_);
176#endif
177 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
178 MultiVecTraitsTpetra::MvTimesMatAddMv(alpha, toTpetra(A), B, beta, toTpetra(mv));
179 return;
180 }
181 }
182
183 static void MvAddMv( Scalar alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, Scalar beta, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
184 {
185 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
186 MultiVecTraitsTpetra::MvAddMv(alpha, toTpetra(A), beta, toTpetra(B), toTpetra(mv));
187 return;
188 }
189 }
190
191 static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha )
192 {
193 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
194 MultiVecTraitsTpetra::MvScale(toTpetra(mv), alpha);
195 return;
196 }
197 }
198
199 static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<BaseScalar>& alphas )
200 {
201 std::vector<Scalar> alphas_mp(alphas.size());
202 const size_t sz = alphas.size();
203 for (size_t i=0; i<sz; ++i)
204 alphas_mp[i] = alphas[i];
205 MvScale (mv, alphas_mp);
206 }
207
208 static void MvScale ( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, const std::vector<Scalar>& alphas )
209 {
210 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
211 MultiVecTraitsTpetra::MvScale(toTpetra(mv), alphas);
212 return;
213 }
214 }
215
216 static void MvTransMv( dot_type alpha, const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, Teuchos::SerialDenseMatrix<int,dot_type>& C)
217 {
218#ifdef HAVE_BELOS_XPETRA_TIMERS
219 Teuchos::TimeMonitor lcltimer(*mvTransMvTimer_);
220#endif
221
222 if (A.getMap()->lib() == Xpetra::UseTpetra) {
223 MultiVecTraitsTpetra::MvTransMv(alpha, toTpetra(A), toTpetra(B), C);
224 return;
225 }
226 }
227
228 static void MvDot( const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const Xpetra::MultiVector<Scalar,LO,GO,Node>& B, std::vector<dot_type> &dots)
229 {
230 if (A.getMap()->lib() == Xpetra::UseTpetra) {
231 MultiVecTraitsTpetra::MvDot(toTpetra(A), toTpetra(B), dots);
232 return;
233 }
234 }
235
236 static void MvNorm(const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::vector<mag_type> &normvec, NormType type=TwoNorm)
237 {
238 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
239 MultiVecTraitsTpetra::MvNorm(toTpetra(mv), normvec, type);
240 return;
241 }
242 }
243
244 static void SetBlock( const Xpetra::MultiVector<Scalar,LO,GO,Node>& A, const std::vector<int>& index, Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
245 {
246 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
247 MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
248 return;
249 }
250 }
251
252 static void
253 SetBlock (const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
254 const Teuchos::Range1D& index,
255 Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
256 {
257 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
258 MultiVecTraitsTpetra::SetBlock(toTpetra(A), index, toTpetra(mv));
259 return;
260 }
261 }
262
263 static void
264 Assign (const Xpetra::MultiVector<Scalar,LO,GO,Node>& A,
265 Xpetra::MultiVector<Scalar,LO,GO,Node>& mv)
266 {
267 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
268 MultiVecTraitsTpetra::Assign(toTpetra(A), toTpetra(mv));
269 return;
270 }
271 }
272
273 static void MvRandom( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv )
274 {
275 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
276 MultiVecTraitsTpetra::MvRandom(toTpetra(mv));
277 return;
278 }
279 }
280
281 static void MvInit( Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, Scalar alpha = Teuchos::ScalarTraits<Scalar>::zero() )
282 {
283 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
284 MultiVecTraitsTpetra::MvInit(toTpetra(mv), alpha);
285 return;
286 }
287 }
288
289 static void MvPrint( const Xpetra::MultiVector<Scalar,LO,GO,Node>& mv, std::ostream& os )
290 {
291 if (mv.getMap()->lib() == Xpetra::UseTpetra) {
292 MultiVecTraitsTpetra::MvPrint(toTpetra(mv), os);
293 return;
294 }
295 }
296
297 };
298
299} // end of Belos namespace
300
301#endif
302
303#endif