Xpetra Version of the Day
Loading...
Searching...
No Matches
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>
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"
89#include "Xpetra_CrsGraph.hpp"
90#include "Xpetra_CrsMatrixWrap.hpp"
92
93#include "Xpetra_Map.hpp"
94#include "Xpetra_StridedMap.hpp"
95#include "Xpetra_StridedMapFactory.hpp"
96#include "Xpetra_MapExtractor.hpp"
98
99#include <Teuchos_MatrixMarket_Raw_Writer.hpp>
100#include <Teuchos_MatrixMarket_Raw_Reader.hpp>
101#include <string>
102
103
104namespace Xpetra {
105
106
107#ifdef HAVE_XPETRA_EPETRA
108 // This non-member templated function exists so that the matrix-matrix multiply will compile if Epetra, Tpetra, and ML are enabled.
109 template<class SC,class LO,class GO,class NO>
110 RCP<Xpetra::CrsMatrixWrap<SC,LO,GO,NO> >
113 "Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap cannot be used with Scalar != double, LocalOrdinal != int, GlobalOrdinal != int");
114 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
115 }
116
117 // specialization for the case of ScalarType=double and LocalOrdinal=GlobalOrdinal=int
118 template<>
119 inline RCP<Xpetra::CrsMatrixWrap<double,int,int,Xpetra::EpetraNode> > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<double,int,int,Xpetra::EpetraNode> (RCP<Epetra_CrsMatrix> &epAB) {
120 typedef double SC;
121 typedef int LO;
122 typedef int GO;
123 typedef Xpetra::EpetraNode NO;
124
125 RCP<Xpetra::EpetraCrsMatrixT<GO,NO> > tmpC1 = rcp(new Xpetra::EpetraCrsMatrixT<GO,NO>(epAB));
126 RCP<Xpetra::CrsMatrix<SC,LO,GO,NO> > tmpC2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<SC,LO,GO,NO> >(tmpC1);
127 RCP<Xpetra::CrsMatrixWrap<SC,LO,GO,NO> > tmpC3 = rcp(new Xpetra::CrsMatrixWrap<SC,LO,GO,NO>(tmpC2));
128
129 return tmpC3;
130 }
131
132
133 template<class SC,class LO,class GO,class NO>
134 RCP<Xpetra::MultiVector<SC,LO,GO,NO> >
137 "Convert_Epetra_MultiVector_ToXpetra_MultiVector cannot be used with Scalar != double, LocalOrdinal != int, GlobalOrdinal != int");
138 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
139 }
140
141 // specialization for the case of ScalarType=double and LocalOrdinal=GlobalOrdinal=int
142 template<>
143 inline RCP<Xpetra::MultiVector<double,int,int,Xpetra::EpetraNode> > Convert_Epetra_MultiVector_ToXpetra_MultiVector<double,int,int,Xpetra::EpetraNode> (RCP<Epetra_MultiVector> &epX) {
144 typedef double SC;
145 typedef int LO;
146 typedef int GO;
147 typedef Xpetra::EpetraNode NO;
148
149 RCP<Xpetra::MultiVector<SC,LO,GO,NO >> tmp = Xpetra::toXpetra<GO,NO>(epX);
150 return tmp;
151 }
152
153#endif
154
159 template <class Scalar,
160 class LocalOrdinal = int,
161 class GlobalOrdinal = LocalOrdinal,
163 class IO {
164
165 private:
166#undef XPETRA_IO_SHORT
168
169 public:
170
171#ifdef HAVE_XPETRA_EPETRA
173 // @{
174 /*static RCP<const Epetra_MultiVector> MV2EpetraMV(RCP<MultiVector> const vec);
175 static RCP< Epetra_MultiVector> MV2NonConstEpetraMV(RCP<MultiVector> vec);
176
177 static const Epetra_MultiVector& MV2EpetraMV(const MultiVector& vec);
178 static Epetra_MultiVector& MV2NonConstEpetraMV(MultiVector& vec);
179
180 static RCP<const Epetra_CrsMatrix> Op2EpetraCrs(RCP<const Matrix> Op);
181 static RCP< Epetra_CrsMatrix> Op2NonConstEpetraCrs(RCP<Matrix> Op);
182
183 static const Epetra_CrsMatrix& Op2EpetraCrs(const Matrix& Op);
184 static Epetra_CrsMatrix& Op2NonConstEpetraCrs(Matrix& Op);*/
185
187 RCP<const Xpetra::EpetraMapT<GlobalOrdinal,Node> > xeMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal,Node> >(Teuchos::rcpFromRef(map));
188 if (xeMap == Teuchos::null)
189 throw Exceptions::BadCast("Utils::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
190 return xeMap->getEpetra_Map();
191 }
192 // @}
193#endif
194
195#ifdef HAVE_XPETRA_TPETRA
197 // @{
198 /*static RCP<const Tpetra::MultiVector<SC,LO,GO,NO> > MV2TpetraMV(RCP<MultiVector> const vec);
199 static RCP< Tpetra::MultiVector<SC,LO,GO,NO> > MV2NonConstTpetraMV(RCP<MultiVector> vec);
200 static RCP< Tpetra::MultiVector<SC,LO,GO,NO> > MV2NonConstTpetraMV2(MultiVector& vec);
201
202 static const Tpetra::MultiVector<SC,LO,GO,NO>& MV2TpetraMV(const MultiVector& vec);
203 static Tpetra::MultiVector<SC,LO,GO,NO>& MV2NonConstTpetraMV(MultiVector& vec);
204
205 static RCP<const Tpetra::CrsMatrix<SC,LO,GO,NO> > Op2TpetraCrs(RCP<const Matrix> Op);
206 static RCP< Tpetra::CrsMatrix<SC,LO,GO,NO> > Op2NonConstTpetraCrs(RCP<Matrix> Op);
207
208 static const Tpetra::CrsMatrix<SC,LO,GO,NO>& Op2TpetraCrs(const Matrix& Op);
209 static Tpetra::CrsMatrix<SC,LO,GO,NO>& Op2NonConstTpetraCrs(Matrix& Op);
210
211 static RCP<const Tpetra::RowMatrix<SC,LO,GO,NO> > Op2TpetraRow(RCP<const Matrix> Op);
212 static RCP< Tpetra::RowMatrix<SC,LO,GO,NO> > Op2NonConstTpetraRow(RCP<Matrix> Op);*/
213
214
216 const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >& tmp_TMap = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
217 if (tmp_TMap == Teuchos::null)
218 throw Exceptions::BadCast("Utils::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
219 return tmp_TMap->getTpetra_Map();
220 }
221#endif
222
223
225
226
227 static void Write(const std::string& fileName, const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> & M) {
229#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
230 const RCP<const Xpetra::EpetraMapT<GlobalOrdinal,Node> >& tmp_EMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal,Node> >(tmp_Map);
231 if (tmp_EMap != Teuchos::null) {
232 int rv = EpetraExt::BlockMapToMatrixMarketFile(fileName.c_str(), tmp_EMap->getEpetra_Map());
233 if (rv != 0)
234 throw Exceptions::RuntimeError("EpetraExt::BlockMapToMatrixMarketFile() return value of " + Teuchos::toString(rv));
235 return;
236 }
237#endif // HAVE_XPETRA_EPETRAEXT
238
239#ifdef HAVE_XPETRA_TPETRA
241 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node> >(tmp_Map);
242 if (tmp_TMap != Teuchos::null) {
243 RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > TMap = tmp_TMap->getTpetra_Map();
244 Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >::writeMapFile(fileName, *TMap);
245 return;
246 }
247#endif // HAVE_XPETRA_TPETRA
248
249 throw Exceptions::BadCast("Could not cast to EpetraMap or TpetraMap in map writing");
250
251 } //Write
252
254 static void Write(const std::string& fileName, const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & vec) {
255 std::string mapfile = "map_" + fileName;
256 Write(mapfile, *(vec.getMap()));
257
259#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
260 const RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >& tmp_EVec = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(tmp_Vec);
261 if (tmp_EVec != Teuchos::null) {
262 int rv = EpetraExt::MultiVectorToMatrixMarketFile(fileName.c_str(), *(tmp_EVec->getEpetra_MultiVector()));
263 if (rv != 0)
264 throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
265 return;
266 }
267#endif // HAVE_XPETRA_EPETRA
268
269#ifdef HAVE_XPETRA_TPETRA
271 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmp_Vec);
272 if (tmp_TVec != Teuchos::null) {
273 RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> > TVec = tmp_TVec->getTpetra_MultiVector();
274 Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >::writeDenseFile(fileName, TVec);
275 return;
276 }
277#endif // HAVE_XPETRA_TPETRA
278
279 throw Exceptions::BadCast("Could not cast to EpetraMultiVector or TpetraMultiVector in multivector writing");
280
281 } //Write
282
283
285 static void Write(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Op, const bool &writeAllMaps = false) {
286
287 Write("rowmap_" + fileName, *(Op.getRowMap()));
288 if ( !Op.getDomainMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps )
289 Write("domainmap_" + fileName, *(Op.getDomainMap()));
290 if ( !Op.getRangeMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps )
291 Write("rangemap_" + fileName, *(Op.getRangeMap()));
292 if ( !Op.getColMap()->isSameAs(*(Op.getDomainMap())) || writeAllMaps )
293 Write("colmap_" + fileName, *(Op.getColMap()));
294
298#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
299 const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >& tmp_ECrsMtx = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(tmp_CrsMtx);
300 if (tmp_ECrsMtx != Teuchos::null) {
301 RCP<const Epetra_CrsMatrix> A = tmp_ECrsMtx->getEpetra_CrsMatrix();
302 int rv = EpetraExt::RowMatrixToMatrixMarketFile(fileName.c_str(), *A);
303 if (rv != 0)
304 throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
305 return;
306 }
307#endif
308
309#ifdef HAVE_XPETRA_TPETRA
311 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmp_CrsMtx);
312 if (tmp_TCrsMtx != Teuchos::null) {
313 RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A = tmp_TCrsMtx->getTpetra_CrsMatrix();
314 Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >::writeSparseFile(fileName, A);
315 return;
316 }
318 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmp_CrsMtx);
319 if(tmp_BlockCrs != Teuchos::null) {
320 std::ofstream outstream (fileName,std::ofstream::out);
321 Teuchos::FancyOStream ofs(Teuchos::rcpFromRef(outstream));
322 tmp_BlockCrs->getTpetra_BlockCrsMatrix()->describe(ofs,Teuchos::VERB_EXTREME);
323 return;
324 }
325
326#endif // HAVE_XPETRA_TPETRA
327
328 throw Exceptions::BadCast("Could not cast to EpetraCrsMatrix or TpetraCrsMatrix in matrix writing");
329 } //Write
330
331
333 static void WriteLocal(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Op) {
337
338 ArrayRCP<const size_t> rowptr_RCP;
339 ArrayRCP<LocalOrdinal> rowptr2_RCP;
341 ArrayRCP<const Scalar> vals_RCP;
342 tmp_CrsMtx->getAllValues(rowptr_RCP, colind_RCP, vals_RCP);
343
344 ArrayView<const size_t> rowptr = rowptr_RCP();
345 ArrayView<const LocalOrdinal> colind = colind_RCP();
346 ArrayView<const Scalar> vals = vals_RCP();
347
348 rowptr2_RCP.resize(rowptr.size());
349 ArrayView<LocalOrdinal> rowptr2 = rowptr2_RCP();
350 for (size_t j = 0; j<rowptr.size(); j++)
351 rowptr2[j] = rowptr[j];
352
354 writer.writeFile(fileName + "." + std::to_string(Op.getRowMap()->getComm()->getSize()) + "." + std::to_string(Op.getRowMap()->getComm()->getRank()),
355 rowptr2,colind,vals,
356 rowptr.size()-1,Op.getColMap()->getLocalNumElements());
357 } //WriteLocal
358
359
374 static void WriteBlockedCrsMatrix(const std::string& fileName, const Xpetra::BlockedCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Op, const bool &writeAllMaps = false) {
376
377 // Write all matrix blocks with their maps
378 for (size_t row = 0; row < Op.Rows(); ++row) {
379 for (size_t col = 0; col < Op.Cols(); ++col) {
380 RCP<const Matrix > m = Op.getMatrix(row,col);
381 if(m != Teuchos::null) { // skip empty blocks
382 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(m) == Teuchos::null, Xpetra::Exceptions::BadCast,
383 "Sub block matrix (" << row << "," << col << ") is not of type CrsMatrixWrap.");
384 XpIO::Write(fileName + toString(row) + toString(col) + ".m", *m, writeAllMaps);
385 }
386 }
387 }
388
389 // write map information of map extractors
390 RCP<const MapExtractor> rangeMapExtractor = Op.getRangeMapExtractor();
391 RCP<const MapExtractor> domainMapExtractor = Op.getDomainMapExtractor();
392
393 for(size_t row = 0; row < rangeMapExtractor->NumMaps(); ++row) {
394 RCP<const Map> map = rangeMapExtractor->getMap(row);
395 XpIO::Write("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", *map);
396 }
397 XpIO::Write("fullRangeMap_" + fileName +".m",*(rangeMapExtractor->getFullMap()));
398
399 for(size_t col = 0; col < domainMapExtractor->NumMaps(); ++col) {
400 RCP<const Map> map = domainMapExtractor->getMap(col);
401 XpIO::Write("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", *map);
402 }
403 XpIO::Write("fullDomainMap_" + fileName+ ".m",*(domainMapExtractor->getFullMap()));
404 } // WriteBlockedCrsMatrix
405
407 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) {
408 if (binary == false) {
409 // Matrix Market file format (ASCII)
410 if (lib == Xpetra::UseEpetra) {
411#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
413 const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(comm);
414 int rv = EpetraExt::MatrixMarketFileToCrsMatrix(fileName.c_str(), *epcomm, eA);
415 if (rv != 0)
416 throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
417
418 RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
419
421 Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
422 return A;
423#else
424 throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
425#endif
426 } else if (lib == Xpetra::UseTpetra) {
427#ifdef HAVE_XPETRA_TPETRA
428 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
429
430 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
431
432 bool callFillComplete = true;
433
434 RCP<sparse_matrix_type> tA = reader_type::readSparseFile(fileName, comm, callFillComplete);
435
436 if (tA.is_null())
437 throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
438
440 RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmpA1);
442
443 return A;
444#else
445 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
446#endif
447 } else {
448 throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
449 }
450 } else {
451 // Custom file format (binary)
452 std::ifstream ifs(fileName.c_str(), std::ios::binary);
453 TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << fileName << "\"");
454 int m, n, nnz;
455 ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
456 ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
457 ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
458
459 int myRank = comm->getRank();
460
461 GO indexBase = 0;
462 RCP<Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, m, (myRank == 0 ? m : 0), indexBase, comm), rangeMap = rowMap;
463 RCP<Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > colMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, n, (myRank == 0 ? n : 0), indexBase, comm), domainMap = colMap;
465
466 if (myRank == 0) {
469 // Scan matrix to determine the exact nnz per row.
470 Teuchos::ArrayRCP<size_t> numEntriesPerRow(m,(size_t)(0));
471 for (int i = 0; i < m; i++) {
472 int row, rownnz;
473 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
474 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
475 numEntriesPerRow[row] = rownnz;
476 for (int j = 0; j < rownnz; j++) {
477 int index;
478 ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
479 }
480 for (int j = 0; j < rownnz; j++) {
481 double value;
482 ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
483 }
484 }
485
487
488 // Now that nnz per row are known, reread and store the matrix.
489 ifs.seekg(0, ifs.beg); //rewind to beginning of file
490 int junk; //skip header info
491 ifs.read(reinterpret_cast<char*>(&junk), sizeof(junk));
492 ifs.read(reinterpret_cast<char*>(&junk), sizeof(junk));
493 ifs.read(reinterpret_cast<char*>(&junk), sizeof(junk));
494 for (int i = 0; i < m; i++) {
495 int row, rownnz;
496 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
497 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
498 inds.resize(rownnz);
499 vals.resize(rownnz);
500 for (int j = 0; j < rownnz; j++) {
501 int index;
502 ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
503 inds[j] = Teuchos::as<GlobalOrdinal>(index);
504 }
505 for (int j = 0; j < rownnz; j++) {
506 double value;
507 ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
508 vals[j] = Teuchos::as<Scalar>(value);
509 }
510 A->insertGlobalValues(row, inds, vals);
511 }
512 } //if (myRank == 0)
513 else {
514 Teuchos::ArrayRCP<size_t> numEntriesPerRow(0,(size_t)(0));
516 }
517
518 A->fillComplete(domainMap, rangeMap);
519
520 return A;
521 } //if (binary == false) ... else
522
523 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
524
525 } //Read()
526
527
533 Read(const std::string& filename,
535 RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > colMap = Teuchos::null,
536 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > domainMap = Teuchos::null,
537 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > rangeMap = Teuchos::null,
538 const bool callFillComplete = true,
539 const bool binary = false,
540 const bool tolerant = false,
541 const bool debug = false) {
542 TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::Read() : rowMap cannot be null");
543
544 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domain = (domainMap.is_null() ? rowMap : domainMap);
545 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > range = (rangeMap .is_null() ? rowMap : rangeMap);
546
547 const Xpetra::UnderlyingLib lib = rowMap->lib();
548 if (binary == false) {
549 if (lib == Xpetra::UseEpetra) {
550#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
552 const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(rowMap->getComm());
554 const Epetra_Map& epetraDomainMap = (domainMap.is_null() ? epetraRowMap : Xpetra::IO<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Map2EpetraMap(*domainMap));
555 const Epetra_Map& epetraRangeMap = (rangeMap .is_null() ? epetraRowMap : Xpetra::IO<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Map2EpetraMap(*rangeMap));
556 int rv;
557 if (colMap.is_null()) {
558 rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraRangeMap, epetraDomainMap, eA);
559
560 } else {
561 const Epetra_Map& epetraColMap = Map2EpetraMap(*colMap);
562 rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraColMap, epetraRangeMap, epetraDomainMap, eA);
563 }
564
565 if (rv != 0)
566 throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
567
568 RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
570 Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
571
572 return A;
573#else
574 throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
575#endif
576 } else if (lib == Xpetra::UseTpetra) {
577#ifdef HAVE_XPETRA_TPETRA
578 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
579 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
580 typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
581
582 const RCP<const map_type> tpetraRowMap = Map2TpetraMap(*rowMap);
583 RCP<const map_type> tpetraColMap = (colMap.is_null() ? Teuchos::null : Map2TpetraMap(*colMap));
584 const RCP<const map_type> tpetraRangeMap = (rangeMap.is_null() ? tpetraRowMap : Map2TpetraMap(*rangeMap));
585 const RCP<const map_type> tpetraDomainMap = (domainMap.is_null() ? tpetraRowMap : Map2TpetraMap(*domainMap));
586
587 RCP<sparse_matrix_type> tA = reader_type::readSparseFile(filename, tpetraRowMap, tpetraColMap, tpetraDomainMap, tpetraRangeMap,
588 callFillComplete, tolerant, debug);
589 if (tA.is_null())
590 throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
591
593 RCP<Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(tmpA1);
595
596 return A;
597#else
598 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
599#endif
600 } else {
601 throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
602 }
603 } else {
604
605 // Read in on rank 0.
606 auto tempA = Read(filename, lib, rowMap->getComm(), binary);
607
609 auto importer = Xpetra::ImportFactory<LocalOrdinal,GlobalOrdinal,Node>::Build(tempA->getRowMap(), rowMap);
610 A->doImport(*tempA, *importer, Xpetra::INSERT);
611 if (callFillComplete)
612 A->fillComplete(domainMap, rangeMap);
613
614 return A;
615 }
616
617 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
618 }
619
631 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > domainMap = Teuchos::null,
632 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > rangeMap = Teuchos::null,
633 const bool callFillComplete = true,
634 const bool binary = false,
635 const bool tolerant = false,
636 const bool debug = false) {
637 TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : rowMap cannot be null");
638 TEUCHOS_TEST_FOR_EXCEPTION(colMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : colMap cannot be null");
639
643
644 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domain = (domainMap.is_null() ? rowMap : domainMap);
645 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > range = (rangeMap .is_null() ? rowMap : rangeMap);
646
647 std::string rankFilename = filename + "." + std::to_string(rowMap->getComm()->getSize()) + "." + std::to_string(rowMap->getComm()->getRank());
648 RCP<matrix_type> A = rcp(new crs_wrap_type(rowMap, colMap, 0));
649
650 if (binary == false) {
651
653 params->set("Parse tolerantly", tolerant);
654 params->set("Debug mode", debug);
655
656 LocalOrdinal numRows = rowMap->getLocalNumElements();
657 LocalOrdinal numCols = colMap->getLocalNumElements();
658
659 ArrayRCP<LocalOrdinal> rowptr2_RCP;
660 ArrayRCP<LocalOrdinal> colind2_RCP;
661 ArrayRCP<Scalar> vals2_RCP;
662
664 reader.readFile(rowptr2_RCP,colind2_RCP,vals2_RCP,
665 numRows,numCols,
666 rankFilename);
667
668
669 RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
670
671 ArrayRCP<size_t> rowptr_RCP;
672 ArrayRCP<LocalOrdinal> colind_RCP;
673 ArrayRCP<Scalar> vals_RCP;
674 ACrs->allocateAllValues(colind2_RCP.size(), rowptr_RCP, colind_RCP, vals_RCP);
675
676 rowptr_RCP.assign(rowptr2_RCP.begin(), rowptr2_RCP.end());
677 colind_RCP = colind2_RCP;
678 vals_RCP = vals2_RCP;
679
680 ACrs->setAllValues(rowptr_RCP, colind_RCP, vals_RCP);
681 } else {
682 // Custom file format (binary)
683 std::ifstream ifs = std::ifstream(rankFilename.c_str(), std::ios::binary);
684 TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << filename << "\"");
685
686 int m, n, nnz;
687 ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
688 ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
689 ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
690
691 TEUCHOS_ASSERT_EQUALITY(Teuchos::as<int>(rowMap->getLocalNumElements()), m);
692
696
697 RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
698
699 ACrs->allocateAllValues(nnz, rowptrRCP, indicesRCP, valuesRCP);
700
701 Teuchos::ArrayView<size_t> rowptr = rowptrRCP();
702 Teuchos::ArrayView<LocalOrdinal> indices = indicesRCP();
703 Teuchos::ArrayView<Scalar> values = valuesRCP();
704
705 bool sorted = true;
706
707 // Read in rowptr
708 for (int i = 0; i < m; i++) {
709 int row, rownnz;
710 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
711 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
712
713 rowptr[row+1] += rownnz;
714 ifs.seekg(sizeof(int)*rownnz + sizeof(double)*rownnz, ifs.cur);
715 }
716 for (int i = 0; i < m; i++)
717 rowptr[i+1] += rowptr[i];
718 TEUCHOS_ASSERT(Teuchos::as<int>(rowptr[m]) == nnz);
719
720 // reset to where the data starts
721 ifs.seekg(sizeof(int)*3, ifs.beg);
722
723 // read in entries
724 for (int i = 0; i < m; i++) {
725 int row, rownnz;
726 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
727 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
728 size_t ptr = rowptr[row];
729 for (int j = 0; j < rownnz; j++) {
730 int index;
731 ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
732 indices[ptr] = Teuchos::as<LocalOrdinal>(index);
733 if (j>0)
734 sorted = sorted & (indices[ptr-1] < indices[ptr]);
735 ++ptr;
736 }
737 ptr = rowptr[row];
738 for (int j = 0; j < rownnz; j++) {
739 double value;
740 ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
741 values[ptr] = Teuchos::as<Scalar>(value);
742 ++ptr;
743 }
744 rowptr[row] += rownnz;
745 }
746 for (int i = m; i > 0; i--)
747 rowptr[i] = rowptr[i-1];
748 rowptr[0] = 0;
749
750#ifdef HAVE_XPETRA_TPETRA
751 if (!sorted) {
752 for (LocalOrdinal lclRow = 0; lclRow < m; lclRow++) {
753 size_t rowBegin = rowptr[lclRow];
754 size_t rowEnd = rowptr[lclRow+1];
755 Tpetra::sort2(&indices[rowBegin], &indices[rowEnd], &values[rowBegin]);
756 }
757 }
758#else
759 TEUCHOS_ASSERT(sorted);
760#endif
761
762 ACrs->setAllValues(rowptrRCP, indicesRCP, valuesRCP);
763
764 }
765
766 if (callFillComplete)
767 A->fillComplete(domainMap, rangeMap);
768 return A;
769
770 }
772
773
774 static RCP<MultiVector> ReadMultiVector (const std::string& fileName,
775 const RCP<const Map>& map,
776 const bool binary=false) {
777 Xpetra::UnderlyingLib lib = map->lib();
778
779 if (lib == Xpetra::UseEpetra) {
780 TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
781
782 } else if (lib == Xpetra::UseTpetra) {
783#ifdef HAVE_XPETRA_TPETRA
784 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
785 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
786 typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
787 typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> multivector_type;
788
789 RCP<const map_type> temp = toTpetra(map);
790 RCP<multivector_type> TMV = reader_type::readDenseFile(fileName,map->getComm(),temp,false,false,binary);
792 return rmv;
793#else
794 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
795#endif
796 } else {
797 throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
798 }
799
800 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
801 }
802
803 static RCP<const Map> ReadMap (const std::string& fileName,
805 const RCP<const Teuchos::Comm<int> >& comm,
806 const bool binary=false) {
807 if (lib == Xpetra::UseEpetra) {
808 TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
809 } else if (lib == Xpetra::UseTpetra) {
810#ifdef HAVE_XPETRA_TPETRA
811 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
812 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
813
814 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > tMap = reader_type::readMapFile(fileName, comm, false, false, binary);
815 if (tMap.is_null())
816 throw Exceptions::RuntimeError("The Tpetra::Map returned from readSparseFile() is null.");
817
818 return Xpetra::toXpetra(tMap);
819#else
820 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
821#endif
822 } else {
823 throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
824 }
825
826 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
827
828 }
829
845
846 size_t numBlocks = 2; // TODO user parameter?
847
848 std::vector<RCP<const Map> > rangeMapVec;
849 for(size_t row = 0; row < numBlocks; ++row) {
850 RCP<const Map> map = XpIO::ReadMap("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", lib, comm);
851 rangeMapVec.push_back(map);
852 }
853 RCP<const Map> fullRangeMap = XpIO::ReadMap("fullRangeMap_" + fileName + ".m", lib, comm);
854
855 std::vector<RCP<const Map> > domainMapVec;
856 for(size_t col = 0; col < numBlocks; ++col) {
857 RCP<const Map> map = XpIO::ReadMap("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", lib, comm);
858 domainMapVec.push_back(map);
859 }
860 RCP<const Map> fullDomainMap = XpIO::ReadMap("fullDomainMap_" + fileName + ".m", lib, comm);
861
862 /*std::vector<RCP<const XpMap> > testRgMapVec;
863 for(size_t r = 0; r < numBlocks; ++r) {
864 RCP<const XpMap> map = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(r) + "0.m", lib, comm);
865 testRgMapVec.push_back(map);
866 }
867 std::vector<RCP<const XpMap> > testDoMapVec;
868 for(size_t c = 0; c < numBlocks; ++c) {
869 RCP<const XpMap> map = XpIO::ReadMap("domainmap_" + fileName + "0" + XpIO::toString<size_t>(c) + ".m", lib, comm);
870 testDoMapVec.push_back(map);
871 }*/
872
873 // create map extractors
874
875 // range map extractor
876 bool bRangeUseThyraStyleNumbering = false;
877 /*GlobalOrdinal gMinGids = 0;
878 for(size_t v = 0; v < testRgMapVec.size(); ++v) {
879 gMinGids += testRgMapVec[v]->getMinAllGlobalIndex();
880 }
881 if ( gMinGids==0 && testRgMapVec.size() > 1 ) bRangeUseThyraStyleNumbering = true;
882 */
883 RCP<const MapExtractor> rangeMapExtractor =
884 Teuchos::rcp(new MapExtractor(fullRangeMap, rangeMapVec, bRangeUseThyraStyleNumbering));
885
886
887 // domain map extractor
888 bool bDomainUseThyraStyleNumbering = false;
889 /*gMinGids = 0;
890 for(size_t v = 0; v < testDoMapVec.size(); ++v) {
891 gMinGids += testDoMapVec[v]->getMinAllGlobalIndex();
892 }
893 if ( gMinGids==0 && testDoMapVec.size() > 1 ) bDomainUseThyraStyleNumbering = true;
894 */
895 RCP<const MapExtractor> domainMapExtractor =
896 Teuchos::rcp(new MapExtractor(fullDomainMap, domainMapVec, bDomainUseThyraStyleNumbering));
897
898 RCP<BlockedCrsMatrix> bOp = Teuchos::rcp(new BlockedCrsMatrix(rangeMapExtractor, domainMapExtractor,33));
899
900 // Read all sub-matrices with their maps and place into blocked operator
901 for (size_t row = 0; row < numBlocks; ++row) {
902 for (size_t col = 0; col < numBlocks; ++col) {
903 RCP<const Map> rowSubMap = XpIO::ReadMap("rowmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
904 RCP<const Map> colSubMap = XpIO::ReadMap("colmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
905 RCP<const Map> domSubMap = XpIO::ReadMap("domainmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
906 RCP<const Map> ranSubMap = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
907 RCP<Matrix> mat = XpIO::Read(fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", rowSubMap, colSubMap, domSubMap, ranSubMap);
908 bOp->setMatrix(row, col, mat);
909 }
910 }
911
912 bOp->fillComplete();
913
914 return bOp;
915 } // ReadBlockedCrsMatrix
916
917
919 template<class T>
920 static std::string toString(const T& what) {
921 std::ostringstream buf;
922 buf << what;
923 return buf.str();
924 }
925 };
926
927
928#ifdef HAVE_XPETRA_EPETRA
938 template <class Scalar>
939 class IO<Scalar,int,int,EpetraNode> {
940 public:
941 typedef int LocalOrdinal;
942 typedef int GlobalOrdinal;
944
945#ifdef HAVE_XPETRA_EPETRA
947 // @{
949 RCP<const Xpetra::EpetraMapT<GlobalOrdinal,Node> > xeMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal,Node> >(Teuchos::rcpFromRef(map));
950 if (xeMap == Teuchos::null)
951 throw Exceptions::BadCast("IO::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
952 return xeMap->getEpetra_Map();
953 }
954 // @}
955#endif
956
957#ifdef HAVE_XPETRA_TPETRA
959 // @{
961 const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >& tmp_TMap = Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
962 if (tmp_TMap == Teuchos::null)
963 throw Exceptions::BadCast("IO::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
964 return tmp_TMap->getTpetra_Map();
965 }
966#endif
967
968
970
971
972 static void Write(const std::string& fileName, const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> & M) {
974#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
975 const RCP<const Xpetra::EpetraMapT<GlobalOrdinal,Node> >& tmp_EMap = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal,Node> >(tmp_Map);
976 if (tmp_EMap != Teuchos::null) {
977 int rv = EpetraExt::BlockMapToMatrixMarketFile(fileName.c_str(), tmp_EMap->getEpetra_Map());
978 if (rv != 0)
979 throw Exceptions::RuntimeError("EpetraExt::BlockMapToMatrixMarketFile() return value of " + Teuchos::toString(rv));
980 return;
981 }
982#endif // HAVE_XPETRA_EPETRA
983
984#ifdef HAVE_XPETRA_TPETRA
985# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
986 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
987 // do nothing
988# else
990 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal, GlobalOrdinal, Node> >(tmp_Map);
991 if (tmp_TMap != Teuchos::null) {
992 RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > TMap = tmp_TMap->getTpetra_Map();
993 Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >::writeMapFile(fileName, *TMap);
994 return;
995 }
996# endif
997#endif // HAVE_XPETRA_TPETRA
998 throw Exceptions::BadCast("Could not cast to EpetraMap or TpetraMap in map writing");
999 }
1000
1002 static void Write(const std::string& fileName, const Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & vec) {
1003 std::string mapfile = "map_" + fileName;
1004 Write(mapfile, *(vec.getMap()));
1005
1007#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1008 const RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >& tmp_EVec = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(tmp_Vec);
1009 if (tmp_EVec != Teuchos::null) {
1010 int rv = EpetraExt::MultiVectorToMatrixMarketFile(fileName.c_str(), *(tmp_EVec->getEpetra_MultiVector()));
1011 if (rv != 0)
1012 throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
1013 return;
1014 }
1015#endif // HAVE_XPETRA_EPETRAEXT
1016
1017#ifdef HAVE_XPETRA_TPETRA
1018# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1019 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1020 // do nothin
1021# else
1023 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmp_Vec);
1024 if (tmp_TVec != Teuchos::null) {
1025 RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> > TVec = tmp_TVec->getTpetra_MultiVector();
1026 Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >::writeDenseFile(fileName, TVec);
1027 return;
1028 }
1029# endif
1030#endif // HAVE_XPETRA_TPETRA
1031
1032 throw Exceptions::BadCast("Could not cast to EpetraMultiVector or TpetraMultiVector in multivector writing");
1033
1034 } //Write
1035
1036
1037
1039 static void Write(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Op, const bool &writeAllMaps = false) {
1040
1041 Write("rowmap_" + fileName, *(Op.getRowMap()));
1042 if ( !Op.getDomainMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps )
1043 Write("domainmap_" + fileName, *(Op.getDomainMap()));
1044 if ( !Op.getRangeMap()->isSameAs(*(Op.getRowMap())) || writeAllMaps )
1045 Write("rangemap_" + fileName, *(Op.getRangeMap()));
1046 if ( !Op.getColMap()->isSameAs(*(Op.getDomainMap())) || writeAllMaps )
1047 Write("colmap_" + fileName, *(Op.getColMap()));
1048
1052#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1053 const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >& tmp_ECrsMtx = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(tmp_CrsMtx);
1054 if (tmp_ECrsMtx != Teuchos::null) {
1055 RCP<const Epetra_CrsMatrix> A = tmp_ECrsMtx->getEpetra_CrsMatrix();
1056 int rv = EpetraExt::RowMatrixToMatrixMarketFile(fileName.c_str(), *A);
1057 if (rv != 0)
1058 throw Exceptions::RuntimeError("EpetraExt::RowMatrixToMatrixMarketFile return value of " + Teuchos::toString(rv));
1059 return;
1060 }
1061#endif // endif HAVE_XPETRA_EPETRA
1062
1063#ifdef HAVE_XPETRA_TPETRA
1064# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1065 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1066 // do nothin
1067# else
1069 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmp_CrsMtx);
1070 if (tmp_TCrsMtx != Teuchos::null) {
1071 RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A = tmp_TCrsMtx->getTpetra_CrsMatrix();
1072 Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >::writeSparseFile(fileName, A);
1073 return;
1074 }
1076 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmp_CrsMtx);
1077 if(tmp_BlockCrs != Teuchos::null) {
1078 std::ofstream outstream (fileName,std::ofstream::out);
1079 Teuchos::FancyOStream ofs(Teuchos::rcpFromRef(outstream));
1080 tmp_BlockCrs->getTpetra_BlockCrsMatrix()->describe(ofs,Teuchos::VERB_EXTREME);
1081 return;
1082 }
1083
1084# endif
1085#endif // HAVE_XPETRA_TPETRA
1086
1087 throw Exceptions::BadCast("Could not cast to EpetraCrsMatrix or TpetraCrsMatrix in matrix writing");
1088 } //Write
1089
1090
1092 static void Write(const std::string& fileName, const Xpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> & graph, const bool &writeAllMaps = false) {
1093
1094 Write("rowmap_" + fileName, *(graph.getRowMap()));
1095 if ( !graph.getDomainMap()->isSameAs(*(graph.getRowMap())) || writeAllMaps )
1096 Write("domainmap_" + fileName, *(graph.getDomainMap()));
1097 if ( !graph.getRangeMap()->isSameAs(*(graph.getRowMap())) || writeAllMaps )
1098 Write("rangemap_" + fileName, *(graph.getRangeMap()));
1099 if ( !graph.getColMap()->isSameAs(*(graph.getDomainMap())) || writeAllMaps )
1100 Write("colmap_" + fileName, *(graph.getColMap()));
1101
1103
1104#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1105 const RCP<const Xpetra::EpetraCrsGraphT<GlobalOrdinal,Node> >& tmp_ECrsGraph = Teuchos::rcp_dynamic_cast<const Xpetra::EpetraCrsGraphT<GlobalOrdinal,Node> >(tmp_Graph);
1106 if (tmp_ECrsGraph != Teuchos::null) {
1107 throw Exceptions::BadCast("Writing not implemented for EpetraCrsGraphT");
1108 }
1109#endif // endif HAVE_XPETRA_EPETRA
1110
1111#ifdef HAVE_XPETRA_TPETRA
1112# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1113 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1114 // do nothin
1115# else
1117 Teuchos::rcp_dynamic_cast<const Xpetra::TpetraCrsGraph<LocalOrdinal, GlobalOrdinal, Node> >(tmp_Graph);
1118 if (tmp_TCrsGraph != Teuchos::null) {
1119 RCP<const Tpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node> > G = tmp_TCrsGraph->getTpetra_CrsGraph();
1120 Tpetra::MatrixMarket::Writer<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >::writeSparseGraphFile(fileName, G);
1121 return;
1122 }
1123# endif
1124#endif // HAVE_XPETRA_TPETRA
1125
1126 throw Exceptions::BadCast("Could not cast to EpetraCrsMatrix or TpetraCrsMatrix in matrix writing");
1127 } //Write
1128
1129
1131 static void WriteLocal(const std::string& fileName, const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Op) {
1135
1136 ArrayRCP<const size_t> rowptr_RCP;
1137 ArrayRCP<LocalOrdinal> rowptr2_RCP;
1139 ArrayRCP<const Scalar> vals_RCP;
1140 tmp_CrsMtx->getAllValues(rowptr_RCP, colind_RCP, vals_RCP);
1141
1142 ArrayView<const size_t> rowptr = rowptr_RCP();
1143 ArrayView<const LocalOrdinal> colind = colind_RCP();
1144 ArrayView<const Scalar> vals = vals_RCP();
1145
1146 rowptr2_RCP.resize(rowptr.size());
1147 ArrayView<LocalOrdinal> rowptr2 = rowptr2_RCP();
1148 for (size_t j = 0; j<Teuchos::as<size_t>(rowptr.size()); j++)
1149 rowptr2[j] = rowptr[j];
1150
1152 writer.writeFile(fileName + "." + std::to_string(Op.getRowMap()->getComm()->getSize()) + "." + std::to_string(Op.getRowMap()->getComm()->getRank()),
1153 rowptr2,colind,vals,
1154 rowptr.size()-1,Op.getColMap()->getLocalNumElements());
1155 } //WriteLocal
1156
1173 static void WriteBlockedCrsMatrix(const std::string& fileName, const Xpetra::BlockedCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Op, const bool &writeAllMaps = false) {
1174#include "Xpetra_UseShortNames.hpp"
1176
1177 // write all matrices with their maps
1178 for (size_t row = 0; row < Op.Rows(); ++row) {
1179 for (size_t col = 0; col < Op.Cols(); ++col) {
1180 RCP<const Matrix> m = Op.getMatrix(row, col);
1181 if(m != Teuchos::null) { // skip empty blocks
1182 TEUCHOS_TEST_FOR_EXCEPTION(Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(m) == Teuchos::null, Xpetra::Exceptions::BadCast,
1183 "Sub block matrix (" << row << "," << col << ") is not of type CrsMatrixWrap.");
1184 XpIO::Write(fileName + toString(row) + toString(col) + ".m", *m, writeAllMaps);
1185 }
1186 }
1187 }
1188
1189 // write map information of map extractors
1190 RCP<const MapExtractor> rangeMapExtractor = Op.getRangeMapExtractor();
1191 RCP<const MapExtractor> domainMapExtractor = Op.getDomainMapExtractor();
1192
1193 for(size_t row = 0; row < rangeMapExtractor->NumMaps(); ++row) {
1194 RCP<const Map> map = rangeMapExtractor->getMap(row);
1195 XpIO::Write("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", *map);
1196 }
1197 XpIO::Write("fullRangeMap_" + fileName +".m", *(rangeMapExtractor->getFullMap()));
1198
1199 for(size_t col = 0; col < domainMapExtractor->NumMaps(); ++col) {
1200 RCP<const Map> map = domainMapExtractor->getMap(col);
1201 XpIO::Write("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", *map);
1202 }
1203 XpIO::Write("fullDomainMap_" + fileName+ ".m", *(domainMapExtractor->getFullMap()));
1204 } // WriteBlockedCrsMatrix
1205
1207 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) {
1208 if (binary == false) {
1209 // Matrix Market file format (ASCII)
1210 if (lib == Xpetra::UseEpetra) {
1211#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1212 Epetra_CrsMatrix *eA;
1213 const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(comm);
1214 int rv = EpetraExt::MatrixMarketFileToCrsMatrix(fileName.c_str(), *epcomm, eA);
1215 if (rv != 0)
1216 throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
1217
1218 RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
1219
1221 Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
1222 return A;
1223#else
1224 throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1225#endif
1226 } else if (lib == Xpetra::UseTpetra) {
1227#ifdef HAVE_XPETRA_TPETRA
1228# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1229 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1230 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int enabled.");
1231# else
1232 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1233
1234 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1235
1236 bool callFillComplete = true;
1237
1238 RCP<sparse_matrix_type> tA = reader_type::readSparseFile(fileName, comm, callFillComplete);
1239
1240 if (tA.is_null())
1241 throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
1242
1244 RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(tmpA1);
1246
1247 return A;
1248# endif
1249#else
1250 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1251#endif
1252 } else {
1253 throw Exceptions::RuntimeError("Xpetra:IO: you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1254 }
1255 } else {
1256 // Custom file format (binary)
1257 std::ifstream ifs(fileName.c_str(), std::ios::binary);
1258 TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << fileName << "\"");
1259 int m, n, nnz;
1260 ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
1261 ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
1262 ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
1263
1264 int myRank = comm->getRank();
1265
1266 GlobalOrdinal indexBase = 0;
1267 RCP<Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, m, (myRank == 0 ? m : 0), indexBase, comm), rangeMap = rowMap;
1268 RCP<Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > colMap = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>::Build(lib, n, (myRank == 0 ? n : 0), indexBase, comm), domainMap = colMap;
1269
1271
1272 if (myRank == 0) {
1275 // Scan matrix to determine the exact nnz per row.
1276 Teuchos::ArrayRCP<size_t> numEntriesPerRow(m);
1277 for (int i = 0; i < m; i++) {
1278 int row, rownnz;
1279 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1280 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1281 numEntriesPerRow[i] = rownnz;
1282 for (int j = 0; j < rownnz; j++) {
1283 int index;
1284 ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
1285 }
1286 for (int j = 0; j < rownnz; j++) {
1287 double value;
1288 ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
1289 }
1290 }
1291
1293
1294 // Now that nnz per row are known, reread and store the matrix.
1295 ifs.seekg(0, ifs.beg); //rewind to beginning of file
1296 int junk; //skip header info
1297 ifs.read(reinterpret_cast<char*>(&m), sizeof(junk));
1298 ifs.read(reinterpret_cast<char*>(&n), sizeof(junk));
1299 ifs.read(reinterpret_cast<char*>(&nnz), sizeof(junk));
1300 for (int i = 0; i < m; i++) {
1301 int row, rownnz;
1302 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1303 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1304 inds.resize(rownnz);
1305 vals.resize(rownnz);
1306 for (int j = 0; j < rownnz; j++) {
1307 int index;
1308 ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
1309 inds[j] = Teuchos::as<GlobalOrdinal>(index);
1310 }
1311 for (int j = 0; j < rownnz; j++) {
1312 double value;
1313 ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
1314 vals[j] = Teuchos::as<Scalar>(value);
1315 }
1316 A->insertGlobalValues(row, inds, vals);
1317 }
1318 } //if (myRank == 0)
1319 else {
1320 Teuchos::ArrayRCP<size_t> numEntriesPerRow(0);
1322 }
1323
1324 A->fillComplete(domainMap, rangeMap);
1325
1326 return A;
1327 }
1328
1329 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1330
1331 } //Read()
1332
1333
1340 RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > colMap = Teuchos::null,
1341 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > domainMap = Teuchos::null,
1342 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > rangeMap = Teuchos::null,
1343 const bool callFillComplete = true,
1344 const bool binary = false,
1345 const bool tolerant = false,
1346 const bool debug = false) {
1347 TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::Read() : rowMap cannot be null");
1348
1349 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domain = (domainMap.is_null() ? rowMap : domainMap);
1350 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > range = (rangeMap .is_null() ? rowMap : rangeMap);
1351
1352 const Xpetra::UnderlyingLib lib = rowMap->lib();
1353 if (binary == false) {
1354 if (lib == Xpetra::UseEpetra) {
1355#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1356 Epetra_CrsMatrix *eA;
1357 const RCP<const Epetra_Comm> epcomm = Xpetra::toEpetra(rowMap->getComm());
1359 const Epetra_Map& epetraDomainMap = (domainMap.is_null() ? epetraRowMap : Xpetra::IO<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Map2EpetraMap(*domainMap));
1360 const Epetra_Map& epetraRangeMap = (rangeMap .is_null() ? epetraRowMap : Xpetra::IO<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Map2EpetraMap(*rangeMap));
1361 int rv;
1362 if (colMap.is_null()) {
1363 rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraRangeMap, epetraDomainMap, eA);
1364
1365 } else {
1366 const Epetra_Map& epetraColMap = Map2EpetraMap(*colMap);
1367 rv = EpetraExt::MatrixMarketFileToCrsMatrix(filename.c_str(), epetraRowMap, epetraColMap, epetraRangeMap, epetraDomainMap, eA);
1368 }
1369
1370 if (rv != 0)
1371 throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToCrsMatrix return value of " + Teuchos::toString(rv));
1372
1373 RCP<Epetra_CrsMatrix> tmpA = rcp(eA);
1375 Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(tmpA);
1376
1377 return A;
1378#else
1379 throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1380#endif
1381 } else if (lib == Xpetra::UseTpetra) {
1382#ifdef HAVE_XPETRA_TPETRA
1383# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1384 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1385 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int support.");
1386# else
1387 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1388 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1389 typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
1390
1391 const RCP<const map_type> tpetraRowMap = Map2TpetraMap(*rowMap);
1392 RCP<const map_type> tpetraColMap = (colMap.is_null() ? Teuchos::null : Map2TpetraMap(*colMap));
1393 const RCP<const map_type> tpetraRangeMap = (rangeMap.is_null() ? tpetraRowMap : Map2TpetraMap(*rangeMap));
1394 const RCP<const map_type> tpetraDomainMap = (domainMap.is_null() ? tpetraRowMap : Map2TpetraMap(*domainMap));
1395
1396 RCP<sparse_matrix_type> tA = reader_type::readSparseFile(filename, tpetraRowMap, tpetraColMap, tpetraDomainMap, tpetraRangeMap,
1397 callFillComplete, tolerant, debug);
1398 if (tA.is_null())
1399 throw Exceptions::RuntimeError("The Tpetra::CrsMatrix returned from readSparseFile() is null.");
1400
1402 RCP<Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpA2 = Teuchos::rcp_implicit_cast<Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(tmpA1);
1404
1405 return A;
1406# endif
1407#else
1408 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1409#endif
1410 } else {
1411 throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1412 }
1413 } else {
1414
1415 // Read in on rank 0.
1416 auto tempA = Read(filename, lib, rowMap->getComm(), binary);
1417
1419 auto importer = Xpetra::ImportFactory<LocalOrdinal,GlobalOrdinal,Node>::Build(tempA->getRowMap(), rowMap);
1420 A->doImport(*tempA, *importer, Xpetra::INSERT);
1421 if (callFillComplete)
1422 A->fillComplete(domainMap, rangeMap);
1423
1424 return A;
1425 }
1426
1427 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1428 }
1429
1441 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > domainMap = Teuchos::null,
1442 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > rangeMap = Teuchos::null,
1443 const bool callFillComplete = true,
1444 const bool binary = false,
1445 const bool tolerant = false,
1446 const bool debug = false) {
1447 TEUCHOS_TEST_FOR_EXCEPTION(rowMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : rowMap cannot be null");
1448 TEUCHOS_TEST_FOR_EXCEPTION(colMap.is_null(), Exceptions::RuntimeError, "Utils::ReadLocal() : colMap cannot be null");
1449
1453
1454 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domain = (domainMap.is_null() ? rowMap : domainMap);
1455 RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > range = (rangeMap .is_null() ? rowMap : rangeMap);
1456
1457 std::string rankFilename = filename + "." + std::to_string(rowMap->getComm()->getSize()) + "." + std::to_string(rowMap->getComm()->getRank());
1458 RCP<matrix_type> A = rcp(new crs_wrap_type(rowMap, colMap, 0));
1459
1460 if (binary == false) {
1461
1463 params->set("Parse tolerantly", tolerant);
1464 params->set("Debug mode", debug);
1465
1466 LocalOrdinal numRows = rowMap->getLocalNumElements();
1467 LocalOrdinal numCols = colMap->getLocalNumElements();
1468
1469 ArrayRCP<LocalOrdinal> rowptr2_RCP;
1470 ArrayRCP<LocalOrdinal> colind2_RCP;
1471 ArrayRCP<Scalar> vals2_RCP;
1472
1474 reader.readFile(rowptr2_RCP,colind2_RCP,vals2_RCP,
1475 numRows,numCols,
1476 rankFilename);
1477
1478
1479 RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
1480
1481 ArrayRCP<size_t> rowptr_RCP;
1482 ArrayRCP<LocalOrdinal> colind_RCP;
1483 ArrayRCP<Scalar> vals_RCP;
1484 ACrs->allocateAllValues(colind2_RCP.size(), rowptr_RCP, colind_RCP, vals_RCP);
1485
1486 rowptr_RCP.assign(rowptr2_RCP.begin(), rowptr2_RCP.end());
1487 colind_RCP = colind2_RCP;
1488 vals_RCP = vals2_RCP;
1489
1490 ACrs->setAllValues(rowptr_RCP, colind_RCP, vals_RCP);
1491 } else {
1492 // Custom file format (binary)
1493 std::ifstream ifs = std::ifstream(rankFilename.c_str(), std::ios::binary);
1494 TEUCHOS_TEST_FOR_EXCEPTION(!ifs.good(), Exceptions::RuntimeError, "Can not read \"" << filename << "\"");
1495
1496 int m, n, nnz;
1497 ifs.read(reinterpret_cast<char*>(&m), sizeof(m));
1498 ifs.read(reinterpret_cast<char*>(&n), sizeof(n));
1499 ifs.read(reinterpret_cast<char*>(&nnz), sizeof(nnz));
1500
1501 TEUCHOS_ASSERT_EQUALITY(Teuchos::as<int>(rowMap->getLocalNumElements()), m);
1502
1503 Teuchos::ArrayRCP<size_t> rowptrRCP;
1505 Teuchos::ArrayRCP<Scalar> valuesRCP;
1506
1507 RCP<crs_type> ACrs = Teuchos::rcp_dynamic_cast<crs_wrap_type>(A)->getCrsMatrix();
1508
1509 ACrs->allocateAllValues(nnz, rowptrRCP, indicesRCP, valuesRCP);
1510
1511 Teuchos::ArrayView<size_t> rowptr = rowptrRCP();
1512 Teuchos::ArrayView<LocalOrdinal> indices = indicesRCP();
1513 Teuchos::ArrayView<Scalar> values = valuesRCP();
1514
1515 bool sorted = true;
1516
1517 // Read in rowptr
1518 for (int i = 0; i < m; i++) {
1519 int row, rownnz;
1520 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1521 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1522
1523 rowptr[row+1] += rownnz;
1524 ifs.seekg(sizeof(int)*rownnz + sizeof(double)*rownnz, ifs.cur);
1525 }
1526 for (int i = 0; i < m; i++)
1527 rowptr[i+1] += rowptr[i];
1528 TEUCHOS_ASSERT(Teuchos::as<int>(rowptr[m]) == nnz);
1529
1530 // reset to where the data starts
1531 ifs.seekg(sizeof(int)*3, ifs.beg);
1532
1533 // read in entries
1534 for (int i = 0; i < m; i++) {
1535 int row, rownnz;
1536 ifs.read(reinterpret_cast<char*>(&row), sizeof(row));
1537 ifs.read(reinterpret_cast<char*>(&rownnz), sizeof(rownnz));
1538 size_t ptr = rowptr[row];
1539 for (int j = 0; j < rownnz; j++) {
1540 int index;
1541 ifs.read(reinterpret_cast<char*>(&index), sizeof(index));
1542 indices[ptr] = Teuchos::as<LocalOrdinal>(index);
1543 if (j>0)
1544 sorted = sorted & (indices[ptr-1] < indices[ptr]);
1545 ++ptr;
1546 }
1547 ptr = rowptr[row];
1548 for (int j = 0; j < rownnz; j++) {
1549 double value;
1550 ifs.read(reinterpret_cast<char*>(&value), sizeof(value));
1551 values[ptr] = Teuchos::as<Scalar>(value);
1552 ++ptr;
1553 }
1554 rowptr[row] += rownnz;
1555 }
1556 for (int i = m; i > 0; i--)
1557 rowptr[i] = rowptr[i-1];
1558 rowptr[0] = 0;
1559
1560#ifdef HAVE_XPETRA_TPETRA
1561 if (!sorted) {
1562 for (LocalOrdinal lclRow = 0; lclRow < m; lclRow++) {
1563 size_t rowBegin = rowptr[lclRow];
1564 size_t rowEnd = rowptr[lclRow+1];
1565 Tpetra::sort2(&indices[rowBegin], &indices[rowEnd], &values[rowBegin]);
1566 }
1567 }
1568#else
1569 TEUCHOS_ASSERT(sorted);
1570#endif
1571
1572 ACrs->setAllValues(rowptrRCP, indicesRCP, valuesRCP);
1573
1574 }
1575
1576 if (callFillComplete)
1577 A->fillComplete(domainMap, rangeMap);
1578 return A;
1579
1580 }
1582
1583
1586 const bool binary=false) {
1587 Xpetra::UnderlyingLib lib = map->lib();
1588
1589 if (lib == Xpetra::UseEpetra) {
1590 // taw: Oct 9 2015: do we need a specialization for <double,int,int>??
1591 //TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
1592#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1593 TEUCHOS_ASSERT(!binary);
1594 Epetra_MultiVector * MV;
1595 int rv = EpetraExt::MatrixMarketFileToMultiVector(fileName.c_str(), toEpetra(map), MV);
1596 if(rv != 0) throw Exceptions::RuntimeError("EpetraExt::MatrixMarketFileToMultiVector failed");
1597 RCP<Epetra_MultiVector> MVrcp = rcp(MV);
1598 return Convert_Epetra_MultiVector_ToXpetra_MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>(MVrcp);
1599#else
1600 throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1601#endif
1602 } else if (lib == Xpetra::UseTpetra) {
1603#ifdef HAVE_XPETRA_TPETRA
1604# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1605 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1606 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int support.");
1607# else
1608 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1609 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1610 typedef Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> map_type;
1611 typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> multivector_type;
1612
1613 RCP<const map_type> temp = toTpetra(map);
1614 RCP<multivector_type> TMV = reader_type::readDenseFile(fileName,map->getComm(),temp,false,false,binary);
1616 return rmv;
1617# endif
1618#else
1619 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1620#endif
1621 } else {
1622 throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1623 }
1624
1625 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1626
1627 }
1628
1629
1632 const RCP<const Teuchos::Comm<int> >& comm,
1633 const bool binary=false) {
1634 if (lib == Xpetra::UseEpetra) {
1635 // do we need another specialization for <double,int,int> ??
1636 //TEUCHOS_TEST_FOR_EXCEPTION(true, ::Xpetra::Exceptions::BadCast, "Epetra can only be used with Scalar=double and Ordinal=int");
1637#if defined(HAVE_XPETRA_EPETRA) && defined(HAVE_XPETRA_EPETRAEXT)
1638 TEUCHOS_ASSERT(!binary);
1639 Epetra_Map *eMap;
1640 int rv = EpetraExt::MatrixMarketFileToMap(fileName.c_str(), *(Xpetra::toEpetra(comm)), eMap);
1641 if (rv != 0)
1642 throw Exceptions::RuntimeError("Error reading map from file " + fileName + " with EpetraExt::MatrixMarketToMap (returned " + Teuchos::toString(rv) + ")");
1643
1644 RCP<Epetra_Map> eMap1 = rcp(new Epetra_Map(*eMap));
1645 return Xpetra::toXpetra<int,Node>(*eMap1);
1646#else
1647 throw Exceptions::RuntimeError("Xpetra has not been compiled with Epetra and EpetraExt support.");
1648#endif
1649 } else if (lib == Xpetra::UseTpetra) {
1650#ifdef HAVE_XPETRA_TPETRA
1651# if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
1652 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
1653 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra GO=int support.");
1654# else
1655 typedef Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> sparse_matrix_type;
1656 typedef Tpetra::MatrixMarket::Reader<sparse_matrix_type> reader_type;
1657
1658 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > tMap = reader_type::readMapFile(fileName, comm, false, false, binary);
1659 if (tMap.is_null())
1660 throw Exceptions::RuntimeError("The Tpetra::Map returned from readSparseFile() is null.");
1661
1662 return Xpetra::toXpetra(tMap);
1663# endif
1664#else
1665 throw Exceptions::RuntimeError("Xpetra has not been compiled with Tpetra support.");
1666#endif
1667 } else {
1668 throw Exceptions::RuntimeError("Utils::Read : you must specify Xpetra::UseEpetra or Xpetra::UseTpetra.");
1669 }
1670
1671 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
1672
1673 }
1674
1691#include "Xpetra_UseShortNames.hpp"
1693
1694 size_t numBlocks = 2; // TODO user parameter?
1695
1696 std::vector<RCP<const Map> > rangeMapVec;
1697 for(size_t row = 0; row < numBlocks; ++row) {
1698 RCP<const Map> map = XpIO::ReadMap("subRangeMap_" + fileName + XpIO::toString<size_t>(row) + ".m", lib, comm);
1699 rangeMapVec.push_back(map);
1700 }
1701 RCP<const Map> fullRangeMap = XpIO::ReadMap("fullRangeMap_" + fileName + ".m", lib, comm);
1702
1703 std::vector<RCP<const Map> > domainMapVec;
1704 for(size_t col = 0; col < numBlocks; ++col) {
1705 RCP<const Map> map = XpIO::ReadMap("subDomainMap_" + fileName + XpIO::toString<size_t>(col) + ".m", lib, comm);
1706 domainMapVec.push_back(map);
1707 }
1708 RCP<const Map> fullDomainMap = XpIO::ReadMap("fullDomainMap_" + fileName + ".m", lib, comm);
1709
1710 /*std::vector<RCP<const XpMap> > testRgMapVec;
1711 for(size_t r = 0; r < numBlocks; ++r) {
1712 RCP<const XpMap> map = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(r) + "0.m", lib, comm);
1713 testRgMapVec.push_back(map);
1714 }
1715 std::vector<RCP<const XpMap> > testDoMapVec;
1716 for(size_t c = 0; c < numBlocks; ++c) {
1717 RCP<const XpMap> map = XpIO::ReadMap("domainmap_" + fileName + "0" + XpIO::toString<size_t>(c) + ".m", lib, comm);
1718 testDoMapVec.push_back(map);
1719 }*/
1720
1721 // create map extractors
1722
1723 // range map extractor
1724 bool bRangeUseThyraStyleNumbering = false;
1725 /*
1726 GlobalOrdinal gMinGids = 0;
1727 for(size_t v = 0; v < testRgMapVec.size(); ++v) {
1728 gMinGids += testRgMapVec[v]->getMinAllGlobalIndex();
1729 }
1730 if ( gMinGids==0 && testRgMapVec.size() > 1 ) bRangeUseThyraStyleNumbering = true;*/
1731 RCP<const MapExtractor> rangeMapExtractor =
1732 rcp(new MapExtractor(fullRangeMap, rangeMapVec, bRangeUseThyraStyleNumbering));
1733
1734 // domain map extractor
1735 bool bDomainUseThyraStyleNumbering = false;
1736 /*gMinGids = 0;
1737 for(size_t v = 0; v < testDoMapVec.size(); ++v) {
1738 gMinGids += testDoMapVec[v]->getMinAllGlobalIndex();
1739 }
1740 if ( gMinGids==0 && testDoMapVec.size() > 1) bDomainUseThyraStyleNumbering = true;*/
1741 RCP<const MapExtractor> domainMapExtractor =
1742 rcp(new MapExtractor(fullDomainMap, domainMapVec, bDomainUseThyraStyleNumbering));
1743
1744 RCP<BlockedCrsMatrix> bOp = Teuchos::rcp(new BlockedCrsMatrix(rangeMapExtractor, domainMapExtractor, 33));
1745
1746 // Read all matrices with their maps and create the BlockedCrsMatrix
1747 for (size_t row = 0; row < numBlocks; ++row) {
1748 for (size_t col = 0; col < numBlocks; ++col) {
1749 RCP<const Map> rowSubMap = XpIO::ReadMap("rowmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1750 RCP<const Map> colSubMap = XpIO::ReadMap("colmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1751 RCP<const Map> domSubMap = XpIO::ReadMap("domainmap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1752 RCP<const Map> ranSubMap = XpIO::ReadMap("rangemap_" + fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", lib, comm);
1753 RCP<Matrix> mat = XpIO::Read(fileName + XpIO::toString<size_t>(row) + XpIO::toString<size_t>(col) + ".m", rowSubMap, colSubMap, domSubMap, ranSubMap);
1754 bOp->setMatrix(row, col, mat);
1755 }
1756 }
1757
1758 bOp->fillComplete();
1759
1760 return bOp;
1761 } // ReadBlockedCrsMatrix
1762
1764 template<class T>
1765 static std::string toString(const T& what) {
1766 std::ostringstream buf;
1767 buf << what;
1768 return buf.str();
1769 }
1770 };
1771#endif // HAVE_XPETRA_EPETRA
1772
1773
1774} // end namespace Xpetra
1775
1776#define XPETRA_IO_SHORT
1777
1778#endif /* PACKAGES_XPETRA_SUP_UTILS_XPETRA_IO_HPP_ */
size_type size() const
void resize(const size_type n, const T &val=T())
void assign(size_type n, const T &val)
iterator begin() const
iterator end() const
size_type size() const
void resize(size_type new_size, const value_type &x=value_type())
bool readFile(ArrayRCP< Ordinal > &rowptr, ArrayRCP< Ordinal > &colind, ArrayRCP< Scalar > &values, Ordinal &numRows, Ordinal &numCols, const std::string &filename)
void writeFile(const std::string &filename, const ArrayView< const OrdinalType > &rowptr, const ArrayView< const OrdinalType > &colind, const ArrayView< const ScalarType > &values, const OrdinalType numRows, const OrdinalType numCols)
bool is_null() const
Teuchos::RCP< Matrix > getMatrix(size_t r, size_t c) const
return block (r,c)
virtual size_t Rows() const
number of row blocks
RCP< const MapExtractor > getRangeMapExtractor() const
Returns map extractor class for range map.
RCP< const MapExtractor > getDomainMapExtractor() const
Returns map extractor for domain map.
virtual size_t Cols() const
number of column blocks
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const =0
Returns the Map associated with the domain of this graph.
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const =0
Returns the Map that describes the row distribution in this graph.
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const =0
Returns the Map associated with the domain of this graph.
virtual RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const =0
Returns the Map that describes the column distribution in this graph.
Concrete implementation of Xpetra::Matrix.
RCP< CrsMatrix > getCrsMatrix() const
virtual Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const =0
The Map describing the parallel distribution of this object.
Exception indicating invalid cast attempted.
Exception throws to report errors in the internal logical of the program.
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.
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)
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.
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.
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.
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.
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.
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.
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.
static std::string toString(const T &what)
Little helper function to convert non-string types to strings.
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.
static void Write(const std::string &fileName, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &M)
Read/Write methods.
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.
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)
static void Write(const std::string &fileName, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
Save vector to file in Matrix Market format.
Xpetra utility class containing IO routines to read/write vectors, matrices etc...
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.
static RCP< MultiVector > ReadMultiVector(const std::string &fileName, const RCP< const Map > &map, const bool binary=false)
static RCP< const Map > ReadMap(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int > > &comm, const bool binary=false)
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.
static std::string toString(const T &what)
Little helper function to convert non-string types to strings.
static void Write(const std::string &fileName, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
Save vector to file in Matrix Market format.
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.
static void Write(const std::string &fileName, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &M)
Read/Write methods.
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.
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.
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.
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.
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.
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.
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.
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.
static RCP< Matrix > Build(const RCP< const Map > &rowMap)
Xpetra-specific matrix class.
virtual const RCP< const Map > & getColMap() const
Returns the Map that describes the column distribution in this matrix. This might be null until fillC...
virtual const RCP< const Map > & getRowMap() const
Returns the Map that describes the row distribution in this matrix.
virtual Teuchos::RCP< const Map > getRangeMap() const =0
The Map associated with the range of this operator, which must be compatible with Y....
virtual Teuchos::RCP< const Map > getDomainMap() const =0
The Map associated with the domain of this operator, which must be compatible with X....
#define TEUCHOS_ASSERT(assertion_test)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)
#define TEUCHOS_UNREACHABLE_RETURN(dummyReturnVal)
std::string toString(const T &t)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Xpetra namespace
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > Convert_Epetra_MultiVector_ToXpetra_MultiVector(RCP< Epetra_MultiVector > &epX)
RCP< Xpetra::MultiVector< double, int, int, Xpetra::EpetraNode > > Convert_Epetra_MultiVector_ToXpetra_MultiVector< double, int, int, Xpetra::EpetraNode >(RCP< Epetra_MultiVector > &epX)
RCP< Xpetra::CrsMatrixWrap< double, int, int, Xpetra::EpetraNode > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap< double, int, int, Xpetra::EpetraNode >(RCP< Epetra_CrsMatrix > &epAB)
const Epetra_CrsGraph & toEpetra(const RCP< const CrsGraph< int, GlobalOrdinal, Node > > &graph)
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > toTpetra(const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph)
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &)