Thyra Version of the Day
Loading...
Searching...
No Matches
Thyra_MultiVectorStdOpsTester_def.hpp
1// @HEADER
2// *****************************************************************************
3// Thyra: Interfaces and Support for Abstract Numerical Algorithms
4//
5// Copyright 2004 NTESS and the Thyra contributors.
6// SPDX-License-Identifier: BSD-3-Clause
7// *****************************************************************************
8// @HEADER
9
10#ifndef THYRA_MULTI_VECTOR_STD_OPS_TESTER_HPP
11#define THYRA_MULTI_VECTOR_STD_OPS_TESTER_HPP
12
13#include "Thyra_MultiVectorStdOpsTester_decl.hpp"
14#include "Thyra_MultiVectorStdOps.hpp"
15#include "Thyra_VectorStdOps.hpp"
16#include "Thyra_TestingTools.hpp"
17
18namespace Thyra {
19
20// MultiVectorStdOpsTester
21
22template <class Scalar>
24 const ScalarMag &warning_tol_in
25 ,const ScalarMag &error_tol_in
26 ,const int num_mv_cols_in
27 )
28 :warning_tol_(warning_tol_in)
29 ,error_tol_(error_tol_in)
30 ,num_mv_cols_(num_mv_cols_in)
31{}
32
33template <class Scalar>
35 const VectorSpaceBase<Scalar> &vecSpc
36 ,std::ostream *out
37 ,const bool &/* dumpAll */
38 )
39{
40 using Teuchos::as;
41 using Teuchos::ptr;
42 using Teuchos::tuple;
44
45 if(out)
46 *out << "\n*** Entering MultiVectorStdOpsTester<"<<ST::name()<<">::checkStdOps(...) ...\n"
47 << "using a \'" << vecSpc.description() << "\' object ...\n";
48
49 bool success = true;
50 if(out) *out << "\nvecSpc.dim() = " << vecSpc.dim() << std::endl;
51
52 const Ordinal n = vecSpc.dim();
53
54 const Scalar
55 two = as<Scalar>(2.0),
56 three = as<Scalar>(3.0),
57 four = as<Scalar>(4.0);
58
59 int tc = 0;
60
61 if(out) *out << "\nCreating MultiVectorBase objects V1, V2, and V3 ...\n";
63 V1 = createMembers(vecSpc,num_mv_cols()),
64 V2 = createMembers(vecSpc,num_mv_cols()),
65 V3 = createMembers(vecSpc,num_mv_cols());
66
67 if(out) *out << "\nassign(V1.ptr(),-2.0);\n";
68 assign(V1.ptr(),Scalar(-two));
69
70 Teuchos::Array<Scalar> scalars1(num_mv_cols());
71 Teuchos::Array<Scalar> scalars2(num_mv_cols());
72 Teuchos::Array<ScalarMag> mags1(num_mv_cols());
73 Teuchos::Array<ScalarMag> mags2(num_mv_cols());
74
75 // sums
76 if(out) *out << "\n"<<tc<<") sums(*V1);\n";
77 ++tc;
78 {
79 sums(*V1, scalars1());
80 for (Ordinal i = 0; i < scalars2.size(); ++i)
81 scalars2[i] = -two*as<Scalar>(vecSpc.dim());
82 if (!testRelErrors<Scalar, Scalar, ScalarMag>(
83 "sums(*V1)", scalars1(),
84 "-2.0*n", scalars2(),
85 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
86 )
87 ) success = false;
88 }
89
90 // norms_1
91 if(out) *out << "\n"<<tc<<") norms_1(*V1);\n";
92 ++tc;
93 {
94 norms_1(*V1, mags1());
95 for (Ordinal i = 0; i < mags2.size(); ++i)
96 mags2[i] = ST::magnitude(two)*as<ScalarMag>(vecSpc.dim());
97 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
98 "norms_1(*V1)", mags1(),
99 "2.0*n", mags2(),
100 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
101 )
102 ) success = false;
103 }
104
105 // norms_2
106 if(out) *out << "\n"<<tc<<") norms_2(*V1);\n";
107 ++tc;
108 {
109 norms_2(*V1, mags1());
110 for (Ordinal i = 0; i < mags2.size(); ++i)
111 mags2[i] = ST::magnitude(two * ST::squareroot(as<Scalar>(n)));
112 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
113 "norms_2(*V1)", mags1(),
114 "2.0*sqrt(n)", mags2(),
115 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
116 )
117 ) success = false;
118 }
119
120 // norms_inf
121 if(out) *out << "\n"<<tc<<") norms_inf(*V1);\n";
122 ++tc;
123 {
124 norms_inf(*V1, mags1());
125 for (Ordinal i = 0; i < mags2.size(); ++i)
126 mags2[i] = ST::magnitude(two);
127 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
128 "norms_inf(*V1)", mags1(),
129 "2.0", mags2(),
130 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
131 )
132 ) success = false;
133 }
134
135 // assign(scalar)
136 if(out) *out << "\n"<<tc<<") assign(V2.ptr(), alpha);\n";
137 ++tc;
138 {
139 assign(V2.ptr(), three);
140 norms_2(*V2, mags1());
141 for (Ordinal i = 0; i < mags2.size(); ++i)
142 mags2[i] = ST::magnitude(three * ST::squareroot(as<Scalar>(n)));
143 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
144 "norms_2(*V2)", mags1(),
145 "3.0*sqrt(n)", mags2(),
146 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
147 )
148 ) success = false;
149 }
150
151 // assign(MV)
152 if(out) *out << "\n"<<tc<<") assign(V2.ptr(), *V1);\n";
153 ++tc;
154 assign(V2.ptr(), *V1);
155 norms_2(*V1, mags1());
156 norms_2(*V2, mags2());
157 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
158 "norms_2(*V1)", mags1(),
159 "norms_2(*V2)", mags2(),
160 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
161 )
162 ) success = false;
163
164 // scale
165 if(out) *out << "\n"<<tc<<") scale(alpha,V2.ptr());\n";
166 ++tc;
167 {
168 Scalar alpha = as<Scalar>(1.2345);
169 assign(V2.ptr(), *V1);
170 scale(alpha, V2.ptr());
171 norms_2(*V1, mags1());
172 for (Ordinal i = 0; i < mags1.size(); ++i)
173 mags1[i] *= ST::magnitude(alpha);
174 norms_2(*V2, mags2());
175 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
176 "norms_2(alpha*V1)", mags1(),
177 "alpha*norms_2(V1)", mags2(),
178 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
179 )
180 ) success = false;
181 }
182
183 // scaleUpdate
184 if(out) *out << "\n"<<tc<<") scaleUpdate(a,V1,V2.ptr());\n";
185 ++tc;
186 {
187 Teuchos::RCP<VectorBase<Scalar> > a = createMember(vecSpc);
188 assign(a.ptr(), two);
189 assign(V2.ptr(), four);
190 scaleUpdate(*a, *V1, V2.ptr()); // V2(i,j) = 2.0*(-2.0) + 4.0
191 norms_2(*V2, mags1());
192 if (!testMaxErrors<Scalar>(
193 "norms_2(*V2)", mags1(),
194 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
195 )
196 ) success = false;
197 }
198
199 // update(alpha, U, V.ptr())
200 if(out) *out << "\n"<<tc<<") update(a,V1,V2.ptr());\n";
201 ++tc;
202 {
203 Scalar alpha = as<Scalar>(1.2345);
204 assign(V2.ptr(), three);
205 assign(V3.ptr(), *V1);
206 scale(alpha, V3.ptr());
207 Vp_V(V3.ptr(), *V2);
208 update(alpha, *V1, V2.ptr());
209 norms_2(*V2, mags1());
210 norms_2(*V3, mags2());
211 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
212 "norms_2(*V2)", mags1(),
213 "norms_2(*V3)", mags2(),
214 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
215 )
216 ) success = false;
217 }
218
219 // update(alpha, beta, U, V.ptr())
220 if(out) *out << "\n"<<tc<<") update(alpha,beta,*V1,V2.ptr());\n";
221 ++tc;
222 {
223 Teuchos::Array<Scalar> alpha(num_mv_cols());
224 for (Ordinal i = 0; i < alpha.size(); ++i)
225 alpha[i] = as<Scalar>(i+1);
226 Scalar beta = as<Scalar>(1.2345);
227 assign(V2.ptr(), three);
228 assign(V3.ptr(), *V1);
229 scale(beta, V3.ptr());
230 for (Ordinal i = 0; i < alpha.size(); ++i)
231 scale(alpha[i], V3->col(i).ptr());
232 Vp_V(V3.ptr(), *V2);
233 Teuchos::ArrayView<const Scalar> alphaView = alpha();
234 update(alphaView, beta, *V1, V2.ptr());
235 norms_2(*V2, mags1());
236 norms_2(*V3, mags2());
237 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
238 "norms_2(*V2)", mags1(),
239 "norms_2(*V3)", mags2(),
240 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
241 )
242 ) success = false;
243 }
244
245 // update(U, alpha, beta, V.ptr())
246 if(out) *out << "\n"<<tc<<") update(*V1,alpha,beta,V2.ptr());\n";
247 ++tc;
248 {
249 Teuchos::Array<Scalar> alpha(num_mv_cols());
250 for (Ordinal i = 0; i < alpha.size(); ++i)
251 alpha[i] = as<Scalar>(i+1);
252 Scalar beta = as<Scalar>(1.2345);
253 assign(V2.ptr(), three);
254 assign(V3.ptr(), *V2);
255 scale(beta, V3.ptr());
256 for (Ordinal i = 0; i < alpha.size(); ++i)
257 scale(alpha[i], V3->col(i).ptr());
258 Vp_V(V3.ptr(), *V1);
259 Teuchos::ArrayView<const Scalar> alphaView = alpha();
260 update(*V1, alphaView, beta, V2.ptr());
261 norms_2(*V2, mags1());
262 norms_2(*V3, mags2());
263 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
264 "norms_2(*V2)", mags1(),
265 "norms_2(*V3)", mags2(),
266 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
267 )
268 ) success = false;
269 }
270
271 // linear_combination
272 if(out) *out << "\n"<<tc<<") linear_combination({alpha,beta,gamma},{V1.ptr(),V2.ptr(),V3.ptr()},0.0,V4.ptr());\n";
273 ++tc;
274 {
275 Scalar alpha = two, beta = -three, gamma = three;
276 Teuchos::RCP<MultiVectorBase<Scalar> > V4 = createMembers(vecSpc,num_mv_cols());
277 assign(V2.ptr(), two);
278 assign(V3.ptr(), four);
279 linear_combination<Scalar>(
280 tuple<Scalar>(alpha, beta, gamma),
281 tuple<Ptr<const MultiVectorBase<Scalar> > >(V1.ptr(), V2.ptr(), V3.ptr()),
282 as<Scalar>(0.0),
283 V4.ptr()); // V4(i,j) = 2.0(-2.0) + (-3.0)(2.0) + 3.0(4.0)
284 norms_2(*V4, mags1());
285 for (Ordinal i = 0; i < mags2.size(); ++i)
286 mags2[i] = ST::magnitude(two * ST::squareroot(as<Scalar>(n)));
287 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
288 "norms_2(*V4)", mags1(),
289 "2.0*sqrt(n)", mags2(),
290 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
291 )
292 ) success = false;
293 }
294
295 if(out) *out << "\n"<<tc<<") linear_combination({alpha,beta,gamma},{V1.ptr(),V2.ptr(),V3.ptr()},0.5,V4.ptr());\n";
296 ++tc;
297 {
298 Scalar alpha = two, beta = -three, gamma = three;
299 Teuchos::RCP<MultiVectorBase<Scalar> > V4 = createMembers(vecSpc,num_mv_cols());
300 assign(V2.ptr(), two);
301 assign(V3.ptr(), four);
302 assign(V4.ptr(), -four);
303 linear_combination<Scalar>(
304 tuple<Scalar>(alpha, beta, gamma),
305 tuple<Ptr<const MultiVectorBase<Scalar> > >(V1.ptr(), V2.ptr(), V3.ptr()),
306 as<Scalar>(0.5),
307 V4.ptr()); // V4(i,j) = 0.5(-4.0) + 2.0(-2.0) + (-3.0)(2.0) + 3.0(4.0)
308 norms_2(*V4, mags1());
309 if (!testMaxErrors<Scalar>(
310 "norms_2(*V4)", mags1(),
311 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
312 )
313 ) success = false;
314 }
315
316 // Vt_S
317 if(out) *out << "\n"<<tc<<") Vt_S(V1.ptr(),alpha);\n";
318 ++tc;
319 {
320 Scalar alpha = as<Scalar>(1.2345);
321 assign(V2.ptr(), *V1);
322 Vt_S(V2.ptr(), alpha);
323 norms_2(*V1, mags1());
324 for (Ordinal i = 0; i < mags1.size(); ++i)
325 mags1[i] *= ST::magnitude(alpha);
326 norms_2(*V2, mags2());
327 if (!testRelErrors<ScalarMag, ScalarMag, ScalarMag>(
328 "norms_2(alpha*V1)", mags1(),
329 "alpha*norms_2(V1)", mags2(),
330 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
331 )
332 ) success = false;
333 }
334
335 // Vp_S
336 if(out) *out << "\n"<<tc<<") Vp_S(V2.ptr(),alpha);\n";
337 ++tc;
338 {
339 assign(V2.ptr(), *V1);
340 Vp_S(V2.ptr(), two); // V2(i,j) = -2.0 + 2.0
341 norms_2(*V2, mags1());
342 if (!testMaxErrors<Scalar>(
343 "norms_2(V2)", mags1(),
344 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
345 )
346 ) success = false;
347 }
348
349 // Vp_V
350 if(out) *out << "\n"<<tc<<") Vp_V(V2.ptr(),*V1);\n";
351 ++tc;
352 {
353 assign(V2.ptr(), two);
354 Vp_V(V2.ptr(), *V1); // V2(i,j) = -2.0 + 2.0
355 norms_2(*V2, mags1());
356 if (!testMaxErrors<Scalar>(
357 "norms_2(V2)", mags1(),
358 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
359 )
360 ) success = false;
361 }
362
363
364 // V_VpV
365 if(out) *out << "\n"<<tc<<") V_VpV(V3.ptr(),*V1,*V2);\n";
366 ++tc;
367 {
368 assign(V2.ptr(), two);
369 V_VpV(V3.ptr(), *V1, *V2); // V3(i,j) = -2.0 + 2.0
370 norms_2(*V3, mags1());
371 if (!testMaxErrors<Scalar>(
372 "norms_2(V3)", mags1(),
373 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
374 )
375 ) success = false;
376 }
377
378
379 // V_VmV
380 if(out) *out << "\n"<<tc<<") V_VmV(V3.ptr(),*V1,*V2);\n";
381 ++tc;
382 {
383 assign(V2.ptr(), -two);
384 V_VmV(V3.ptr(), *V1, *V2); // V3(i,j) = -2.0 - (-2.0)
385 norms_2(*V3, mags1());
386 if (!testMaxErrors<Scalar>(
387 "norms_2(V3)", mags1(),
388 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
389 )
390 ) success = false;
391 }
392
393
394 // V_StVpV
395 if(out) *out << "\n"<<tc<<") V_StVpV(V3.ptr(),alpha,*V1,*V2);\n";
396 ++tc;
397 {
398 Scalar alpha = as<Scalar>(1.2345);
399 assign(V2.ptr(), three);
400 V_StVpV(V3.ptr(), alpha, *V1, *V2);
401 scale(alpha, V1.ptr());
402 Vp_V(V2.ptr(), *V1);
403 V_VmV(V3.ptr(), *V2, *V3);
404 norms_2(*V3, mags1());
405 if (!testMaxErrors<Scalar>(
406 "norms_2(V3)", mags1(),
407 "error_tol", error_tol(), "warning_tol", warning_tol(), ptr(out)
408 )
409 ) success = false;
410 }
411
412
413 if(out) *out
414 << "\n*** Leaving MultiVectorStdOpsTester<"<<ST::name()<<">::checkStdOps(...) ...\n";
415
416 return success;
417
418}
419
420} // namespace Thyra
421
422#endif // THYRA_MULTI_VECTOR_STD_OPS_TESTER_HPP
size_type size() const
virtual std::string description() const
Ptr< T > ptr() const
Interface for a collection of column vectors called a multi-vector.
bool checkStdOps(const VectorSpaceBase< Scalar > &vecSpc, std::ostream *out=0, const bool &dumpAll=false)
Run the tests using a vector space.
MultiVectorStdOpsTester(const ScalarMag &warning_tol=0, const ScalarMag &error_tol=0, const int num_mv_cols=4)
Teuchos::ScalarTraits< Scalar >::magnitudeType ScalarMag
Abstract interface for objects that represent a space for vectors.
virtual Ordinal dim() const =0
Return the dimension of the vector space.
Teuchos::Ordinal Ordinal
Type for the dimension of a vector space. `*.
TypeTo as(const TypeFrom &t)