46 #include "Thyra_DefaultScaledAdjointLinearOp.hpp" 47 #include "Thyra_DefaultMultipliedLinearOp.hpp" 48 #include "Thyra_DefaultAddedLinearOp.hpp" 49 #include "Thyra_DefaultDiagonalLinearOp.hpp" 50 #include "Thyra_VectorStdOps.hpp" 51 #include "Thyra_TestingTools.hpp" 52 #include "Thyra_LinearOpTester.hpp" 56 #include "Thyra_DefaultIdentityLinearOp.hpp" 57 #include "EpetraExt_readEpetraLinearSystem.h" 58 #include "Teuchos_GlobalMPISession.hpp" 59 #include "Teuchos_VerboseObject.hpp" 60 #include "Teuchos_XMLParameterListHelpers.hpp" 61 #include "Teuchos_CommandLineProcessor.hpp" 62 #include "Teuchos_StandardCatchMacros.hpp" 65 # include "Epetra_MpiComm.h" 67 # include "Epetra_SerialComm.h" 70 #include "Teuchos_UnitTestHarness.hpp" 79 using Thyra::epetraExtAddTransformer;
80 using Thyra::VectorBase;
81 using Thyra::LinearOpBase;
82 using Thyra::createMember;
83 using Thyra::LinearOpTester;
85 using Thyra::multiply;
86 using Thyra::diagonal;
89 std::string matrixFile =
"";
90 std::string matrixFile2 =
"";
95 Teuchos::UnitTestRepository::getCLP().setOption(
96 "matrix-file", &matrixFile,
97 "Defines the Epetra_CrsMatrix to read in." );
98 Teuchos::UnitTestRepository::getCLP().setOption(
99 "matrix-file-2", &matrixFile2,
100 "Defines the Epetra_CrsMatrix to read in." );
103 const Teuchos::RCP<const Thyra::LinearOpBase<double> >
104 buildAddOperator(
int scenario,
const Teuchos::RCP<
const Thyra::LinearOpBase<double> > & A,
105 const Teuchos::RCP<
const Thyra::LinearOpBase<double> > & B)
108 RCP<const Thyra::LinearOpBase<double> > M;
112 M = Thyra::add(A,B,
"A+B");
115 M = Thyra::add(A,B,
"A+adj(B)");
118 M = Thyra::add(A,B,
"adj(A)+B");
121 M = Thyra::add(A,B,
"adb(A)+adb(B)");
124 TEUCHOS_ASSERT(
false);
138 out <<
"\nReading linear system in Epetra format from the file \'"<<matrixFile<<
"\' ...\n";
141 Epetra_MpiComm comm(MPI_COMM_WORLD);
143 Epetra_SerialComm comm;
145 RCP<Epetra_CrsMatrix> crsMat;
146 EpetraExt::readEpetraLinearSystem( matrixFile, comm, &crsMat, NULL, NULL, NULL );
155 const RCP<const Thyra::LinearOpBase<double> > A = Thyra::scale<double>(scaleMat,Thyra::epetraLinearOp(crsMat));
156 RCP<VectorBase<double> > d = createMember(A->domain(),
"d");
158 const RCP<const Thyra::LinearOpBase<double> > B = diagonal(d);
160 out <<
"\nA = " << *A;
161 out <<
"\nB = " << *B;
163 for(
int scenario=0;scenario<4;scenario++) {
168 const RCP<const Thyra::LinearOpBase<double> > M = buildAddOperator(scenario,A,B);
174 const RCP<EpetraExtAddTransformer> ApB_transformer = epetraExtAddTransformer();
176 TEST_ASSERT(ApB_transformer != null);
178 const RCP<LinearOpBase<double> > M_explicit = ApB_transformer->createOutputOp();
179 ApB_transformer->transform( *M, M_explicit.ptr() );
181 out <<
"\nM_explicit = " << *M_explicit;
187 LinearOpTester<double> M_explicit_tester;
188 M_explicit_tester.show_all_tests(
true);;
190 const bool result = M_explicit_tester.compare( *M, *M_explicit, Teuchos::inOutArg(out) );
191 if (!result) success =
false;
194 for(
int scenario=0;scenario<4;scenario++) {
199 const RCP<const Thyra::LinearOpBase<double> > M = buildAddOperator(scenario,B,A);
205 const RCP<EpetraExtAddTransformer> ApB_transformer = epetraExtAddTransformer();
207 TEST_ASSERT(ApB_transformer != null);
209 const RCP<LinearOpBase<double> > M_explicit = ApB_transformer->createOutputOp();
210 ApB_transformer->transform( *M, M_explicit.ptr() );
212 out <<
"\nM_explicit = " << *M_explicit;
218 LinearOpTester<double> M_explicit_tester;
219 M_explicit_tester.show_all_tests(
true);;
221 const bool result = M_explicit_tester.compare( *M, *M_explicit, Teuchos::inOutArg(out) );
222 if (!result) success =
false;
233 out <<
"\nReading linear system in Epetra format from the file \'"<<matrixFile<<
"\' ...\n";
236 Epetra_MpiComm comm(MPI_COMM_WORLD);
238 Epetra_SerialComm comm;
240 RCP<Epetra_CrsMatrix> crsMat;
241 EpetraExt::readEpetraLinearSystem( matrixFile, comm, &crsMat, NULL, NULL, NULL );
250 const RCP<const Thyra::LinearOpBase<double> > A = Thyra::scale<double>(scaleMat,Thyra::epetraLinearOp(crsMat));
251 const RCP<const Thyra::LinearOpBase<double> > B = identity(A->domain(),
"id");
253 out <<
"\nA = " << *A;
254 out <<
"\nB = " << *B;
256 for(
int scenario=0;scenario<4;scenario++) {
261 const RCP<const Thyra::LinearOpBase<double> > M = buildAddOperator(scenario,A,B);
267 const RCP<EpetraExtAddTransformer> ApB_transformer = epetraExtAddTransformer();
269 TEST_ASSERT(ApB_transformer != null);
271 const RCP<LinearOpBase<double> > M_explicit = ApB_transformer->createOutputOp();
272 ApB_transformer->transform( *M, M_explicit.ptr() );
274 out <<
"\nM_explicit = " << *M_explicit;
280 LinearOpTester<double> M_explicit_tester;
281 M_explicit_tester.show_all_tests(
true);;
283 const bool result = M_explicit_tester.compare( *M, *M_explicit, Teuchos::inOutArg(out) );
284 if (!result) success =
false;
287 for(
int scenario=0;scenario<4;scenario++) {
292 const RCP<const Thyra::LinearOpBase<double> > M = buildAddOperator(scenario,B,A);
298 const RCP<EpetraExtAddTransformer> ApB_transformer = epetraExtAddTransformer();
300 TEST_ASSERT(ApB_transformer != null);
302 const RCP<LinearOpBase<double> > M_explicit = ApB_transformer->createOutputOp();
303 ApB_transformer->transform( *M, M_explicit.ptr() );
305 out <<
"\nM_explicit = " << *M_explicit;
311 LinearOpTester<double> M_explicit_tester;
312 M_explicit_tester.show_all_tests(
true);;
314 const bool result = M_explicit_tester.compare( *M, *M_explicit, Teuchos::inOutArg(out) );
315 if (!result) success =
false;
326 out <<
"\nReading linear system in Epetra format from the file \'"<<matrixFile<<
"\' ...\n";
327 out <<
"\nReading linear system in Epetra format from the file \'"<<matrixFile2<<
"\' ...\n";
330 Epetra_MpiComm comm(MPI_COMM_WORLD);
332 Epetra_SerialComm comm;
334 RCP<Epetra_CrsMatrix> epetra_A;
335 RCP<Epetra_CrsMatrix> epetra_B;
336 EpetraExt::readEpetraLinearSystem( matrixFile, comm, &epetra_A, NULL, NULL, NULL );
337 EpetraExt::readEpetraLinearSystem( matrixFile2, comm, &epetra_B, NULL, NULL, NULL );
345 const RCP<const Thyra::LinearOpBase<double> > A = Thyra::scale<double>(scaleA,Thyra::epetraLinearOp(epetra_B));
346 const RCP<const Thyra::LinearOpBase<double> > B = Thyra::scale<double>(scaleB,Thyra::epetraLinearOp(epetra_B));
348 out <<
"\nA = " << *A;
349 out <<
"\nB = " << *B;
351 for(
int scenario=0;scenario<4;scenario++) {
356 const RCP<const Thyra::LinearOpBase<double> > M = buildAddOperator(scenario,A,B);
362 const RCP<EpetraExtAddTransformer> ApB_transformer = epetraExtAddTransformer();
364 TEST_ASSERT(ApB_transformer != null);
366 const RCP<LinearOpBase<double> > M_explicit = ApB_transformer->createOutputOp();
367 ApB_transformer->transform( *M, M_explicit.ptr() );
369 out <<
"\nM_explicit = " << *M_explicit;
375 LinearOpTester<double> M_explicit_tester;
376 M_explicit_tester.show_all_tests(
true);;
378 const bool result = M_explicit_tester.compare( *M, *M_explicit, Teuchos::inOutArg(out) );
379 if (!result) success =
false;
390 out <<
"\nReading linear system in Epetra format from the file \'"<<matrixFile<<
"\' ...\n";
391 out <<
"\nReading linear system in Epetra format from the file \'"<<matrixFile2<<
"\' ...\n";
394 Epetra_MpiComm comm(MPI_COMM_WORLD);
396 Epetra_SerialComm comm;
398 RCP<Epetra_CrsMatrix> epetra_A;
399 RCP<Epetra_CrsMatrix> epetra_B;
400 EpetraExt::readEpetraLinearSystem( matrixFile, comm, &epetra_A, NULL, NULL, NULL );
401 EpetraExt::readEpetraLinearSystem( matrixFile2, comm, &epetra_B, NULL, NULL, NULL );
409 const RCP<const Thyra::LinearOpBase<double> > A = Thyra::scale<double>(scaleA,Thyra::epetraLinearOp(epetra_B));
410 const RCP<const Thyra::LinearOpBase<double> > B = Thyra::scale<double>(scaleB,Thyra::epetraLinearOp(epetra_B));
412 out <<
"\nA = " << *A;
413 out <<
"\nB = " << *B;
415 for(
int scenario=0;scenario<4;scenario++) {
420 const RCP<const Thyra::LinearOpBase<double> > M = buildAddOperator(scenario,A,B);
426 const RCP<EpetraExtAddTransformer> ApB_transformer = epetraExtAddTransformer();
428 TEST_ASSERT(ApB_transformer != null);
430 const RCP<LinearOpBase<double> > M_explicit = ApB_transformer->createOutputOp();
431 const RCP<LinearOpBase<double> > M_explicit_orig = M_explicit;
432 ApB_transformer->transform( *M, M_explicit.ptr() );
435 double * view;
int numEntries;
436 epetra_B->Scale(3.2);
437 TEUCHOS_ASSERT(epetra_B->ExtractMyRowView(3,numEntries,view)==0);
438 for(
int i=0;i<numEntries;i++) view[i] += view[i]*
double(i+1);
441 ApB_transformer->transform( *M, M_explicit.ptr() );
443 out <<
"\nM_explicit = " << *M_explicit;
452 LinearOpTester<double> M_explicit_tester;
453 M_explicit_tester.show_all_tests(
true);;
455 const bool result = M_explicit_tester.compare( *M, *M_explicit, Teuchos::inoutArg(out) );
456 if (!result) success =
false;
TEUCHOS_UNIT_TEST(EpetraOperatorWrapper, basic)
Teuchos::RCP< Epetra_Operator > get_Epetra_Operator(LinearOpBase< double > &op)
Full specialization for Scalar=double.