46 #include "Thyra_VectorSpaceTester.hpp" 47 #include "Thyra_VectorStdOpsTester.hpp" 48 #include "Thyra_VectorStdOps.hpp" 49 #include "Thyra_MultiVectorStdOps.hpp" 50 #include "Thyra_LinearOpTester.hpp" 51 #include "Thyra_DefaultProductVector.hpp" 52 #include "Thyra_TestingTools.hpp" 53 #include "Thyra_ScaledLinearOpBase.hpp" 54 #include "Thyra_RowStatLinearOpBase.hpp" 55 #include "Thyra_VectorStdOps.hpp" 56 #include "Tpetra_CrsMatrix.hpp" 57 #include "Teuchos_UnitTestHarness.hpp" 58 #include "Teuchos_DefaultComm.hpp" 59 #include "Teuchos_Tuple.hpp" 74 using Teuchos::ArrayView;
75 using Teuchos::rcp_dynamic_cast;
76 using Teuchos::inOutArg;
87 RCP<const TpetraMap_t>
90 typedef Teuchos::OrdinalTraits<Tpetra::global_size_t> OT;
91 return Teuchos::rcp(
new TpetraMap_t(OT::invalid(), localDim, 0,
92 Teuchos::DefaultComm<int>::getComm()));
99 template<
class Scalar>
100 RCP<const VectorSpaceBase<Scalar> >
107 template<
class Scalar>
108 RCP<Tpetra::Operator<Scalar> >
111 typedef Tpetra::global_size_t global_size_t;
112 typedef Tpetra::Map<>::global_ordinal_type GO;
116 const size_t numMyElements = map->getNodeNumElements();
117 const global_size_t numGlobalElements = map->getGlobalNumElements();
119 ArrayView<const GO> myGlobalElements = map->getNodeElementList();
125 Teuchos::ArrayRCP<size_t> numNz = Teuchos::arcp<size_t>(numMyElements);
130 for (
size_t i=0; i < numMyElements; ++i) {
131 if (myGlobalElements[i] == 0 || static_cast<global_size_t>(myGlobalElements[i]) == numGlobalElements-1) {
141 RCP< Tpetra::CrsMatrix<Scalar> > A =
142 Teuchos::rcp(
new Tpetra::CrsMatrix<Scalar>(map, numNz, Tpetra::StaticProfile) );
145 numNz = Teuchos::null;
149 const Scalar two =
static_cast<Scalar
>( 2.0);
150 const Scalar posOne =
static_cast<Scalar
>(+1.0);
151 const Scalar negOne =
static_cast<Scalar
>(-1.0);
152 for (
size_t i = 0; i < numMyElements; i++) {
153 if (myGlobalElements[i] == 0) {
154 A->insertGlobalValues( myGlobalElements[i],
155 tuple<GO>(myGlobalElements[i], myGlobalElements[i]+1)(),
156 tuple<Scalar> (two, posOne)()
159 else if (static_cast<global_size_t>(myGlobalElements[i]) == numGlobalElements-1) {
160 A->insertGlobalValues( myGlobalElements[i],
161 tuple<GO>(myGlobalElements[i]-1, myGlobalElements[i])(),
162 tuple<Scalar> (negOne, two)()
166 A->insertGlobalValues( myGlobalElements[i],
167 tuple<GO>(myGlobalElements[i]-1, myGlobalElements[i], myGlobalElements[i]+1)(),
168 tuple<Scalar> (negOne, two, posOne)()
188 Teuchos::UnitTestRepository::getCLP().setOption(
189 "show-all-tests",
"no-show-all-tests", &
showAllTests,
"Show all tests or not" );
190 Teuchos::UnitTestRepository::getCLP().setOption(
191 "dump-all",
"no-dump-all", &
dumpAll,
"Dump all objects being tested" );
192 Teuchos::UnitTestRepository::getCLP().setOption(
209 RCP<const Comm<int> > tpetraComm = Teuchos::DefaultComm<int>::getComm();
211 TEST_ASSERT(nonnull(thyraComm));
223 const RCP<const VectorSpaceBase<Scalar> > vs =
224 Thyra::createVectorSpace<Scalar>(tpetraMap);
225 TEST_ASSERT(nonnull(vs));
226 out <<
"vs = " << *vs;
227 const RCP<const SpmdVectorSpaceBase<Scalar> > vs_spmd =
228 rcp_dynamic_cast<
const SpmdVectorSpaceBase<Scalar> >(vs,
true);
229 TEST_EQUALITY(vs_spmd->localSubDim(),
g_localDim);
230 TEST_EQUALITY(vs->dim(), as<Ordinal>(tpetraMap->getGlobalNumElements()));
245 const RCP<const VectorSpaceBase<Scalar> > vs =
246 Thyra::createVectorSpace<Scalar>(tpetraMap);
248 const RCP<Tpetra::Vector<Scalar> > tpetraVector =
249 rcp(
new Tpetra::Vector<Scalar>(tpetraMap));
252 const RCP<VectorBase<Scalar> > thyraVector =
createVector(tpetraVector, vs);
253 TEST_EQUALITY(thyraVector->space(), vs);
254 const RCP<Tpetra::Vector<Scalar> > tpetraVector2 =
255 ConverterT::getTpetraVector(thyraVector);
256 TEST_EQUALITY(tpetraVector2, tpetraVector);
261 TEST_INEQUALITY(thyraVector->space(), vs);
262 TEST_ASSERT(thyraVector->space()->isCompatible(*vs));
263 const RCP<Tpetra::Vector<Scalar> > tpetraVector2 =
264 ConverterT::getTpetraVector(thyraVector);
265 TEST_EQUALITY(tpetraVector2, tpetraVector);
282 const RCP<const VectorSpaceBase<Scalar> > vs =
283 Thyra::createVectorSpace<Scalar>(tpetraMap);
285 const RCP<const Tpetra::Vector<Scalar> > tpetraVector =
286 rcp(
new Tpetra::Vector<Scalar>(tpetraMap));
289 const RCP<const VectorBase<Scalar> > thyraVector =
291 TEST_EQUALITY(thyraVector->space(), vs);
292 const RCP<const Tpetra::Vector<Scalar> > tpetraVector2 =
293 ConverterT::getConstTpetraVector(thyraVector);
294 TEST_EQUALITY(tpetraVector2, tpetraVector);
298 const RCP<const VectorBase<Scalar> > thyraVector =
300 TEST_INEQUALITY(thyraVector->space(), vs);
301 TEST_ASSERT(thyraVector->space()->isCompatible(*vs));
302 const RCP<const Tpetra::Vector<Scalar> > tpetraVector2 =
303 ConverterT::getConstTpetraVector(thyraVector);
304 TEST_EQUALITY(tpetraVector2, tpetraVector);
317 typedef Tpetra::Map<>::local_ordinal_type LO;
318 typedef Tpetra::Map<>::global_ordinal_type GO;
321 const int numCols = 3;
324 const RCP<const VectorSpaceBase<Scalar> > rangeVs =
325 Thyra::createVectorSpace<Scalar>(tpetraMap);
327 const RCP<const TpetraMap_t> tpetraLocRepMap =
328 Tpetra::createLocalMapWithNode<LO, GO>(
329 numCols, tpetraMap->getComm(), tpetraMap->getNode());
330 const RCP<const VectorSpaceBase<Scalar> > domainVs =
331 Thyra::createVectorSpace<Scalar>(tpetraLocRepMap);
333 const RCP<Tpetra::MultiVector<Scalar> > tpetraMultiVector =
334 rcp(
new Tpetra::MultiVector<Scalar>(tpetraMap, numCols));
337 const RCP<MultiVectorBase<Scalar> > thyraMultiVector =
339 TEST_EQUALITY(thyraMultiVector->range(), rangeVs);
340 TEST_EQUALITY(thyraMultiVector->domain(), domainVs);
341 const RCP<Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
342 ConverterT::getTpetraMultiVector(thyraMultiVector);
343 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
347 const RCP<MultiVectorBase<Scalar> > thyraMultiVector =
349 TEST_INEQUALITY(thyraMultiVector->range(), rangeVs);
350 TEST_INEQUALITY(thyraMultiVector->domain(), domainVs);
351 TEST_ASSERT(thyraMultiVector->range()->isCompatible(*rangeVs));
352 TEST_ASSERT(thyraMultiVector->domain()->isCompatible(*domainVs));
353 const RCP<Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
354 ConverterT::getTpetraMultiVector(thyraMultiVector);
355 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
368 typedef Tpetra::Map<>::local_ordinal_type LO;
369 typedef Tpetra::Map<>::global_ordinal_type GO;
372 const int numCols = 3;
375 const RCP<const VectorSpaceBase<Scalar> > rangeVs =
376 Thyra::createVectorSpace<Scalar>(tpetraMap);
378 const RCP<const TpetraMap_t> tpetraLocRepMap =
379 Tpetra::createLocalMapWithNode<LO,GO>(
380 numCols, tpetraMap->getComm(), tpetraMap->getNode());
381 const RCP<const VectorSpaceBase<Scalar> > domainVs =
382 Thyra::createVectorSpace<Scalar>(tpetraLocRepMap);
384 const RCP<const Tpetra::MultiVector<Scalar> > tpetraMultiVector =
385 rcp(
new Tpetra::MultiVector<Scalar>(tpetraMap, numCols));
388 const RCP<const MultiVectorBase<Scalar> > thyraMultiVector =
390 TEST_EQUALITY(thyraMultiVector->range(), rangeVs);
391 TEST_EQUALITY(thyraMultiVector->domain(), domainVs);
392 const RCP<const Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
393 ConverterT::getConstTpetraMultiVector(thyraMultiVector);
394 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
398 const RCP<const MultiVectorBase<Scalar> > thyraMultiVector =
400 TEST_INEQUALITY(thyraMultiVector->range(), rangeVs);
401 TEST_INEQUALITY(thyraMultiVector->domain(), domainVs);
402 TEST_ASSERT(thyraMultiVector->range()->isCompatible(*rangeVs));
403 TEST_ASSERT(thyraMultiVector->domain()->isCompatible(*domainVs));
404 const RCP<const Tpetra::MultiVector<Scalar> > tpetraMultiVector2 =
405 ConverterT::getConstTpetraMultiVector(thyraMultiVector);
406 TEST_EQUALITY(tpetraMultiVector2, tpetraMultiVector);
420 const RCP<const VectorSpaceBase<Scalar> > vs =
422 const RCP<VectorBase<Scalar> > v = createMember(vs);
423 TEST_ASSERT(nonnull(v));
424 TEST_EQUALITY(v->space(), vs);
435 const RCP<const VectorSpaceBase<Scalar> > vs
436 = createTpetraVectorSpace<Scalar>(
g_localDim);
437 Thyra::VectorSpaceTester<Scalar> vectorSpaceTester;
439 vectorSpaceTester.dump_all(
dumpAll);
440 TEST_ASSERT(vectorSpaceTester.check(*vs, &out));
474 const int numCols = 3;
475 const RCP<const VectorSpaceBase<Scalar> > vs
476 = createTpetraVectorSpace<Scalar>(
g_localDim);
479 const RCP<MultiVectorBase<Scalar> > mv = createMembers(vs, numCols);
480 const RCP<Tpetra::MultiVector<Scalar> > tmv =
481 ConverterT::getTpetraMultiVector(mv);
482 TEST_ASSERT(nonnull(tmv));
483 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
487 const RCP<VectorBase<Scalar> > v = createMember(vs);
488 const RCP<Tpetra::MultiVector<Scalar> > tmv =
489 ConverterT::getTpetraMultiVector(v);
490 TEST_ASSERT(nonnull(tmv));
491 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
495 const RCP<VectorBase<Scalar> > pv = Thyra::defaultProductVector<Scalar>();
496 TEST_THROW(ConverterT::getTpetraMultiVector(pv), std::logic_error);
511 const int numCols = 3;
512 const RCP<const VectorSpaceBase<Scalar> > vs
513 = createTpetraVectorSpace<Scalar>(
g_localDim);
516 const RCP<const MultiVectorBase<Scalar> > mv = createMembers(vs, numCols);
517 const RCP<const Tpetra::MultiVector<Scalar> > tmv =
518 ConverterT::getConstTpetraMultiVector(mv);
519 TEST_ASSERT(nonnull(tmv));
520 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
524 const RCP<const VectorBase<Scalar> > v = createMember(vs);
525 const RCP<const Tpetra::MultiVector<Scalar> > tmv =
526 ConverterT::getConstTpetraMultiVector(v);
527 TEST_ASSERT(nonnull(tmv));
528 TEST_EQUALITY(as<Ordinal>(tmv->getMap()->getGlobalNumElements()), vs->dim());
532 const RCP<const VectorBase<Scalar> > pv = Thyra::defaultProductVector<Scalar>();
533 TEST_THROW(ConverterT::getConstTpetraMultiVector(pv), std::logic_error);
547 typedef Teuchos::ScalarTraits<Scalar> ST;
550 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
551 createTriDiagonalTpetraOperator<Scalar>(
g_localDim);
552 out <<
"tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
553 TEST_ASSERT(nonnull(tpetraOp));
555 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
556 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
557 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
558 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
559 const RCP<const LinearOpBase<Scalar> > thyraLinearOp =
560 Thyra::tpetraLinearOp(rangeSpace, domainSpace, tpetraOp);
561 TEST_ASSERT(nonnull(thyraLinearOp));
563 out <<
"\nCheck that operator returns the right thing ...\n";
564 const RCP<VectorBase<Scalar> > x = createMember(thyraLinearOp->domain());
565 Thyra::V_S(x.ptr(), ST::one());
566 const RCP<VectorBase<Scalar> > y = createMember(thyraLinearOp->range());
567 Thyra::apply<Scalar>(*thyraLinearOp, Thyra::NOTRANS, *x, y.ptr());
568 const Scalar sum_y = sum(*y);
569 TEST_FLOATING_EQUALITY( sum_y, as<Scalar>(3+1+2*(y->space()->dim()-2)),
572 out <<
"\nCheck the general LinearOp interface ...\n";
573 Thyra::LinearOpTester<Scalar> linearOpTester;
575 linearOpTester.dump_all(
dumpAll);
577 TEST_ASSERT(linearOpTester.check(*thyraLinearOp, Teuchos::inOutArg(out)));
593 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
594 createTriDiagonalTpetraOperator<Scalar>(
g_localDim);
595 out <<
"tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
597 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
598 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
600 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
601 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
604 const RCP<LinearOpBase<Scalar> > thyraOp =
606 TEST_EQUALITY(thyraOp->range(), rangeSpace);
607 TEST_EQUALITY(thyraOp->domain(), domainSpace);
608 const RCP<Tpetra::Operator<Scalar> > tpetraOp2 =
609 ConverterT::getTpetraOperator(thyraOp);
610 TEST_EQUALITY(tpetraOp2, tpetraOp);
614 const RCP<LinearOpBase<Scalar> > thyraOp =
616 TEST_INEQUALITY(thyraOp->range(), rangeSpace);
617 TEST_INEQUALITY(thyraOp->domain(), domainSpace);
618 TEST_ASSERT(thyraOp->range()->isCompatible(*rangeSpace));
619 TEST_ASSERT(thyraOp->domain()->isCompatible(*domainSpace));
620 const RCP<Tpetra::Operator<Scalar> > tpetraOp2 =
621 ConverterT::getTpetraOperator(thyraOp);
622 TEST_EQUALITY(tpetraOp2, tpetraOp);
638 const RCP<const Tpetra::Operator<Scalar> > tpetraOp =
639 createTriDiagonalTpetraOperator<Scalar>(
g_localDim);
640 out <<
"tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
642 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
643 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
645 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
646 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
649 const RCP<const LinearOpBase<Scalar> > thyraOp =
651 TEST_EQUALITY(thyraOp->range(), rangeSpace);
652 TEST_EQUALITY(thyraOp->domain(), domainSpace);
653 const RCP<const Tpetra::Operator<Scalar> > tpetraOp2 =
654 ConverterT::getConstTpetraOperator(thyraOp);
655 TEST_EQUALITY(tpetraOp2, tpetraOp);
659 const RCP<const LinearOpBase<Scalar> > thyraOp =
661 TEST_INEQUALITY(thyraOp->range(), rangeSpace);
662 TEST_INEQUALITY(thyraOp->domain(), domainSpace);
663 TEST_ASSERT(thyraOp->range()->isCompatible(*rangeSpace));
664 TEST_ASSERT(thyraOp->domain()->isCompatible(*domainSpace));
665 const RCP<const Tpetra::Operator<Scalar> > tpetraOp2 =
666 ConverterT::getConstTpetraOperator(thyraOp);
667 TEST_EQUALITY(tpetraOp2, tpetraOp);
680 Teuchos::RCP<Teuchos::Time> timer = Teuchos::TimeMonitor::lookupCounter(label);
681 TEUCHOS_TEST_FOR_EXCEPTION(timer == null,
683 "lookupAndAssertTimer(): timer \"" << label <<
"\" was not present in Teuchos::TimeMonitor." 684 " Unit test not valid.");
689 #define CHECK_TPETRA_FUNC_CALL_INCREMENT( timerStr, tpetraCode, thyraCode ) \ 691 out << "\nTesting that Thyra calls down to " << timerStr << "\n"; \ 693 const RCP<const Time> timer = lookupAndAssertTimer(timerStr); \ 694 const int countBefore = timer->numCalls(); \ 696 const int countAfter = timer->numCalls(); \ 697 TEST_EQUALITY( countAfter, countBefore+1 ); \ 705 typedef Teuchos::ScalarTraits<Scalar> ST;
706 typedef typename ST::magnitudeType Magnitude;
707 typedef VectorSpaceBase<Scalar> VectorSpace;
708 typedef MultiVectorBase<Scalar> MultiVec;
710 typedef Tpetra::MultiVector<Scalar> TpetraMultiVec;
713 const int numCols = 3;
715 const RCP<const VectorSpace> vs =
718 A = createMembers(vs, numCols),
719 B = createMembers(vs, numCols);
720 const RCP<TpetraMultiVec>
721 tA = TOVE::getTpetraMultiVector(A),
722 tB = TOVE::getTpetraMultiVector(B);
723 Array<Scalar> C(numCols*numCols,ST::zero());
725 Teuchos::Array<Magnitude> avMag(numCols);
726 Teuchos::Array<Scalar> avScal(numCols);
729 "Tpetra::MultiVector::putScalar()",
730 tA->putScalar(ST::zero()),
731 Thyra::assign(A.ptr(), ST::zero())
735 "Tpetra::MultiVector::dot()",
736 tA->dot(*tB, avScal() ),
737 Thyra::norms( *A, avMag() )
741 "Tpetra::MultiVector::dot()",
742 tA->dot(*tB, avScal() ),
743 A->range()->scalarProds(*A, *B, avScal() )
809 #ifdef TPETRA_TEUCHOS_TIME_MONITOR 810 # define TPETRA_TIMER_TESTS(SCALAR) \ 811 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, UseTpetraImplementations, SCALAR ) 813 # define TPETRA_TIMER_TESTS(SCALAR) 822 #ifdef HAVE_THYRA_TPETRA_EPETRA 830 using Teuchos::outArg;
831 using Teuchos::rcp_dynamic_cast;
832 using Teuchos::Array;
833 typedef Teuchos::ScalarTraits<Scalar> ST;
835 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
836 createTriDiagonalTpetraOperator<Scalar>(
g_localDim);
838 const RCP<LinearOpBase<Scalar> > thyraOp =
841 const RCP<Thyra::TpetraLinearOp<Scalar> > thyraTpetraOp =
844 RCP<const Epetra_Operator> epetraOp;
845 Thyra::EOpTransp epetraOpTransp;
849 thyraTpetraOp->getEpetraOpView( outArg(epetraOp), outArg(epetraOpTransp),
850 outArg(epetraOpApplyAs), outArg(epetraOpAdjointSupport) );
852 if (
typeid(Scalar) ==
typeid(
double)) {
853 TEST_ASSERT(nonnull(epetraOp));
854 const RCP<const Epetra_RowMatrix> epetraRowMatrix =
855 rcp_dynamic_cast<
const Epetra_RowMatrix>(epetraOp,
true);
856 int numRowEntries = -1;
857 epetraRowMatrix->NumMyRowEntries(1, numRowEntries);
858 TEST_EQUALITY_CONST(numRowEntries, 3);
859 Array<double> row_values(numRowEntries);
860 Array<int> row_indices(numRowEntries);
861 epetraRowMatrix->ExtractMyRowCopy(1, numRowEntries, numRowEntries,
862 row_values.getRawPtr(), row_indices.getRawPtr());
863 TEST_EQUALITY_CONST(row_values[0], -1.0);
864 TEST_EQUALITY_CONST(row_values[1], 2.0);
865 TEST_EQUALITY_CONST(row_values[2], 1.0);
869 TEST_ASSERT(is_null(epetraOp));
886 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
887 createTriDiagonalTpetraOperator<Scalar>(
g_localDim);
888 out <<
"tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
889 TEST_ASSERT(nonnull(tpetraOp));
891 const RCP<Tpetra::CrsMatrix<Scalar> > tpetraCrsMatrix =
892 Teuchos::rcp_dynamic_cast<Tpetra::CrsMatrix<Scalar> >(tpetraOp,
true);
894 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
895 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
896 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
897 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
898 const RCP<LinearOpBase<Scalar> > thyraLinearOp =
899 Thyra::tpetraLinearOp(rangeSpace, domainSpace, tpetraOp);
900 TEST_ASSERT(nonnull(thyraLinearOp));
902 const Teuchos::RCP<Thyra::RowStatLinearOpBase<Scalar> > rowStatOp =
903 Teuchos::rcp_dynamic_cast<Thyra::RowStatLinearOpBase<Scalar> >(thyraLinearOp,
true);
907 const RCP<VectorBase<Scalar> > inv_row_sums =
908 createMember<Scalar>(thyraLinearOp->range());
909 const RCP<VectorBase<Scalar> > row_sums =
910 createMember<Scalar>(thyraLinearOp->range());
912 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_INV_ROW_SUM,
914 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
917 out <<
"inv_row_sums = " << *inv_row_sums;
918 out <<
"row_sums = " << *row_sums;
920 TEST_FLOATING_EQUALITY(
921 Thyra::sum<Scalar>(*row_sums),
922 Teuchos::as<Scalar>(4.0 * thyraLinearOp->domain()->dim() - 2.0),
923 Teuchos::as<Scalar>(10.0 * Teuchos::ScalarTraits<Scalar>::eps())
926 TEST_FLOATING_EQUALITY(
927 Thyra::sum<Scalar>(*inv_row_sums),
928 Teuchos::as<Scalar>( 1.0 / 4.0 * (thyraLinearOp->domain()->dim() - 2) + 2.0 / 3.0 ),
929 Teuchos::as<Scalar>(10.0 * Teuchos::ScalarTraits<Scalar>::eps())
936 const RCP<Tpetra::Operator<Scalar> > tpetraOp =
937 createTriDiagonalTpetraOperator<Scalar>(
g_localDim);
938 out <<
"tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
939 TEST_ASSERT(nonnull(tpetraOp));
941 const RCP<Tpetra::CrsMatrix<Scalar> > tpetraCrsMatrix =
942 Teuchos::rcp_dynamic_cast<Tpetra::CrsMatrix<Scalar> >(tpetraOp,
true);
944 const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
945 Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
946 const RCP<const VectorSpaceBase<Scalar> > domainSpace =
947 Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
948 const RCP<LinearOpBase<Scalar> > thyraLinearOp =
949 Thyra::tpetraLinearOp(rangeSpace, domainSpace, tpetraOp);
950 TEST_ASSERT(nonnull(thyraLinearOp));
952 const Teuchos::RCP<Thyra::RowStatLinearOpBase<Scalar> > rowStatOp =
953 Teuchos::rcp_dynamic_cast<Thyra::RowStatLinearOpBase<Scalar> >(thyraLinearOp,
true);
957 const RCP<VectorBase<Scalar> > inv_row_sums =
958 createMember<Scalar>(thyraLinearOp->range());
959 const RCP<VectorBase<Scalar> > row_sums =
960 createMember<Scalar>(thyraLinearOp->range());
962 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_INV_ROW_SUM,
964 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
967 out <<
"inv_row_sums = " << *inv_row_sums;
968 out <<
"row_sums = " << *row_sums;
970 const Teuchos::RCP<Thyra::ScaledLinearOpBase<Scalar> > scaledOp =
971 Teuchos::rcp_dynamic_cast<Thyra::ScaledLinearOpBase<Scalar> >(thyraLinearOp,
true);
973 TEUCHOS_ASSERT(scaledOp->supportsScaleLeft());
975 scaledOp->scaleLeft(*inv_row_sums);
977 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
980 out <<
"row_sums after left scaling by inv_row_sum = " << *row_sums;
983 TEST_FLOATING_EQUALITY(
984 Scalar(row_sums->space()->dim()),
985 Thyra::sum<Scalar>(*row_sums),
986 as<Scalar>(10.0 * Teuchos::ScalarTraits<Scalar>::eps())
992 TEUCHOS_ASSERT(scaledOp->supportsScaleRight());
993 scaledOp->scaleRight(*inv_row_sums);
994 rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,row_sums.ptr());
995 out <<
"row_sums after right scaling by inv_row_sum = " << *row_sums;
998 #endif // HAVE_THYRA_TPETRA_EPETRA 1005 #define THYRA_TPETRA_THYRA_WRAPPERS_INSTANT(SCALAR) \ 1007 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1008 convertTpetraToThyraComm, SCALAR ) \ 1010 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1011 createVectorSpace, SCALAR ) \ 1013 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1014 createVector, SCALAR ) \ 1016 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1017 createConstVector, SCALAR ) \ 1019 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1020 createMultiVector, SCALAR ) \ 1022 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1023 createConstMultiVector, SCALAR ) \ 1025 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1026 TeptraVectorSpace, SCALAR ) \ 1028 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1029 vectorSpaceTester, SCALAR ) \ 1031 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1032 getTpetraMultiVector, SCALAR ) \ 1034 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1035 getConstTpetraMultiVector, SCALAR ) \ 1037 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1038 TpetraLinearOp, SCALAR ) \ 1040 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1041 createLinearOp, SCALAR ) \ 1043 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1044 createConstLinearOp, SCALAR ) \ 1046 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1047 TpetraLinearOp_EpetraRowMatrix, SCALAR ) \ 1049 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1050 TpetraLinearOp_RowStatLinearOpBase, SCALAR ) \ 1052 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \ 1053 TpetraLinearOp_ScaledLinearOpBase, SCALAR ) \ 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)
#define THYRA_TPETRA_THYRA_WRAPPERS_INSTANT(SCALAR)
RCP< Tpetra::Operator< Scalar > > createTriDiagonalTpetraOperator(const int numLocalRows)
#define CHECK_TPETRA_FUNC_CALL_INCREMENT(timerStr, tpetraCode, thyraCode)
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 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< 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< VectorBase< Scalar > > createVector(const RCP< Tpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &tpetraVector, const RCP< const VectorSpaceBase< Scalar > > space=Teuchos::null)
EAdjointEpetraOp
Determine if adjoints are supported on Epetra_Opeator or not.
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)
Teuchos::RCP< Teuchos::Time > lookupAndAssertTimer(const std::string &label)
TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL(TpetraThyraWrappers, convertTpetraToThyraComm, Scalar)
EApplyEpetraOpAs
Determine how the apply an Epetra_Operator as a linear operator.
RCP< const VectorSpaceBase< Scalar > > createVectorSpace(const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > &tpetraMap)
Create a Thyra::VectorSpaceBase object given a Tpetra::Map.
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...
Concrete Thyra::LinearOpBase subclass for Tpetra::Operator.
RCP< const TpetraMap_t > createTpetraMap(const int localDim)
RCP< const VectorSpaceBase< Scalar > > createTpetraVectorSpace(const int localDim)