Xpetra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Xpetra_IO.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Xpetra: A linear algebra interface package
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact
39 // Jonathan Hu (jhu@sandia.gov)
40 // Andrey Prokopenko (aprokop@sandia.gov)
41 // Ray Tuminaro (rstumin@sandia.gov)
42 // Tobias Wiesner (tawiesn@sandia.gov)
43 //
44 // ***********************************************************************
45 //
46 // @HEADER
47 
48 #ifndef PACKAGES_XPETRA_SUP_UTILS_XPETRA_IO_HPP_
49 #define PACKAGES_XPETRA_SUP_UTILS_XPETRA_IO_HPP_
50 
51 #include <fstream>
52 #include "Xpetra_ConfigDefs.hpp"
53 
54 #ifdef HAVE_XPETRA_EPETRA
55 #ifdef HAVE_MPI
56 #include "Epetra_MpiComm.h"
57 #endif
58 #endif
59 
60 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
61 #include <EpetraExt_MatrixMatrix.h>
62 #include <EpetraExt_RowMatrixOut.h>
63 #include <EpetraExt_MultiVectorOut.h>
64 #include <EpetraExt_CrsMatrixIn.h>
65 #include <EpetraExt_MultiVectorIn.h>
66 #include <EpetraExt_BlockMapIn.h>
67 #include <Xpetra_EpetraUtils.hpp>
68 #include <Xpetra_EpetraMultiVector.hpp>
69 #include <EpetraExt_BlockMapOut.h>
70 #endif
71 
72 #ifdef HAVE_XPETRA_TPETRA
73 #include <MatrixMarket_Tpetra.hpp>
74 #include <Tpetra_RowMatrixTransposer.hpp>
75 #include <TpetraExt_MatrixMatrix.hpp>
76 #include <Xpetra_TpetraMultiVector.hpp>
77 #include <Xpetra_TpetraCrsGraph.hpp>
78 #include <Xpetra_TpetraCrsMatrix.hpp>
79 #include <Xpetra_TpetraBlockCrsMatrix.hpp>
80 #include "Tpetra_Util.hpp"
81 #endif
82 
83 #ifdef HAVE_XPETRA_EPETRA
84 #include <Xpetra_EpetraMap.hpp>
85 #endif
86 
87 #include "Xpetra_Matrix.hpp"
88 #include "Xpetra_MatrixMatrix.hpp"
89 #include "Xpetra_CrsGraph.hpp"
90 #include "Xpetra_CrsMatrixWrap.hpp"
91 #include "Xpetra_BlockedCrsMatrix.hpp"
92 
93 #include "Xpetra_Map.hpp"
94 #include "Xpetra_StridedMap.hpp"
95 #include "Xpetra_StridedMapFactory.hpp"
96 #include "Xpetra_MapExtractor.hpp"
97 #include "Xpetra_MatrixFactory.hpp"
98 
99 #include <Teuchos_MatrixMarket_Raw_Writer.hpp>
100 #include <Teuchos_MatrixMarket_Raw_Reader.hpp>
101 #include <string>
102 
103 namespace Xpetra {
104 
105 #ifdef HAVE_XPETRA_EPETRA
106 // This non-member templated function exists so that the matrix-matrix multiply will compile if Epetra, Tpetra, and ML are enabled.
107 template <class SC, class LO, class GO, class NO>
108 RCP<Xpetra::CrsMatrixWrap<SC, LO, GO, NO>>
109 Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP<Epetra_CrsMatrix>& /* epAB */) {
110  TEUCHOS_TEST_FOR_EXCEPTION(true, Exceptions::RuntimeError,
111  "Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap cannot be used with Scalar != double, LocalOrdinal != int, GlobalOrdinal != int");
112  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
113 }
114 
115 // specialization for the case of ScalarType=double and LocalOrdinal=GlobalOrdinal=int
116 template <>
117 inline RCP<Xpetra::CrsMatrixWrap<double, int, int, Xpetra::EpetraNode>> Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<double, int, int, Xpetra::EpetraNode>(RCP<Epetra_CrsMatrix>& epAB) {
118  typedef double SC;
119  typedef int LO;
120  typedef int GO;
121  typedef Xpetra::EpetraNode NO;
122 
123  RCP<Xpetra::EpetraCrsMatrixT<GO, NO>> tmpC1 = rcp(new Xpetra::EpetraCrsMatrixT<GO, NO>(epAB));
124  RCP<Xpetra::CrsMatrix<SC, LO, GO, NO>> tmpC2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<SC, LO, GO, NO>>(tmpC1);
125  RCP<Xpetra::CrsMatrixWrap<SC, LO, GO, NO>> tmpC3 = rcp(new Xpetra::CrsMatrixWrap<SC, LO, GO, NO>(tmpC2));
126 
127  return tmpC3;
128 }
129 
130 template <class SC, class LO, class GO, class NO>
131 RCP<Xpetra::MultiVector<SC, LO, GO, NO>>
132 Convert_Epetra_MultiVector_ToXpetra_MultiVector(RCP<Epetra_MultiVector>& epX) {
133  TEUCHOS_TEST_FOR_EXCEPTION(true, Exceptions::RuntimeError,
134  "Convert_Epetra_MultiVector_ToXpetra_MultiVector cannot be used with Scalar != double, LocalOrdinal != int, GlobalOrdinal != int");
135  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
136 }
137 
138 // specialization for the case of ScalarType=double and LocalOrdinal=GlobalOrdinal=int
139 template <>
140 inline RCP<Xpetra::MultiVector<double, int, int, Xpetra::EpetraNode>> Convert_Epetra_MultiVector_ToXpetra_MultiVector<double, int, int, Xpetra::EpetraNode>(RCP<Epetra_MultiVector>& epX) {
141  typedef double SC;
142  typedef int LO;
143  typedef int GO;
144  typedef Xpetra::EpetraNode NO;
145 
146  RCP<Xpetra::MultiVector<SC, LO, GO, NO>> tmp = Xpetra::toXpetra<GO, NO>(epX);
147  return tmp;
148 }
149 
150 #endif
151 
156 template <class Scalar,
157  class LocalOrdinal = int,
158  class GlobalOrdinal = LocalOrdinal,
159  class Node = Tpetra::KokkosClassic::DefaultNode::DefaultNodeType>
160 class IO {
161  private:
162 #undef XPETRA_IO_SHORT
163 #include "Xpetra_UseShortNames.hpp"
164 
165  public:
166 #ifdef HAVE_XPETRA_EPETRA
167  // @{
169  /*static RCP<const Epetra_MultiVector> MV2EpetraMV(RCP<MultiVector> const vec);
170  static RCP< Epetra_MultiVector> MV2NonConstEpetraMV(RCP<MultiVector> vec);
171 
172  static const Epetra_MultiVector& MV2EpetraMV(const MultiVector& vec);
173  static Epetra_MultiVector& MV2NonConstEpetraMV(MultiVector& vec);
174 
175  static RCP<const Epetra_CrsMatrix> Op2EpetraCrs(RCP<const Matrix> Op);
176  static RCP< Epetra_CrsMatrix> Op2NonConstEpetraCrs(RCP<Matrix> Op);
177 
178  static const Epetra_CrsMatrix& Op2EpetraCrs(const Matrix& Op);
179  static Epetra_CrsMatrix& Op2NonConstEpetraCrs(Matrix& Op);*/
180 
181  static const Epetra_Map& Map2EpetraMap(const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>& map) {
182  RCP<const Xpetra::EpetraMapT<GlobalOrdinal, Node>> xeMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal, Node>>(Teuchos::rcpFromRef(map));
183  if (xeMap == Teuchos::null)
184  throw Exceptions::BadCast("Utils::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
185  return xeMap->getEpetra_Map();
186  }
187  // @}
188 #endif
189 
190 #ifdef HAVE_XPETRA_TPETRA
191  // @{
193  /*static RCP<const Tpetra::MultiVector<SC,LO,GO,NO> > MV2TpetraMV(RCP<MultiVector> const vec);
194  static RCP< Tpetra::MultiVector<SC,LO,GO,NO> > MV2NonConstTpetraMV(RCP<MultiVector> vec);
195  static RCP< Tpetra::MultiVector<SC,LO,GO,NO> > MV2NonConstTpetraMV2(MultiVector& vec);
196 
197  static const Tpetra::MultiVector<SC,LO,GO,NO>& MV2TpetraMV(const MultiVector& vec);
198  static Tpetra::MultiVector<SC,LO,GO,NO>& MV2NonConstTpetraMV(MultiVector& vec);
199 
200  static RCP<const Tpetra::CrsMatrix<SC,LO,GO,NO> > Op2TpetraCrs(RCP<const Matrix> Op);
201  static RCP< Tpetra::CrsMatrix<SC,LO,GO,NO> > Op2NonConstTpetraCrs(RCP<Matrix> Op);
202 
203  static const Tpetra::CrsMatrix<SC,LO,GO,NO>& Op2TpetraCrs(const Matrix& Op);
204  static Tpetra::CrsMatrix<SC,LO,GO,NO>& Op2NonConstTpetraCrs(Matrix& Op);
205 
206  static RCP<const Tpetra::RowMatrix<SC,LO,GO,NO> > Op2TpetraRow(RCP<const Matrix> Op);
207  static RCP< Tpetra::RowMatrix<SC,LO,GO,NO> > Op2NonConstTpetraRow(RCP<Matrix> Op);*/
208 
209  static const RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> Map2TpetraMap(const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>& map) {
210  const RCP<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>& tmp_TMap = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>(rcpFromRef(map));
211  if (tmp_TMap == Teuchos::null)
212  throw Exceptions::BadCast("Utils::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
213  return tmp_TMap->getTpetra_Map();
214  }
215 #endif
216 
218 
219 
220  static void Write(const std::string& fileName, const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>& M) {
221  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> tmp_Map = rcpFromRef(M);
222 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
223  const RCP<const Xpetra::EpetraMapT<GlobalOrdinal, Node>>& tmp_EMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal, Node>>(tmp_Map);
224  if (tmp_EMap != Teuchos::null) {
225  int rv = EpetraExt::BlockMapToMatrixMarketFile(fileName.c_str(), tmp_EMap->getEpetra_Map());
226  if (rv != 0)
227  throw Exceptions::RuntimeError("EpetraExt::BlockMapToMatrixMarketFile() return value of " + Teuchos::toString(rv));
228  return;
229  }
230 #endif // HAVE_XPETRA_EPETRAEXT
231 
232 #ifdef HAVE_XPETRA_TPETRA
233  const RCP<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>& tmp_TMap =
234  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>(tmp_Map);
235  if (tmp_TMap != Teuchos::null) {
236  RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> TMap = tmp_TMap->getTpetra_Map();
237  Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>::writeMapFile(fileName, *TMap);
238  return;
239  }
240 #endif // HAVE_XPETRA_TPETRA
241 
242  throw Exceptions::BadCast("Could not cast to EpetraMap or TpetraMap in map writing");
243 
244  } // Write
245 
247  static void Write(const std::string& fileName, const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& vec) {
248  std::string mapfile = "map_" + fileName;
249  Write(mapfile, *(vec.getMap()));
250 
251  RCP<const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmp_Vec = Teuchos::rcpFromRef(vec);
252 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
253  const RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal, Node>>& tmp_EVec = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal, Node>>(tmp_Vec);
254  if (tmp_EVec != Teuchos::null) {
255  int rv = EpetraExt::MultiVectorToMatrixMarketFile(fileName.c_str(), *(tmp_EVec->getEpetra_MultiVector()));
256  if (rv != 0)
257  throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
258  return;
259  }
260 #endif // HAVE_XPETRA_EPETRA
261 
262 #ifdef HAVE_XPETRA_TPETRA
263  const RCP<const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>& tmp_TVec =
264  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmp_Vec);
265  if (tmp_TVec != Teuchos::null) {
266  RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> TVec = tmp_TVec->getTpetra_MultiVector();
267  Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>::writeDenseFile(fileName, TVec);
268  return;
269  }
270 #endif // HAVE_XPETRA_TPETRA
271 
272  throw Exceptions::BadCast("Could not cast to EpetraMultiVector or TpetraMultiVector in multivector writing");
273 
274  } // Write
275 
277  static void Write(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Op, const bool& writeAllMaps = false) {
278  Write("rowmap_" + fileName, *(Op.getRowMap()));
279  if (!Op.getDomainMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps)
280  Write("domainmap_" + fileName, *(Op.getDomainMap()));
281  if (!Op.getRangeMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps)
282  Write("rangemap_" + fileName, *(Op.getRangeMap()));
283  if (!Op.getColMap()->isSameAs(*(Op.getDomainMap())) || writeAllMaps)
284  Write("colmap_" + fileName, *(Op.getColMap()));
285 
286  const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>& crsOp =
287  dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>&>(Op);
288  RCP<const Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmp_CrsMtx = crsOp.getCrsMatrix();
289 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
290  const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal, Node>>& tmp_ECrsMtx = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal, Node>>(tmp_CrsMtx);
291  if (tmp_ECrsMtx != Teuchos::null) {
292  RCP<const Epetra_CrsMatrix> A = tmp_ECrsMtx->getEpetra_CrsMatrix();
293  int rv = EpetraExt::RowMatrixToMatrixMarketFile(fileName.c_str(), *A);
294  if (rv != 0)
295  throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
296  return;
297  }
298 #endif
299 
300 #ifdef HAVE_XPETRA_TPETRA
301  const RCP<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>& tmp_TCrsMtx =
302  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmp_CrsMtx);
303  if (tmp_TCrsMtx != Teuchos::null) {
304  RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A = tmp_TCrsMtx->getTpetra_CrsMatrix();
305  Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>::writeSparseFile(fileName, A);
306  return;
307  }
308  const RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>& tmp_BlockCrs =
309  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmp_CrsMtx);
310  if (tmp_BlockCrs != Teuchos::null) {
311  std::ofstream outstream(fileName, std::ofstream::out);
312  Teuchos::FancyOStream ofs(Teuchos::rcpFromRef(outstream));
313  tmp_BlockCrs->getTpetra_BlockCrsMatrix()->describe(ofs, Teuchos::VERB_EXTREME);
314  return;
315  }
316 
317 #endif // HAVE_XPETRA_TPETRA
318 
319  throw Exceptions::BadCast("Could not cast to EpetraCrsMatrix or TpetraCrsMatrix in matrix writing");
320  } // Write
321 
323  static void WriteLocal(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Op) {
324  const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>& crsOp =
325  dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>&>(Op);
326  RCP<const Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmp_CrsMtx = crsOp.getCrsMatrix();
327 
328  ArrayRCP<const size_t> rowptr_RCP;
329  ArrayRCP<LocalOrdinal> rowptr2_RCP;
330  ArrayRCP<const LocalOrdinal> colind_RCP;
331  ArrayRCP<const Scalar> vals_RCP;
332  tmp_CrsMtx->getAllValues(rowptr_RCP, colind_RCP, vals_RCP);
333 
334  ArrayView<const size_t> rowptr = rowptr_RCP();
335  ArrayView<const LocalOrdinal> colind = colind_RCP();
336  ArrayView<const Scalar> vals = vals_RCP();
337 
338  rowptr2_RCP.resize(rowptr.size());
339  ArrayView<LocalOrdinal> rowptr2 = rowptr2_RCP();
340  for (size_t j = 0; j < rowptr.size(); j++)
341  rowptr2[j] = rowptr[j];
342 
343  Teuchos::MatrixMarket::Raw::Writer<Scalar, LocalOrdinal> writer;
344  writer.writeFile(fileName + "." + std::to_string(Op.getRowMap()->getComm()->getSize()) + "." + std::to_string(Op.getRowMap()->getComm()->getRank()),
345  rowptr2, colind, vals,
346  rowptr.size() - 1, Op.getColMap()->getLocalNumElements());
347  } // WriteLocal
348 
363  static void WriteBlockedCrsMatrix(const std::string& fileName, const Xpetra::BlockedCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Op, const bool& writeAllMaps = false) {
364  using XpIO = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
365 
366  // Write all matrix blocks with their maps
367  for (size_t row = 0; row < Op.Rows(); ++row) {
368  for (size_t col = 0; col < Op.Cols(); ++col) {
369  RCP<const Matrix> m = Op.getMatrix(row, col);
370  if (m != Teuchos::null) { // skip empty blocks
371  TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(m) == Teuchos::null, Xpetra::Exceptions::BadCast,
372  "Sub block matrix (" << row << "," << col << ") is not of type CrsMatrixWrap.");
373  XpIO::Write(fileName + toString(row) + toString(col) + ".m", *m, writeAllMaps);
374  }
375  }
376  }
377 
378  // write map information of map extractors
379  RCP<const MapExtractor> rangeMapExtractor = Op.getRangeMapExtractor();
380  RCP<const MapExtractor> domainMapExtractor = Op.getDomainMapExtractor();
381 
382  for (size_t row = 0; row < rangeMapExtractor->NumMaps(); ++row) {
383  RCP<const Map> map = rangeMapExtractor->getMap(row);
384  XpIO::Write("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", *map);
385  }
386  XpIO::Write("fullRangeMap_" + fileName + ".m", *(rangeMapExtractor->getFullMap()));
387 
388  for (size_t col = 0; col < domainMapExtractor->NumMaps(); ++col) {
389  RCP<const Map> map = domainMapExtractor->getMap(col);
390  XpIO::Write("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", *map);
391  }
392  XpIO::Write("fullDomainMap_" + fileName + ".m", *(domainMapExtractor->getFullMap()));
393  } // WriteBlockedCrsMatrix
394 
396  static Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> Read(const std::string& fileName, Xpetra::UnderlyingLib lib, const RCP<const Teuchos::Comm<int>>& comm, bool binary = false) {
397  if (binary == false) {
398  // Matrix Market file format (ASCII)
399  if (lib == Xpetra::UseEpetra) {
400 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
401  Epetra_CrsMatrix* eA;
402  const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(comm);
403  int rv = EpetraExt::MatrixMarketFileToCrsMatrix(fileName.c_str(), *epcomm, eA);
404  if (rv != 0)
405  throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
406 
407  RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
408 
409  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A =
410  Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
411  return A;
412 #else
413  throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
414 #endif
415  } else if (lib == Xpetra::UseTpetra) {
416 #ifdef HAVE_XPETRA_TPETRA
417  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
418 
419  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
420 
421  bool callFillComplete = true;
422 
423  RCP<sparse_matrix_type> tA = reader_type::readSparseFile(fileName, comm, callFillComplete);
424 
425  if (tA.is_null())
426  throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
427 
428  RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA1 = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tA));
429  RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmpA1);
430  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A = rcp(new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA2));
431 
432  return A;
433 #else
434  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
435 #endif
436  } else {
437  throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
438  }
439  } else {
440  // Custom file format (binary)
441  std::ifstream ifs(fileName.c_str(), std::ios::binary);
442  TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << fileName << "\"");
443  int m, n, nnz;
444  ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
445  ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
446  ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
447 
448  int myRank = comm->getRank();
449 
450  GO indexBase = 0;
451  RCP<Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rowMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, m, (myRank == 0 ? m : 0), indexBase, comm), rangeMap = rowMap;
452  RCP<Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> colMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, n, (myRank == 0 ? n : 0), indexBase, comm), domainMap = colMap;
453  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A;
454 
455  if (myRank == 0) {
456  Teuchos::Array<GlobalOrdinal> inds;
457  Teuchos::Array<Scalar> vals;
458  // Scan matrix to determine the exact nnz per row.
459  Teuchos::ArrayRCP<size_t> numEntriesPerRow(m, (size_t)(0));
460  for (int i = 0; i < m; i++) {
461  int row, rownnz;
462  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
463  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
464  numEntriesPerRow[row] = rownnz;
465  for (int j = 0; j < rownnz; j++) {
466  int index;
467  ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
468  }
469  for (int j = 0; j < rownnz; j++) {
470  double value;
471  ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
472  }
473  }
474 
475  A = Xpetra::MatrixFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Build(rowMap, colMap, numEntriesPerRow);
476 
477  // Now that nnz per row are known, reread and store the matrix.
478  ifs.seekg(0, ifs.beg); // rewind to beginning of file
479  int junk; // skip header info
480  ifs.read(reinterpret_cast<char*>(&junk), sizeof(junk));
481  ifs.read(reinterpret_cast<char*>(&junk), sizeof(junk));
482  ifs.read(reinterpret_cast<char*>(&junk), sizeof(junk));
483  for (int i = 0; i < m; i++) {
484  int row, rownnz;
485  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
486  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
487  inds.resize(rownnz);
488  vals.resize(rownnz);
489  for (int j = 0; j < rownnz; j++) {
490  int index;
491  ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
492  inds[j] = Teuchos::as<GlobalOrdinal>(index);
493  }
494  for (int j = 0; j < rownnz; j++) {
495  double value;
496  ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
497  vals[j] = Teuchos::as<Scalar>(value);
498  }
499  A->insertGlobalValues(row, inds, vals);
500  }
501  } // if (myRank == 0)
502  else {
503  Teuchos::ArrayRCP<size_t> numEntriesPerRow(0, (size_t)(0));
504  A = Xpetra::MatrixFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Build(rowMap, colMap, numEntriesPerRow);
505  }
506 
507  A->fillComplete(domainMap, rangeMap);
508 
509  return A;
510  } // if (binary == false) ... else
511 
512  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
513 
514  } // Read()
515 
520  static Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>
521  Read(const std::string& filename,
522  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rowMap,
523  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> colMap = Teuchos::null,
524  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domainMap = Teuchos::null,
525  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rangeMap = Teuchos::null,
526  const bool callFillComplete = true,
527  const bool binary = false,
528  const bool tolerant = false,
529  const bool debug = false) {
530  TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::Read() : rowMap cannot be null");
531 
532  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domain = (domainMap.is_null() ? rowMap : domainMap);
533  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> range = (rangeMap.is_null() ? rowMap : rangeMap);
534 
535  const Xpetra::UnderlyingLib lib = rowMap->lib();
536  if (binary == false) {
537  if (lib == Xpetra::UseEpetra) {
538 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
539  Epetra_CrsMatrix* eA;
540  const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(rowMap->getComm());
541  const Epetra_Map& epetraRowMap = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2EpetraMap(*rowMap);
542  const Epetra_Map& epetraDomainMap = (domainMap.is_null() ? epetraRowMap : Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2EpetraMap(*domainMap));
543  const Epetra_Map& epetraRangeMap = (rangeMap.is_null() ? epetraRowMap : Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2EpetraMap(*rangeMap));
544  int rv;
545  if (colMap.is_null()) {
546  rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraRangeMap, epetraDomainMap, eA);
547 
548  } else {
549  const Epetra_Map& epetraColMap = Map2EpetraMap(*colMap);
550  rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraColMap, epetraRangeMap, epetraDomainMap, eA);
551  }
552 
553  if (rv != 0)
554  throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
555 
556  RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
557  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A =
558  Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
559 
560  return A;
561 #else
562  throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
563 #endif
564  } else if (lib == Xpetra::UseTpetra) {
565 #ifdef HAVE_XPETRA_TPETRA
566  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
567  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
568  typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
569 
570  const RCP<const map_type> tpetraRowMap = Map2TpetraMap(*rowMap);
571  RCP<const map_type> tpetraColMap = (colMap.is_null() ? Teuchos::null : Map2TpetraMap(*colMap));
572  const RCP<const map_type> tpetraRangeMap = (rangeMap.is_null() ? tpetraRowMap : Map2TpetraMap(*rangeMap));
573  const RCP<const map_type> tpetraDomainMap = (domainMap.is_null() ? tpetraRowMap : Map2TpetraMap(*domainMap));
574 
575  RCP<sparse_matrix_type> tA = reader_type::readSparseFile(filename, tpetraRowMap, tpetraColMap, tpetraDomainMap, tpetraRangeMap,
576  callFillComplete, tolerant, debug);
577  if (tA.is_null())
578  throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
579 
580  RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA1 = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tA));
581  RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmpA1);
582  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A = rcp(new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA2));
583 
584  return A;
585 #else
586  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
587 #endif
588  } else {
589  throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
590  }
591  } else {
592  // Read in on rank 0.
593  auto tempA = Read(filename, lib, rowMap->getComm(), binary);
594 
595  auto A = Xpetra::MatrixFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Build(rowMap, colMap, 0);
596  auto importer = Xpetra::ImportFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(tempA->getRowMap(), rowMap);
597  A->doImport(*tempA, *importer, Xpetra::INSERT);
598  if (callFillComplete)
599  A->fillComplete(domainMap, rangeMap);
600 
601  return A;
602  }
603 
604  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
605  }
606 
615  static Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> ReadLocal(const std::string& filename,
616  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rowMap,
617  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> colMap,
618  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domainMap = Teuchos::null,
619  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rangeMap = Teuchos::null,
620  const bool callFillComplete = true,
621  const bool binary = false,
622  const bool tolerant = false,
623  const bool debug = false) {
624  TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : rowMap cannot be null");
625  TEUCHOS_TEST_FOR_EXCEPTION(colMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : colMap cannot be null");
626 
627  using matrix_type = Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
628  using crs_wrap_type = Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
629  using crs_type = Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
630 
631  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domain = (domainMap.is_null() ? rowMap : domainMap);
632  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> range = (rangeMap.is_null() ? rowMap : rangeMap);
633 
634  std::string rankFilename = filename + "." + std::to_string(rowMap->getComm()->getSize()) + "." + std::to_string(rowMap->getComm()->getRank());
635  RCP<matrix_type> A = rcp(new crs_wrap_type(rowMap, colMap, 0));
636 
637  if (binary == false) {
638  RCP<Teuchos::ParameterList> params = rcp(new Teuchos::ParameterList());
639  params->set("Parse tolerantly", tolerant);
640  params->set("Debug mode", debug);
641 
642  LocalOrdinal numRows = rowMap->getLocalNumElements();
643  LocalOrdinal numCols = colMap->getLocalNumElements();
644 
645  ArrayRCP<LocalOrdinal> rowptr2_RCP;
646  ArrayRCP<LocalOrdinal> colind2_RCP;
647  ArrayRCP<Scalar> vals2_RCP;
648 
649  Teuchos::MatrixMarket::Raw::Reader<Scalar, LocalOrdinal> reader;
650  reader.readFile(rowptr2_RCP, colind2_RCP, vals2_RCP,
651  numRows, numCols,
652  rankFilename);
653 
654  RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
655 
656  ArrayRCP<size_t> rowptr_RCP;
657  ArrayRCP<LocalOrdinal> colind_RCP;
658  ArrayRCP<Scalar> vals_RCP;
659  ACrs->allocateAllValues(colind2_RCP.size(), rowptr_RCP, colind_RCP, vals_RCP);
660 
661  rowptr_RCP.assign(rowptr2_RCP.begin(), rowptr2_RCP.end());
662  colind_RCP = colind2_RCP;
663  vals_RCP = vals2_RCP;
664 
665  ACrs->setAllValues(rowptr_RCP, colind_RCP, vals_RCP);
666  } else {
667  // Custom file format (binary)
668  std::ifstream ifs = std::ifstream(rankFilename.c_str(), std::ios::binary);
669  TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << filename << "\"");
670 
671  int m, n, nnz;
672  ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
673  ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
674  ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
675 
676  TEUCHOS_ASSERT_EQUALITY(Teuchos::as<int>(rowMap->getLocalNumElements()), m);
677 
678  Teuchos::ArrayRCP<size_t> rowptrRCP;
679  Teuchos::ArrayRCP<LocalOrdinal> indicesRCP;
680  Teuchos::ArrayRCP<Scalar> valuesRCP;
681 
682  RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
683 
684  ACrs->allocateAllValues(nnz, rowptrRCP, indicesRCP, valuesRCP);
685 
686  Teuchos::ArrayView<size_t> rowptr = rowptrRCP();
687  Teuchos::ArrayView<LocalOrdinal> indices = indicesRCP();
688  Teuchos::ArrayView<Scalar> values = valuesRCP();
689 
690  bool sorted = true;
691 
692  // Read in rowptr
693  for (int i = 0; i < m; i++) {
694  int row, rownnz;
695  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
696  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
697 
698  rowptr[row + 1] += rownnz;
699  ifs.seekg(sizeof(int) * rownnz + sizeof(double) * rownnz, ifs.cur);
700  }
701  for (int i = 0; i < m; i++)
702  rowptr[i + 1] += rowptr[i];
703  TEUCHOS_ASSERT(Teuchos::as<int>(rowptr[m]) == nnz);
704 
705  // reset to where the data starts
706  ifs.seekg(sizeof(int) * 3, ifs.beg);
707 
708  // read in entries
709  for (int i = 0; i < m; i++) {
710  int row, rownnz;
711  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
712  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
713  size_t ptr = rowptr[row];
714  for (int j = 0; j < rownnz; j++) {
715  int index;
716  ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
717  indices[ptr] = Teuchos::as<LocalOrdinal>(index);
718  if (j > 0)
719  sorted = sorted & (indices[ptr - 1] < indices[ptr]);
720  ++ptr;
721  }
722  ptr = rowptr[row];
723  for (int j = 0; j < rownnz; j++) {
724  double value;
725  ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
726  values[ptr] = Teuchos::as<Scalar>(value);
727  ++ptr;
728  }
729  rowptr[row] += rownnz;
730  }
731  for (int i = m; i > 0; i--)
732  rowptr[i] = rowptr[i - 1];
733  rowptr[0] = 0;
734 
735 #ifdef HAVE_XPETRA_TPETRA
736  if (!sorted) {
737  for (LocalOrdinal lclRow = 0; lclRow < m; lclRow++) {
738  size_t rowBegin = rowptr[lclRow];
739  size_t rowEnd = rowptr[lclRow + 1];
740  Tpetra::sort2(&indices[rowBegin], &indices[rowEnd], &values[rowBegin]);
741  }
742  }
743 #else
744  TEUCHOS_ASSERT(sorted);
745 #endif
746 
747  ACrs->setAllValues(rowptrRCP, indicesRCP, valuesRCP);
748  }
749 
750  if (callFillComplete)
751  A->fillComplete(domainMap, rangeMap);
752  return A;
753  }
755 
756  static RCP<MultiVector> ReadMultiVector(const std::string& fileName,
757  const RCP<const Map>& map,
758  const bool binary = false) {
759  Xpetra::UnderlyingLib lib = map->lib();
760 
761  if (lib == Xpetra::UseEpetra) {
762  TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
763 
764  } else if (lib == Xpetra::UseTpetra) {
765 #ifdef HAVE_XPETRA_TPETRA
766  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
767  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
768  typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
769  typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> multivector_type;
770 
771  RCP<const map_type> temp = toTpetra(map);
772  RCP<multivector_type> TMV = reader_type::readDenseFile(fileName, map->getComm(), temp, false, false, binary);
773  RCP<MultiVector> rmv = Xpetra::toXpetra(TMV);
774  return rmv;
775 #else
776  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
777 #endif
778  } else {
779  throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
780  }
781 
782  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
783  }
784 
785  static RCP<const Map> ReadMap(const std::string& fileName,
786  Xpetra::UnderlyingLib lib,
787  const RCP<const Teuchos::Comm<int>>& comm,
788  const bool binary = false) {
789  if (lib == Xpetra::UseEpetra) {
790  TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
791  } else if (lib == Xpetra::UseTpetra) {
792 #ifdef HAVE_XPETRA_TPETRA
793  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
794  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
795 
796  RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> tMap = reader_type::readMapFile(fileName, comm, false, false, binary);
797  if (tMap.is_null())
798  throw Exceptions::RuntimeError("The Tpetra::Map returned from readSparseFile() is null.");
799 
800  return Xpetra::toXpetra(tMap);
801 #else
802  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
803 #endif
804  } else {
805  throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
806  }
807 
808  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
809  }
810 
824  static RCP<const Xpetra::BlockedCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> ReadBlockedCrsMatrix(const std::string& fileName, Xpetra::UnderlyingLib lib, const RCP<const Teuchos::Comm<int>>& comm) {
825  using XpIO = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
826 
827  size_t numBlocks = 2; // TODO user parameter?
828 
829  std::vector<RCP<const Map>> rangeMapVec;
830  for (size_t row = 0; row < numBlocks; ++row) {
831  RCP<const Map> map = XpIO::ReadMap("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", lib, comm);
832  rangeMapVec.push_back(map);
833  }
834  RCP<const Map> fullRangeMap = XpIO::ReadMap("fullRangeMap_" + fileName + ".m", lib, comm);
835 
836  std::vector<RCP<const Map>> domainMapVec;
837  for (size_t col = 0; col < numBlocks; ++col) {
838  RCP<const Map> map = XpIO::ReadMap("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", lib, comm);
839  domainMapVec.push_back(map);
840  }
841  RCP<const Map> fullDomainMap = XpIO::ReadMap("fullDomainMap_" + fileName + ".m", lib, comm);
842 
843  /*std::vector<RCP<const XpMap> > testRgMapVec;
844  for(size_t r = 0; r < numBlocks; ++r) {
845  RCP<const XpMap> map = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(r) + "0.m", lib, comm);
846  testRgMapVec.push_back(map);
847  }
848  std::vector<RCP<const XpMap> > testDoMapVec;
849  for(size_t c = 0; c < numBlocks; ++c) {
850  RCP<const XpMap> map = XpIO::ReadMap("domainmap_" + fileName + "0" + XpIO::toString<size_t>(c) + ".m", lib, comm);
851  testDoMapVec.push_back(map);
852  }*/
853 
854  // create map extractors
855 
856  // range map extractor
857  bool bRangeUseThyraStyleNumbering = false;
858  /*GlobalOrdinal gMinGids = 0;
859  for(size_t v = 0; v < testRgMapVec.size(); ++v) {
860  gMinGids += testRgMapVec[v]->getMinAllGlobalIndex();
861  }
862  if ( gMinGids==0 && testRgMapVec.size() > 1 ) bRangeUseThyraStyleNumbering = true;
863  */
864  RCP<const MapExtractor> rangeMapExtractor =
865  Teuchos::rcp(new MapExtractor(fullRangeMap, rangeMapVec, bRangeUseThyraStyleNumbering));
866 
867  // domain map extractor
868  bool bDomainUseThyraStyleNumbering = false;
869  /*gMinGids = 0;
870  for(size_t v = 0; v < testDoMapVec.size(); ++v) {
871  gMinGids += testDoMapVec[v]->getMinAllGlobalIndex();
872  }
873  if ( gMinGids==0 && testDoMapVec.size() > 1 ) bDomainUseThyraStyleNumbering = true;
874  */
875  RCP<const MapExtractor> domainMapExtractor =
876  Teuchos::rcp(new MapExtractor(fullDomainMap, domainMapVec, bDomainUseThyraStyleNumbering));
877 
878  RCP<BlockedCrsMatrix> bOp = Teuchos::rcp(new BlockedCrsMatrix(rangeMapExtractor, domainMapExtractor, 33));
879 
880  // Read all sub-matrices with their maps and place into blocked operator
881  for (size_t row = 0; row < numBlocks; ++row) {
882  for (size_t col = 0; col < numBlocks; ++col) {
883  RCP<const Map> rowSubMap = XpIO::ReadMap("rowmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
884  RCP<const Map> colSubMap = XpIO::ReadMap("colmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
885  RCP<const Map> domSubMap = XpIO::ReadMap("domainmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
886  RCP<const Map> ranSubMap = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
887  RCP<Matrix> mat = XpIO::Read(fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", rowSubMap, colSubMap, domSubMap, ranSubMap);
888  bOp->setMatrix(row, col, mat);
889  }
890  }
891 
892  bOp->fillComplete();
893 
894  return bOp;
895  } // ReadBlockedCrsMatrix
896 
898  template <class T>
899  static std::string toString(const T& what) {
900  std::ostringstream buf;
901  buf << what;
902  return buf.str();
903  }
904 };
905 
906 #ifdef HAVE_XPETRA_EPETRA
907 
916 template <class Scalar>
917 class IO<Scalar, int, int, EpetraNode> {
918  public:
919  typedef int LocalOrdinal;
920  typedef int GlobalOrdinal;
921  typedef EpetraNode Node;
922 
923 #ifdef HAVE_XPETRA_EPETRA
924  // @{
926  static const Epetra_Map& Map2EpetraMap(const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>& map) {
927  RCP<const Xpetra::EpetraMapT<GlobalOrdinal, Node>> xeMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal, Node>>(Teuchos::rcpFromRef(map));
928  if (xeMap == Teuchos::null)
929  throw Exceptions::BadCast("IO::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
930  return xeMap->getEpetra_Map();
931  }
932  // @}
933 #endif
934 
935 #ifdef HAVE_XPETRA_TPETRA
936  // @{
938  static const RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> Map2TpetraMap(const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>& map) {
939  const RCP<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>& tmp_TMap = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>(rcpFromRef(map));
940  if (tmp_TMap == Teuchos::null)
941  throw Exceptions::BadCast("IO::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
942  return tmp_TMap->getTpetra_Map();
943  }
944 #endif
945 
947 
948 
949  static void Write(const std::string& fileName, const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>& M) {
950  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> tmp_Map = rcpFromRef(M);
951 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
952  const RCP<const Xpetra::EpetraMapT<GlobalOrdinal, Node>>& tmp_EMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal, Node>>(tmp_Map);
953  if (tmp_EMap != Teuchos::null) {
954  int rv = EpetraExt::BlockMapToMatrixMarketFile(fileName.c_str(), tmp_EMap->getEpetra_Map());
955  if (rv != 0)
956  throw Exceptions::RuntimeError("EpetraExt::BlockMapToMatrixMarketFile() return value of " + Teuchos::toString(rv));
957  return;
958  }
959 #endif // HAVE_XPETRA_EPETRA
960 
961 #ifdef HAVE_XPETRA_TPETRA
962 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
963  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
964  // do nothing
965 #else
966  const RCP<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>& tmp_TMap =
967  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node>>(tmp_Map);
968  if (tmp_TMap != Teuchos::null) {
969  RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> TMap = tmp_TMap->getTpetra_Map();
970  Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>::writeMapFile(fileName, *TMap);
971  return;
972  }
973 #endif
974 #endif // HAVE_XPETRA_TPETRA
975  throw Exceptions::BadCast("Could not cast to EpetraMap or TpetraMap in map writing");
976  }
977 
979  static void Write(const std::string& fileName, const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& vec) {
980  std::string mapfile = "map_" + fileName;
981  Write(mapfile, *(vec.getMap()));
982 
983  RCP<const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmp_Vec = Teuchos::rcpFromRef(vec);
984 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
985  const RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal, Node>>& tmp_EVec = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal, Node>>(tmp_Vec);
986  if (tmp_EVec != Teuchos::null) {
987  int rv = EpetraExt::MultiVectorToMatrixMarketFile(fileName.c_str(), *(tmp_EVec->getEpetra_MultiVector()));
988  if (rv != 0)
989  throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
990  return;
991  }
992 #endif // HAVE_XPETRA_EPETRAEXT
993 
994 #ifdef HAVE_XPETRA_TPETRA
995 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
996  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
997  // do nothin
998 #else
999  const RCP<const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>& tmp_TVec =
1000  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmp_Vec);
1001  if (tmp_TVec != Teuchos::null) {
1002  RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> TVec = tmp_TVec->getTpetra_MultiVector();
1003  Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>::writeDenseFile(fileName, TVec);
1004  return;
1005  }
1006 #endif
1007 #endif // HAVE_XPETRA_TPETRA
1008 
1009  throw Exceptions::BadCast("Could not cast to EpetraMultiVector or TpetraMultiVector in multivector writing");
1010 
1011  } // Write
1012 
1014  static void Write(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Op, const bool& writeAllMaps = false) {
1015  Write("rowmap_" + fileName, *(Op.getRowMap()));
1016  if (!Op.getDomainMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps)
1017  Write("domainmap_" + fileName, *(Op.getDomainMap()));
1018  if (!Op.getRangeMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps)
1019  Write("rangemap_" + fileName, *(Op.getRangeMap()));
1020  if (!Op.getColMap()->isSameAs(*(Op.getDomainMap())) || writeAllMaps)
1021  Write("colmap_" + fileName, *(Op.getColMap()));
1022 
1023  const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>& crsOp =
1024  dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>&>(Op);
1025  RCP<const Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmp_CrsMtx = crsOp.getCrsMatrix();
1026 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1027  const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal, Node>>& tmp_ECrsMtx = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal, Node>>(tmp_CrsMtx);
1028  if (tmp_ECrsMtx != Teuchos::null) {
1029  RCP<const Epetra_CrsMatrix> A = tmp_ECrsMtx->getEpetra_CrsMatrix();
1030  int rv = EpetraExt::RowMatrixToMatrixMarketFile(fileName.c_str(), *A);
1031  if (rv != 0)
1032  throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
1033  return;
1034  }
1035 #endif // endif HAVE_XPETRA_EPETRA
1036 
1037 #ifdef HAVE_XPETRA_TPETRA
1038 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1039  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1040  // do nothin
1041 #else
1042  const RCP<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>& tmp_TCrsMtx =
1043  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmp_CrsMtx);
1044  if (tmp_TCrsMtx != Teuchos::null) {
1045  RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A = tmp_TCrsMtx->getTpetra_CrsMatrix();
1046  Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>::writeSparseFile(fileName, A);
1047  return;
1048  }
1049  const RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>& tmp_BlockCrs =
1050  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmp_CrsMtx);
1051  if (tmp_BlockCrs != Teuchos::null) {
1052  std::ofstream outstream(fileName, std::ofstream::out);
1053  Teuchos::FancyOStream ofs(Teuchos::rcpFromRef(outstream));
1054  tmp_BlockCrs->getTpetra_BlockCrsMatrix()->describe(ofs, Teuchos::VERB_EXTREME);
1055  return;
1056  }
1057 
1058 #endif
1059 #endif // HAVE_XPETRA_TPETRA
1060 
1061  throw Exceptions::BadCast("Could not cast to EpetraCrsMatrix or TpetraCrsMatrix in matrix writing");
1062  } // Write
1063 
1065  static void Write(const std::string& fileName, const Xpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node>& graph, const bool& writeAllMaps = false) {
1066  Write("rowmap_" + fileName, *(graph.getRowMap()));
1067  if (!graph.getDomainMap()->isSameAs(*(graph.getRowMap())) || writeAllMaps)
1068  Write("domainmap_" + fileName, *(graph.getDomainMap()));
1069  if (!graph.getRangeMap()->isSameAs(*(graph.getRowMap())) || writeAllMaps)
1070  Write("rangemap_" + fileName, *(graph.getRangeMap()));
1071  if (!graph.getColMap()->isSameAs(*(graph.getDomainMap())) || writeAllMaps)
1072  Write("colmap_" + fileName, *(graph.getColMap()));
1073 
1074  RCP<const Xpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node>> tmp_Graph = rcpFromRef(graph);
1075 
1076 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1077  const RCP<const Xpetra::EpetraCrsGraphT<GlobalOrdinal, Node>>& tmp_ECrsGraph = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraCrsGraphT<GlobalOrdinal, Node>>(tmp_Graph);
1078  if (tmp_ECrsGraph != Teuchos::null) {
1079  throw Exceptions::BadCast("Writing not implemented for EpetraCrsGraphT");
1080  }
1081 #endif // endif HAVE_XPETRA_EPETRA
1082 
1083 #ifdef HAVE_XPETRA_TPETRA
1084 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1085  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1086  // do nothin
1087 #else
1088  RCP<const Xpetra::TpetraCrsGraph<LocalOrdinal, GlobalOrdinal, Node>> tmp_TCrsGraph =
1089  Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsGraph<LocalOrdinal, GlobalOrdinal, Node>>(tmp_Graph);
1090  if (tmp_TCrsGraph != Teuchos::null) {
1091  RCP<const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node>> G = tmp_TCrsGraph->getTpetra_CrsGraph();
1092  Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>::writeSparseGraphFile(fileName, G);
1093  return;
1094  }
1095 #endif
1096 #endif // HAVE_XPETRA_TPETRA
1097 
1098  throw Exceptions::BadCast("Could not cast to EpetraCrsMatrix or TpetraCrsMatrix in matrix writing");
1099  } // Write
1100 
1102  static void WriteLocal(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Op) {
1103  const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>& crsOp =
1104  dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>&>(Op);
1105  RCP<const Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmp_CrsMtx = crsOp.getCrsMatrix();
1106 
1107  ArrayRCP<const size_t> rowptr_RCP;
1108  ArrayRCP<LocalOrdinal> rowptr2_RCP;
1109  ArrayRCP<const LocalOrdinal> colind_RCP;
1110  ArrayRCP<const Scalar> vals_RCP;
1111  tmp_CrsMtx->getAllValues(rowptr_RCP, colind_RCP, vals_RCP);
1112 
1113  ArrayView<const size_t> rowptr = rowptr_RCP();
1114  ArrayView<const LocalOrdinal> colind = colind_RCP();
1115  ArrayView<const Scalar> vals = vals_RCP();
1116 
1117  rowptr2_RCP.resize(rowptr.size());
1118  ArrayView<LocalOrdinal> rowptr2 = rowptr2_RCP();
1119  for (size_t j = 0; j < Teuchos::as<size_t>(rowptr.size()); j++)
1120  rowptr2[j] = rowptr[j];
1121 
1122  Teuchos::MatrixMarket::Raw::Writer<Scalar, LocalOrdinal> writer;
1123  writer.writeFile(fileName + "." + std::to_string(Op.getRowMap()->getComm()->getSize()) + "." + std::to_string(Op.getRowMap()->getComm()->getRank()),
1124  rowptr2, colind, vals,
1125  rowptr.size() - 1, Op.getColMap()->getLocalNumElements());
1126  } // WriteLocal
1127 
1144  static void WriteBlockedCrsMatrix(const std::string& fileName, const Xpetra::BlockedCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Op, const bool& writeAllMaps = false) {
1145 #include "Xpetra_UseShortNames.hpp"
1146  using XpIO = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
1147 
1148  // write all matrices with their maps
1149  for (size_t row = 0; row < Op.Rows(); ++row) {
1150  for (size_t col = 0; col < Op.Cols(); ++col) {
1151  RCP<const Matrix> m = Op.getMatrix(row, col);
1152  if (m != Teuchos::null) { // skip empty blocks
1153  TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(m) == Teuchos::null, Xpetra::Exceptions::BadCast,
1154  "Sub block matrix (" << row << "," << col << ") is not of type CrsMatrixWrap.");
1155  XpIO::Write(fileName + toString(row) + toString(col) + ".m", *m, writeAllMaps);
1156  }
1157  }
1158  }
1159 
1160  // write map information of map extractors
1161  RCP<const MapExtractor> rangeMapExtractor = Op.getRangeMapExtractor();
1162  RCP<const MapExtractor> domainMapExtractor = Op.getDomainMapExtractor();
1163 
1164  for (size_t row = 0; row < rangeMapExtractor->NumMaps(); ++row) {
1165  RCP<const Map> map = rangeMapExtractor->getMap(row);
1166  XpIO::Write("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", *map);
1167  }
1168  XpIO::Write("fullRangeMap_" + fileName + ".m", *(rangeMapExtractor->getFullMap()));
1169 
1170  for (size_t col = 0; col < domainMapExtractor->NumMaps(); ++col) {
1171  RCP<const Map> map = domainMapExtractor->getMap(col);
1172  XpIO::Write("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", *map);
1173  }
1174  XpIO::Write("fullDomainMap_" + fileName + ".m", *(domainMapExtractor->getFullMap()));
1175  } // WriteBlockedCrsMatrix
1176 
1178  static Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> Read(const std::string& fileName, Xpetra::UnderlyingLib lib, const RCP<const Teuchos::Comm<int>>& comm, bool binary = false) {
1179  if (binary == false) {
1180  // Matrix Market file format (ASCII)
1181  if (lib == Xpetra::UseEpetra) {
1182 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1183  Epetra_CrsMatrix* eA;
1184  const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(comm);
1185  int rv = EpetraExt::MatrixMarketFileToCrsMatrix(fileName.c_str(), *epcomm, eA);
1186  if (rv != 0)
1187  throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
1188 
1189  RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
1190 
1191  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A =
1192  Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
1193  return A;
1194 #else
1195  throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1196 #endif
1197  } else if (lib == Xpetra::UseTpetra) {
1198 #ifdef HAVE_XPETRA_TPETRA
1199 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1200  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1201  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int enabled.");
1202 #else
1203  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1204 
1205  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1206 
1207  bool callFillComplete = true;
1208 
1209  RCP<sparse_matrix_type> tA = reader_type::readSparseFile(fileName, comm, callFillComplete);
1210 
1211  if (tA.is_null())
1212  throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
1213 
1214  RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA1 = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tA));
1215  RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmpA1);
1216  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A = rcp(new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA2));
1217 
1218  return A;
1219 #endif
1220 #else
1221  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1222 #endif
1223  } else {
1224  throw Exceptions::RuntimeError("Xpetra:IO: you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1225  }
1226  } else {
1227  // Custom file format (binary)
1228  std::ifstream ifs(fileName.c_str(), std::ios::binary);
1229  TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << fileName << "\"");
1230  int m, n, nnz;
1231  ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
1232  ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
1233  ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
1234 
1235  int myRank = comm->getRank();
1236 
1237  GlobalOrdinal indexBase = 0;
1238  RCP<Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rowMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, m, (myRank == 0 ? m : 0), indexBase, comm), rangeMap = rowMap;
1239  RCP<Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> colMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, n, (myRank == 0 ? n : 0), indexBase, comm), domainMap = colMap;
1240 
1241  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A;
1242 
1243  if (myRank == 0) {
1244  Teuchos::Array<GlobalOrdinal> inds;
1245  Teuchos::Array<Scalar> vals;
1246  // Scan matrix to determine the exact nnz per row.
1247  Teuchos::ArrayRCP<size_t> numEntriesPerRow(m);
1248  for (int i = 0; i < m; i++) {
1249  int row, rownnz;
1250  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1251  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1252  numEntriesPerRow[i] = rownnz;
1253  for (int j = 0; j < rownnz; j++) {
1254  int index;
1255  ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
1256  }
1257  for (int j = 0; j < rownnz; j++) {
1258  double value;
1259  ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
1260  }
1261  }
1262 
1263  A = Xpetra::MatrixFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Build(rowMap, colMap, numEntriesPerRow);
1264 
1265  // Now that nnz per row are known, reread and store the matrix.
1266  ifs.seekg(0, ifs.beg); // rewind to beginning of file
1267  int junk; // skip header info
1268  ifs.read(reinterpret_cast<char*>(&m), sizeof(junk));
1269  ifs.read(reinterpret_cast<char*>(&n), sizeof(junk));
1270  ifs.read(reinterpret_cast<char*>(&nnz), sizeof(junk));
1271  for (int i = 0; i < m; i++) {
1272  int row, rownnz;
1273  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1274  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1275  inds.resize(rownnz);
1276  vals.resize(rownnz);
1277  for (int j = 0; j < rownnz; j++) {
1278  int index;
1279  ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
1280  inds[j] = Teuchos::as<GlobalOrdinal>(index);
1281  }
1282  for (int j = 0; j < rownnz; j++) {
1283  double value;
1284  ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
1285  vals[j] = Teuchos::as<Scalar>(value);
1286  }
1287  A->insertGlobalValues(row, inds, vals);
1288  }
1289  } // if (myRank == 0)
1290  else {
1291  Teuchos::ArrayRCP<size_t> numEntriesPerRow(0);
1292  A = Xpetra::MatrixFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Build(rowMap, colMap, numEntriesPerRow);
1293  }
1294 
1295  A->fillComplete(domainMap, rangeMap);
1296 
1297  return A;
1298  }
1299 
1300  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1301 
1302  } // Read()
1303 
1308  static Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> Read(const std::string& filename,
1309  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rowMap,
1310  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> colMap = Teuchos::null,
1311  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domainMap = Teuchos::null,
1312  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rangeMap = Teuchos::null,
1313  const bool callFillComplete = true,
1314  const bool binary = false,
1315  const bool tolerant = false,
1316  const bool debug = false) {
1317  TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::Read() : rowMap cannot be null");
1318 
1319  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domain = (domainMap.is_null() ? rowMap : domainMap);
1320  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> range = (rangeMap.is_null() ? rowMap : rangeMap);
1321 
1322  const Xpetra::UnderlyingLib lib = rowMap->lib();
1323  if (binary == false) {
1324  if (lib == Xpetra::UseEpetra) {
1325 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1326  Epetra_CrsMatrix* eA;
1327  const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(rowMap->getComm());
1328  const Epetra_Map& epetraRowMap = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2EpetraMap(*rowMap);
1329  const Epetra_Map& epetraDomainMap = (domainMap.is_null() ? epetraRowMap : Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2EpetraMap(*domainMap));
1330  const Epetra_Map& epetraRangeMap = (rangeMap.is_null() ? epetraRowMap : Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2EpetraMap(*rangeMap));
1331  int rv;
1332  if (colMap.is_null()) {
1333  rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraRangeMap, epetraDomainMap, eA);
1334 
1335  } else {
1336  const Epetra_Map& epetraColMap = Map2EpetraMap(*colMap);
1337  rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraColMap, epetraRangeMap, epetraDomainMap, eA);
1338  }
1339 
1340  if (rv != 0)
1341  throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
1342 
1343  RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
1344  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A =
1345  Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
1346 
1347  return A;
1348 #else
1349  throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1350 #endif
1351  } else if (lib == Xpetra::UseTpetra) {
1352 #ifdef HAVE_XPETRA_TPETRA
1353 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1354  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1355  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int support.");
1356 #else
1357  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1358  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1359  typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
1360 
1361  const RCP<const map_type> tpetraRowMap = Map2TpetraMap(*rowMap);
1362  RCP<const map_type> tpetraColMap = (colMap.is_null() ? Teuchos::null : Map2TpetraMap(*colMap));
1363  const RCP<const map_type> tpetraRangeMap = (rangeMap.is_null() ? tpetraRowMap : Map2TpetraMap(*rangeMap));
1364  const RCP<const map_type> tpetraDomainMap = (domainMap.is_null() ? tpetraRowMap : Map2TpetraMap(*domainMap));
1365 
1366  RCP<sparse_matrix_type> tA = reader_type::readSparseFile(filename, tpetraRowMap, tpetraColMap, tpetraDomainMap, tpetraRangeMap,
1367  callFillComplete, tolerant, debug);
1368  if (tA.is_null())
1369  throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
1370 
1371  RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA1 = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tA));
1372  RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>>(tmpA1);
1373  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> A = rcp(new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA2));
1374 
1375  return A;
1376 #endif
1377 #else
1378  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1379 #endif
1380  } else {
1381  throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1382  }
1383  } else {
1384  // Read in on rank 0.
1385  auto tempA = Read(filename, lib, rowMap->getComm(), binary);
1386 
1387  auto A = Xpetra::MatrixFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Build(rowMap, colMap, 0);
1388  auto importer = Xpetra::ImportFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(tempA->getRowMap(), rowMap);
1389  A->doImport(*tempA, *importer, Xpetra::INSERT);
1390  if (callFillComplete)
1391  A->fillComplete(domainMap, rangeMap);
1392 
1393  return A;
1394  }
1395 
1396  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1397  }
1398 
1407  static Teuchos::RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> ReadLocal(const std::string& filename,
1408  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rowMap,
1409  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> colMap,
1410  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domainMap = Teuchos::null,
1411  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> rangeMap = Teuchos::null,
1412  const bool callFillComplete = true,
1413  const bool binary = false,
1414  const bool tolerant = false,
1415  const bool debug = false) {
1416  TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : rowMap cannot be null");
1417  TEUCHOS_TEST_FOR_EXCEPTION(colMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : colMap cannot be null");
1418 
1419  using matrix_type = Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
1420  using crs_wrap_type = Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
1421  using crs_type = Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
1422 
1423  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> domain = (domainMap.is_null() ? rowMap : domainMap);
1424  RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> range = (rangeMap.is_null() ? rowMap : rangeMap);
1425 
1426  std::string rankFilename = filename + "." + std::to_string(rowMap->getComm()->getSize()) + "." + std::to_string(rowMap->getComm()->getRank());
1427  RCP<matrix_type> A = rcp(new crs_wrap_type(rowMap, colMap, 0));
1428 
1429  if (binary == false) {
1430  RCP<Teuchos::ParameterList> params = rcp(new Teuchos::ParameterList());
1431  params->set("Parse tolerantly", tolerant);
1432  params->set("Debug mode", debug);
1433 
1434  LocalOrdinal numRows = rowMap->getLocalNumElements();
1435  LocalOrdinal numCols = colMap->getLocalNumElements();
1436 
1437  ArrayRCP<LocalOrdinal> rowptr2_RCP;
1438  ArrayRCP<LocalOrdinal> colind2_RCP;
1439  ArrayRCP<Scalar> vals2_RCP;
1440 
1441  Teuchos::MatrixMarket::Raw::Reader<Scalar, LocalOrdinal> reader;
1442  reader.readFile(rowptr2_RCP, colind2_RCP, vals2_RCP,
1443  numRows, numCols,
1444  rankFilename);
1445 
1446  RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
1447 
1448  ArrayRCP<size_t> rowptr_RCP;
1449  ArrayRCP<LocalOrdinal> colind_RCP;
1450  ArrayRCP<Scalar> vals_RCP;
1451  ACrs->allocateAllValues(colind2_RCP.size(), rowptr_RCP, colind_RCP, vals_RCP);
1452 
1453  rowptr_RCP.assign(rowptr2_RCP.begin(), rowptr2_RCP.end());
1454  colind_RCP = colind2_RCP;
1455  vals_RCP = vals2_RCP;
1456 
1457  ACrs->setAllValues(rowptr_RCP, colind_RCP, vals_RCP);
1458  } else {
1459  // Custom file format (binary)
1460  std::ifstream ifs = std::ifstream(rankFilename.c_str(), std::ios::binary);
1461  TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << filename << "\"");
1462 
1463  int m, n, nnz;
1464  ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
1465  ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
1466  ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
1467 
1468  TEUCHOS_ASSERT_EQUALITY(Teuchos::as<int>(rowMap->getLocalNumElements()), m);
1469 
1470  Teuchos::ArrayRCP<size_t> rowptrRCP;
1471  Teuchos::ArrayRCP<LocalOrdinal> indicesRCP;
1472  Teuchos::ArrayRCP<Scalar> valuesRCP;
1473 
1474  RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
1475 
1476  ACrs->allocateAllValues(nnz, rowptrRCP, indicesRCP, valuesRCP);
1477 
1478  Teuchos::ArrayView<size_t> rowptr = rowptrRCP();
1479  Teuchos::ArrayView<LocalOrdinal> indices = indicesRCP();
1480  Teuchos::ArrayView<Scalar> values = valuesRCP();
1481 
1482  bool sorted = true;
1483 
1484  // Read in rowptr
1485  for (int i = 0; i < m; i++) {
1486  int row, rownnz;
1487  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1488  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1489 
1490  rowptr[row + 1] += rownnz;
1491  ifs.seekg(sizeof(int) * rownnz + sizeof(double) * rownnz, ifs.cur);
1492  }
1493  for (int i = 0; i < m; i++)
1494  rowptr[i + 1] += rowptr[i];
1495  TEUCHOS_ASSERT(Teuchos::as<int>(rowptr[m]) == nnz);
1496 
1497  // reset to where the data starts
1498  ifs.seekg(sizeof(int) * 3, ifs.beg);
1499 
1500  // read in entries
1501  for (int i = 0; i < m; i++) {
1502  int row, rownnz;
1503  ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1504  ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1505  size_t ptr = rowptr[row];
1506  for (int j = 0; j < rownnz; j++) {
1507  int index;
1508  ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
1509  indices[ptr] = Teuchos::as<LocalOrdinal>(index);
1510  if (j > 0)
1511  sorted = sorted & (indices[ptr - 1] < indices[ptr]);
1512  ++ptr;
1513  }
1514  ptr = rowptr[row];
1515  for (int j = 0; j < rownnz; j++) {
1516  double value;
1517  ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
1518  values[ptr] = Teuchos::as<Scalar>(value);
1519  ++ptr;
1520  }
1521  rowptr[row] += rownnz;
1522  }
1523  for (int i = m; i > 0; i--)
1524  rowptr[i] = rowptr[i - 1];
1525  rowptr[0] = 0;
1526 
1527 #ifdef HAVE_XPETRA_TPETRA
1528  if (!sorted) {
1529  for (LocalOrdinal lclRow = 0; lclRow < m; lclRow++) {
1530  size_t rowBegin = rowptr[lclRow];
1531  size_t rowEnd = rowptr[lclRow + 1];
1532  Tpetra::sort2(&indices[rowBegin], &indices[rowEnd], &values[rowBegin]);
1533  }
1534  }
1535 #else
1536  TEUCHOS_ASSERT(sorted);
1537 #endif
1538 
1539  ACrs->setAllValues(rowptrRCP, indicesRCP, valuesRCP);
1540  }
1541 
1542  if (callFillComplete)
1543  A->fillComplete(domainMap, rangeMap);
1544  return A;
1545  }
1547 
1548  static RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> ReadMultiVector(const std::string& fileName,
1549  const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>>& map,
1550  const bool binary = false) {
1551  Xpetra::UnderlyingLib lib = map->lib();
1552 
1553  if (lib == Xpetra::UseEpetra) {
1554  // taw: Oct 9 2015: do we need a specialization for <double,int,int>??
1555  // TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
1556 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1557  TEUCHOS_ASSERT(!binary);
1558  Epetra_MultiVector* MV;
1559  int rv = EpetraExt::MatrixMarketFileToMultiVector(fileName.c_str(), toEpetra(map), MV);
1560  if (rv != 0) throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToMultiVector failed");
1561  RCP<Epetra_MultiVector> MVrcp = rcp(MV);
1562  return Convert_Epetra_MultiVector_ToXpetra_MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(MVrcp);
1563 #else
1564  throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1565 #endif
1566  } else if (lib == Xpetra::UseTpetra) {
1567 #ifdef HAVE_XPETRA_TPETRA
1568 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1569  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1570  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int support.");
1571 #else
1572  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1573  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1574  typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
1575  typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> multivector_type;
1576 
1577  RCP<const map_type> temp = toTpetra(map);
1578  RCP<multivector_type> TMV = reader_type::readDenseFile(fileName, map->getComm(), temp, false, false, binary);
1579  RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>> rmv = Xpetra::toXpetra(TMV);
1580  return rmv;
1581 #endif
1582 #else
1583  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1584 #endif
1585  } else {
1586  throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1587  }
1588 
1589  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1590  }
1591 
1592  static RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> ReadMap(const std::string& fileName,
1593  Xpetra::UnderlyingLib lib,
1594  const RCP<const Teuchos::Comm<int>>& comm,
1595  const bool binary = false) {
1596  if (lib == Xpetra::UseEpetra) {
1597  // do we need another specialization for <double,int,int> ??
1598  // TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
1599 #if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1600  TEUCHOS_ASSERT(!binary);
1601  Epetra_Map* eMap;
1602  int rv = EpetraExt::MatrixMarketFileToMap(fileName.c_str(), *(Xpetra::toEpetra(comm)), eMap);
1603  if (rv != 0)
1604  throw Exceptions::RuntimeError("Error reading map from file " + fileName + " with EpetraExt::MatrixMarketToMap (returned " + Teuchos::toString(rv) + ")");
1605 
1606  RCP<Epetra_Map> eMap1 = rcp(new Epetra_Map(*eMap));
1607  return Xpetra::toXpetra<int, Node>(*eMap1);
1608 #else
1609  throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1610 #endif
1611  } else if (lib == Xpetra::UseTpetra) {
1612 #ifdef HAVE_XPETRA_TPETRA
1613 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1614  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1615  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int support.");
1616 #else
1617  typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1618  typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1619 
1620  RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> tMap = reader_type::readMapFile(fileName, comm, false, false, binary);
1621  if (tMap.is_null())
1622  throw Exceptions::RuntimeError("The Tpetra::Map returned from readSparseFile() is null.");
1623 
1624  return Xpetra::toXpetra(tMap);
1625 #endif
1626 #else
1627  throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1628 #endif
1629  } else {
1630  throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1631  }
1632 
1633  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1634  }
1635 
1651  static RCP<const Xpetra::BlockedCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> ReadBlockedCrsMatrix(const std::string& fileName, Xpetra::UnderlyingLib lib, const RCP<const Teuchos::Comm<int>>& comm) {
1652 #include "Xpetra_UseShortNames.hpp"
1653  using XpIO = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
1654 
1655  size_t numBlocks = 2; // TODO user parameter?
1656 
1657  std::vector<RCP<const Map>> rangeMapVec;
1658  for (size_t row = 0; row < numBlocks; ++row) {
1659  RCP<const Map> map = XpIO::ReadMap("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", lib, comm);
1660  rangeMapVec.push_back(map);
1661  }
1662  RCP<const Map> fullRangeMap = XpIO::ReadMap("fullRangeMap_" + fileName + ".m", lib, comm);
1663 
1664  std::vector<RCP<const Map>> domainMapVec;
1665  for (size_t col = 0; col < numBlocks; ++col) {
1666  RCP<const Map> map = XpIO::ReadMap("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", lib, comm);
1667  domainMapVec.push_back(map);
1668  }
1669  RCP<const Map> fullDomainMap = XpIO::ReadMap("fullDomainMap_" + fileName + ".m", lib, comm);
1670 
1671  /*std::vector<RCP<const XpMap> > testRgMapVec;
1672  for(size_t r = 0; r < numBlocks; ++r) {
1673  RCP<const XpMap> map = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(r) + "0.m", lib, comm);
1674  testRgMapVec.push_back(map);
1675  }
1676  std::vector<RCP<const XpMap> > testDoMapVec;
1677  for(size_t c = 0; c < numBlocks; ++c) {
1678  RCP<const XpMap> map = XpIO::ReadMap("domainmap_" + fileName + "0" + XpIO::toString<size_t>(c) + ".m", lib, comm);
1679  testDoMapVec.push_back(map);
1680  }*/
1681 
1682  // create map extractors
1683 
1684  // range map extractor
1685  bool bRangeUseThyraStyleNumbering = false;
1686  /*
1687  GlobalOrdinal gMinGids = 0;
1688  for(size_t v = 0; v < testRgMapVec.size(); ++v) {
1689  gMinGids += testRgMapVec[v]->getMinAllGlobalIndex();
1690  }
1691  if ( gMinGids==0 && testRgMapVec.size() > 1 ) bRangeUseThyraStyleNumbering = true;*/
1692  RCP<const MapExtractor> rangeMapExtractor =
1693  rcp(new MapExtractor(fullRangeMap, rangeMapVec, bRangeUseThyraStyleNumbering));
1694 
1695  // domain map extractor
1696  bool bDomainUseThyraStyleNumbering = false;
1697  /*gMinGids = 0;
1698  for(size_t v = 0; v < testDoMapVec.size(); ++v) {
1699  gMinGids += testDoMapVec[v]->getMinAllGlobalIndex();
1700  }
1701  if ( gMinGids==0 && testDoMapVec.size() > 1) bDomainUseThyraStyleNumbering = true;*/
1702  RCP<const MapExtractor> domainMapExtractor =
1703  rcp(new MapExtractor(fullDomainMap, domainMapVec, bDomainUseThyraStyleNumbering));
1704 
1705  RCP<BlockedCrsMatrix> bOp = Teuchos::rcp(new BlockedCrsMatrix(rangeMapExtractor, domainMapExtractor, 33));
1706 
1707  // Read all matrices with their maps and create the BlockedCrsMatrix
1708  for (size_t row = 0; row < numBlocks; ++row) {
1709  for (size_t col = 0; col < numBlocks; ++col) {
1710  RCP<const Map> rowSubMap = XpIO::ReadMap("rowmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1711  RCP<const Map> colSubMap = XpIO::ReadMap("colmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1712  RCP<const Map> domSubMap = XpIO::ReadMap("domainmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1713  RCP<const Map> ranSubMap = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1714  RCP<Matrix> mat = XpIO::Read(fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", rowSubMap, colSubMap, domSubMap, ranSubMap);
1715  bOp->setMatrix(row, col, mat);
1716  }
1717  }
1718 
1719  bOp->fillComplete();
1720 
1721  return bOp;
1722  } // ReadBlockedCrsMatrix
1723 
1725  template <class T>
1726  static std::string toString(const T& what) {
1727  std::ostringstream buf;
1728  buf << what;
1729  return buf.str();
1730  }
1731 };
1732 #endif // HAVE_XPETRA_EPETRA
1733 
1734 } // end namespace Xpetra
1735 
1736 #define XPETRA_IO_SHORT
1737 
1738 #endif /* PACKAGES_XPETRA_SUP_UTILS_XPETRA_IO_HPP_ */
Xpetra::CrsMatrixWrap::getCrsMatrix
RCP< CrsMatrix > getCrsMatrix() const
Definition: Xpetra_CrsMatrixWrap_def.hpp:470
Xpetra::IO::Read
static Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Read(const std::string &filename, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rowMap, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> colMap=Teuchos::null, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> domainMap=Teuchos::null, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rangeMap=Teuchos::null, const bool callFillComplete=true, const bool binary=false, const bool tolerant=false, const bool debug=false)
Read matrix from file in Matrix Market or binary format.
Definition: Xpetra_IO.hpp:521
Xpetra::toString
std::string toString(Xpetra::UnderlyingLib lib)
Convert a Xpetra::UnderlyingLib to a std::string.
Definition: Xpetra_Utils.cpp:51
Xpetra::IO< Scalar, int, int, EpetraNode >::ReadLocal
static Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > ReadLocal(const std::string &filename, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rowMap, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> colMap, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> domainMap=Teuchos::null, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rangeMap=Teuchos::null, const bool callFillComplete=true, const bool binary=false, const bool tolerant=false, const bool debug=false)
Read matrix from local files in Matrix Market or binary format.
Definition: Xpetra_IO.hpp:1407
Xpetra::IO< Scalar, int, int, EpetraNode >::Write
static void Write(const std::string &fileName, const Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > &graph, const bool &writeAllMaps=false)
Save CrsGraph to file in Matrix Market format.
Definition: Xpetra_IO.hpp:1065
Xpetra::IO::ReadBlockedCrsMatrix
static RCP< const Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > ReadBlockedCrsMatrix(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int >> &comm)
Read block matrix from one file per block in Matrix Market format.
Definition: Xpetra_IO.hpp:824
Xpetra::CrsGraph::getColMap
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const =0
Returns the Map that describes the column distribution in this graph.
Xpetra_UseShortNames.hpp
Xpetra::IO::Write
static void Write(const std::string &fileName, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &M)
Read/Write methods.
Definition: Xpetra_IO.hpp:220
Xpetra::MatrixFactory::Build
static RCP< Matrix > Build(const RCP< const Map > &rowMap)
Definition: Xpetra_MatrixFactory.hpp:235
Xpetra_EpetraUtils.hpp
Xpetra::IO::ReadMultiVector
static RCP< MultiVector > ReadMultiVector(const std::string &fileName, const RCP< const Map > &map, const bool binary=false)
Definition: Xpetra_IO.hpp:756
Xpetra::IO< Scalar, int, int, EpetraNode >::WriteLocal
static void WriteLocal(const std::string &fileName, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
Save local parts of matrix to files in Matrix Market format.
Definition: Xpetra_IO.hpp:1102
Xpetra::Exceptions::RuntimeError
Exception throws to report errors in the internal logical of the program.
Definition: Xpetra_Exceptions.hpp:101
Xpetra_Matrix.hpp
Xpetra_MatrixMatrix.hpp
Xpetra::TpetraCrsGraph::getTpetra_CrsGraph
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsGraph() const
Get the underlying Tpetra graph.
Definition: Xpetra_TpetraCrsGraph_def.hpp:469
Xpetra_CrsGraph.hpp
Xpetra::Matrix::getColMap
virtual const RCP< const Map > & getColMap() const
Returns the Map that describes the column distribution in this matrix. This might be null until fillC...
Definition: Xpetra_Matrix.hpp:325
Xpetra::Exceptions::BadCast
Exception indicating invalid cast attempted.
Definition: Xpetra_Exceptions.hpp:86
Xpetra::IO::ReadLocal
static Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > ReadLocal(const std::string &filename, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rowMap, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> colMap, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> domainMap=Teuchos::null, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rangeMap=Teuchos::null, const bool callFillComplete=true, const bool binary=false, const bool tolerant=false, const bool debug=false)
Read matrix from local files in Matrix Market or binary format.
Definition: Xpetra_IO.hpp:615
Xpetra::DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node >::getMap
virtual Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const =0
The Map describing the parallel distribution of this object.
Xpetra::MapFactory::Build
static Teuchos::RCP< Map< LocalOrdinal, GlobalOrdinal, Node > > Build(UnderlyingLib lib, global_size_t numGlobalElements, GlobalOrdinal indexBase, const Teuchos::RCP< const Teuchos::Comm< int >> &comm, LocalGlobal lg=Xpetra::GloballyDistributed)
Map constructor with Xpetra-defined contiguous uniform distribution.
Definition: Xpetra_MapFactory_def.hpp:65
Xpetra::BlockedCrsMatrix::getMatrix
Teuchos::RCP< Matrix > getMatrix(size_t r, size_t c) const
return block (r,c)
Definition: Xpetra_BlockedCrsMatrix.hpp:1380
Xpetra::BlockedCrsMatrix::Cols
virtual size_t Cols() const
number of column blocks
Definition: Xpetra_BlockedCrsMatrix.hpp:1344
Xpetra::Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &)
Definition: Xpetra_IO.hpp:109
Xpetra::toTpetra
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > toTpetra(const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph)
Definition: Xpetra_TpetraCrsGraph_decl.hpp:408
Xpetra::IO< Scalar, int, int, EpetraNode >::ReadBlockedCrsMatrix
static RCP< const Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > ReadBlockedCrsMatrix(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int >> &comm)
Read block matrix from one file per block in Matrix Market format.
Definition: Xpetra_IO.hpp:1651
Xpetra::IO< Scalar, int, int, EpetraNode >::WriteBlockedCrsMatrix
static void WriteBlockedCrsMatrix(const std::string &fileName, const Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const bool &writeAllMaps=false)
Save block matrix to one file per block in Matrix Market format.
Definition: Xpetra_IO.hpp:1144
Xpetra::EpetraCrsMatrixT::getEpetra_CrsMatrix
RCP< const Epetra_CrsMatrix > getEpetra_CrsMatrix() const
Definition: Xpetra_EpetraCrsMatrix.hpp:241
Xpetra::Operator::getRangeMap
virtual const Teuchos::RCP< const map_type > getRangeMap() const =0
The Map associated with the range of this operator, which must be compatible with Y...
Xpetra::CrsGraph::getDomainMap
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const =0
Returns the Map associated with the domain of this graph.
Xpetra::IO< Scalar, int, int, EpetraNode >::ReadMap
static RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > ReadMap(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int >> &comm, const bool binary=false)
Definition: Xpetra_IO.hpp:1592
Xpetra::IO::WriteBlockedCrsMatrix
static void WriteBlockedCrsMatrix(const std::string &fileName, const Xpetra::BlockedCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const bool &writeAllMaps=false)
Save block matrix to one file per block in Matrix Market format.
Definition: Xpetra_IO.hpp:363
Xpetra::IO< Scalar, int, int, EpetraNode >::toString
static std::string toString(const T &what)
Little helper function to convert non-string types to strings.
Definition: Xpetra_IO.hpp:1726
Xpetra::IO< Scalar, int, int, EpetraNode >::Read
static Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Read(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int >> &comm, bool binary=false)
Read matrix from file in Matrix Market or binary format.
Definition: Xpetra_IO.hpp:1178
Xpetra::TpetraMap::getTpetra_Map
RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getTpetra_Map() const
Get the underlying Tpetra map.
Definition: Xpetra_TpetraMap_def.hpp:274
Xpetra_EpetraMap.hpp
Xpetra_BlockedCrsMatrix.hpp
Xpetra_MatrixFactory.hpp
Xpetra::BlockedCrsMatrix::getDomainMapExtractor
RCP< const MapExtractor > getDomainMapExtractor() const
Returns map extractor for domain map.
Definition: Xpetra_BlockedCrsMatrix.hpp:1123
Xpetra::EpetraNode
Tpetra::KokkosCompat::KokkosSerialWrapperNode EpetraNode
Definition: Xpetra_Map_decl.hpp:71
Xpetra::IO::Write
static void Write(const std::string &fileName, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const bool &writeAllMaps=false)
Save matrix to file in Matrix Market format.
Definition: Xpetra_IO.hpp:277
Xpetra::IO::WriteLocal
static void WriteLocal(const std::string &fileName, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
Save local parts of matrix to files in Matrix Market format.
Definition: Xpetra_IO.hpp:323
Xpetra::IO< Scalar, int, int, EpetraNode >::ReadMultiVector
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > ReadMultiVector(const std::string &fileName, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> &map, const bool binary=false)
Definition: Xpetra_IO.hpp:1548
Xpetra::Convert_Epetra_MultiVector_ToXpetra_MultiVector
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > Convert_Epetra_MultiVector_ToXpetra_MultiVector(RCP< Epetra_MultiVector > &epX)
Definition: Xpetra_IO.hpp:132
Xpetra::TpetraCrsMatrix::getTpetra_CrsMatrix
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
Definition: Xpetra_TpetraCrsMatrix_def.hpp:647
Xpetra_EpetraMultiVector.hpp
Xpetra_ConfigDefs.hpp
Xpetra::IO::Map2EpetraMap
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
Helper utility to pull out the underlying Epetra objects from an Xpetra object.
Definition: Xpetra_IO.hpp:181
Xpetra::CrsMatrixWrap
Concrete implementation of Xpetra::Matrix.
Definition: Xpetra_CrsMatrixWrap_fwd.hpp:51
Xpetra::IO::Map2TpetraMap
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
Definition: Xpetra_IO.hpp:209
Xpetra::BlockedCrsMatrix::Rows
virtual size_t Rows() const
number of row blocks
Definition: Xpetra_BlockedCrsMatrix.hpp:1338
Xpetra::Operator::getDomainMap
virtual const Teuchos::RCP< const map_type > getDomainMap() const =0
The Map associated with the domain of this operator, which must be compatible with X...
Xpetra::TpetraMultiVector::getTpetra_MultiVector
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_MultiVector() const
Get the underlying Tpetra multivector.
Definition: Xpetra_TpetraMultiVector_def.hpp:555
Xpetra::Matrix::getRowMap
virtual const RCP< const Map > & getRowMap() const
Returns the Map that describes the row distribution in this matrix.
Definition: Xpetra_Matrix.hpp:315
Xpetra::CrsGraph::getRangeMap
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const =0
Returns the Map associated with the domain of this graph.
Xpetra::IO< Scalar, int, int, EpetraNode >::Map2EpetraMap
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
Helper utility to pull out the underlying Epetra objects from an Xpetra object.
Definition: Xpetra_IO.hpp:926
Xpetra::IO< Scalar, int, int, EpetraNode >::Map2TpetraMap
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
Definition: Xpetra_IO.hpp:938
Xpetra::IO
Xpetra utility class containing IO routines to read/write vectors, matrices etc...
Definition: Xpetra_IO.hpp:160
Xpetra::IO< Scalar, int, int, EpetraNode >::Write
static void Write(const std::string &fileName, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &M)
Read/Write methods.
Definition: Xpetra_IO.hpp:949
Xpetra::IO::Write
static void Write(const std::string &fileName, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
Save vector to file in Matrix Market format.
Definition: Xpetra_IO.hpp:247
Xpetra::IO< Scalar, int, int, EpetraNode >::Write
static void Write(const std::string &fileName, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const bool &writeAllMaps=false)
Save matrix to file in Matrix Market format.
Definition: Xpetra_IO.hpp:1014
Xpetra::IO< Scalar, int, int, EpetraNode >::Write
static void Write(const std::string &fileName, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
Save vector to file in Matrix Market format.
Definition: Xpetra_IO.hpp:979
Xpetra::IO< Scalar, int, int, EpetraNode >::Read
static Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Read(const std::string &filename, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rowMap, RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> colMap=Teuchos::null, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> domainMap=Teuchos::null, const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node >> rangeMap=Teuchos::null, const bool callFillComplete=true, const bool binary=false, const bool tolerant=false, const bool debug=false)
Read matrix from file in Matrix Market or binary format.
Definition: Xpetra_IO.hpp:1308
Xpetra::ImportFactory::Build
static RCP< Import< LocalOrdinal, GlobalOrdinal, Node > > Build(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &source, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &target, const Teuchos::RCP< Teuchos::ParameterList > &plist=Teuchos::null)
Constructor specifying the number of non-zeros for all rows.
Definition: Xpetra_ImportFactory.hpp:74
Xpetra::toEpetra
const Epetra_CrsGraph & toEpetra(const RCP< const CrsGraph< int, GlobalOrdinal, Node > > &graph)
Definition: Xpetra_EpetraCrsGraph.cpp:57
Xpetra::IO::Read
static Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Read(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int >> &comm, bool binary=false)
Read matrix from file in Matrix Market or binary format.
Definition: Xpetra_IO.hpp:396
Xpetra::Matrix
Xpetra-specific matrix class.
Definition: Xpetra_Matrix_fwd.hpp:51
Xpetra::CrsGraph::getRowMap
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const =0
Returns the Map that describes the row distribution in this graph.
Xpetra::IO::ReadMap
static RCP< const Map > ReadMap(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int >> &comm, const bool binary=false)
Definition: Xpetra_IO.hpp:785
Xpetra::toXpetra
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
Definition: Xpetra_EpetraCrsGraph.cpp:65
Xpetra::BlockedCrsMatrix::getRangeMapExtractor
RCP< const MapExtractor > getRangeMapExtractor() const
Returns map extractor class for range map.
Definition: Xpetra_BlockedCrsMatrix.hpp:1117
Xpetra::IO::toString
static std::string toString(const T &what)
Little helper function to convert non-string types to strings.
Definition: Xpetra_IO.hpp:899