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"
53 #include "Xpetra_Matrix.hpp"
54 #include "Xpetra_CrsMatrixWrap.hpp"
55 #include "Xpetra_Map.hpp"
56 #include "Xpetra_Vector.hpp"
57 #include "Xpetra_Exceptions.hpp"
58 
59 namespace Xpetra {
60  template <class Scalar = Matrix<>::scalar_type,
61  class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type,
62  class GlobalOrdinal =
64  class Node =
66  class MatrixFactory {
67 #undef XPETRA_MATRIXFACTORY_SHORT
68 #include "Xpetra_UseShortNames.hpp"
69 
70  private:
73 
74  public:
75 
77  static RCP<Matrix> Build(const RCP<const Map>& rowMap, size_t maxNumEntriesPerRow, Xpetra::ProfileType pftype = Xpetra::DynamicProfile) {
78  return rcp(new CrsMatrixWrap(rowMap, maxNumEntriesPerRow, pftype));
79  }
80 
82  static RCP<Matrix> Build(const RCP<const Map>& rowMap, const RCP<const Map>& colMap, size_t maxNumEntriesPerRow, Xpetra::ProfileType pftype = Xpetra::DynamicProfile) {
83  return rcp(new CrsMatrixWrap(rowMap, colMap, maxNumEntriesPerRow, pftype));
84  }
85 
87  static RCP<Matrix> Build(const RCP<const Map>& rowMap, const RCP<const Map>& colMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, Xpetra::ProfileType pftype = Xpetra::DynamicProfile) {
88  return rcp(new CrsMatrixWrap(rowMap, colMap, NumEntriesPerRowToAlloc, pftype));
89  }
90 
91 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
92  static RCP<Matrix> Build (
94  const Teuchos::RCP<const Map>& rowMap,
95  const Teuchos::RCP<const Map>& colMap,
97  const Teuchos::RCP<Teuchos::ParameterList>& params = null) {
98  XPETRA_MONITOR("MatrixFactory::Build");
99  return rcp(new CrsMatrixWrap(rowMap, colMap, lclMatrix, params));
100  }
101 #endif
102 
104  static RCP<Matrix> Build(const RCP<const Map> &rowMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc, ProfileType pftype = Xpetra::DynamicProfile) {
105  return rcp( new CrsMatrixWrap(rowMap, NumEntriesPerRowToAlloc, pftype) );
106  }
107 
109  static RCP<Matrix> Build(const RCP<const CrsGraph>& graph, const RCP<ParameterList>& paramList = Teuchos::null) {
110  return rcp(new CrsMatrixWrap(graph, paramList));
111  }
112 
114  static RCP<Matrix> Build(const RCP<const Vector>& diagonal) {
115  Teuchos::ArrayRCP<const Scalar> vals = diagonal->getData(0);
116  LocalOrdinal NumMyElements = diagonal->getMap()->getNodeNumElements();
117  Teuchos::ArrayView<const GlobalOrdinal> MyGlobalElements = diagonal->getMap()->getNodeElementList();
118 
120 
121  for (LocalOrdinal i = 0; i < NumMyElements; ++i) {
122  mtx->insertGlobalValues(MyGlobalElements[i],
123  Teuchos::tuple<GlobalOrdinal>(MyGlobalElements[i]),
124  Teuchos::tuple<Scalar>(vals[i]) );
125  }
126  mtx->fillComplete();
127  return mtx;
128  }
129 
131  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) {
132  RCP<const CrsMatrixWrap> crsOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(sourceMatrix);
133  if (crsOp == Teuchos::null)
134  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
135 
136  RCP<CrsMatrix> originalCrs = crsOp->getCrsMatrix();
137  RCP<CrsMatrix> newCrs = CrsMatrixFactory::Build(originalCrs, importer, domainMap, rangeMap, params);
138  if (newCrs->hasMatrix())
139  return rcp(new CrsMatrixWrap(newCrs));
140  else
141  return Teuchos::null;
142  }
143 
145  static RCP<Matrix> Build(const RCP<const Matrix> & sourceMatrix, const Export &exporter, const RCP<const Map> & domainMap = Teuchos::null, const RCP<const Map> & rangeMap = Teuchos::null,const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null) {
146  RCP<const CrsMatrixWrap> crsOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(sourceMatrix);
147  if (crsOp == Teuchos::null)
148  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
149 
150  RCP<CrsMatrix> originalCrs = crsOp->getCrsMatrix();
151  return rcp(new CrsMatrixWrap(CrsMatrixFactory::Build(originalCrs, exporter, domainMap, rangeMap, params)));
152  }
153  };
154 #define XPETRA_MATRIXFACTORY_SHORT
155 
156 
157  /*template <class Scalar = Matrix<>::scalar_type,
158  class LocalOrdinal = typename Matrix<Scalar>::local_ordinal_type,
159  class GlobalOrdinal =
160  typename Matrix<Scalar, LocalOrdinal>::global_ordinal_type,
161  class Node =
162  typename Matrix<Scalar, LocalOrdinal, GlobalOrdinal>::node_type>*/
163  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
164  class MatrixFactory2 {
165 #undef XPETRA_MATRIXFACTORY2_SHORT
166 #include "Xpetra_UseShortNames.hpp"
167 
168  public:
170  RCP<const CrsMatrixWrap> oldOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(A);
171  if (oldOp == Teuchos::null)
172  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
173 
174  RCP<const CrsMatrix> oldCrsOp = oldOp->getCrsMatrix();
175 
176  UnderlyingLib lib = A->getRowMap()->lib();
177 
179  "Not Epetra or Tpetra matrix");
180 
181 #ifdef HAVE_XPETRA_EPETRA
182  if (lib == UseEpetra) {
183  // NOTE: The proper Epetra conversion in Xpetra_MatrixFactory.cpp
184  throw Exceptions::RuntimeError("Xpetra::BuildCopy(): matrix templates are incompatible with Epetra");
185  }
186 #endif
187 
188 #ifdef HAVE_XPETRA_TPETRA
189  if (lib == UseTpetra) {
190  // Underlying matrix is Tpetra
191  RCP<const TpetraCrsMatrix> oldTCrsOp = Teuchos::rcp_dynamic_cast<const TpetraCrsMatrix>(oldCrsOp);
192 
193  if (oldTCrsOp != Teuchos::null) {
194  RCP<TpetraCrsMatrix> newTCrsOp(new TpetraCrsMatrix(*oldTCrsOp));
196 
197  return newOp;
198  } else {
199  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::TpetraCrsMatrix failed");
200  }
201  }
202 #endif
203 
204  return Teuchos::null;
205  }
206  };
207 #define XPETRA_MATRIXFACTORY2_SHORT
208 
209  //template<>
210  //class MatrixFactory2<double,int,int,typename Xpetra::Matrix<double, int, int>::node_type> {
211  template<class Node>
212  class MatrixFactory2<double,int,int,Node> {
213  typedef double Scalar;
214  typedef int LocalOrdinal;
215  typedef int GlobalOrdinal;
216  //typedef Matrix<double, int, GlobalOrdinal>::node_type Node;
217 #undef XPETRA_MATRIXFACTORY2_SHORT
218 #include "Xpetra_UseShortNames.hpp"
219  public:
221  RCP<const CrsMatrixWrap> oldOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(A);
222  if (oldOp == Teuchos::null)
223  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
224 
225  RCP<const CrsMatrix> oldCrsOp = oldOp->getCrsMatrix();
226 
227  #ifdef HAVE_XPETRA_EPETRA
228  #ifndef XPETRA_EPETRA_NO_32BIT_GLOBAL_INDICES
229  RCP<const EpetraCrsMatrixT<GlobalOrdinal,Node> > oldECrsOp = Teuchos::rcp_dynamic_cast<const EpetraCrsMatrixT<GlobalOrdinal,Node> >(oldCrsOp);
230  if (oldECrsOp != Teuchos::null) {
231  // Underlying matrix is Epetra
232  RCP<CrsMatrix> newECrsOp(new EpetraCrsMatrixT<GlobalOrdinal,Node>(*oldECrsOp));
233  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap (newECrsOp));
234 
235  return newOp;
236  }
237  #endif
238  #endif
239 
240  #ifdef HAVE_XPETRA_TPETRA
241  // Underlying matrix is Tpetra
242  RCP<const TpetraCrsMatrix> oldTCrsOp = Teuchos::rcp_dynamic_cast<const TpetraCrsMatrix>(oldCrsOp);
243  if (oldTCrsOp != Teuchos::null) {
244  RCP<CrsMatrix> newTCrsOp(new TpetraCrsMatrix(*oldTCrsOp));
245  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap(newTCrsOp));
246 
247  return newOp;
248  }
249  #else
250  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::EpetraCrsMatrix or Xpetra::TpetraCrsMatrix failed");
251  #endif
252 
253  return Teuchos::null; // make compiler happy
254  }
255  };
256 
257 #define XPETRA_MATRIXFACTORY2_SHORT
258 
259 #ifdef HAVE_XPETRA_INT_LONG_LONG
260  //template<>
261  //class MatrixFactory2<double,int,long long,typename Xpetra::Matrix<double, int, long long>::node_type> {
262  template<class Node>
263  class MatrixFactory2<double, int, long long, Node> {
264  typedef double Scalar;
265  typedef int LocalOrdinal;
266  typedef long long GlobalOrdinal;
267  //typedef Matrix<double, int, GlobalOrdinal>::node_type Node;
268 #undef XPETRA_MATRIXFACTORY2_SHORT
269 #include "Xpetra_UseShortNames.hpp"
270  public:
272  RCP<const CrsMatrixWrap> oldOp = Teuchos::rcp_dynamic_cast<const CrsMatrixWrap>(A);
273  if (oldOp == Teuchos::null)
274  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
275 
276  RCP<const CrsMatrix> oldCrsOp = oldOp->getCrsMatrix();
277 
278  #ifdef HAVE_XPETRA_EPETRA
279  #ifndef XPETRA_EPETRA_NO_64BIT_GLOBAL_INDICES
280  RCP<const EpetraCrsMatrixT<GlobalOrdinal,Node> > oldECrsOp = Teuchos::rcp_dynamic_cast<const EpetraCrsMatrixT<GlobalOrdinal,Node> >(oldCrsOp);
281  if (oldECrsOp != Teuchos::null) {
282  // Underlying matrix is Epetra
283  RCP<CrsMatrix> newECrsOp(new EpetraCrsMatrixT<GlobalOrdinal,Node>(*oldECrsOp));
284  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap (newECrsOp));
285 
286  return newOp;
287  }
288  #endif
289  #endif
290 
291  #ifdef HAVE_XPETRA_TPETRA
292  // Underlying matrix is Tpetra
293  RCP<const TpetraCrsMatrix> oldTCrsOp = Teuchos::rcp_dynamic_cast<const TpetraCrsMatrix>(oldCrsOp);
294  if (oldTCrsOp != Teuchos::null) {
295  RCP<CrsMatrix> newTCrsOp(new TpetraCrsMatrix(*oldTCrsOp));
296  RCP<CrsMatrixWrap> newOp (new CrsMatrixWrap(newTCrsOp));
297 
298  return newOp;
299  }
300  #else
301  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::EpetraCrsMatrix or Xpetra::TpetraCrsMatrix failed");
302  #endif
303 
304  return Teuchos::null; // make compiler happy
305  }
306  };
307 #endif // HAVE_XPETRA_INT_LONG_LONG
308 
309 #define XPETRA_MATRIXFACTORY2_SHORT
310 
311 }
312 
313 #define XPETRA_MATRIXFACTORY_SHORT
314 #define XPETRA_MATRIXFACTORY2_SHORT
315 #endif
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< 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.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
static RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Build(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, Xpetra::ProfileType pftype=Xpetra::DynamicProfile, const Teuchos::RCP< Teuchos::ParameterList > &plist=Teuchos::null)
Constructor specifying fixed number of entries for each row.
Xpetra namespace
LocalOrdinal local_ordinal_type
Exception throws to report errors in the internal logical of the program.
static RCP< Matrix > Build(const RCP< const Map > &rowMap, const RCP< const Map > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, Xpetra::ProfileType pftype=Xpetra::DynamicProfile)
Constructor specifying the (possibly different) number of entries per row and providing column map...
static RCP< Matrix > Build(const RCP< const Map > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, ProfileType pftype=Xpetra::DynamicProfile)
Constructor specifying (possibly different) number of entries in each row.
Exception indicating invalid cast attempted.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
static RCP< Matrix > Build(const RCP< const Map > &rowMap, const RCP< const Map > &colMap, size_t maxNumEntriesPerRow, Xpetra::ProfileType pftype=Xpetra::DynamicProfile)
Constructor specifying the max number of non-zeros 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)
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > BuildCopy(const RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > A)
TypeTo as(const TypeFrom &t)
Concrete implementation of Xpetra::Matrix.
#define XPETRA_MONITOR(funcName)
static RCP< Matrix > Build(const RCP< const Map > &rowMap, size_t maxNumEntriesPerRow, Xpetra::ProfileType pftype=Xpetra::DynamicProfile)
Constructor specifying the number of non-zeros for all rows.
Xpetra-specific matrix class.
MatrixFactory()
Private constructor. This is a static class.
static RCP< Matrix > Build(const RCP< const Matrix > &sourceMatrix, const Export &exporter, 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...
GlobalOrdinal global_ordinal_type