40#ifndef TPETRA_BLOCKCRSMATRIX_HELPERS_DEF_HPP
41#define TPETRA_BLOCKCRSMATRIX_HELPERS_DEF_HPP
45#include "Tpetra_BlockCrsMatrix.hpp"
46#include "Tpetra_CrsMatrix.hpp"
47#include "Tpetra_HashTable.hpp"
48#include "Tpetra_Import.hpp"
49#include "Tpetra_Map.hpp"
50#include "Tpetra_MultiVector.hpp"
51#include "Teuchos_ParameterList.hpp"
52#include "Teuchos_ScalarTraits.hpp"
58 template<
class Scalar,
class LO,
class GO,
class Node>
60 Teuchos::ParameterList pl;
62 out.open(fileName.c_str());
66 template<
class Scalar,
class LO,
class GO,
class Node>
69 out.open(fileName.c_str());
73 template<
class Scalar,
class LO,
class GO,
class Node>
75 Teuchos::ParameterList pl;
79 template<
class Scalar,
class LO,
class GO,
class Node>
85 typedef Teuchos::OrdinalTraits<GO> TOT;
92 RCP<const map_type>
const rowMap = A.
getRowMap();
93 RCP<const Teuchos::Comm<int> > comm = rowMap->getComm();
94 const int myRank = comm->getRank();
95 const size_t numProcs = comm->getSize();
98 bool alwaysUseParallelAlgorithm =
false;
99 if (params.isParameter(
"always use parallel algorithm"))
100 alwaysUseParallelAlgorithm = params.get<
bool>(
"always use parallel algorithm");
101 bool printMatrixMarketHeader =
true;
102 if (params.isParameter(
"print MatrixMarket header"))
103 printMatrixMarketHeader = params.get<
bool>(
"print MatrixMarket header");
105 if (printMatrixMarketHeader && myRank==0) {
106 std::time_t now = std::time(NULL);
108 const std::string dataTypeStr =
109 Teuchos::ScalarTraits<Scalar>::isComplex ?
"complex" :
"real";
119 os <<
"%%MatrixMarket matrix coordinate " << dataTypeStr <<
" general" << std::endl;
120 os <<
"% time stamp: " << ctime(&now);
121 os <<
"% written from " << numProcs <<
" processes" << std::endl;
122 os <<
"% point representation of Tpetra::BlockCrsMatrix" << std::endl;
125 os <<
"% " << numRows <<
" block rows, " << numCols <<
" block columns" << std::endl;
127 os <<
"% block size " << blockSize << std::endl;
128 os << numRows*blockSize <<
" " << numCols*blockSize <<
" " << A.
getGlobalNumEntries()*blockSize*blockSize << std::endl;
131 if (numProcs==1 && !alwaysUseParallelAlgorithm) {
134 size_t numRows = rowMap->getLocalNumElements();
137 RCP<const map_type> allMeshGidsMap = rcp(
new map_type(TOT::invalid(), numRows, A.
getIndexBase(), comm));
140 mv_type allMeshGids(allMeshGidsMap,1);
141 Teuchos::ArrayRCP<GO> allMeshGidsData = allMeshGids.getDataNonConst(0);
143 for (
size_t i=0; i<numRows; i++)
144 allMeshGidsData[i] = rowMap->getGlobalElement(i);
145 allMeshGidsData = Teuchos::null;
148 size_t stripSize = allMeshGids.getGlobalLength() / numProcs;
149 size_t remainder = allMeshGids.getGlobalLength() % numProcs;
151 size_t curStripSize = 0;
152 Teuchos::Array<GO> importMeshGidList;
153 for (
size_t i=0; i<numProcs; i++) {
155 curStripSize = stripSize;
156 if (i<remainder) curStripSize++;
157 importMeshGidList.resize(curStripSize);
158 for (
size_t j=0; j<curStripSize; j++) importMeshGidList[j] = j + curStart + A.
getIndexBase();
159 curStart += curStripSize;
162 TEUCHOS_TEST_FOR_EXCEPTION(myRank>0 && curStripSize!=0,
163 std::runtime_error,
"Tpetra::blockCrsMatrixWriter: (pid "
164 << myRank <<
") map size should be zero, but is " << curStripSize);
165 RCP<map_type> importMeshGidMap = rcp(
new map_type(TOT::invalid(), importMeshGidList(), A.
getIndexBase(), comm));
166 import_type gidImporter(allMeshGidsMap, importMeshGidMap);
167 mv_type importMeshGids(importMeshGidMap, 1);
168 importMeshGids.doImport(allMeshGids, gidImporter,
INSERT);
174 Teuchos::ArrayRCP<const GO> importMeshGidsData = importMeshGids.getData(0);
175 Teuchos::Array<GO> importMeshGidsGO;
176 importMeshGidsGO.reserve(importMeshGidsData.size());
177 for (
typename Teuchos::ArrayRCP<const GO>::size_type j=0; j<importMeshGidsData.size(); ++j)
178 importMeshGidsGO.push_back(importMeshGidsData[j]);
179 RCP<const map_type> importMap = rcp(
new map_type(TOT::invalid(), importMeshGidsGO(), rowMap->getIndexBase(), comm) );
181 import_type importer(rowMap,importMap );
183 RCP<crs_graph_type> graph =
createCrsGraph(importMap,numEntriesPerRow);
184 RCP<const map_type> domainMap = A.getCrsGraph().
getDomainMap();
185 graph->doImport(A.getCrsGraph(), importer,
INSERT);
186 graph->fillComplete(domainMap, importMap);
188 block_crs_matrix_type importA(*graph, A.
getBlockSize());
189 importA.doImport(A, importer,
INSERT);
197 template<
class Scalar,
class LO,
class GO,
class Node>
202 using bcrs_local_inds_host_view_type =
typename bcrs_type::local_inds_host_view_type;
203 using bcrs_values_host_view_type =
typename bcrs_type::values_host_view_type;
204 using impl_scalar_type =
typename bcrs_type::impl_scalar_type;
207 RCP<const map_type> rowMap = A.
getRowMap();
208 RCP<const map_type> colMap = A.
getColMap();
209 RCP<const Teuchos::Comm<int> > comm = rowMap->getComm();
210 const int myRank = comm->getRank();
212 const size_t meshRowOffset = rowMap->getIndexBase();
213 const size_t meshColOffset = colMap->getIndexBase();
214 TEUCHOS_TEST_FOR_EXCEPTION(meshRowOffset != meshColOffset,
215 std::runtime_error,
"Tpetra::writeMatrixStrip: "
216 "mesh row index base != mesh column index base");
221 std::runtime_error,
"Tpetra::writeMatrixStrip: pid "
224 std::runtime_error,
"Tpetra::writeMatrixStrip: pid "
225 << myRank <<
" should have 0 columns but has " << A.
getLocalNumCols());
230 std::runtime_error,
"Tpetra::writeMatrixStrip: "
231 "number of rows on pid 0 does not match global number of rows");
237 bool precisionChanged=
false;
238 int oldPrecision = 0;
239 if (params.isParameter(
"precision")) {
240 oldPrecision = os.precision(params.get<
int>(
"precision"));
241 precisionChanged=
true;
244 if (params.isParameter(
"zero-based indexing")) {
245 if (params.get<
bool>(
"zero-based indexing") ==
true)
250 for (localRowInd = 0; localRowInd < numLocalRows; ++localRowInd) {
253 bcrs_local_inds_host_view_type localColInds;
254 bcrs_values_host_view_type vals;
256 A.
getLocalRowView (localRowInd, localColInds, vals); numEntries = localColInds.extent(0);
257 GO globalMeshRowID = rowMap->getGlobalElement(localRowInd) - meshRowOffset;
259 for (LO k = 0; k < numEntries; ++k) {
260 GO globalMeshColID = colMap->getGlobalElement(localColInds[k]) - meshColOffset;
261 const impl_scalar_type* curBlock = vals.data() + blockSize * blockSize * k;
263 for (LO j = 0; j < blockSize; ++j) {
264 GO globalPointRowID = globalMeshRowID * blockSize + j + pointOffset;
265 for (LO i = 0; i < blockSize; ++i) {
266 GO globalPointColID = globalMeshColID * blockSize + i + pointOffset;
267 const impl_scalar_type curVal = curBlock[i + j * blockSize];
269 os << globalPointRowID <<
" " << globalPointColID <<
" ";
270 if (Teuchos::ScalarTraits<impl_scalar_type>::isComplex) {
273 os << Teuchos::ScalarTraits<impl_scalar_type>::real (curVal) <<
" "
274 << Teuchos::ScalarTraits<impl_scalar_type>::imag (curVal);
284 if (precisionChanged)
285 os.precision(oldPrecision);
286 TEUCHOS_TEST_FOR_EXCEPTION(err != 0,
287 std::runtime_error,
"Tpetra::writeMatrixStrip: "
288 "error getting view of local row " << localRowInd);
294 template<
class Scalar,
class LO,
class GO,
class Node>
295 Teuchos::RCP<BlockCrsMatrix<Scalar, LO, GO, Node> >
307 using Teuchos::Array;
308 using Teuchos::ArrayView;
316 const map_type &pointRowMap = *(pointMatrix.
getRowMap());
317 RCP<const map_type> meshRowMap = createMeshMap<LO,GO,Node>(blockSize, pointRowMap);
319 const map_type &pointColMap = *(pointMatrix.
getColMap());
320 RCP<const map_type> meshColMap = createMeshMap<LO,GO,Node>(blockSize, pointColMap);
321 if(meshColMap.is_null())
throw std::runtime_error(
"ERROR: Cannot create mesh colmap");
323 const map_type &pointDomainMap = *(pointMatrix.
getDomainMap());
324 RCP<const map_type> meshDomainMap = createMeshMap<LO,GO,Node>(blockSize, pointDomainMap);
326 const map_type &pointRangeMap = *(pointMatrix.
getRangeMap());
327 RCP<const map_type> meshRangeMap = createMeshMap<LO,GO,Node>(blockSize, pointRangeMap);
332 RCP<crs_graph_type> meshCrsGraph = rcp(
new crs_graph_type(meshRowMap, meshColMap,
338 typename crs_matrix_type::local_inds_host_view_type pointColInds;
339 typename crs_matrix_type::values_host_view_type pointVals;
340 Array<GO> meshColGids;
345 GO indexBase = pointColMap.getIndexBase();
347 for (
int j=0; j<blockSize; ++j) {
348 LO rowLid = i*blockSize+j;
351 for (
size_t k=0; k<pointColInds.size(); ++k) {
352 GO meshColInd = (pointColMap.getGlobalElement(pointColInds[k]) - indexBase) / blockSize + indexBase;
353 if (meshColMap->getLocalElement(meshColInd) == Teuchos::OrdinalTraits<GO>::invalid()) {
354 std::ostringstream oss;
355 oss<<
"["<<i<<
"] ERROR: meshColId "<< meshColInd <<
" is not in the column map. Correspnds to pointColId = "<<pointColInds[k]<<std::endl;
356 throw std::runtime_error(oss.str());
359 meshColGids.push_back(meshColInd);
364 std::sort(meshColGids.begin(), meshColGids.end());
365 meshColGids.erase( std::unique(meshColGids.begin(), meshColGids.end()), meshColGids.end() );
366 meshCrsGraph->insertGlobalIndices(meshRowMap->getGlobalElement(i), meshColGids());
369 meshCrsGraph->fillComplete(meshDomainMap,meshRangeMap);
372 RCP<block_crs_matrix_type> blockMatrix = rcp(
new block_crs_matrix_type(*meshCrsGraph, blockSize));
375 Array<Scalar> tmpBlock(blockSize*blockSize);
378 int maxBlockEntries = blockMatrix->getLocalMaxNumRowEntries();
379 Array<Array<Scalar>> blocks(maxBlockEntries);
380 for (
int i=0; i<maxBlockEntries; ++i)
381 blocks[i].reserve(blockSize*blockSize);
382 std::map<int,int> bcol2bentry;
383 std::map<int,int>::iterator iter;
393 for (
int j=0; j<blockSize; ++j) {
394 LO rowLid = i*blockSize+j;
396 for (
size_t k=0; k<pointColInds.size(); ++k) {
398 LO meshColInd = pointColInds[k] / blockSize;
399 iter = bcol2bentry.find(meshColInd);
400 if (iter == bcol2bentry.end()) {
402 bcol2bentry[meshColInd] = blkCnt;
403 blocks[blkCnt].push_back(pointVals[k]);
407 int littleBlock = iter->second;
408 blocks[littleBlock].push_back(pointVals[k]);
414 for (iter=bcol2bentry.begin(); iter != bcol2bentry.end(); ++iter) {
415 LO localBlockCol = iter->first;
416 Scalar *vals = (blocks[iter->second]).getRawPtr();
417 if (std::is_same<typename block_crs_matrix_type::little_block_type::array_layout,Kokkos::LayoutLeft>::value) {
419 for (LO ii=0;ii<blockSize;++ii)
420 for (LO jj=0;jj<blockSize;++jj)
421 tmpBlock[ii+jj*blockSize] = vals[ii*blockSize+jj];
422 Scalar *tmp_vals = tmpBlock.getRawPtr();
423 blockMatrix->replaceLocalValues(i, &localBlockCol, tmp_vals, 1);
426 blockMatrix->replaceLocalValues(i, &localBlockCol, vals, 1);
431 for (
int j=0; j<maxBlockEntries; ++j)
443 template<
class LO,
class GO,
class Node>
444 Teuchos::RCP<const Tpetra::Map<LO,GO,Node> >
447 typedef Teuchos::OrdinalTraits<Tpetra::global_size_t> TOT;
452 Teuchos::Array<GO> meshGids;
458 meshGids.reserve(pointGids.size());
460 for (
int i=0; i<pointGids.size(); ++i) {
461 GO meshGid = (pointGids[i]-indexBase) / blockSize + indexBase;
462 if (hashTable.
get(meshGid) == -1) {
463 hashTable.
add(meshGid,1);
464 meshGids.push_back(meshGid);
468 Teuchos::RCP<const map_type> meshMap = Teuchos::rcp(
new map_type(TOT::invalid(), meshGids(), 0, pointMap.
getComm()) );
474 template<
class LO,
class GO,
class Node>
475 Teuchos::RCP<const Tpetra::Map<LO,GO,Node> >
478 typedef Teuchos::OrdinalTraits<Tpetra::global_size_t> TOT;
483 Teuchos::Array<GO> pointGids(blockGids.size() * blockSize);
486 for(LO i=0, ct=0; i<(LO)blockGids.size(); i++) {
487 GO base = (blockGids[i] - indexBase)* blockSize + indexBase;
488 for(LO j=0; j<blockSize; j++, ct++)
489 pointGids[i*blockSize+j] = base+j;
492 Teuchos::RCP<const map_type> pointMap = Teuchos::rcp(
new map_type(TOT::invalid(), pointGids(), 0, blockMap.
getComm()) );
498 template<
class Scalar,
class LO,
class GO,
class Node>
499 Teuchos::RCP<CrsMatrix<Scalar, LO, GO, Node> >
502 using Teuchos::Array;
503 using Teuchos::ArrayView;
510 using crs_graph_type =
typename block_crs_matrix_type::crs_graph_type;
511 using local_graph_device_type =
typename crs_matrix_type::local_graph_device_type;
512 using local_matrix_device_type =
typename crs_matrix_type::local_matrix_device_type;
513 using row_map_type =
typename local_graph_device_type::row_map_type::non_const_type;
514 using entries_type =
typename local_graph_device_type::entries_type::non_const_type;
515 using values_type =
typename local_matrix_device_type::values_type::non_const_type;
517 using row_map_type_const =
typename local_graph_device_type::row_map_type;
518 using entries_type_const =
typename local_graph_device_type::entries_type;
520 using little_block_type =
typename block_crs_matrix_type::const_little_block_type;
521 using offset_type =
typename row_map_type::non_const_value_type;
522 using column_type =
typename entries_type::non_const_value_type;
524 using execution_space =
typename Node::execution_space;
525 using range_type = Kokkos::RangePolicy<execution_space, LO>;
529 const offset_type bs2 = blocksize * blocksize;
531 size_t point_nnz = block_nnz * bs2;
534 RCP<const map_type> pointDomainMap = blockMatrix.
getDomainMap();
535 RCP<const map_type> pointRangeMap = blockMatrix.
getRangeMap();
538 RCP<const map_type> blockRowMap = blockMatrix.
getRowMap();
539 RCP<const map_type> pointRowMap = createPointMap<LO,GO,Node>(blocksize, *blockRowMap);
541 RCP<const map_type> blockColMap = blockMatrix.
getColMap();
542 RCP<const map_type> pointColMap = createPointMap<LO,GO,Node>(blocksize, *blockColMap);
546 const crs_graph_type & blockGraph = blockMatrix.getCrsGraph();
547 LO point_rows = (LO) pointRowMap->getLocalNumElements();
548 LO block_rows = (LO) blockRowMap->getLocalNumElements();
549 auto blockValues = blockMatrix.getValuesDevice();
550 auto blockLocalGraph = blockGraph.getLocalGraphDevice();
551 row_map_type_const blockRowptr = blockLocalGraph.row_map;
552 entries_type_const blockColind = blockLocalGraph.entries;
555 row_map_type rowptr(
"row_map", point_rows+1);
556 entries_type colind(
"entries", point_nnz);
557 values_type values(
"values", point_nnz);
560 Kokkos::parallel_for(
"fillRowPtr",range_type(0,block_rows),KOKKOS_LAMBDA(
const LO i) {
561 offset_type base = blockRowptr[i];
562 offset_type diff = blockRowptr[i+1] - base;
563 for(LO j=0; j<blocksize; j++) {
564 rowptr[i*blocksize +j] = base*bs2 + j*diff*blocksize;
568 if(i==block_rows-1) {
569 rowptr[block_rows*blocksize] = blockRowptr[i+1] * bs2;
574 Kokkos::parallel_for(
"copyEntriesAndValues",range_type(0,block_rows),KOKKOS_LAMBDA(
const LO i) {
575 const offset_type blkBeg = blockRowptr[i];
576 const offset_type numBlocks = blockRowptr[i+1] - blkBeg;
579 for (offset_type block=0; block < numBlocks; block++) {
580 column_type point_col_base = blockColind[blkBeg + block] * blocksize;
581 little_block_type my_block(blockValues.data () + (blkBeg+block) * bs2, blocksize, blocksize);
584 for(LO little_row=0; little_row<blocksize; little_row++) {
585 offset_type point_row_offset = rowptr[i*blocksize + little_row];
586 for(LO little_col=0; little_col<blocksize; little_col++) {
587 colind[point_row_offset + block*blocksize + little_col] = point_col_base + little_col;
588 values[point_row_offset + block*blocksize + little_col] = my_block(little_row,little_col);
596 RCP<crs_matrix_type> pointCrsMatrix = rcp(
new crs_matrix_type(pointRowMap, pointColMap, 0));
597 pointCrsMatrix->setAllValues(rowptr,colind,values);
600 pointCrsMatrix->expertStaticFillComplete(pointDomainMap,pointRangeMap);
602 return pointCrsMatrix;
614#define TPETRA_BLOCKCRSMATRIX_HELPERS_INSTANT(S,LO,GO,NODE) \
615 template void blockCrsMatrixWriter(BlockCrsMatrix<S,LO,GO,NODE> const &A, std::string const &fileName); \
616 template void blockCrsMatrixWriter(BlockCrsMatrix<S,LO,GO,NODE> const &A, std::string const &fileName, Teuchos::ParameterList const ¶ms); \
617 template void blockCrsMatrixWriter(BlockCrsMatrix<S,LO,GO,NODE> const &A, std::ostream &os); \
618 template void blockCrsMatrixWriter(BlockCrsMatrix<S,LO,GO,NODE> const &A, std::ostream &os, Teuchos::ParameterList const ¶ms); \
619 template void writeMatrixStrip(BlockCrsMatrix<S,LO,GO,NODE> const &A, std::ostream &os, Teuchos::ParameterList const ¶ms); \
620 template Teuchos::RCP<BlockCrsMatrix<S, LO, GO, NODE> > convertToBlockCrsMatrix(const CrsMatrix<S, LO, GO, NODE>& pointMatrix, const LO &blockSize); \
621 template Teuchos::RCP<CrsMatrix<S, LO, GO, NODE> > convertToCrsMatrix(const BlockCrsMatrix<S, LO, GO, NODE>& blockMatrix);
626#define TPETRA_CREATEMESHMAP_INSTANT(LO,GO,NODE) \
627 template Teuchos::RCP<const Map<LO,GO,NODE> > createMeshMap (const LO& blockSize, const Map<LO,GO,NODE>& pointMap); \
628 template Teuchos::RCP<const Map<LO,GO,NODE> > createPointMap (const LO& blockSize, const Map<LO,GO,NODE>& blockMap);
Sparse matrix whose entries are small dense square blocks, all of the same dimensions.
virtual size_t getLocalNumCols() const override
The number of columns needed to apply the forward operator on this node.
Teuchos::RCP< const map_type > getRangeMap() const override
Get the (point) range Map of this matrix.
size_t getLocalNumRows() const override
get the local number of block rows
virtual size_t getLocalNumEntries() const override
The local number of stored (structurally nonzero) entries.
Teuchos::RCP< const map_type > getColMap() const override
get the (mesh) map for the columns of this block matrix.
virtual global_size_t getGlobalNumCols() const override
The global number of columns of this matrix.
Teuchos::RCP< const map_type > getRowMap() const override
get the (mesh) map for the rows of this block matrix.
void getLocalRowView(LO LocalRow, local_inds_host_view_type &indices, values_host_view_type &values) const override
Get a view of the (mesh, i.e., block) row, using local (mesh, i.e., block) indices.
virtual GO getIndexBase() const override
The index base for global indices in this matrix.
virtual global_size_t getGlobalNumEntries() const override
The global number of stored (structurally nonzero) entries.
global_size_t getGlobalNumRows() const override
get the global number of block rows
Teuchos::RCP< const map_type > getDomainMap() const override
Get the (point) domain Map of this matrix.
LO getBlockSize() const
The number of degrees of freedom per mesh point.
A distributed graph accessed by rows (adjacency lists) and stored sparsely.
Teuchos::RCP< const map_type > getDomainMap() const override
Returns the Map associated with the domain of this graph.
size_t getGlobalMaxNumRowEntries() const override
Maximum number of entries in any row of the graph, over all processes in the graph's communicator.
Sparse matrix that presents a row-oriented interface that lets users read or modify entries.
size_t getGlobalMaxNumRowEntries() const override
Maximum number of entries in any row of the matrix, over all processes in the matrix's communicator.
Teuchos::RCP< const map_type > getDomainMap() const override
The domain Map of this matrix.
Teuchos::RCP< const map_type > getRangeMap() const override
The range Map of this matrix.
void getLocalRowView(LocalOrdinal LocalRow, local_inds_host_view_type &indices, values_host_view_type &values) const override
Get a constant view of a row of this matrix, using local row and column indices.
Teuchos::RCP< const map_type > getColMap() const override
The Map that describes the column distribution in this matrix.
size_t getLocalNumRows() const override
The number of matrix rows owned by the calling process.
Teuchos::RCP< const map_type > getRowMap() const override
The Map that describes the row distribution in this matrix.
ValueType get(const KeyType key)
Get the value corresponding to the given key.
void add(const KeyType key, const ValueType value)
Add a key and its value to the hash table.
Communication plan for data redistribution from a uniquely-owned to a (possibly) multiply-owned distr...
A parallel distribution of indices over processes.
Teuchos::ArrayView< const global_ordinal_type > getLocalElementList() const
Return a NONOWNING view of the global indices owned by this process.
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
Accessors for the Teuchos::Comm and Kokkos Node objects.
global_ordinal_type getIndexBase() const
The index base for this Map.
One or more distributed dense vectors.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
void writeMatrixStrip(BlockCrsMatrix< Scalar, LO, GO, Node > const &A, std::ostream &os, Teuchos::ParameterList const ¶ms)
Helper function called by blockCrsMatrixWriter.
Teuchos::RCP< const Tpetra::Map< LO, GO, Node > > createPointMap(LO const &blockSize, const Tpetra::Map< LO, GO, Node > &blockMap)
Helper function to generate a point map from a block map (with a given block size) GIDs associated wi...
Teuchos::RCP< CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > createCrsGraph(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t maxNumEntriesPerRow=0, const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Nonmember function to create an empty CrsGraph given a row Map and the max number of entries allowed ...
Teuchos::RCP< CrsMatrix< Scalar, LO, GO, Node > > convertToCrsMatrix(const Tpetra::BlockCrsMatrix< Scalar, LO, GO, Node > &blockMatrix)
Non-member constructor that creates a point CrsMatrix from an existing BlockCrsMatrix.
Teuchos::RCP< BlockCrsMatrix< Scalar, LO, GO, Node > > convertToBlockCrsMatrix(const Tpetra::CrsMatrix< Scalar, LO, GO, Node > &pointMatrix, const LO &blockSize)
Non-member constructor that creates a BlockCrsMatrix from an existing point CrsMatrix.
void blockCrsMatrixWriter(BlockCrsMatrix< Scalar, LO, GO, Node > const &A, std::string const &fileName)
Helper function to write a BlockCrsMatrix. Calls the 3-argument version.
@ INSERT
Insert new values that don't currently exist.
Teuchos::RCP< const Tpetra::Map< LO, GO, Node > > createMeshMap(LO const &blockSize, const Tpetra::Map< LO, GO, Node > &pointMap)
Helper function to generate a mesh map from a point map. Important! It's assumed that point GIDs asso...