Import redistributes from a uniquely owned (one-to-one) Map to a possibly not uniquely owned Map. Export redistributes from a possibly not uniquely owned to a uniquely owned Map. We distinguish between these cases both for historical reasons and for performance reasons.
In both cases, Import and Export let the user specify how to combine incoming new data with existing data that has the same global index. For example, one may replace old data with new data or sum them together.
#include <Teuchos_TimeMonitor.hpp>
#include <Tpetra_CrsMatrix.hpp>
#include <iostream>
Teuchos::RCP<Teuchos::Time> exportTimer;
template<class CrsMatrixType>
Teuchos::RCP<const CrsMatrixType>
createMatrix (const Teuchos::RCP<const typename CrsMatrixType::map_type>& map)
{
using Teuchos::arcp;
using Teuchos::ArrayRCP;
using Teuchos::ArrayView;
using Teuchos::RCP;
using Teuchos::rcp;
using Teuchos::Time;
using Teuchos::TimeMonitor;
using Teuchos::tuple;
typedef typename CrsMatrixType::scalar_type scalar_type;
typedef typename CrsMatrixType::local_ordinal_type LO;
typedef typename CrsMatrixType::global_ordinal_type GO;
RCP<Time> timer = TimeMonitor::getNewCounter ("Sparse matrix creation");
TimeMonitor monitor (*timer);
RCP<CrsMatrixType> A (new CrsMatrixType (map, 3));
const scalar_type two = static_cast<scalar_type> ( 2.0);
const scalar_type negOne = static_cast<scalar_type> (-1.0);
const GST numGlobalIndices = map->getGlobalNumElements ();
ArrayView<const GO> myGlobalElements = map->getLocalElementList ();
typedef typename ArrayView<const GO>::const_iterator iter_type;
for (iter_type it = myGlobalElements.begin(); it != myGlobalElements.end(); ++it) {
const LO i_local = *it;
const GO i_global = map->getGlobalElement (i_local);
if (i_global == 0) {
A->insertGlobalValues (i_global,
tuple (i_global, i_global+1),
tuple (two, negOne));
} else if (static_cast<GST> (i_global) == numGlobalIndices - 1) {
A->insertGlobalValues (i_global,
tuple (i_global-1, i_global),
tuple (negOne, two));
} else {
A->insertGlobalValues (i_global,
tuple (i_global-1, i_global, i_global+1),
tuple (negOne, two, negOne));
}
}
A->fillComplete ();
return A;
}
void
example (const Teuchos::RCP<const Teuchos::Comm<int> >& comm,
std::ostream& out)
{
using std::endl;
using Teuchos::ParameterList;
using Teuchos::RCP;
using Teuchos::rcp;
using Teuchos::Time;
using Teuchos::TimeMonitor;
const int myRank = comm->getRank ();
const GST numGlobalIndices = 10 * comm->getSize ();
const global_ordinal_type indexBase = 0;
if (myRank == 0) {
out << "Construct Process 0 Map" << endl;
}
RCP<const map_type> procZeroMap;
{
const size_t numLocalIndices = (myRank == 0) ? numGlobalIndices : 0;
procZeroMap = rcp (new map_type (numGlobalIndices, numLocalIndices,
indexBase, comm));
}
if (myRank == 0) {
out << "Construct global Map" << endl;
}
RCP<const map_type> globalMap =
rcp (new map_type (numGlobalIndices, indexBase, comm,
Tpetra::GloballyDistributed));
if (myRank == 0) {
out << "Create sparse matrix using Process 0 Map" << endl;
}
RCP<const crs_matrix_type> A = createMatrix<crs_matrix_type> (procZeroMap);
RCP<crs_matrix_type> B;
if (myRank == 0) {
out << "Redistribute sparse matrix" << endl;
}
{
TimeMonitor monitor (*exportTimer);
export_type exporter (procZeroMap, globalMap);
B = rcp (new crs_matrix_type (globalMap, 0));
}
B->fillComplete ();
}
int
main (int argc, char *argv[])
{
using Teuchos::TimeMonitor;
{
const int myRank = comm->getRank ();
exportTimer =
TimeMonitor::getNewCounter ("Sparse matrix redistribution");
example (comm, std::cout);
TimeMonitor::summarize (std::cout);
exportTimer = Teuchos::null;
if (myRank == 0) {
std::cout << "End Result: TEST PASSED" << std::endl;
}
}
return 0;
}
Functions for initializing and finalizing Tpetra.
Struct that holds views of the contents of a CrsMatrix.
Scope guard whose destructor automatically calls Tpetra::finalize for you.
Teuchos::RCP< const Teuchos::Comm< int > > getDefaultComm()
Get Tpetra's default communicator.
@ INSERT
Insert new values that don't currently exist.