85 const ParameterList& pL = GetParameterList();
86 const size_t row = Teuchos::as<size_t>(pL.get<
int>(
"block row"));
87 const size_t col = Teuchos::as<size_t>(pL.get<
int>(
"block col"));
90 RCP<Matrix> Ain = currentLevel.
Get<RCP<Matrix>>(
"A", this->GetFactory(
"A").get());
91 RCP<BlockedCrsMatrix> A = rcp_dynamic_cast<BlockedCrsMatrix>(Ain);
94 TEUCHOS_TEST_FOR_EXCEPTION(A.is_null(),
Exceptions::BadCast,
"Input matrix A is not a BlockedCrsMatrix.");
95 TEUCHOS_TEST_FOR_EXCEPTION(row > A->Rows(),
Exceptions::RuntimeError,
"row [" << row <<
"] > A.Rows() [" << A->Rows() <<
"].");
96 TEUCHOS_TEST_FOR_EXCEPTION(col > A->Cols(),
Exceptions::RuntimeError,
"col [" << col <<
"] > A.Cols() [" << A->Cols() <<
"].");
99 RCP<const StridedMap> stridedRangeMap = Teuchos::null;
100 RCP<const StridedMap> stridedDomainMap = Teuchos::null;
103 RCP<const MapExtractor> rangeMapExtractor = A->getRangeMapExtractor();
104 RCP<const MapExtractor> domainMapExtractor = A->getDomainMapExtractor();
109 if (stridedRangeMap.is_null()) {
111 "Range map extractor contains non-strided maps in block row " << row <<
". This should not be.");
112 stridedRangeMap = rcp_dynamic_cast<const StridedMap>(rangeMapExtractor->getMap(row));
115 if (stridedDomainMap.is_null()) {
117 "Domain map extractor contains non-strided maps in block row " << row <<
". This should not be.");
118 stridedDomainMap = rcp_dynamic_cast<const StridedMap>(domainMapExtractor->getMap(col));
121 TEUCHOS_TEST_FOR_EXCEPTION(stridedRangeMap.is_null(),
Exceptions::BadCast,
"rangeMap " << row <<
" is not a strided map.");
122 TEUCHOS_TEST_FOR_EXCEPTION(stridedDomainMap.is_null(),
Exceptions::BadCast,
"domainMap " << col <<
" is not a strided map.");
124 RCP<Matrix> Op = MatrixFactory::Build(stridedRangeMap, stridedDomainMap,
static_cast<size_t>(0));
125 TEUCHOS_ASSERT(!Op.is_null());
127 Op->fillComplete(stridedDomainMap, stridedRangeMap);
128 TEUCHOS_ASSERT(Op->isFillComplete());
130 GetOStream(
Statistics1) <<
"A(" << row <<
"," << col <<
") is a single block and has strided maps:"
131 <<
"\n range map fixed block size = " << stridedRangeMap ->getFixedBlockSize() <<
", strided block id = " << stridedRangeMap ->getStridedBlockId()
132 <<
"\n domain map fixed block size = " << stridedDomainMap->getFixedBlockSize() <<
", strided block id = " << stridedDomainMap->getStridedBlockId() << std::endl;
133 GetOStream(
Statistics2) <<
"A(" << row <<
"," << col <<
") has " << Op->getGlobalNumRows() <<
"x" << Op->getGlobalNumCols() <<
" rows and columns." << std::endl;
135 if (Op->IsView(
"stridedMaps") ==
true)
136 Op->RemoveView(
"stridedMaps");
137 Op->CreateView(
"stridedMaps", stridedRangeMap, stridedDomainMap);
139 currentLevel.
Set(
"A", Op,
this);
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access)....