Amesos2 using "amesos2_cholmod" solver generates wrong results
Created by: davoodansari
Hi everyone Spent some time compiling metis, chomod and on top of that trilinos (amesos2 with cholmod enabled). Then I modified one of the examples to see if the "amesos2_cholmod" solver is actually accessible. To my surpride, the solver generates wrong results, i.e. [-2.88889; 2.25; -5.66667; 3.4; 3.77778; 4.66667] instead of [1 2 3 4 5 6]. Below is the code and the actual output. Any ideas what can be wrong? Any advice is highly appreciated Thanks
output:
Kokkos::OpenMP::initialize WARNING: OMP_PROC_BIND environment variable not set In general, for best performance with OpenMP 4.0 or better set OMP_PROC_BIND=spread and OMP_PLACES=threads For best performance with OpenMP 3.1 set OMP_PROC_BIND=true For unit testing set OMP_PROC_BIND=false Amesos2 using 3P solver: amesos2_cholmod Amesos2 version 0.3d - 07/28/2011
code:
Solution : "Tpetra::MultiVector": Template parameters: Scalar: double LocalOrdinal: int GlobalOrdinal: int Node: Serial/Wrapper Global number of rows: 6 Number of columns: 1 Process 0 of 1: Local number of rows: 6 Column stride: 6 Values: [-2.88889; 2.25; -5.66667; 3.4; 3.77778; 4.66667]
int main(int argc, char argv[]) { Tpetra::ScopeGuard tpetraScope(&argc,&argv); typedef double Scalar; typedef Tpetra::Map<>::local_ordinal_type LO; typedef Tpetra::Map<>::global_ordinal_type GO; typedef Tpetra::CrsMatrix<Scalar,LO,GO> MAT; typedef Tpetra::MultiVector<Scalar,LO,GO> MV; using Tpetra::global_size_t; using Teuchos::tuple; using Teuchos::RCP; using Teuchos::rcp; const int s=9; const char solverNames[10] = {"shylubasker", "basker", "klu2", "superlu_dist", "superlu_mt", "superlu", "pardiso_mkl", "lapack", "mumps", "amesos2_cholmod"};
// Before we do anything, check that SuperLU is enabled if( !Amesos2::query(solverNames[s]) ) { std::cerr << solverNames[s] << " not enabled. Exiting..." << std::endl; return EXIT_SUCCESS; // Otherwise CTest will pick it up as failure, which it isn't really } else std::cout << "Amesos2 using 3P solver: " << solverNames[s] << std::endl;
Teuchos::RCP<const Teuchos::Comm > comm = Tpetra::getDefaultComm(); size_t myRank = comm->getRank(); std::ostream &out = std::cout; out << Amesos2::version() << std::endl << std::endl; const size_t numVectors = 1; // create a Map global_size_t nrows = 6; RCP<Tpetra::Map<LO,GO> > map = rcp( new Tpetra::Map<LO,GO>(nrows,0,comm) ); RCP A = rcp( new MAT(map,3) ); // max of three entries in a row /*
- We will solve a system with a known solution, for which we will be using
- the following matrix:
- [ [ 7, 0, -3, 0, -1, 0 ]
- [ 2, 8, 0, 0, 0, 0 ]
- [ 0, 0, 1, 0, 0, 0 ]
- [ -3, 0, 0, 5, 0, 0 ]
- [ 0, -1, 0, 0, 4, 0 ]
- [ 0, 0, 0, -2, 0, 6 ] ]
/ // Construct matrix if( myRank == 0 ) { A->insertGlobalValues(0,tuple(0,2,4),tuple(7,-3,-1)); A->insertGlobalValues(1,tuple(0,1),tuple(2,8)); A->insertGlobalValues(2,tuple(2),tuple(1)); A->insertGlobalValues(3,tuple(0,3),tuple(-3,5)); A->insertGlobalValues(4,tuple(1,4),tuple(-1,4)); A->insertGlobalValues(5,tuple(3,5),tuple(-2,6)); } A->fillComplete(); // Create random X RCP X = rcp(new MV(map,numVectors)); //X->randomize(); / Create B *
- Use RHS:
- [[-7]
- [18]
- [ 3]
- [17]
- [18]
- [28]] / RCP B = rcp(new MV(map,numVectors)); int data[6] = {-7,18,3,17,18,28}; for( int i = 0; i < 6; ++i ) { if( B->getMap()->isNodeGlobalElement(i) ) { B->replaceGlobalValue(i,0,data[i]); } } // Create solver interface to Superlu with Amesos2 factory method RCP<Amesos2::Solver<MAT,MV> > solver = Amesos2::create<MAT,MV>(solverNames[s], A, X, B); solver->symbolicFactorization().numericFactorization().solve(); / Print the solution
- Should be:
- [[1]
- [2]
- [3]
- [4]
- [5]
- [6]] */ RCPTeuchos::FancyOStream fos = Teuchos::fancyOStream(Teuchos::rcpFromRef(out)); *fos << "Solution :" << std::endl; X->describe(*fos,Teuchos::VERB_EXTREME); *fos << std::endl; // We are done. return 0; }