55#ifdef HAVE_TEUCHOS_COMPLEX
56#define STYPE std::complex<double>
63#ifdef HAVE_TEUCHOS_COMPLEX
64#define SCALARMAX STYPE(10,0)
66#define SCALARMAX STYPE(10)
69template<
typename TYPE>
78template<
typename TYPE>
90std::complex<T>
GetRandom( std::complex<T>, std::complex<T> );
103int main(
int argc,
char* argv[])
112 if (argc>1)
if (argv[1][0]==
'-' && argv[1][1]==
'v') verbose =
true;
117 int numberFailedTests = 0;
119 std::string testName =
"", testType =
"";
121#ifdef HAVE_TEUCHOS_COMPLEX
122 testType =
"COMPLEX";
127 if (verbose) std::cout<<std::endl<<
"********** CHECKING TEUCHOS SERIAL DENSE SOLVER - " << testType <<
"-VALUED **********"<<std::endl<<std::endl;
136 testName =
"Generating right-hand side vector using A*x, where x is a random vector:";
137 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
140 testName =
"Generating right-hand side vector using A^T*x, where x is a random vector:";
141 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
143#ifdef HAVE_TEUCHOS_COMPLEX
146 testName =
"Generating right-hand side vector using A^H*x, where x is a random vector:";
147 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
161 returnCode = solver1.
factor();
162 testName =
"Simple solve: factor() random A:";
163 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
166 returnCode = solver1.
solve();
167 testName =
"Simple solve: solve() random A (NO_TRANS):";
169 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
175 returnCode = solver1.
solve();
176 testName =
"Simple solve: solve() random A (TRANS):";
178 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
180#ifdef HAVE_TEUCHOS_COMPLEX
185 returnCode = solver1.
solve();
186 testName =
"Simple solve: solve() random A (CONJ_TRANS):";
188 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
192 returnCode = solver1.
invert();
193 testName =
"Simple solve: invert() random A:";
194 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
198 testName =
"Computing solution using inverted random A (NO_TRANS):";
200 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
203 testName =
"Computing solution using inverted random A (TRANS):";
205 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
207#ifdef HAVE_TEUCHOS_COMPLEX
209 testName =
"Computing solution using inverted random A (CONJ_TRANS):";
211 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
215#ifdef HAVE_TEUCHOSNUMERICS_EIGEN
226#ifdef HAVE_TEUCHOS_COMPLEX
239 returnCode = solver2.
factor();
240 testName =
"Solve with iterative refinement: factor() random A:";
241 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
244 returnCode = solver2.
solve();
245 testName =
"Solve with iterative refinement: solve() random A (NO_TRANS):";
247 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
253 returnCode = solver2.
solve();
254 testName =
"Solve with iterative refinement: solve() random A (TRANS):";
256 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
258#ifdef HAVE_TEUCHOS_COMPLEX
263 returnCode = solver2.
solve();
264 testName =
"Solve with iterative refinement: solve() random A (CONJ_TRANS):";
266 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
280#ifdef HAVE_TEUCHOS_COMPLEX
297 returnCode = solver3.
factor();
298 testName =
"Solve with matrix equilibration: factor() random A:";
299 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
302 returnCode = solver3.
solve();
303 testName =
"Solve with matrix equilibration: solve() random A (NO_TRANS):";
305 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
311 returnCode = solver3.
solve();
312 testName =
"Solve with matrix equilibration: solve() random A (TRANS):";
314 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
316#ifdef HAVE_TEUCHOS_COMPLEX
321 returnCode = solver3.
solve();
322 testName =
"Solve with matrix equilibration: solve() random A (CONJ_TRANS):";
324 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
333 returnCode = solver3.
solve();
334 testName =
"Solve with matrix equilibration: solve() without factor() random A (NO_TRANS):";
336 numberFailedTests +=
ReturnCodeCheck(testName, returnCode, 0, verbose);
341 if(numberFailedTests > 0)
344 std::cout <<
"Number of failed tests: " << numberFailedTests << std::endl;
345 std::cout <<
"End Result: TEST FAILED" << std::endl;
349 if(numberFailedTests == 0)
350 std::cout <<
"End Result: TEST PASSED" << std::endl;
355template<
typename TYPE>
356int PrintTestResults(std::string testName, TYPE calculatedResult, TYPE expectedResult,
bool verbose)
359 if(calculatedResult == expectedResult)
361 if(verbose) std::cout << testName <<
" successful." << std::endl;
366 if(verbose) std::cout << testName <<
" unsuccessful." << std::endl;
372int ReturnCodeCheck(std::string testName,
int returnCode,
int expectedResult,
bool verbose)
375 if(expectedResult == 0)
379 if(verbose) std::cout << testName <<
" test successful." << std::endl;
384 if(verbose) std::cout << testName <<
" test unsuccessful. Return code was " << returnCode <<
"." << std::endl;
392 if(verbose) std::cout << testName <<
" test successful -- failed as expected." << std::endl;
397 if(verbose) std::cout << testName <<
" test unsuccessful -- did not fail as expected. Return code was " << returnCode <<
"." << std::endl;
404template<
typename TYPE>
411std::complex<T>
GetRandom( std::complex<T> Low, std::complex<T> High)
413 T lowMag = Low.real();
414 T highMag = High.real();
417 return std::complex<T>( real, imag );
437 for (
int i=0; i<m; i++)
438 for (
int j=0; j<n; j++)
449 for (
int i=0; i<n; i++)
470 MagnitudeType temp = norm_diff;
474 if (temp > Tolerance)
477 std::cout <<
"COMPARISON FAILED : ";
Reference-counted pointer class and non-member templated function implementations.
Defines basic traits for the scalar field type.
Non-member helper functions on the templated serial, dense matrix/vector classes.
Templated serial dense matrix class.
Templated class for solving dense linear problems.
Templated serial dense vector class.
Smart reference counting pointer class for automatic garbage collection.
This class creates and provides basic support for dense rectangular matrix of templated type.
int multiply(ETransp transa, ETransp transb, ScalarType alpha, const SerialDenseMatrix< OrdinalType, ScalarType > &A, const SerialDenseMatrix< OrdinalType, ScalarType > &B, ScalarType beta)
Multiply A * B and add them to this; this = beta * this + alpha*A*B.
ScalarTraits< ScalarType >::magnitudeType normFrobenius() const
Returns the Frobenius-norm of the matrix.
A class for solving dense linear problems.
void factorWithEquilibration(bool flag)
Causes equilibration to be called just before the matrix factorization as part of the call to factor.
int setMatrix(const RCP< SerialDenseMatrix< OrdinalType, ScalarType > > &A)
Sets the pointers for coefficient matrix.
int invert()
Inverts the this matrix.
int factor()
Computes the in-place LU factorization of the matrix using the LAPACK routine _GETRF.
void solveWithTransposeFlag(Teuchos::ETransp trans)
All subsequent function calls will work with the transpose-type set by this method (Teuchos::NO_TRANS...
void solveToRefinedSolution(bool flag)
Causes all solves to compute solution to best ability using iterative refinement.
int setVectors(const RCP< SerialDenseMatrix< OrdinalType, ScalarType > > &X, const RCP< SerialDenseMatrix< OrdinalType, ScalarType > > &B)
Sets the pointers for left and right hand side vector(s).
int solve()
Computes the solution X to AX = B for the this matrix and the B provided to SetVectors()....
This class creates and provides basic support for dense vectors of templated type as a specialization...
TYPE GetRandom(TYPE, TYPE)
Teuchos::RCP< DMatrix > GetRandomMatrix(int m, int n)
SerialDenseVector< OTYPE, STYPE > DVector
int PrintTestResults(std::string, TYPE, TYPE, bool)
int CompareVectors(const SerialDenseVector< OTYPE, STYPE > &Vector1, const SerialDenseVector< OTYPE, STYPE > &Vector2, ScalarTraits< STYPE >::magnitudeType Tolerance, bool verbose)
SerialDenseMatrix< OTYPE, STYPE > DMatrix
int ReturnCodeCheck(std::string, int, int, bool)
Teuchos::RCP< DVector > GetRandomVector(int n)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Deprecated.
std::string Teuchos_Version()
This structure defines some basic traits for a scalar field type.
static T one()
Returns representation of one for this scalar type.