Xpetra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Xpetra_MatrixFactory.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 //
43 // ***********************************************************************
44 //
45 // @HEADER
46 
47 // WARNING: This code is experimental. Backwards compatibility should not be expected.
48 
49 #ifndef XPETRA_MATRIXFACTORY_HPP
50 #define XPETRA_MATRIXFACTORY_HPP
51 
52 #include "Xpetra_ConfigDefs.hpp"
54 #include "Xpetra_Matrix.hpp"
55 #include "Xpetra_CrsMatrixWrap.hpp"
57 #include "Xpetra_Map.hpp"
58 #include "Xpetra_Vector.hpp"
59 #include "Xpetra_Exceptions.hpp"
60 
61 namespace Xpetra {
62 
63  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node = KokkosClassic::DefaultNode::DefaultNodeType>
64  class MatrixFactory2 {
65 #undef XPETRA_MATRIXFACTORY2_SHORT
66 #include "Xpetra_UseShortNames.hpp"
67 
68  public:
69  static RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > BuildCopy(const RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A, bool setFixedBlockSize = true) {
70  RCP<const CrsMatrixWrap> oldOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(A);
71  if (oldOp == Teuchos::null)
72  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
73 
74  RCP<const CrsMatrix> oldCrsOp = oldOp->getCrsMatrix();
75 
76  UnderlyingLib lib = A->getRowMap()->lib();
77 
78  TEUCHOS_TEST_FOR_EXCEPTION(lib != UseEpetra && lib != UseTpetra, Exceptions::RuntimeError,
79  "Not Epetra or Tpetra matrix");
80 
81 #ifdef HAVE_XPETRA_EPETRA
82  if (lib == UseEpetra) {
83  // NOTE: The proper Epetra conversion in Xpetra_MatrixFactory.cpp
84  throw Exceptions::RuntimeError("Xpetra::BuildCopy(): matrix templates are incompatible with Epetra");
85  }
86 #endif
87 
88 #ifdef HAVE_XPETRA_TPETRA
89  if (lib == UseTpetra) {
90  // Underlying matrix is Tpetra
91  RCP<const TpetraCrsMatrix> oldTCrsOp = Teuchos::rcp_dynamic_cast<const TpetraCrsMatrix>(oldCrsOp);
92 
93  if (oldTCrsOp != Teuchos::null) {
94  RCP<TpetraCrsMatrix> newTCrsOp(new TpetraCrsMatrix(*oldTCrsOp));
95  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap(Teuchos::as<RCP<CrsMatrix> >(newTCrsOp)));
96  if (setFixedBlockSize)
97  newOp->SetFixedBlockSize(A->GetFixedBlockSize());
98 
99  return newOp;
100  } else {
101  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::TpetraCrsMatrix failed");
102  }
103  }
104 #endif
105 
106  return Teuchos::null;
107  }
108  };
109 #define XPETRA_MATRIXFACTORY2_SHORT
110 
111  //template<>
112  //class MatrixFactory2<double,int,int,typename Xpetra::Matrix<double, int, int>::node_type> {
113  template<class Node>
114  class MatrixFactory2<double,int,int,Node> {
115  typedef double Scalar;
116  typedef int LocalOrdinal;
117  typedef int GlobalOrdinal;
118  //typedef Matrix<double, int, GlobalOrdinal>::node_type Node;
119 #undef XPETRA_MATRIXFACTORY2_SHORT
120 #include "Xpetra_UseShortNames.hpp"
121  public:
122  static RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > BuildCopy(const RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A, bool setFixedBlockSize = true) {
123  RCP<const CrsMatrixWrap> oldOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(A);
124  if (oldOp == Teuchos::null)
125  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
126 
127  RCP<const CrsMatrix> oldCrsOp = oldOp->getCrsMatrix();
128 
129 #ifdef HAVE_XPETRA_EPETRA
130 #ifndef XPETRA_EPETRA_NO_32BIT_GLOBAL_INDICES
131  RCP<const EpetraCrsMatrixT<GlobalOrdinal,Node> > oldECrsOp = Teuchos::rcp_dynamic_cast<const EpetraCrsMatrixT<GlobalOrdinal,Node> >(oldCrsOp);
132  if (oldECrsOp != Teuchos::null) {
133  // Underlying matrix is Epetra
134  RCP<CrsMatrix> newECrsOp(new EpetraCrsMatrixT<GlobalOrdinal,Node>(*oldECrsOp));
135  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap (newECrsOp));
136  if (setFixedBlockSize)
137  newOp->SetFixedBlockSize(A->GetFixedBlockSize());
138  return newOp;
139  }
140 #endif
141 #endif
142 
143 #ifdef HAVE_XPETRA_TPETRA
144  // Underlying matrix is Tpetra
145  RCP<const TpetraCrsMatrix> oldTCrsOp = Teuchos::rcp_dynamic_cast<const TpetraCrsMatrix>(oldCrsOp);
146  if (oldTCrsOp != Teuchos::null) {
147  RCP<CrsMatrix> newTCrsOp(new TpetraCrsMatrix(*oldTCrsOp));
148  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap(newTCrsOp));
149  if (setFixedBlockSize)
150  newOp->SetFixedBlockSize(A->GetFixedBlockSize());
151  return newOp;
152  }
153  return Teuchos::null;
154 #else
155  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::EpetraCrsMatrix or Xpetra::TpetraCrsMatrix failed");
156  TEUCHOS_UNREACHABLE_RETURN(Teuchos::null); // make compiler happy
157 #endif
158 
159  } //BuildCopy
160  };
161 
162 #define XPETRA_MATRIXFACTORY2_SHORT
163 
164 #ifdef HAVE_XPETRA_INT_LONG_LONG
165  //template<>
166  //class MatrixFactory2<double,int,long long,typename Xpetra::Matrix<double, int, long long>::node_type> {
167  template<class Node>
168  class MatrixFactory2<double, int, long long, Node> {
169  typedef double Scalar;
170  typedef int LocalOrdinal;
171  typedef long long GlobalOrdinal;
172  //typedef Matrix<double, int, GlobalOrdinal>::node_type Node;
173 #undef XPETRA_MATRIXFACTORY2_SHORT
174 #include "Xpetra_UseShortNames.hpp"
175  public:
176  static RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > BuildCopy(const RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A, bool setFixedBlockSize = true) {
177  RCP<const CrsMatrixWrap> oldOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(A);
178  if (oldOp == Teuchos::null)
179  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
180 
181  RCP<const CrsMatrix> oldCrsOp = oldOp->getCrsMatrix();
182 
183 #ifdef HAVE_XPETRA_EPETRA
184 #ifndef XPETRA_EPETRA_NO_64BIT_GLOBAL_INDICES
185  RCP<const EpetraCrsMatrixT<GlobalOrdinal,Node> > oldECrsOp = Teuchos::rcp_dynamic_cast<const EpetraCrsMatrixT<GlobalOrdinal,Node> >(oldCrsOp);
186  if (oldECrsOp != Teuchos::null) {
187  // Underlying matrix is Epetra
188  RCP<CrsMatrix> newECrsOp(new EpetraCrsMatrixT<GlobalOrdinal,Node>(*oldECrsOp));
189  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap (newECrsOp));
190  if (setFixedBlockSize)
191  newOp->SetFixedBlockSize(A->GetFixedBlockSize());
192  return newOp;
193  }
194 #endif
195 #endif
196 
197 #ifdef HAVE_XPETRA_TPETRA
198  // Underlying matrix is Tpetra
199  RCP<const TpetraCrsMatrix> oldTCrsOp = Teuchos::rcp_dynamic_cast<const TpetraCrsMatrix>(oldCrsOp);
200  if (oldTCrsOp != Teuchos::null) {
201  RCP<CrsMatrix> newTCrsOp(new TpetraCrsMatrix(*oldTCrsOp));
202  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap(newTCrsOp));
203  if (setFixedBlockSize)
204  newOp->SetFixedBlockSize(A->GetFixedBlockSize());
205  return newOp;
206  }
207 #else
208  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::EpetraCrsMatrix or Xpetra::TpetraCrsMatrix failed");
209 #endif
210 
211  return Teuchos::null; // make compiler happy
212  }
213  };
214 #endif // HAVE_XPETRA_INT_LONG_LONG
215 
216 #define XPETRA_MATRIXFACTORY2_SHORT
217 
218 
219  template <class Scalar,
220  class LocalOrdinal,
221  class GlobalOrdinal,
222  class Node>
223  class MatrixFactory {
224 #undef XPETRA_MATRIXFACTORY_SHORT
225 #include "Xpetra_UseShortNames.hpp"
226 
227  private:
230 
231  public:
234  static RCP<Matrix> Build(const RCP<const Map>& rowMap) {
235  return rcp(new CrsMatrixWrap(rowMap));
236  }
237 
239  static RCP<Matrix> Build(const RCP<const Map>& rowMap, size_t maxNumEntriesPerRow) {
240  return rcp(new CrsMatrixWrap(rowMap, maxNumEntriesPerRow));
241  }
242 
244  static RCP<Matrix> Build(const RCP<const Map>& rowMap, const RCP<const Map>& colMap, size_t maxNumEntriesPerRow) {
245  return rcp(new CrsMatrixWrap(rowMap, colMap, maxNumEntriesPerRow));
246  }
247 
249  static RCP<Matrix> Build(const RCP<const Map>& rowMap, const RCP<const Map>& colMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc) {
250  return rcp(new CrsMatrixWrap(rowMap, colMap, NumEntriesPerRowToAlloc));
251  }
252 
253 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
254  static RCP<Matrix> Build (
256  const Teuchos::RCP<const Map>& rowMap,
257  const Teuchos::RCP<const Map>& colMap,
259  const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
260  XPETRA_MONITOR("MatrixFactory::Build");
261  return rcp(new CrsMatrixWrap(rowMap, colMap, lclMatrix, params));
262  }
264  static RCP<Matrix> Build (
266  const Teuchos::RCP<const Map>& rowMap,
267  const Teuchos::RCP<const Map>& colMap,
268  const Teuchos::RCP<const Map>& domainMap = Teuchos::null,
269  const Teuchos::RCP<const Map>& rangeMap = Teuchos::null,
270  const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
271  XPETRA_MONITOR("MatrixFactory::Build");
272  return rcp(new CrsMatrixWrap(lclMatrix, rowMap, colMap, domainMap, rangeMap, params));
273  }
274 #endif
275 
277  static RCP<Matrix> Build(const RCP<const Map> &rowMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc) {
278  return rcp( new CrsMatrixWrap(rowMap, NumEntriesPerRowToAlloc) );
279  }
280 
282  static RCP<Matrix> Build(const RCP<const CrsGraph>& graph, const RCP<ParameterList>& paramList = Teuchos::null) {
283  return rcp(new CrsMatrixWrap(graph, paramList));
284  }
285 
287  static RCP<Matrix> Build(const RCP<const Vector>& diagonal) {
288  Teuchos::ArrayRCP<const Scalar> vals = diagonal->getData(0);
289  LocalOrdinal NumMyElements = diagonal->getMap()->getNodeNumElements();
290  Teuchos::ArrayView<const GlobalOrdinal> MyGlobalElements = diagonal->getMap()->getNodeElementList();
291 
292  Teuchos::RCP<CrsMatrixWrap> mtx = Teuchos::rcp(new CrsMatrixWrap(diagonal->getMap(), 1));
293 
294  for (LocalOrdinal i = 0; i < NumMyElements; ++i) {
295  mtx->insertGlobalValues(MyGlobalElements[i],
296  Teuchos::tuple<GlobalOrdinal>(MyGlobalElements[i]),
297  Teuchos::tuple<Scalar>(vals[i]) );
298  }
299  mtx->fillComplete();
300  return mtx;
301  }
302 
304  static RCP<Matrix> Build(const RCP<const Matrix>& sourceMatrix, const Import& importer, const RCP<const Map>& domainMap = Teuchos::null, const RCP<const Map>& rangeMap = Teuchos::null, const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
305  RCP<const CrsMatrixWrap> crsOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(sourceMatrix);
306  if (crsOp == Teuchos::null)
307  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
308 
309  RCP<CrsMatrix> originalCrs = crsOp->getCrsMatrix();
310  RCP<CrsMatrix> newCrs = CrsMatrixFactory::Build(originalCrs, importer, domainMap, rangeMap, params);
311  if (newCrs->hasMatrix())
312  return rcp(new CrsMatrixWrap(newCrs));
313  else
314  return Teuchos::null;
315  }
316 
318  static RCP<Matrix> Build(const RCP<const Matrix> & sourceMatrix, const Export &exporter, const RCP<const Map> & domainMap, const RCP<const Map> & rangeMap,const Teuchos::RCP<Teuchos::ParameterList>& params) {
319  RCP<const CrsMatrixWrap> crsOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(sourceMatrix);
320  if (crsOp == Teuchos::null)
321  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
322 
323  RCP<CrsMatrix> originalCrs = crsOp->getCrsMatrix();
324  return rcp(new CrsMatrixWrap(CrsMatrixFactory::Build(originalCrs, exporter, domainMap, rangeMap, params)));
325  }
326 
328  static RCP<Matrix> Build(const RCP<const Matrix>& sourceMatrix, const Import& RowImporter, const Import& DomainImporter, const RCP<const Map>& domainMap, const RCP<const Map>& rangeMap, const Teuchos::RCP<Teuchos::ParameterList>& params) {
329  RCP<const CrsMatrixWrap> crsOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(sourceMatrix);
330  if (crsOp == Teuchos::null)
331  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
332 
333  RCP<CrsMatrix> originalCrs = crsOp->getCrsMatrix();
334  RCP<CrsMatrix> newCrs = CrsMatrixFactory::Build(originalCrs, RowImporter, Teuchos::rcpFromRef(DomainImporter), domainMap, rangeMap, params);
335  if (newCrs->hasMatrix())
336  return rcp(new CrsMatrixWrap(newCrs));
337  else
338  return Teuchos::null;
339  }
340 
342  static RCP<Matrix> Build(const RCP<const Matrix> & sourceMatrix, const Export &RowExporter, const Export &DomainExporter, const RCP<const Map> & domainMap = Teuchos::null, const RCP<const Map> & rangeMap = Teuchos::null,const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
343  RCP<const CrsMatrixWrap> crsOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(sourceMatrix);
344  if (crsOp == Teuchos::null)
345  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
346 
347  RCP<CrsMatrix> originalCrs = crsOp->getCrsMatrix();
348  RCP<CrsMatrix> newCrs = CrsMatrixFactory::Build(originalCrs, RowExporter, Teuchos::rcpFromRef(DomainExporter), domainMap, rangeMap, params);
349  if (newCrs->hasMatrix())
350  return rcp(new CrsMatrixWrap(newCrs));
351  else
352  return Teuchos::null;
353  }
354 
355 
358  static RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > BuildCopy(const RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A, bool setFixedBlockSize = true) {
359  RCP<const BlockedCrsMatrix> input = Teuchos::rcp_dynamic_cast<const BlockedCrsMatrix>(A);
360  if(input == Teuchos::null)
362 
363  // deep copy of MapExtractors (and underlying maps)
364  RCP<const MapExtractor> rgMapExt = Teuchos::rcp(new MapExtractor(*(input->getRangeMapExtractor())));
365  RCP<const MapExtractor> doMapExt = Teuchos::rcp(new MapExtractor(*(input->getDomainMapExtractor())));
366 
367  // create new BlockedCrsMatrix object
368  RCP<BlockedCrsMatrix> bop = Teuchos::rcp(new BlockedCrsMatrix(rgMapExt, doMapExt, input->getNodeMaxNumRowEntries()));
369 
370  for (size_t r = 0; r < input->Rows(); ++r) {
371  for (size_t c = 0; c < input->Cols(); ++c)
372  if(input->getMatrix(r,c) != Teuchos::null) {
373  // make a deep copy of the matrix
374  // This is a recursive call to this function
375  RCP<Matrix> mat =
376  Xpetra::MatrixFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::BuildCopy(input->getMatrix(r,c),setFixedBlockSize);
377  bop->setMatrix(r,c,mat);
378  }
379  }
380 
381  if(input->isFillComplete())
382  bop->fillComplete();
383  return bop;
384  }
385  };
386 #define XPETRA_MATRIXFACTORY_SHORT
387 
388 }
389 
390 #define XPETRA_MATRIXFACTORY_SHORT
391 #define XPETRA_MATRIXFACTORY2_SHORT
392 #endif
static RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Build(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap)
Constructor for empty matrix (intended use is an import/export target - can&#39;t insert entries directly...
static RCP< Matrix > Build(const RCP< const Vector > &diagonal)
Constructor for creating a diagonal Xpetra::Matrix using the entries of a given vector for the diagon...
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > BuildCopy(const RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > A, bool setFixedBlockSize=true)
static RCP< Matrix > Build(const RCP< const Matrix > &sourceMatrix, const Import &importer, const RCP< const Map > &domainMap=Teuchos::null, const RCP< const Map > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor to create a Matrix using a fusedImport-style construction. The originalMatrix must be a X...
static RCP< Matrix > Build(const RCP< const CrsGraph > &graph, const RCP< ParameterList > &paramList=Teuchos::null)
Constructor specifying graph.
static RCP< Matrix > Build(const RCP< const Map > &rowMap, const RCP< const Map > &colMap, size_t maxNumEntriesPerRow)
Constructor specifying the max number of non-zeros per row and providing column map.
static RCP< Matrix > Build(const RCP< const Map > &rowMap, size_t maxNumEntriesPerRow)
Constructor specifying the number of non-zeros for all rows.
static RCP< Matrix > Build(const RCP< const Map > &rowMap)
Exception throws to report errors in the internal logical of the program.
Exception indicating invalid cast attempted.
static RCP< Matrix > Build(const RCP< const Map > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc)
Constructor specifying (possibly different) number of entries in each row.
static RCP< Matrix > Build(const RCP< const Map > &rowMap, const RCP< const Map > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc)
Constructor specifying the (possibly different) number of entries per row and providing column map...
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > BuildCopy(const RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > A, bool setFixedBlockSize=true)
static RCP< Matrix > Build(const RCP< const Matrix > &sourceMatrix, const Export &exporter, const RCP< const Map > &domainMap, const RCP< const Map > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor to create a Matrix using a fusedExport-style construction. The originalMatrix must be a X...
Concrete implementation of Xpetra::Matrix.
#define XPETRA_MONITOR(funcName)
static RCP< Matrix > Build(const RCP< const Matrix > &sourceMatrix, const Export &RowExporter, const Export &DomainExporter, const RCP< const Map > &domainMap=Teuchos::null, const RCP< const Map > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor to create a Matrix using a fusedExport-style construction. The originalMatrix must be a X...
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > BuildCopy(const RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > A, bool setFixedBlockSize=true)
Xpetra-specific matrix class.
MatrixFactory()
Private constructor. This is a static class.
static RCP< Matrix > Build(const RCP< const Matrix > &sourceMatrix, const Import &RowImporter, const Import &DomainImporter, const RCP< const Map > &domainMap, const RCP< const Map > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor to create a Matrix using a fusedImport-style construction. The originalMatrix must be a X...