Skip to content

Commit

Permalink
#215: finish cxx_main.cpp
Browse files Browse the repository at this point in the history
  • Loading branch information
cwschilly committed Sep 26, 2023
1 parent 875f269 commit 09ec914
Show file tree
Hide file tree
Showing 3 changed files with 39 additions and 534 deletions.
Original file line number Diff line number Diff line change
@@ -1,13 +1,4 @@

# TRIBITS_ADD_EXECUTABLE_AND_TEST(
# Tpetra_GeneralizedDavidson_solvertest
# SOURCES cxx_main_solvertest.cpp
# ARGS
# "--verbose"
# "--debug"
# COMM serial mpi
# )

ASSERT_DEFINED(${PACKAGE_NAME}_ENABLE_Triutils)
IF (${PACKAGE_NAME}_ENABLE_Triutils)
TRIBITS_ADD_EXECUTABLE_AND_TEST(
Expand Down
74 changes: 39 additions & 35 deletions packages/anasazi/tpetra/test/GeneralizedDavidson/cxx_main.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -60,46 +60,49 @@
#include "AnasaziTpetraAdapter.hpp"
#include "AnasaziBasicEigenproblem.hpp"
#include "AnasaziGeneralizedDavidsonSolMgr.hpp"

#include <Teuchos_CommandLineProcessor.hpp>

#include <Tpetra_Map.hpp>
#include <Tpetra_Core.hpp>
#include <Tpetra_Operator.hpp>
#include <Tpetra_CrsMatrix.hpp>
#include <Tpetra_MultiVector.hpp>

// I/O for Harwell-Boeing files
#include <Trilinos_Util_iohb.h>

int main(int argc, char *argv[])
{
// #ifndef HAVE_TPETRA_COMPLEX_DOUBLE
// # error "Anasazi: This test requires Scalar = std::complex<double> to be enabled in Tpetra."
// #else
template <typename ScalarType>
int run(int argc, char *argv[]) {

using ST = typename Tpetra::MultiVector<ScalarType>::scalar_type;
using LO = typename Tpetra::MultiVector<>::local_ordinal_type;
using GO = typename Tpetra::MultiVector<>::global_ordinal_type;
using NT = typename Tpetra::MultiVector<>::node_type;

using SCT = typename Teuchos::ScalarTraits<ScalarType>;
using MT = typename SCT::magnitudeType;

using OP = Tpetra::Operator<ST,LO,GO,NT>;
using MV = Tpetra::MultiVector<ST,LO,GO,NT>;
using OPT = Anasazi::OperatorTraits<ST,MV,OP>;
using MVT = Anasazi::MultiVecTraits<ST,MV>;

using tmap_t = Tpetra::Map<LO,GO,NT>;
using tcrsmatrix_t = Tpetra::CrsMatrix<ST,LO,GO,NT>;

using namespace Teuchos;
using Tpetra::CrsMatrix;
using Tpetra::Map;
using Tpetra::MultiVector;
using Tpetra::Operator;
using std::vector;
using std::cout;
using std::endl;

typedef double ST;
typedef ScalarTraits<ST> SCT;
typedef SCT::magnitudeType MT;
typedef MultiVector<ST> MV;
typedef MultiVector<ST>::global_ordinal_type GO;
typedef Operator<ST> OP;
typedef Anasazi::MultiVecTraits<ST,MV> MVT;
typedef Anasazi::OperatorTraits<ST,MV,OP> OPT;

Tpetra::ScopeGuard tpetraScope (&argc, &argv);

const ST ONE = SCT::one();
int info = 0;
int MyPID = 0;

RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm ();

MyPID = rank(*comm);
const auto comm = Tpetra::getDefaultComm ();
const int MyPID = comm->getRank();

bool testFailed;
bool verbose = false;
Expand Down Expand Up @@ -156,20 +159,20 @@ int main(int argc, char *argv[])
return -1;
}
// create map
RCP<const Map<> > map = rcp (new Map<> (dim, 0, comm));
RCP<CrsMatrix<ST> > K = rcp(new CrsMatrix<ST> (map, rnnzmax));
RCP<const tmap_t> map = rcp (new tmap_t (dim, 0, comm));
RCP<tcrsmatrix_t> K = rcp(new tcrsmatrix_t (map, rnnzmax));
if (MyPID == 0) {
// Convert interleaved doubles to complex values
// HB format is compressed column. CrsMatrix is compressed row.
const double *dptr = dvals;
const ST *dptr = dvals;
const int *rptr = rowind;
for (int c=0; c<dim; ++c) {
for (int colnnz=0; colnnz < colptr[c+1]-colptr[c]; ++colnnz) {
K->insertGlobalValues (*rptr++ - 1, tuple<GO> (c), tuple (ST (dptr[0], dptr[1])));
dptr += 2;
for (int c = 0; c < dim; ++c) {
for (int colnnz = 0; colnnz < colptr[c + 1] - colptr[c]; ++colnnz) {
ST value = dptr[0];
K->insertGlobalValues(*rptr++ - 1, tuple<GO>(c), tuple(value));
dptr++;
}
}
}

if (MyPID == 0) {
// Clean up.
free( dvals );
Expand All @@ -184,9 +187,7 @@ int main(int argc, char *argv[])

// Create eigenproblem
RCP<Anasazi::BasicEigenproblem<ST,MV,OP> > problem =
rcp( new Anasazi::BasicEigenproblem<ST,MV,OP>() );
problem->setA(K);
problem->setInitVec(ivec);
rcp( new Anasazi::BasicEigenproblem<ST,MV,OP>(K,ivec) );
//
// Inform the eigenproblem that the operator K is symmetric
problem->setHermitian(true);
Expand Down Expand Up @@ -302,6 +303,9 @@ int main(int argc, char *argv[])
cout << "End Result: TEST PASSED" << endl;
}
return 0;
}

// #endif // HAVE_TPETRA_COMPLEX_DOUBLE
int main(int argc, char *argv[]) {
return run<double>(argc,argv);
// run<float>(argc,argv);
}
Loading

0 comments on commit 09ec914

Please sign in to comment.