Thyra Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
TpetraThyraWrappers_UnitTests.cpp
Go to the documentation of this file.
1/*
2// @HEADER
3// ***********************************************************************
4//
5// Thyra: Interfaces and Support for Abstract Numerical Algorithms
6// Copyright (2004) Sandia Corporation
7//
8// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9// license for use of this work by or on behalf of the U.S. Government.
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 Roscoe A. Bartlett (bartlettra@ornl.gov)
39//
40// ***********************************************************************
41// @HEADER
42*/
43
44
46#include "Thyra_VectorSpaceTester.hpp"
47#include "Thyra_VectorStdOpsTester.hpp"
48#include "Thyra_MultiVectorStdOpsTester.hpp"
49#include "Thyra_VectorStdOps.hpp"
50#include "Thyra_MultiVectorStdOps.hpp"
51#include "Thyra_LinearOpTester.hpp"
52#include "Thyra_DefaultProductVector.hpp"
53#include "Thyra_TestingTools.hpp"
54#include "Thyra_ScaledLinearOpBase.hpp"
55#include "Thyra_RowStatLinearOpBase.hpp"
56#include "Thyra_VectorStdOps.hpp"
57#include "Tpetra_CrsMatrix.hpp"
60#include "Teuchos_Tuple.hpp"
61
62
63namespace Thyra {
64
65
66//
67// Helper code and declarations
68//
69
70
71using Teuchos::as;
72using Teuchos::null;
73using Teuchos::RCP;
74using Teuchos::rcp;
76using Teuchos::rcp_dynamic_cast;
77using Teuchos::inOutArg;
78using Teuchos::Comm;
79using Teuchos::tuple;
80
81
82const int g_localDim = 4; // ToDo: Make variable!
83
84
85typedef Tpetra::Map<> TpetraMap_t;
86
87
89createTpetraMap(const int localDim)
90{
92 return Teuchos::rcp(new TpetraMap_t(OT::invalid(), localDim, 0,
94 // ToDo: Pass in the comm?
95}
96
97// ToDo: Above: Vary the LocalOrdinal and GlobalOrdinal types?
98
99
100template<class Scalar>
101RCP<const VectorSpaceBase<Scalar> >
102createTpetraVectorSpace(const int localDim)
103{
104 return Thyra::createVectorSpace<Scalar>(createTpetraMap(g_localDim));
105}
106
107
108template<class Scalar>
109RCP<Tpetra::Operator<Scalar> >
110createTriDiagonalTpetraOperator(const int numLocalRows)
111{
112 typedef Tpetra::global_size_t global_size_t;
113 typedef Tpetra::Map<>::global_ordinal_type GO;
114
115 RCP<const Tpetra::Map<> > map = createTpetraMap(numLocalRows);
116
117 const size_t numMyElements = map->getLocalNumElements();
118 const global_size_t numGlobalElements = map->getGlobalNumElements();
119
120 ArrayView<const GO> myGlobalElements = map->getLocalElementList();
121
122 // Create an OTeger vector numNz that is used to build the Petra Matrix.
123 // numNz[i] is the Number of OFF-DIAGONAL term for the ith global equation
124 // on this processor
125
126 Teuchos::Array<size_t> numNz (numMyElements);
127
128 // We are building a tridiagonal matrix where each row has (-1 2 -1)
129 // So we need 2 off-diagonal terms (except for the first and last equation)
130
131 for (size_t i=0; i < numMyElements; ++i) {
132 if (myGlobalElements[i] == 0 || static_cast<global_size_t>(myGlobalElements[i]) == numGlobalElements-1) {
133 // boundary
134 numNz[i] = 2;
135 }
136 else {
137 numNz[i] = 3;
138 }
139 }
140
141 // Create a Tpetra::Matrix using the Map, with a static allocation dictated by numNz
143 Teuchos::rcp( new Tpetra::CrsMatrix<Scalar>(map, numNz ()) );
144
145 // We are done with NumNZ
146 {
148 swap (empty, numNz); // classic idiom for freeing a container
149 }
150
151 // Add rows one-at-a-time
152 // Off diagonal values will always be -1
153 const Scalar two = static_cast<Scalar>( 2.0);
154 const Scalar posOne = static_cast<Scalar>(+1.0);
155 const Scalar negOne = static_cast<Scalar>(-1.0);
156 for (size_t i = 0; i < numMyElements; i++) {
157 if (myGlobalElements[i] == 0) {
158 A->insertGlobalValues( myGlobalElements[i],
159 tuple<GO>(myGlobalElements[i], myGlobalElements[i]+1)(),
160 tuple<Scalar> (two, posOne)()
161 );
162 }
163 else if (static_cast<global_size_t>(myGlobalElements[i]) == numGlobalElements-1) {
164 A->insertGlobalValues( myGlobalElements[i],
165 tuple<GO>(myGlobalElements[i]-1, myGlobalElements[i])(),
166 tuple<Scalar> (negOne, two)()
167 );
168 }
169 else {
170 A->insertGlobalValues( myGlobalElements[i],
171 tuple<GO>(myGlobalElements[i]-1, myGlobalElements[i], myGlobalElements[i]+1)(),
172 tuple<Scalar> (negOne, two, posOne)()
173 );
174 }
175 }
176
177 // Finish up
178 A->fillComplete();
179
180 return A;
181
182}
183
184
185bool showAllTests = false;
186bool dumpAll = false;
188
189
191{
192 Teuchos::UnitTestRepository::getCLP().setOption(
193 "show-all-tests", "no-show-all-tests", &showAllTests, "Show all tests or not" );
194 Teuchos::UnitTestRepository::getCLP().setOption(
195 "dump-all", "no-dump-all", &dumpAll, "Dump all objects being tested" );
196 Teuchos::UnitTestRepository::getCLP().setOption(
197 "run-linear-op-tester", "no-run-linear-op-tester", &runLinearOpTester, "..." );
198}
199
200
201//
202// Unit Tests
203//
204
205
206//
207// convertTpetraToThyraComm
208//
209
211 Scalar )
212{
215 TEST_ASSERT(nonnull(thyraComm));
216}
217
218
219//
220// createVectorSpace
221//
222
224 Scalar )
225{
228 Thyra::createVectorSpace<Scalar>(tpetraMap);
229 TEST_ASSERT(nonnull(vs));
230 out << "vs = " << *vs;
232 rcp_dynamic_cast<const SpmdVectorSpaceBase<Scalar> >(vs, true);
233 TEST_EQUALITY(vs_spmd->localSubDim(), g_localDim);
234 TEST_EQUALITY(vs->dim(), as<Ordinal>(tpetraMap->getGlobalNumElements()));
235
237 RCP<const TpetraMap_t> tpetraMap2 = ConverterT::getTpetraMap(vs);
238 TEST_EQUALITY(tpetraMap2, tpetraMap);
239}
240
241
242//
243// createVector
244//
245
247 Scalar )
248{
249
251
254 Thyra::createVectorSpace<Scalar>(tpetraMap);
255
256 const RCP<Tpetra::Vector<Scalar> > tpetraVector =
257 rcp(new Tpetra::Vector<Scalar>(tpetraMap));
258
259 {
260 const RCP<VectorBase<Scalar> > thyraVector = createVector(tpetraVector, vs);
261 TEST_EQUALITY(thyraVector->space(), vs);
262 const RCP<Tpetra::Vector<Scalar> > tpetraVector2 =
263 ConverterT::getTpetraVector(thyraVector);
264 TEST_EQUALITY(tpetraVector2, tpetraVector);
265 }
266
267 {
268 const RCP<VectorBase<Scalar> > thyraVector = Thyra::createVector(tpetraVector);
269 TEST_INEQUALITY(thyraVector->space(), vs);
270 TEST_ASSERT(thyraVector->space()->isCompatible(*vs));
271 const RCP<Tpetra::Vector<Scalar> > tpetraVector2 =
272 ConverterT::getTpetraVector(thyraVector);
273 TEST_EQUALITY(tpetraVector2, tpetraVector);
274 }
275
276}
277
278
279//
280// createConstVector
281//
282
284 Scalar )
285{
286
288
291 Thyra::createVectorSpace<Scalar>(tpetraMap);
292
293 const RCP<const Tpetra::Vector<Scalar> > tpetraVector =
294 rcp(new Tpetra::Vector<Scalar>(tpetraMap));
295
296 {
297 const RCP<const VectorBase<Scalar> > thyraVector =
298 createConstVector(tpetraVector, vs);
299 TEST_EQUALITY(thyraVector->space(), vs);
300 const RCP<const Tpetra::Vector<Scalar> > tpetraVector2 =
301 ConverterT::getConstTpetraVector(thyraVector);
302 TEST_EQUALITY(tpetraVector2, tpetraVector);
303 }
304
305 {
306 const RCP<const VectorBase<Scalar> > thyraVector =
307 Thyra::createConstVector(tpetraVector);
308 TEST_INEQUALITY(thyraVector->space(), vs);
309 TEST_ASSERT(thyraVector->space()->isCompatible(*vs));
310 const RCP<const Tpetra::Vector<Scalar> > tpetraVector2 =
311 ConverterT::getConstTpetraVector(thyraVector);
312 TEST_EQUALITY(tpetraVector2, tpetraVector);
313 }
314
315}
316
317
318//
319// createMultiVector
320//
321
323 Scalar )
324{
325 typedef Tpetra::Map<>::local_ordinal_type LO;
326 typedef Tpetra::Map<>::global_ordinal_type GO;
327 typedef Tpetra::Map<>::node_type NODE;
329
330 const int numCols = 3;
331
333 const RCP<const VectorSpaceBase<Scalar> > rangeVs =
334 Thyra::createVectorSpace<Scalar>(tpetraMap);
335
336 const RCP<const TpetraMap_t> tpetraLocRepMap =
337 Tpetra::createLocalMapWithNode<LO,GO,NODE>(
338 numCols, tpetraMap->getComm());
339 const RCP<const VectorSpaceBase<Scalar> > domainVs =
340 Thyra::createVectorSpace<Scalar>(tpetraLocRepMap);
341
342 const RCP<Tpetra::MultiVector<Scalar> > tpetraMultiVector =
343 rcp(new Tpetra::MultiVector<Scalar>(tpetraMap, numCols));
344
345 {
346 const RCP<MultiVectorBase<Scalar> > thyraMultiVector =
347 createMultiVector(tpetraMultiVector, rangeVs, domainVs);
348 TEST_EQUALITY(thyraMultiVector->range(), rangeVs);
349 TEST_EQUALITY(thyraMultiVector->domain(), domainVs);
350 const RCP<Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
351 ConverterT::getTpetraMultiVector(thyraMultiVector);
352 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
353 }
354
355 {
356 const RCP<MultiVectorBase<Scalar> > thyraMultiVector =
357 Thyra::createMultiVector(tpetraMultiVector);
358 TEST_INEQUALITY(thyraMultiVector->range(), rangeVs);
359 TEST_INEQUALITY(thyraMultiVector->domain(), domainVs);
360 TEST_ASSERT(thyraMultiVector->range()->isCompatible(*rangeVs));
361 TEST_ASSERT(thyraMultiVector->domain()->isCompatible(*domainVs));
362 const RCP<Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
363 ConverterT::getTpetraMultiVector(thyraMultiVector);
364 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
365 }
366
367}
368
369
370//
371// createConstMultiVector
372//
373
375 Scalar )
376{
377 typedef Tpetra::Map<>::local_ordinal_type LO;
378 typedef Tpetra::Map<>::global_ordinal_type GO;
379 typedef Tpetra::Map<>::node_type NODE;
381
382 const int numCols = 3;
383
385 const RCP<const VectorSpaceBase<Scalar> > rangeVs =
386 Thyra::createVectorSpace<Scalar>(tpetraMap);
387
388 const RCP<const TpetraMap_t> tpetraLocRepMap =
389 Tpetra::createLocalMapWithNode<LO,GO,NODE>(
390 numCols, tpetraMap->getComm());
391 const RCP<const VectorSpaceBase<Scalar> > domainVs =
392 Thyra::createVectorSpace<Scalar>(tpetraLocRepMap);
393
394 const RCP<const Tpetra::MultiVector<Scalar> > tpetraMultiVector =
395 rcp(new Tpetra::MultiVector<Scalar>(tpetraMap, numCols));
396
397 {
398 const RCP<const MultiVectorBase<Scalar> > thyraMultiVector =
399 createConstMultiVector(tpetraMultiVector, rangeVs, domainVs);
400 TEST_EQUALITY(thyraMultiVector->range(), rangeVs);
401 TEST_EQUALITY(thyraMultiVector->domain(), domainVs);
402 const RCP<const Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
403 ConverterT::getConstTpetraMultiVector(thyraMultiVector);
404 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
405 }
406
407 {
408 const RCP<const MultiVectorBase<Scalar> > thyraMultiVector =
409 Thyra::createConstMultiVector(tpetraMultiVector);
410 TEST_INEQUALITY(thyraMultiVector->range(), rangeVs);
411 TEST_INEQUALITY(thyraMultiVector->domain(), domainVs);
412 TEST_ASSERT(thyraMultiVector->range()->isCompatible(*rangeVs));
413 TEST_ASSERT(thyraMultiVector->domain()->isCompatible(*domainVs));
414 const RCP<const Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
415 ConverterT::getConstTpetraMultiVector(thyraMultiVector);
416 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
417 }
418
419}
420
421
422//
423// TeptraVectorSpace
424//
425
426
427TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TeptraVectorSpace,
428 Scalar )
429{
431 Thyra::createVectorSpace<Scalar>(createTpetraMap(g_localDim));
432 const RCP<VectorBase<Scalar> > v = createMember(vs);
434 TEST_EQUALITY(v->space(), vs);
435}
436
437
438//
439// vectorSpaceTester
440//
441
442TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, vectorSpaceTester,
443 Scalar )
444{
446 = createTpetraVectorSpace<Scalar>(g_localDim);
447 Thyra::VectorSpaceTester<Scalar> vectorSpaceTester;
448 vectorSpaceTester.show_all_tests(showAllTests);
449 vectorSpaceTester.dump_all(dumpAll);
450 TEST_ASSERT(vectorSpaceTester.check(*vs, &out));
451}
452
453
454//
455// vectorStdOpsTester
456//
457
458
459TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, vectorStdOpsTester,
460 Scalar )
461{
463 Thyra::createVectorSpace<Scalar>(createTpetraMap(g_localDim));
464 Thyra::VectorStdOpsTester<Scalar> vectorStdOpsTester;
465 vectorStdOpsTester.warning_tol(5.0e-13);
466 vectorStdOpsTester.error_tol(5.0e-14);
467 TEST_ASSERT(vectorStdOpsTester.checkStdOps(*vs, &out));
468}
469
470
471//
472// multiVectorStdOpsTester
473//
474
475
476TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, multiVectorStdOpsTester,
477 Scalar )
478{
480 Thyra::createVectorSpace<Scalar>(createTpetraMap(g_localDim));
481 Thyra::MultiVectorStdOpsTester<Scalar> mvStdOpsTester;
482 mvStdOpsTester.warning_tol(5.0e-13);
483 mvStdOpsTester.error_tol(5.0e-14);
484 TEST_ASSERT(mvStdOpsTester.checkStdOps(*vs, &out));
485}
486
487
488//
489// getTpetraMultiVector
490//
491
492TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, getTpetraMultiVector,
493 Scalar )
494{
496
497 const int numCols = 3;
499 = createTpetraVectorSpace<Scalar>(g_localDim);
500
501 {
502 const RCP<MultiVectorBase<Scalar> > mv = createMembers(vs, numCols);
504 ConverterT::getTpetraMultiVector(mv);
505 TEST_ASSERT(nonnull(tmv));
506 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
507 }
508
509 {
510 const RCP<VectorBase<Scalar> > v = createMember(vs);
512 ConverterT::getTpetraMultiVector(v);
513 TEST_ASSERT(nonnull(tmv));
514 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
515 }
516
517#ifdef THYRA_DEBUG
518 const RCP<VectorBase<Scalar> > pv = Thyra::defaultProductVector<Scalar>();
519 TEST_THROW(ConverterT::getTpetraMultiVector(pv), std::logic_error);
520#endif
521
522}
523
524
525//
526// getConstTpetraMultiVector
527//
528
529TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, getConstTpetraMultiVector,
530 Scalar )
531{
533
534 const int numCols = 3;
536 = createTpetraVectorSpace<Scalar>(g_localDim);
537
538 {
539 const RCP<const MultiVectorBase<Scalar> > mv = createMembers(vs, numCols);
541 ConverterT::getConstTpetraMultiVector(mv);
542 TEST_ASSERT(nonnull(tmv));
543 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
544 }
545
546 {
547 const RCP<const VectorBase<Scalar> > v = createMember(vs);
549 ConverterT::getConstTpetraMultiVector(v);
550 TEST_ASSERT(nonnull(tmv));
551 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
552 }
553
554#ifdef THYRA_DEBUG
555 const RCP<const VectorBase<Scalar> > pv = Thyra::defaultProductVector<Scalar>();
556 TEST_THROW(ConverterT::getConstTpetraMultiVector(pv), std::logic_error);
557#endif
558
559}
560
561
562//
563// TpetraLinearOp
564//
565
567 Scalar )
568{
569
571 using Teuchos::as;
572
573 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
574 createTriDiagonalTpetraOperator<Scalar>(g_localDim);
575 out << "tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
576 TEST_ASSERT(nonnull(tpetraOp));
577
578 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
579 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
580 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
581 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
582 const RCP<const LinearOpBase<Scalar> > thyraLinearOp =
583 Thyra::tpetraLinearOp(rangeSpace, domainSpace, tpetraOp);
584 TEST_ASSERT(nonnull(thyraLinearOp));
585
586 out << "\nCheck that operator returns the right thing ...\n";
587 const RCP<VectorBase<Scalar> > x = createMember(thyraLinearOp->domain());
588 Thyra::V_S(x.ptr(), ST::one());
589 const RCP<VectorBase<Scalar> > y = createMember(thyraLinearOp->range());
590 Thyra::apply<Scalar>(*thyraLinearOp, Thyra::NOTRANS, *x, y.ptr());
591 const Scalar sum_y = sum(*y);
592 TEST_FLOATING_EQUALITY( sum_y, as<Scalar>(3+1+2*(y->space()->dim()-2)),
593 100.0 * ST::eps() );
594
595 out << "\nCheck the general LinearOp interface ...\n";
596 Thyra::LinearOpTester<Scalar> linearOpTester;
597 linearOpTester.show_all_tests(showAllTests);
598 linearOpTester.dump_all(dumpAll);
599 if (runLinearOpTester) {
600 TEST_ASSERT(linearOpTester.check(*thyraLinearOp, Teuchos::inOutArg(out)));
601 }
602
603}
604
605
606//
607// createLinearOp
608//
609
611 Scalar )
612{
613
615
616 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
617 createTriDiagonalTpetraOperator<Scalar>(g_localDim);
618 out << "tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
619
620 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
621 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
622
623 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
624 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
625
626 {
627 const RCP<LinearOpBase<Scalar> > thyraOp =
628 createLinearOp(tpetraOp, rangeSpace, domainSpace);
629 TEST_EQUALITY(thyraOp->range(), rangeSpace);
630 TEST_EQUALITY(thyraOp->domain(), domainSpace);
631 const RCP<Tpetra::Operator<Scalar> > tpetraOp2 =
632 ConverterT::getTpetraOperator(thyraOp);
633 TEST_EQUALITY(tpetraOp2, tpetraOp);
634 }
635
636 {
637 const RCP<LinearOpBase<Scalar> > thyraOp =
638 Thyra::createLinearOp(tpetraOp);
639 TEST_INEQUALITY(thyraOp->range(), rangeSpace);
640 TEST_INEQUALITY(thyraOp->domain(), domainSpace);
641 TEST_ASSERT(thyraOp->range()->isCompatible(*rangeSpace));
642 TEST_ASSERT(thyraOp->domain()->isCompatible(*domainSpace));
643 const RCP<Tpetra::Operator<Scalar> > tpetraOp2 =
644 ConverterT::getTpetraOperator(thyraOp);
645 TEST_EQUALITY(tpetraOp2, tpetraOp);
646 }
647
648}
649
650
651//
652// createConstLinearOp
653//
654
656 Scalar )
657{
658
660
661 const RCP<const Tpetra::Operator<Scalar> > tpetraOp =
662 createTriDiagonalTpetraOperator<Scalar>(g_localDim);
663 out << "tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
664
665 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
666 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
667
668 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
669 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
670
671 {
672 const RCP<const LinearOpBase<Scalar> > thyraOp =
673 createConstLinearOp(tpetraOp, rangeSpace, domainSpace);
674 TEST_EQUALITY(thyraOp->range(), rangeSpace);
675 TEST_EQUALITY(thyraOp->domain(), domainSpace);
676 const RCP<const Tpetra::Operator<Scalar> > tpetraOp2 =
677 ConverterT::getConstTpetraOperator(thyraOp);
678 TEST_EQUALITY(tpetraOp2, tpetraOp);
679 }
680
681 {
682 const RCP<const LinearOpBase<Scalar> > thyraOp =
684 TEST_INEQUALITY(thyraOp->range(), rangeSpace);
685 TEST_INEQUALITY(thyraOp->domain(), domainSpace);
686 TEST_ASSERT(thyraOp->range()->isCompatible(*rangeSpace));
687 TEST_ASSERT(thyraOp->domain()->isCompatible(*domainSpace));
688 const RCP<const Tpetra::Operator<Scalar> > tpetraOp2 =
689 ConverterT::getConstTpetraOperator(thyraOp);
690 TEST_EQUALITY(tpetraOp2, tpetraOp);
691 }
692
693}
694
695
696//
697// Tpetra-implemented methods
698//
699
700
702{
704 TEUCHOS_TEST_FOR_EXCEPTION(timer == null,
705 std::runtime_error,
706 "lookupAndAssertTimer(): timer \"" << label << "\" was not present in Teuchos::TimeMonitor."
707 " Unit test not valid.");
708 return timer;
709}
710
711
712#define CHECK_TPETRA_FUNC_CALL_INCREMENT( timerStr, tpetraCode, thyraCode ) \
713{ \
714 out << "\nTesting that Thyra calls down to " << timerStr << "\n"; \
715 ECHO(tpetraCode); \
716 const RCP<const Time> timer = lookupAndAssertTimer(timerStr); \
717 const int countBefore = timer->numCalls(); \
718 ECHO(thyraCode); \
719 const int countAfter = timer->numCalls(); \
720 TEST_EQUALITY( countAfter, countBefore+1 ); \
721}
722
723
724TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, UseTpetraImplementations,
725 Scalar )
726{
727 using Teuchos::Time;
729 typedef typename ST::magnitudeType Magnitude;
730 typedef VectorSpaceBase<Scalar> VectorSpace;
731 typedef MultiVectorBase<Scalar> MultiVec;
732 //typedef Tpetra::Map<int> TpetraMap;
733 typedef Tpetra::MultiVector<Scalar> TpetraMultiVec;
735
736 const int numCols = 3;
737
738 const RCP<const VectorSpace> vs =
739 createTpetraVectorSpace<Scalar>(g_localDim);
740 const RCP<MultiVec>
741 A = createMembers(vs, numCols),
742 B = createMembers(vs, numCols);
744 tA = TOVE::getTpetraMultiVector(A),
745 tB = TOVE::getTpetraMultiVector(B);
746 Array<Scalar> C(numCols*numCols,ST::zero());
747
748 Teuchos::Array<Magnitude> avMag(numCols);
749 Teuchos::Array<Scalar> avScal(numCols);
750
752 "Tpetra::MultiVector::putScalar()",
753 tA->putScalar(ST::zero()),
754 Thyra::assign(A.ptr(), ST::zero())
755 );
756
758 "Tpetra::MultiVector::dot()",
759 tA->dot(*tB, avScal() ), // norms calls scalarProd, which calls Tpetra::dot
760 Thyra::norms( *A, avMag() )
761 );
762
764 "Tpetra::MultiVector::dot()",
765 tA->dot(*tB, avScal() ),
766 A->range()->scalarProds(*A, *B, avScal() )
767 );
768
769 /*
770 CHECK_TPETRA_FUNC_CALL_INCREMENT(
771 "Tpetra::MultiVector::scale(alpha)",
772 tB->scale( ST::zero() ),
773 Thyra::scale( ST::zero(), B.ptr() )
774 );
775 */
776
777 /*
778 CHECK_TPETRA_FUNC_CALL_INCREMENT(
779 "Tpetra::MultiVector::operator=()",
780 (*tA) = *tB,
781 Thyra::assign( A.ptr(), *B )
782 );
783 */
784
785 /*
786 {
787 RCP<MultiVec> Ctmvb = Thyra::createMembersView(
788 A->domain(),
789 RTOpPack::SubMultiVectorView<Scalar>(
790 0, numCols, 0, numCols,
791 Teuchos::arcpFromArrayView(C()), numCols
792 )
793 );
794 CHECK_TPETRA_FUNC_CALL_INCREMENT(
795 "Tpetra::multiplyOntoHost()",
796 Tpetra::multiplyOntoHost(Teuchos::CONJ_TRANS,Teuchos::NO_TRANS,ST::one(),*tA,*tB, C()),
797 A->apply(Thyra::CONJTRANS,*B,Ctmvb.ptr(),ST::one(),ST::zero())
798 );
799 }
800 */
801
802 /*
803 {
804 RCP<const MultiVec> Ctmvb = Thyra::createMembersView(
805 A->domain(),
806 RTOpPack::ConstSubMultiVectorView<Scalar>(
807 0, numCols, 0, numCols,
808 Teuchos::arcpFromArrayView(C()), numCols
809 )
810 );
811 const RCP<const TpetraMultiVec>
812 tC = TOVE::getConstTpetraMultiVector(Ctmvb);
813 CHECK_TPETRA_FUNC_CALL_INCREMENT(
814 "Tpetra::MultiVector::multiply()",
815 tB->multiply(Teuchos::NO_TRANS,Teuchos::NO_TRANS,ST::one(),*tA,*tC,ST::zero()),
816 A->apply(Thyra::NOTRANS,*Ctmvb,B.ptr(),ST::one(),ST::zero())
817 );
818 }
819 */
820
821/*
822 RCP<Time>
823 timerUpdate = lookupAndAssertTimer("Tpetra::MultiVector::update(alpha,A,beta,B,gamma)"),
824 timerUpdate2 = lookupAndAssertTimer("Tpetra::MultiVector::update(alpha,A,beta)"),
825
826 // TODO: test update(two vector)
827 // TODO: test update(three vector)
828*/
829}
830
831
832#ifdef TPETRA_TEUCHOS_TIME_MONITOR
833# define TPETRA_TIMER_TESTS(SCALAR) \
834 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, UseTpetraImplementations, SCALAR )
835#else
836# define TPETRA_TIMER_TESTS(SCALAR)
837#endif
838
839
840//
841// TpetraLinearOp_EpetraRowMatrix
842//
843
844
845#ifdef HAVE_THYRA_TPETRA_EPETRA
846
847
848TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_EpetraRowMatrix,
849 Scalar )
850{
851
852 using Teuchos::as;
853 using Teuchos::outArg;
854 using Teuchos::rcp_dynamic_cast;
855 using Teuchos::Array;
857
858 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
859 createTriDiagonalTpetraOperator<Scalar>(g_localDim);
860
861 const RCP<LinearOpBase<Scalar> > thyraOp =
862 Thyra::createLinearOp(tpetraOp);
863
864 const RCP<Thyra::TpetraLinearOp<Scalar> > thyraTpetraOp =
865 Teuchos::rcp_dynamic_cast<Thyra::TpetraLinearOp<Scalar> >(thyraOp);
866
867 RCP<const Epetra_Operator> epetraOp;
868 Thyra::EOpTransp epetraOpTransp;
869 Thyra::EApplyEpetraOpAs epetraOpApplyAs;
870 Thyra::EAdjointEpetraOp epetraOpAdjointSupport;
871
872 thyraTpetraOp->getEpetraOpView( outArg(epetraOp), outArg(epetraOpTransp),
873 outArg(epetraOpApplyAs), outArg(epetraOpAdjointSupport) );
874
875 if (typeid(Scalar) == typeid(double)) {
876 TEST_ASSERT(nonnull(epetraOp));
877 const RCP<const Epetra_RowMatrix> epetraRowMatrix =
878 rcp_dynamic_cast<const Epetra_RowMatrix>(epetraOp, true);
879 int numRowEntries = -1;
880 epetraRowMatrix->NumMyRowEntries(1, numRowEntries);
881 TEST_EQUALITY_CONST(numRowEntries, 3);
882 Array<double> row_values(numRowEntries);
883 Array<int> row_indices(numRowEntries);
884 epetraRowMatrix->ExtractMyRowCopy(1, numRowEntries, numRowEntries,
885 row_values.getRawPtr(), row_indices.getRawPtr());
886 TEST_EQUALITY_CONST(row_values[0], -1.0);
887 TEST_EQUALITY_CONST(row_values[1], 2.0);
888 TEST_EQUALITY_CONST(row_values[2], 1.0);
889 // ToDo: Test column indices!
890 }
891 else {
892 TEST_ASSERT(is_null(epetraOp));
893 }
894
895}
896
897#else
898
899TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_EpetraRowMatrix,
900 Scalar )
901{
902}
903
904TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_RowStatLinearOpBase,
905 Scalar )
906{
907 using Teuchos::as;
908
909 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
910 createTriDiagonalTpetraOperator<Scalar>(g_localDim);
911 out << "tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
912 TEST_ASSERT(nonnull(tpetraOp));
913
914 const RCP<Tpetra::CrsMatrix<Scalar> > tpetraCrsMatrix =
915 Teuchos::rcp_dynamic_cast<Tpetra::CrsMatrix<Scalar> >(tpetraOp,true);
916
917 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
918 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
919 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
920 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
921 const RCP<LinearOpBase<Scalar> > thyraLinearOp =
922 Thyra::tpetraLinearOp(rangeSpace, domainSpace, tpetraOp);
923 TEST_ASSERT(nonnull(thyraLinearOp));
924
926 Teuchos::rcp_dynamic_cast<Thyra::RowStatLinearOpBase<Scalar> >(thyraLinearOp, true);
927
928 // Get the inverse row sums
929
930 const RCP<VectorBase<Scalar> > inv_row_sums =
931 createMember<Scalar>(thyraLinearOp->range());
932 const RCP<VectorBase<Scalar> > row_sums =
933 createMember<Scalar>(thyraLinearOp->range());
934
935 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_INV_ROW_SUM,
936 inv_row_sums.ptr());
937 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
938 row_sums.ptr());
939
940 out << "inv_row_sums = " << *inv_row_sums;
941 out << "row_sums = " << *row_sums;
942
944 Thyra::sum<Scalar>(*row_sums),
945 Teuchos::as<Scalar>(4.0 * thyraLinearOp->domain()->dim() - 2.0),
946 Teuchos::as<Scalar>(10.0 * Teuchos::ScalarTraits<Scalar>::eps())
947 );
948
950 Thyra::sum<Scalar>(*inv_row_sums),
951 Teuchos::as<Scalar>( 1.0 / 4.0 * (thyraLinearOp->domain()->dim() - 2) + 2.0 / 3.0 ),
952 Teuchos::as<Scalar>(10.0 * Teuchos::ScalarTraits<Scalar>::eps())
953 );
954}
955
956TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_ScaledLinearOpBase,
957 Scalar )
958{
959 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
960 createTriDiagonalTpetraOperator<Scalar>(g_localDim);
961 out << "tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
962 TEST_ASSERT(nonnull(tpetraOp));
963
964 const RCP<Tpetra::CrsMatrix<Scalar> > tpetraCrsMatrix =
965 Teuchos::rcp_dynamic_cast<Tpetra::CrsMatrix<Scalar> >(tpetraOp,true);
966
967 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
968 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
969 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
970 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
971 const RCP<LinearOpBase<Scalar> > thyraLinearOp =
972 Thyra::tpetraLinearOp(rangeSpace, domainSpace, tpetraOp);
973 TEST_ASSERT(nonnull(thyraLinearOp));
974
976 Teuchos::rcp_dynamic_cast<Thyra::RowStatLinearOpBase<Scalar> >(thyraLinearOp, true);
977
978 // Get the inverse row sums
979
980 const RCP<VectorBase<Scalar> > inv_row_sums =
981 createMember<Scalar>(thyraLinearOp->range());
982 const RCP<VectorBase<Scalar> > row_sums =
983 createMember<Scalar>(thyraLinearOp->range());
984
985 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_INV_ROW_SUM,
986 inv_row_sums.ptr());
987 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
988 row_sums.ptr());
989
990 out << "inv_row_sums = " << *inv_row_sums;
991 out << "row_sums = " << *row_sums;
992
994 Teuchos::rcp_dynamic_cast<Thyra::ScaledLinearOpBase<Scalar> >(thyraLinearOp, true);
995
996 TEUCHOS_ASSERT(scaledOp->supportsScaleLeft());
997
998 scaledOp->scaleLeft(*inv_row_sums);
999
1000 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
1001 row_sums.ptr());
1002
1003 out << "row_sums after left scaling by inv_row_sum = " << *row_sums;
1004
1005 // scaled row sums should be one for each entry
1007 Scalar(row_sums->space()->dim()),
1008 Thyra::sum<Scalar>(*row_sums),
1009 as<Scalar>(10.0 * Teuchos::ScalarTraits<Scalar>::eps())
1010 );
1011
1012 // Don't currently check the results of right scaling. Tpetra tests
1013 // already check this. Once column sums are supported in tpetra
1014 // adapters, this can be checked easily.
1015 TEUCHOS_ASSERT(scaledOp->supportsScaleRight());
1016 scaledOp->scaleRight(*inv_row_sums);
1017 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,row_sums.ptr());
1018 out << "row_sums after right scaling by inv_row_sum = " << *row_sums;
1019}
1020
1021#endif // HAVE_THYRA_TPETRA_EPETRA
1022
1023
1024//
1025// Unit test instantiations
1026//
1027
1028#define THYRA_TPETRA_THYRA_WRAPPERS_INSTANT(SCALAR) \
1029 \
1030 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1031 convertTpetraToThyraComm, SCALAR ) \
1032 \
1033 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1034 createVectorSpace, SCALAR ) \
1035 \
1036 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1037 createVector, SCALAR ) \
1038 \
1039 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1040 createConstVector, SCALAR ) \
1041 \
1042 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1043 createMultiVector, SCALAR ) \
1044 \
1045 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1046 createConstMultiVector, SCALAR ) \
1047 \
1048 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1049 TeptraVectorSpace, SCALAR ) \
1050 \
1051 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1052 vectorSpaceTester, SCALAR ) \
1053 \
1054 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1055 vectorStdOpsTester, SCALAR ) \
1056 \
1057 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1058 multiVectorStdOpsTester, SCALAR ) \
1059 \
1060 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1061 getTpetraMultiVector, SCALAR ) \
1062 \
1063 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1064 getConstTpetraMultiVector, SCALAR ) \
1065 \
1066 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1067 TpetraLinearOp, SCALAR ) \
1068 \
1069 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1070 createLinearOp, SCALAR ) \
1071 \
1072 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1073 createConstLinearOp, SCALAR ) \
1074 \
1075 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1076 TpetraLinearOp_EpetraRowMatrix, SCALAR ) \
1077 \
1078 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1079 TpetraLinearOp_RowStatLinearOpBase, SCALAR ) \
1080 \
1081 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1082 TpetraLinearOp_ScaledLinearOpBase, SCALAR ) \
1083
1084
1085// We can currently only explicitly instantiate with double support because
1086// Tpetra only supports explicit instantaition with double. As for implicit
1087// instantation, g++ 3.4.6 on my Linux machine was taking more than 30 minutes
1088// to compile this file when all of the types double, float, complex<double>,
1089// and complex<float> where enabled. Therefore, we will only test double for
1090// now until explicit instantation with other types are supported by Tpetra.
1091
1093
1094
1095} // namespace Thyra
#define TEST_EQUALITY_CONST(v1, v2)
#define TEST_THROW(code, ExceptType)
#define TEST_INEQUALITY(v1, v2)
#define TEST_FLOATING_EQUALITY(v1, v2, tol)
#define THYRA_TPETRA_THYRA_WRAPPERS_INSTANT(SCALAR)
#define CHECK_TPETRA_FUNC_CALL_INCREMENT(timerStr, tpetraCode, thyraCode)
void setOption(const char option_true[], const char option_false[], bool *option_val, const char documentation[]=NULL)
static Teuchos::RCP< const Comm< OrdinalType > > getComm()
static RCP< T > lookupCounter(const std::string &name)
Ptr< T > ptr() const
Concrete Thyra::LinearOpBase subclass for Tpetra::Operator.
Traits class that enables the extraction of Tpetra operator/vector objects wrapped in Thyra operator/...
EApplyEpetraOpAs
Determine how the apply an Epetra_Operator as a linear operator.
EAdjointEpetraOp
Determine if adjoints are supported on Epetra_Opeator or not.
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
bool is_null(const boost::shared_ptr< T > &p)
bool nonnull(const boost::shared_ptr< T > &p)
TypeTo as(const TypeFrom &t)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
TEST_ASSERT(castedDep1->getValuesAndValidators().size()==2)
TEST_EQUALITY(rcp_dynamic_cast< const EnhancedNumberValidator< double > >(castedDep1->getValuesAndValidators().find("val1") ->second, true) ->getMax(), double1Vali->getMax())
RCP< const VectorBase< Scalar > > createConstVector(const RCP< const Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVector, const RCP< const VectorSpaceBase< Scalar > > space=Teuchos::null)
RCP< const TpetraMap_t > createTpetraMap(const int localDim)
RCP< const MultiVectorBase< Scalar > > createConstMultiVector(const RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraMultiVector, const RCP< const VectorSpaceBase< Scalar > > rangeSpace=Teuchos::null, const RCP< const VectorSpaceBase< Scalar > > domainSpace=Teuchos::null)
RCP< const VectorSpaceBase< Scalar > > createVectorSpace(const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &tpetraMap)
Create a Thyra::VectorSpaceBase object given a Tpetra::Map.
RCP< VectorBase< Scalar > > createVector(const RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVector, const RCP< const VectorSpaceBase< Scalar > > space=Teuchos::null)
RCP< Tpetra::Operator< Scalar > > createTriDiagonalTpetraOperator(const int numLocalRows)
Teuchos::RCP< Teuchos::Time > lookupAndAssertTimer(const std::string &label)
RCP< MultiVectorBase< Scalar > > createMultiVector(const RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraMultiVector, const RCP< const VectorSpaceBase< Scalar > > rangeSpace=Teuchos::null, const RCP< const VectorSpaceBase< Scalar > > domainSpace=Teuchos::null)
RCP< LinearOpBase< Scalar > > createLinearOp(const RCP< Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraOperator, const RCP< const VectorSpaceBase< Scalar > > rangeSpace=Teuchos::null, const RCP< const VectorSpaceBase< Scalar > > domainSpace=Teuchos::null)
RCP< const LinearOpBase< Scalar > > createConstLinearOp(const RCP< const Tpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraOperator, const RCP< const VectorSpaceBase< Scalar > > rangeSpace=Teuchos::null, const RCP< const VectorSpaceBase< Scalar > > domainSpace=Teuchos::null)
RCP< const Teuchos::Comm< Ordinal > > convertTpetraToThyraComm(const RCP< const Teuchos::Comm< int > > &tpetraComm)
Given an Tpetra Teuchos::Comm<int> object, return an equivalent Teuchos::Comm<Ordinal> object.
RCP< const VectorSpaceBase< Scalar > > createTpetraVectorSpace(const int localDim)
TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL(Comm, ReduceSum, PacketType)