Thyra Package Browser (Single Doxygen Collection)  Version of the Day
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_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"
60 
61 
62 namespace Thyra {
63 
64 
65 //
66 // Helper code and declarations
67 //
68 
69 
70 using Teuchos::as;
71 using Teuchos::null;
72 using Teuchos::RCP;
73 using Teuchos::rcp;
74 using Teuchos::ArrayView;
75 using Teuchos::rcp_dynamic_cast;
76 using Teuchos::inOutArg;
77 using Teuchos::Comm;
78 using Teuchos::tuple;
79 
80 
81 const int g_localDim = 4; // ToDo: Make variable!
82 
83 
84 typedef Tpetra::Map<> TpetraMap_t;
85 
86 
87 RCP<const TpetraMap_t>
88 createTpetraMap(const int localDim)
89 {
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()));
93  // ToDo: Pass in the comm?
94 }
95 
96 // ToDo: Above: Vary the LocalOrdinal and GlobalOrdinal types?
97 
98 
99 template<class Scalar>
100 RCP<const VectorSpaceBase<Scalar> >
101 createTpetraVectorSpace(const int localDim)
102 {
103  return Thyra::createVectorSpace<Scalar>(createTpetraMap(g_localDim));
104 }
105 
106 
107 template<class Scalar>
108 RCP<Tpetra::Operator<Scalar> >
109 createTriDiagonalTpetraOperator(const int numLocalRows)
110 {
111  typedef Tpetra::global_size_t global_size_t;
112  typedef Tpetra::Map<>::global_ordinal_type GO;
113 
114  RCP<const Tpetra::Map<> > map = createTpetraMap(numLocalRows);
115 
116  const size_t numMyElements = map->getNodeNumElements();
117  const global_size_t numGlobalElements = map->getGlobalNumElements();
118 
119  ArrayView<const GO> myGlobalElements = map->getNodeElementList();
120 
121  // Create an OTeger vector numNz that is used to build the Petra Matrix.
122  // numNz[i] is the Number of OFF-DIAGONAL term for the ith global equation
123  // on this processor
124 
125  Teuchos::ArrayRCP<size_t> numNz = Teuchos::arcp<size_t>(numMyElements);
126 
127  // We are building a tridiagonal matrix where each row has (-1 2 -1)
128  // So we need 2 off-diagonal terms (except for the first and last equation)
129 
130  for (size_t i=0; i < numMyElements; ++i) {
131  if (myGlobalElements[i] == 0 || static_cast<global_size_t>(myGlobalElements[i]) == numGlobalElements-1) {
132  // boundary
133  numNz[i] = 2;
134  }
135  else {
136  numNz[i] = 3;
137  }
138  }
139 
140  // Create a Tpetra::Matrix using the Map, with a static allocation dictated by numNz
141  RCP< Tpetra::CrsMatrix<Scalar> > A =
142  Teuchos::rcp( new Tpetra::CrsMatrix<Scalar>(map, numNz, Tpetra::StaticProfile) );
143 
144  // We are done with NumNZ
145  numNz = Teuchos::null;
146 
147  // Add rows one-at-a-time
148  // Off diagonal values will always be -1
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)()
157  );
158  }
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)()
163  );
164  }
165  else {
166  A->insertGlobalValues( myGlobalElements[i],
167  tuple<GO>(myGlobalElements[i]-1, myGlobalElements[i], myGlobalElements[i]+1)(),
168  tuple<Scalar> (negOne, two, posOne)()
169  );
170  }
171  }
172 
173  // Finish up
174  A->fillComplete();
175 
176  return A;
177 
178 }
179 
180 
181 bool showAllTests = false;
182 bool dumpAll = false;
183 bool runLinearOpTester = true;
184 
185 
187 {
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(
193  "run-linear-op-tester", "no-run-linear-op-tester", &runLinearOpTester, "..." );
194 }
195 
196 
197 //
198 // Unit Tests
199 //
200 
201 
202 //
203 // convertTpetraToThyraComm
204 //
205 
207  Scalar )
208 {
209  RCP<const Comm<int> > tpetraComm = Teuchos::DefaultComm<int>::getComm();
210  RCP<const Comm<Ordinal> > thyraComm = Thyra::convertTpetraToThyraComm(tpetraComm);
211  TEST_ASSERT(nonnull(thyraComm));
212 }
213 
214 
215 //
216 // createVectorSpace
217 //
218 
220  Scalar )
221 {
222  const RCP<const TpetraMap_t> tpetraMap = createTpetraMap(g_localDim);
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()));
231 }
232 
233 
234 //
235 // createVector
236 //
237 
239  Scalar )
240 {
241 
243 
244  const RCP<const TpetraMap_t> tpetraMap = createTpetraMap(g_localDim);
245  const RCP<const VectorSpaceBase<Scalar> > vs =
246  Thyra::createVectorSpace<Scalar>(tpetraMap);
247 
248  const RCP<Tpetra::Vector<Scalar> > tpetraVector =
249  rcp(new Tpetra::Vector<Scalar>(tpetraMap));
250 
251  {
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);
257  }
258 
259  {
260  const RCP<VectorBase<Scalar> > thyraVector = Thyra::createVector(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);
266  }
267 
268 }
269 
270 
271 //
272 // createConstVector
273 //
274 
276  Scalar )
277 {
278 
280 
281  const RCP<const TpetraMap_t> tpetraMap = createTpetraMap(g_localDim);
282  const RCP<const VectorSpaceBase<Scalar> > vs =
283  Thyra::createVectorSpace<Scalar>(tpetraMap);
284 
285  const RCP<const Tpetra::Vector<Scalar> > tpetraVector =
286  rcp(new Tpetra::Vector<Scalar>(tpetraMap));
287 
288  {
289  const RCP<const VectorBase<Scalar> > thyraVector =
290  createConstVector(tpetraVector, vs);
291  TEST_EQUALITY(thyraVector->space(), vs);
292  const RCP<const Tpetra::Vector<Scalar> > tpetraVector2 =
293  ConverterT::getConstTpetraVector(thyraVector);
294  TEST_EQUALITY(tpetraVector2, tpetraVector);
295  }
296 
297  {
298  const RCP<const VectorBase<Scalar> > thyraVector =
299  Thyra::createConstVector(tpetraVector);
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);
305  }
306 
307 }
308 
309 
310 //
311 // createMultiVector
312 //
313 
315  Scalar )
316 {
317  typedef Tpetra::Map<>::local_ordinal_type LO;
318  typedef Tpetra::Map<>::global_ordinal_type GO;
320 
321  const int numCols = 3;
322 
323  const RCP<const TpetraMap_t> tpetraMap = createTpetraMap(g_localDim);
324  const RCP<const VectorSpaceBase<Scalar> > rangeVs =
325  Thyra::createVectorSpace<Scalar>(tpetraMap);
326 
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);
332 
333  const RCP<Tpetra::MultiVector<Scalar> > tpetraMultiVector =
334  rcp(new Tpetra::MultiVector<Scalar>(tpetraMap, numCols));
335 
336  {
337  const RCP<MultiVectorBase<Scalar> > thyraMultiVector =
338  createMultiVector(tpetraMultiVector, rangeVs, domainVs);
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);
344  }
345 
346  {
347  const RCP<MultiVectorBase<Scalar> > thyraMultiVector =
348  Thyra::createMultiVector(tpetraMultiVector);
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);
356  }
357 
358 }
359 
360 
361 //
362 // createConstMultiVector
363 //
364 
366  Scalar )
367 {
368  typedef Tpetra::Map<>::local_ordinal_type LO;
369  typedef Tpetra::Map<>::global_ordinal_type GO;
371 
372  const int numCols = 3;
373 
374  const RCP<const TpetraMap_t> tpetraMap = createTpetraMap(g_localDim);
375  const RCP<const VectorSpaceBase<Scalar> > rangeVs =
376  Thyra::createVectorSpace<Scalar>(tpetraMap);
377 
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);
383 
384  const RCP<const Tpetra::MultiVector<Scalar> > tpetraMultiVector =
385  rcp(new Tpetra::MultiVector<Scalar>(tpetraMap, numCols));
386 
387  {
388  const RCP<const MultiVectorBase<Scalar> > thyraMultiVector =
389  createConstMultiVector(tpetraMultiVector, rangeVs, domainVs);
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);
395  }
396 
397  {
398  const RCP<const MultiVectorBase<Scalar> > thyraMultiVector =
399  Thyra::createConstMultiVector(tpetraMultiVector);
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);
407  }
408 
409 }
410 
411 
412 //
413 // TeptraVectorSpace
414 //
415 
416 
417 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TeptraVectorSpace,
418  Scalar )
419 {
420  const RCP<const VectorSpaceBase<Scalar> > vs =
421  Thyra::createVectorSpace<Scalar>(createTpetraMap(g_localDim));
422  const RCP<VectorBase<Scalar> > v = createMember(vs);
423  TEST_ASSERT(nonnull(v));
424  TEST_EQUALITY(v->space(), vs);
425 }
426 
427 
428 //
429 // vectorSpaceTester
430 //
431 
432 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, vectorSpaceTester,
433  Scalar )
434 {
435  const RCP<const VectorSpaceBase<Scalar> > vs
436  = createTpetraVectorSpace<Scalar>(g_localDim);
437  Thyra::VectorSpaceTester<Scalar> vectorSpaceTester;
438  vectorSpaceTester.show_all_tests(showAllTests);
439  vectorSpaceTester.dump_all(dumpAll);
440  TEST_ASSERT(vectorSpaceTester.check(*vs, &out));
441 }
442 
443 
444 // ToDo: Fix the default tolerances for below
445 
446 /*
447 
448 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraVectorSpace, vectorStdOpsTester,
449  Scalar )
450 {
451  const RCP<const VectorSpaceBase<Scalar> > vs =
452  Thyra::createVectorSpace<Scalar>(createTpetraMap(g_localDim));
453  Thyra::VectorStdOpsTester<Scalar> vectorStdOpsTester;
454  //vectorStdOpsTester.show_all_tests(showAllTests);
455  //vectorStdOpsTester.dump_all(dumpAll);
456  TEST_ASSERT(vectorStdOpsTester.checkStdOps(*vs, &out));
457 }
458 
459 TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT_SCALAR_TYPES( TpetraVectorSpace,
460  vectorStdOpsTester )
461 
462 */
463 
464 
465 //
466 // getTpetraMultiVector
467 //
468 
469 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, getTpetraMultiVector,
470  Scalar )
471 {
473 
474  const int numCols = 3;
475  const RCP<const VectorSpaceBase<Scalar> > vs
476  = createTpetraVectorSpace<Scalar>(g_localDim);
477 
478  {
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());
484  }
485 
486  {
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());
492  }
493 
494 #ifdef THYRA_DEBUG
495  const RCP<VectorBase<Scalar> > pv = Thyra::defaultProductVector<Scalar>();
496  TEST_THROW(ConverterT::getTpetraMultiVector(pv), std::logic_error);
497 #endif
498 
499 }
500 
501 
502 //
503 // getConstTpetraMultiVector
504 //
505 
506 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, getConstTpetraMultiVector,
507  Scalar )
508 {
510 
511  const int numCols = 3;
512  const RCP<const VectorSpaceBase<Scalar> > vs
513  = createTpetraVectorSpace<Scalar>(g_localDim);
514 
515  {
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());
521  }
522 
523  {
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());
529  }
530 
531 #ifdef THYRA_DEBUG
532  const RCP<const VectorBase<Scalar> > pv = Thyra::defaultProductVector<Scalar>();
533  TEST_THROW(ConverterT::getConstTpetraMultiVector(pv), std::logic_error);
534 #endif
535 
536 }
537 
538 
539 //
540 // TpetraLinearOp
541 //
542 
544  Scalar )
545 {
546 
547  typedef Teuchos::ScalarTraits<Scalar> ST;
548  using Teuchos::as;
549 
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));
554 
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));
562 
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)),
570  100.0 * ST::eps() );
571 
572  out << "\nCheck the general LinearOp interface ...\n";
573  Thyra::LinearOpTester<Scalar> linearOpTester;
574  linearOpTester.show_all_tests(showAllTests);
575  linearOpTester.dump_all(dumpAll);
576  if (runLinearOpTester) {
577  TEST_ASSERT(linearOpTester.check(*thyraLinearOp, Teuchos::inOutArg(out)));
578  }
579 
580 }
581 
582 
583 //
584 // createLinearOp
585 //
586 
588  Scalar )
589 {
590 
592 
593  const RCP<Tpetra::Operator<Scalar> > tpetraOp =
594  createTriDiagonalTpetraOperator<Scalar>(g_localDim);
595  out << "tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
596 
597  const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
598  Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
599 
600  const RCP<const VectorSpaceBase<Scalar> > domainSpace =
601  Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
602 
603  {
604  const RCP<LinearOpBase<Scalar> > thyraOp =
605  createLinearOp(tpetraOp, rangeSpace, domainSpace);
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);
611  }
612 
613  {
614  const RCP<LinearOpBase<Scalar> > thyraOp =
615  Thyra::createLinearOp(tpetraOp);
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);
623  }
624 
625 }
626 
627 
628 //
629 // createConstLinearOp
630 //
631 
633  Scalar )
634 {
635 
637 
638  const RCP<const Tpetra::Operator<Scalar> > tpetraOp =
639  createTriDiagonalTpetraOperator<Scalar>(g_localDim);
640  out << "tpetraOp = " << Teuchos::describe(*tpetraOp, Teuchos::VERB_HIGH) << std::endl;
641 
642  const RCP<const VectorSpaceBase<Scalar> > rangeSpace =
643  Thyra::createVectorSpace<Scalar>(tpetraOp->getRangeMap());
644 
645  const RCP<const VectorSpaceBase<Scalar> > domainSpace =
646  Thyra::createVectorSpace<Scalar>(tpetraOp->getDomainMap());
647 
648  {
649  const RCP<const LinearOpBase<Scalar> > thyraOp =
650  createConstLinearOp(tpetraOp, rangeSpace, domainSpace);
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);
656  }
657 
658  {
659  const RCP<const LinearOpBase<Scalar> > thyraOp =
660  Thyra::createConstLinearOp(tpetraOp);
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);
668  }
669 
670 }
671 
672 
673 //
674 // Tpetra-implemented methods
675 //
676 
677 
678 Teuchos::RCP<Teuchos::Time> lookupAndAssertTimer(const std::string &label)
679 {
680  Teuchos::RCP<Teuchos::Time> timer = Teuchos::TimeMonitor::lookupCounter(label);
681  TEUCHOS_TEST_FOR_EXCEPTION(timer == null,
682  std::runtime_error,
683  "lookupAndAssertTimer(): timer \"" << label << "\" was not present in Teuchos::TimeMonitor."
684  " Unit test not valid.");
685  return timer;
686 }
687 
688 
689 #define CHECK_TPETRA_FUNC_CALL_INCREMENT( timerStr, tpetraCode, thyraCode ) \
690 { \
691  out << "\nTesting that Thyra calls down to " << timerStr << "\n"; \
692  ECHO(tpetraCode); \
693  const RCP<const Time> timer = lookupAndAssertTimer(timerStr); \
694  const int countBefore = timer->numCalls(); \
695  ECHO(thyraCode); \
696  const int countAfter = timer->numCalls(); \
697  TEST_EQUALITY( countAfter, countBefore+1 ); \
698 }
699 
700 
701 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, UseTpetraImplementations,
702  Scalar )
703 {
704  using Teuchos::Time;
705  typedef Teuchos::ScalarTraits<Scalar> ST;
706  typedef typename ST::magnitudeType Magnitude;
707  typedef VectorSpaceBase<Scalar> VectorSpace;
708  typedef MultiVectorBase<Scalar> MultiVec;
709  //typedef Tpetra::Map<int> TpetraMap;
710  typedef Tpetra::MultiVector<Scalar> TpetraMultiVec;
712 
713  const int numCols = 3;
714 
715  const RCP<const VectorSpace> vs =
716  createTpetraVectorSpace<Scalar>(g_localDim);
717  const RCP<MultiVec>
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());
724 
725  Teuchos::Array<Magnitude> avMag(numCols);
726  Teuchos::Array<Scalar> avScal(numCols);
727 
729  "Tpetra::MultiVector::putScalar()",
730  tA->putScalar(ST::zero()),
731  Thyra::assign(A.ptr(), ST::zero())
732  );
733 
735  "Tpetra::MultiVector::dot()",
736  tA->dot(*tB, avScal() ), // norms calls scalarProd, which calls Tpetra::dot
737  Thyra::norms( *A, avMag() )
738  );
739 
741  "Tpetra::MultiVector::dot()",
742  tA->dot(*tB, avScal() ),
743  A->range()->scalarProds(*A, *B, avScal() )
744  );
745 
746  /*
747  CHECK_TPETRA_FUNC_CALL_INCREMENT(
748  "Tpetra::MultiVector::scale(alpha)",
749  tB->scale( ST::zero() ),
750  Thyra::scale( ST::zero(), B.ptr() )
751  );
752  */
753 
754  /*
755  CHECK_TPETRA_FUNC_CALL_INCREMENT(
756  "Tpetra::MultiVector::operator=()",
757  (*tA) = *tB,
758  Thyra::assign( A.ptr(), *B )
759  );
760  */
761 
762  /*
763  {
764  RCP<MultiVec> Ctmvb = Thyra::createMembersView(
765  A->domain(),
766  RTOpPack::SubMultiVectorView<Scalar>(
767  0, numCols, 0, numCols,
768  Teuchos::arcpFromArrayView(C()), numCols
769  )
770  );
771  CHECK_TPETRA_FUNC_CALL_INCREMENT(
772  "Tpetra::multiplyOntoHost()",
773  Tpetra::multiplyOntoHost(Teuchos::CONJ_TRANS,Teuchos::NO_TRANS,ST::one(),*tA,*tB, C()),
774  A->apply(Thyra::CONJTRANS,*B,Ctmvb.ptr(),ST::one(),ST::zero())
775  );
776  }
777  */
778 
779  /*
780  {
781  RCP<const MultiVec> Ctmvb = Thyra::createMembersView(
782  A->domain(),
783  RTOpPack::ConstSubMultiVectorView<Scalar>(
784  0, numCols, 0, numCols,
785  Teuchos::arcpFromArrayView(C()), numCols
786  )
787  );
788  const RCP<const TpetraMultiVec>
789  tC = TOVE::getConstTpetraMultiVector(Ctmvb);
790  CHECK_TPETRA_FUNC_CALL_INCREMENT(
791  "Tpetra::MultiVector::multiply()",
792  tB->multiply(Teuchos::NO_TRANS,Teuchos::NO_TRANS,ST::one(),*tA,*tC,ST::zero()),
793  A->apply(Thyra::NOTRANS,*Ctmvb,B.ptr(),ST::one(),ST::zero())
794  );
795  }
796  */
797 
798 /*
799  RCP<Time>
800  timerUpdate = lookupAndAssertTimer("Tpetra::MultiVector::update(alpha,A,beta,B,gamma)"),
801  timerUpdate2 = lookupAndAssertTimer("Tpetra::MultiVector::update(alpha,A,beta)"),
802 
803  // TODO: test update(two vector)
804  // TODO: test update(three vector)
805 */
806 }
807 
808 
809 #ifdef TPETRA_TEUCHOS_TIME_MONITOR
810 # define TPETRA_TIMER_TESTS(SCALAR) \
811  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, UseTpetraImplementations, SCALAR )
812 #else
813 # define TPETRA_TIMER_TESTS(SCALAR)
814 #endif
815 
816 
817 //
818 // TpetraLinearOp_EpetraRowMatrix
819 //
820 
821 
822 #ifdef HAVE_THYRA_TPETRA_EPETRA
823 
824 
825 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_EpetraRowMatrix,
826  Scalar )
827 {
828 
829  using Teuchos::as;
830  using Teuchos::outArg;
831  using Teuchos::rcp_dynamic_cast;
832  using Teuchos::Array;
833  typedef Teuchos::ScalarTraits<Scalar> ST;
834 
835  const RCP<Tpetra::Operator<Scalar> > tpetraOp =
836  createTriDiagonalTpetraOperator<Scalar>(g_localDim);
837 
838  const RCP<LinearOpBase<Scalar> > thyraOp =
839  Thyra::createLinearOp(tpetraOp);
840 
841  const RCP<Thyra::TpetraLinearOp<Scalar> > thyraTpetraOp =
842  Teuchos::rcp_dynamic_cast<Thyra::TpetraLinearOp<Scalar> >(thyraOp);
843 
844  RCP<const Epetra_Operator> epetraOp;
845  Thyra::EOpTransp epetraOpTransp;
846  Thyra::EApplyEpetraOpAs epetraOpApplyAs;
847  Thyra::EAdjointEpetraOp epetraOpAdjointSupport;
848 
849  thyraTpetraOp->getEpetraOpView( outArg(epetraOp), outArg(epetraOpTransp),
850  outArg(epetraOpApplyAs), outArg(epetraOpAdjointSupport) );
851 
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);
866  // ToDo: Test column indices!
867  }
868  else {
869  TEST_ASSERT(is_null(epetraOp));
870  }
871 
872 }
873 
874 #else
875 
876 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_EpetraRowMatrix,
877  Scalar )
878 {
879 }
880 
881 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_RowStatLinearOpBase,
882  Scalar )
883 {
884  using Teuchos::as;
885 
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));
890 
891  const RCP<Tpetra::CrsMatrix<Scalar> > tpetraCrsMatrix =
892  Teuchos::rcp_dynamic_cast<Tpetra::CrsMatrix<Scalar> >(tpetraOp,true);
893 
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));
901 
902  const Teuchos::RCP<Thyra::RowStatLinearOpBase<Scalar> > rowStatOp =
903  Teuchos::rcp_dynamic_cast<Thyra::RowStatLinearOpBase<Scalar> >(thyraLinearOp, true);
904 
905  // Get the inverse row sums
906 
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());
911 
912  rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_INV_ROW_SUM,
913  inv_row_sums.ptr());
914  rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
915  row_sums.ptr());
916 
917  out << "inv_row_sums = " << *inv_row_sums;
918  out << "row_sums = " << *row_sums;
919 
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())
924  );
925 
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())
930  );
931 }
932 
933 TEUCHOS_UNIT_TEST_TEMPLATE_1_DECL( TpetraThyraWrappers, TpetraLinearOp_ScaledLinearOpBase,
934  Scalar )
935 {
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));
940 
941  const RCP<Tpetra::CrsMatrix<Scalar> > tpetraCrsMatrix =
942  Teuchos::rcp_dynamic_cast<Tpetra::CrsMatrix<Scalar> >(tpetraOp,true);
943 
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));
951 
952  const Teuchos::RCP<Thyra::RowStatLinearOpBase<Scalar> > rowStatOp =
953  Teuchos::rcp_dynamic_cast<Thyra::RowStatLinearOpBase<Scalar> >(thyraLinearOp, true);
954 
955  // Get the inverse row sums
956 
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());
961 
962  rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_INV_ROW_SUM,
963  inv_row_sums.ptr());
964  rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
965  row_sums.ptr());
966 
967  out << "inv_row_sums = " << *inv_row_sums;
968  out << "row_sums = " << *row_sums;
969 
970  const Teuchos::RCP<Thyra::ScaledLinearOpBase<Scalar> > scaledOp =
971  Teuchos::rcp_dynamic_cast<Thyra::ScaledLinearOpBase<Scalar> >(thyraLinearOp, true);
972 
973  TEUCHOS_ASSERT(scaledOp->supportsScaleLeft());
974 
975  scaledOp->scaleLeft(*inv_row_sums);
976 
977  rowStatOp->getRowStat(Thyra::RowStatLinearOpBaseUtils::ROW_STAT_ROW_SUM,
978  row_sums.ptr());
979 
980  out << "row_sums after left scaling by inv_row_sum = " << *row_sums;
981 
982  // scaled row sums should be one for each entry
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())
987  );
988 
989  // Don't currently check the results of right scaling. Tpetra tests
990  // already check this. Once column sums are supported in tpetra
991  // adapters, this can be checked easily.
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;
996 }
997 
998 #endif // HAVE_THYRA_TPETRA_EPETRA
999 
1000 
1001 //
1002 // Unit test instantiations
1003 //
1004 
1005 #define THYRA_TPETRA_THYRA_WRAPPERS_INSTANT(SCALAR) \
1006  \
1007  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1008  convertTpetraToThyraComm, SCALAR ) \
1009  \
1010  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1011  createVectorSpace, SCALAR ) \
1012  \
1013  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1014  createVector, SCALAR ) \
1015  \
1016  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1017  createConstVector, SCALAR ) \
1018  \
1019  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1020  createMultiVector, SCALAR ) \
1021  \
1022  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1023  createConstMultiVector, SCALAR ) \
1024  \
1025  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1026  TeptraVectorSpace, SCALAR ) \
1027  \
1028  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1029  vectorSpaceTester, SCALAR ) \
1030  \
1031  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1032  getTpetraMultiVector, SCALAR ) \
1033  \
1034  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1035  getConstTpetraMultiVector, SCALAR ) \
1036  \
1037  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1038  TpetraLinearOp, SCALAR ) \
1039  \
1040  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1041  createLinearOp, SCALAR ) \
1042  \
1043  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1044  createConstLinearOp, SCALAR ) \
1045  \
1046  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1047  TpetraLinearOp_EpetraRowMatrix, SCALAR ) \
1048  \
1049  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1050  TpetraLinearOp_RowStatLinearOpBase, SCALAR ) \
1051  \
1052  TEUCHOS_UNIT_TEST_TEMPLATE_1_INSTANT( TpetraThyraWrappers, \
1053  TpetraLinearOp_ScaledLinearOpBase, SCALAR ) \
1054 
1055 
1056 // We can currently only explicitly instantiate with double support because
1057 // Tpetra only supports explicit instantaition with double. As for implicit
1058 // instantation, g++ 3.4.6 on my Linux machine was taking more than 30 minutes
1059 // to compile this file when all of the types double, float, complex<double>,
1060 // and complex<float> where enabled. Therefore, we will only test double for
1061 // now until explicit instantation with other types are supported by Tpetra.
1062 
1064 
1065 
1066 } // namespace Thyra
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)
Traits class that enables the extraction of Tpetra operator/vector objects wrapped in Thyra operator/...
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)