Xpetra Version of the Day
Loading...
Searching...
No Matches
Xpetra_TpetraCrsMatrix_def.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#ifndef XPETRA_TPETRACRSMATRIX_DEF_HPP
47#define XPETRA_TPETRACRSMATRIX_DEF_HPP
48
49#include <Xpetra_MultiVectorFactory.hpp>
51#include "Tpetra_Details_residual.hpp"
52
53namespace Xpetra {
54
55 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
57 : mtx_ (Teuchos::rcp (new Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> (toTpetra(rowMap), maxNumEntriesPerRow, params))) { }
58
59 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
61 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> (toTpetra(rowMap), NumEntriesPerRowToAlloc(), params))) { }
62
63 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
65 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(toTpetra(rowMap), toTpetra(colMap), maxNumEntriesPerRow, params))) { }
66
67 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
69 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >(toTpetra(rowMap), toTpetra(colMap), NumEntriesPerRowToAlloc(), params))) { }
70
71 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
73 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >(toTpetra(graph), params))) { }
74
75
76
77 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
83 {
84 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
85 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, *sourceMatrix, tSourceMatrix, "Xpetra::TpetraCrsMatrix constructor only accepts Xpetra::TpetraCrsMatrix as the input argument.");//TODO: remove and use toTpetra()
87
88 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ? toTpetra(domainMap) : Teuchos::null;
89 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ? toTpetra(rangeMap) : Teuchos::null;
90 mtx_=Tpetra::importAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),toTpetra(importer),myDomainMap,myRangeMap,params);
91 bool restrictComm=false;
92 if(!params.is_null()) restrictComm = params->get("Restrict Communicator",restrictComm);
93 if(restrictComm && mtx_->getRowMap().is_null()) mtx_=Teuchos::null;
94
95 }
96
97 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
100 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
103 {
104 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
105 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, *sourceMatrix, tSourceMatrix, "Xpetra::TpetraCrsMatrix constructor only accepts Xpetra::TpetraCrsMatrix as the input argument.");//TODO: remove and use toTpetra()
106 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSourceMatrix.getTpetra_CrsMatrix();
107
108 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ? toTpetra(domainMap) : Teuchos::null;
109 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ? toTpetra(rangeMap) : Teuchos::null;
110 mtx_=Tpetra::exportAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),toTpetra(exporter),myDomainMap,myRangeMap,params);
111
112 }
113
114 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
116 const Import<LocalOrdinal,GlobalOrdinal,Node> & RowImporter,
117 const Teuchos::RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > DomainImporter,
118 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
121 {
122 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
123 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, *sourceMatrix, tSourceMatrix, "Xpetra::TpetraCrsMatrix constructor only accepts Xpetra::TpetraCrsMatrix as the input argument.");//TODO: remove and use toTpetra()
124 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSourceMatrix.getTpetra_CrsMatrix();
125
126 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ? toTpetra(domainMap) : Teuchos::null;
127 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ? toTpetra(rangeMap) : Teuchos::null;
128
129 mtx_=Tpetra::importAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),toTpetra(RowImporter),toTpetra(*DomainImporter),myDomainMap,myRangeMap,params);
130 bool restrictComm=false;
131 if(!params.is_null()) restrictComm = params->get("Restrict Communicator",restrictComm);
132 if(restrictComm && mtx_->getRowMap().is_null()) mtx_=Teuchos::null;
133 }
134
135 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
137 const Export<LocalOrdinal,GlobalOrdinal,Node> & RowExporter,
138 const Teuchos::RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > DomainExporter,
139 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
142 {
143 typedef Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> MyTpetraCrsMatrix;
144 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, *sourceMatrix, tSourceMatrix, "Xpetra::TpetraCrsMatrix constructor only accepts Xpetra::TpetraCrsMatrix as the input argument.");//TODO: remove and use toTpetra()
145 RCP< const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > v = tSourceMatrix.getTpetra_CrsMatrix();
146
147 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myDomainMap = domainMap!=Teuchos::null ? toTpetra(domainMap) : Teuchos::null;
148 RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > myRangeMap = rangeMap!=Teuchos::null ? toTpetra(rangeMap) : Teuchos::null;
149
150 mtx_=Tpetra::exportAndFillCompleteCrsMatrix<MyTpetraCrsMatrix>(tSourceMatrix.getTpetra_CrsMatrix(),toTpetra(RowExporter),toTpetra(*DomainExporter),myDomainMap,myRangeMap,params);
151 }
152
154
155
156#ifdef HAVE_XPETRA_TPETRA
157 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
160 const local_matrix_type& lclMatrix,
162 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(toTpetra(rowMap), toTpetra(colMap), lclMatrix, params))) { }
163
164 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
166 const local_matrix_type& lclMatrix,
169 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
172 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(lclMatrix, toTpetra(rowMap), toTpetra(colMap), toTpetra(domainMap), toTpetra(rangeMap), params))) { }
173
174 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
176 const local_matrix_type& lclMatrix,
179 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
184 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(lclMatrix, toTpetra(rowMap), toTpetra(colMap), toTpetra(domainMap), toTpetra(rangeMap), toTpetra(importer), toTpetra(exporter), params))) { }
185#endif
186
187 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
189
190 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
191 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals) { XPETRA_MONITOR("TpetraCrsMatrix::insertGlobalValues"); mtx_->insertGlobalValues(globalRow, cols, vals); }
192
193 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
194 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals) { XPETRA_MONITOR("TpetraCrsMatrix::insertLocalValues"); mtx_->insertLocalValues(localRow, cols, vals); }
195
196 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
197 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals) { XPETRA_MONITOR("TpetraCrsMatrix::replaceGlobalValues"); mtx_->replaceGlobalValues(globalRow, cols, vals); }
198
199 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
200 void
203 const ArrayView<const Scalar> &vals)
204 {
205 XPETRA_MONITOR("TpetraCrsMatrix::replaceLocalValues");
206 typedef typename ArrayView<const LocalOrdinal>::size_type size_type;
207 const LocalOrdinal numValid =
208 mtx_->replaceLocalValues (localRow, cols, vals);
210 static_cast<size_type> (numValid) != cols.size (), std::runtime_error,
211 "replaceLocalValues returned " << numValid << " != cols.size() = " <<
212 cols.size () << ".");
213 }
214
215 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
216 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::setAllToScalar(const Scalar &alpha) { XPETRA_MONITOR("TpetraCrsMatrix::setAllToScalar"); mtx_->setAllToScalar(alpha); }
217
218 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
219 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::scale(const Scalar &alpha) { XPETRA_MONITOR("TpetraCrsMatrix::scale"); mtx_->scale(alpha); }
220
221 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
223 { XPETRA_MONITOR("TpetraCrsMatrix::allocateAllValues"); rowptr.resize(getLocalNumRows()+1); colind.resize(numNonZeros); values.resize(numNonZeros);}
224
225 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
227 { XPETRA_MONITOR("TpetraCrsMatrix::setAllValues"); mtx_->setAllValues(rowptr,colind,values); }
228
229 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
231 { XPETRA_MONITOR("TpetraCrsMatrix::getAllValues");
232 // TODO: Change getAllValues interface to return Kokkos views,
233 // TODO: or add getLocalRowPtrsHost, getLocalIndicesHost, etc., interfaces
234 // TODO: and use them in MueLu
235 rowptr = Kokkos::Compat::persistingView(mtx_->getLocalRowPtrsHost());
236 colind = Kokkos::Compat::persistingView(mtx_->getLocalIndicesHost());
237 values = Teuchos::arcp_reinterpret_cast<const Scalar>(
238 Kokkos::Compat::persistingView(mtx_->getLocalValuesHost(
239 Tpetra::Access::ReadOnly)));
240 }
241
242 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
244 { XPETRA_MONITOR("TpetraCrsMatrix::getAllValues");
245 // TODO: Change getAllValues interface to return Kokkos view,
246 // TODO: or add getLocalValuesDevice() interfaces
247 values = Teuchos::arcp_reinterpret_cast<Scalar>(
248 Kokkos::Compat::persistingView(mtx_->getLocalValuesDevice(
249 Tpetra::Access::ReadWrite)));
250 }
251
252 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
254 { return mtx_->haveGlobalConstants();}
255
256 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
257 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::resumeFill(const RCP< ParameterList > &params) { XPETRA_MONITOR("TpetraCrsMatrix::resumeFill"); mtx_->resumeFill(params); }
258
259 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
260 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params) { XPETRA_MONITOR("TpetraCrsMatrix::fillComplete"); mtx_->fillComplete(toTpetra(domainMap), toTpetra(rangeMap), params); }
261
262 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
263 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::fillComplete(const RCP< ParameterList > &params) { XPETRA_MONITOR("TpetraCrsMatrix::fillComplete"); mtx_->fillComplete(params); }
264
265
266 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
268 XPETRA_MONITOR("TpetraCrsMatrix::replaceDomainMapAndImporter");
269 XPETRA_DYNAMIC_CAST( const TpetraImportClass , *newImporter, tImporter, "Xpetra::TpetraCrsMatrix::replaceDomainMapAndImporter only accepts Xpetra::TpetraImport.");
270 RCP<const Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> > myImport = tImporter.getTpetra_Import();
271 mtx_->replaceDomainMapAndImporter( toTpetra(newDomainMap),myImport);
272 }
273
274 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
276 const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & rangeMap,
277 const RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > &importer,
278 const RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > &exporter,
279 const RCP<ParameterList> &params) {
280 XPETRA_MONITOR("TpetraCrsMatrix::expertStaticFillComplete");
283
284 if(importer!=Teuchos::null) {
285 XPETRA_DYNAMIC_CAST( const TpetraImportClass , *importer, tImporter, "Xpetra::TpetraCrsMatrix::expertStaticFillComplete only accepts Xpetra::TpetraImport.");
286 myImport = tImporter.getTpetra_Import();
287 }
288 if(exporter!=Teuchos::null) {
289 XPETRA_DYNAMIC_CAST( const TpetraExportClass , *exporter, tExporter, "Xpetra::TpetraCrsMatrix::expertStaticFillComplete only accepts Xpetra::TpetraExport.");
290 myExport = tExporter.getTpetra_Export();
291 }
292
293 mtx_->expertStaticFillComplete(toTpetra(domainMap),toTpetra(rangeMap),myImport,myExport,params);
294 }
295
296 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
298
299 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
301
302 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
304
305 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
306 global_size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getGlobalNumRows() const { XPETRA_MONITOR("TpetraCrsMatrix::getGlobalNumRows"); return mtx_->getGlobalNumRows(); }
307
308 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
309 global_size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getGlobalNumCols() const { XPETRA_MONITOR("TpetraCrsMatrix::getGlobalNumCols"); return mtx_->getGlobalNumCols(); }
310
311 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
312 size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getLocalNumRows() const { XPETRA_MONITOR("TpetraCrsMatrix::getLocalNumRows"); return mtx_->getLocalNumRows(); }
313
314 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
315 size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getLocalNumCols() const { XPETRA_MONITOR("TpetraCrsMatrix::getLocalNumCols"); return mtx_->getLocalNumCols(); }
316
317 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
318 global_size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getGlobalNumEntries() const { XPETRA_MONITOR("TpetraCrsMatrix::getGlobalNumEntries"); return mtx_->getGlobalNumEntries(); }
319
320 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
321 size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getLocalNumEntries() const { XPETRA_MONITOR("TpetraCrsMatrix::getLocalNumEntries"); return mtx_->getLocalNumEntries(); }
322
323 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
324 size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getNumEntriesInLocalRow(LocalOrdinal localRow) const { XPETRA_MONITOR("TpetraCrsMatrix::getNumEntriesInLocalRow"); return mtx_->getNumEntriesInLocalRow(localRow); }
325
326 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
327 size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const { XPETRA_MONITOR("TpetraCrsMatrix::getNumEntriesInGlobalRow"); return mtx_->getNumEntriesInGlobalRow(globalRow); }
328
329 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
330 size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getGlobalMaxNumRowEntries() const { XPETRA_MONITOR("TpetraCrsMatrix::getGlobalMaxNumRowEntries"); return mtx_->getGlobalMaxNumRowEntries(); }
331
332 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
333 size_t TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getLocalMaxNumRowEntries() const { XPETRA_MONITOR("TpetraCrsMatrix::getLocalMaxNumRowEntries"); return mtx_->getLocalMaxNumRowEntries(); }
334
335 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
336 bool TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::isLocallyIndexed() const { XPETRA_MONITOR("TpetraCrsMatrix::isLocallyIndexed"); return mtx_->isLocallyIndexed(); }
337
338 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
339 bool TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::isGloballyIndexed() const { XPETRA_MONITOR("TpetraCrsMatrix::isGloballyIndexed"); return mtx_->isGloballyIndexed(); }
340
341 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
342 bool TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::isFillComplete() const { XPETRA_MONITOR("TpetraCrsMatrix::isFillComplete"); return mtx_->isFillComplete(); }
343
344 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
345 bool TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::isFillActive() const { XPETRA_MONITOR("TpetraCrsMatrix::isFillActive"); return mtx_->isFillActive(); }
346
347 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
348 typename ScalarTraits< Scalar >::magnitudeType TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getFrobeniusNorm() const { XPETRA_MONITOR("TpetraCrsMatrix::getFrobeniusNorm"); return mtx_->getFrobeniusNorm(); }
349
350 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
351 bool TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::supportsRowViews() const { XPETRA_MONITOR("TpetraCrsMatrix::supportsRowViews"); return mtx_->supportsRowViews(); }
352
353 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
354 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const {
355 XPETRA_MONITOR("TpetraCrsMatrix::getLocalRowCopy");
356 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::nonconst_local_inds_host_view_type indices("indices",Indices.size());
357 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonconst_values_host_view_type values("values", Values.size());
358
359 mtx_->getLocalRowCopy(LocalRow, indices, values, NumEntries);
360 for (size_t i=0;i<NumEntries; ++i) {
361 Indices[i]=indices(i);
362 Values[i]=values(i);
363 }
364 }
365
366 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
367 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const {
368 XPETRA_MONITOR("TpetraCrsMatrix::getGlobalRowCopy");
369 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::nonconst_global_inds_host_view_type indices("indices",Indices.size());
370 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonconst_values_host_view_type values("values", Values.size());
371
372 mtx_->getGlobalRowCopy(GlobalRow, indices, values, NumEntries);
373 for (size_t i=0;i<NumEntries; ++i) {
374 Indices[i]=indices(i);
375 Values[i]=values(i);
376 }
377
378 }
379
380 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
382 XPETRA_MONITOR("TpetraCrsMatrix::getGlobalRowView");
383 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::global_inds_host_view_type indices;
384 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::values_host_view_type values;
385
386 mtx_->getGlobalRowView(GlobalRow, indices, values);
387 Indices = ArrayView<const GlobalOrdinal> (indices.data(), indices.extent(0));
388 Values = ArrayView<const Scalar> (reinterpret_cast<const Scalar*>(values.data()), values.extent(0));
389 }
390
391 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
393 XPETRA_MONITOR("TpetraCrsMatrix::getLocalRowView");
394 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::local_inds_host_view_type indices;
395 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::values_host_view_type values;
396
397 mtx_->getLocalRowView(LocalRow, indices, values);
398 Indices = ArrayView<const LocalOrdinal> (indices.data(), indices.extent(0));
399 Values = ArrayView<const Scalar> (reinterpret_cast<const Scalar*>(values.data()), values.extent(0));
400 }
401
402 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
403 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta) const { XPETRA_MONITOR("TpetraCrsMatrix::apply"); mtx_->apply(toTpetra(X), toTpetra(Y), mode, alpha, beta); }
404
405 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
406 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode, Scalar alpha, Scalar beta, bool sumInterfaceValues, const RCP<Import<LocalOrdinal, GlobalOrdinal, Node> >& regionInterfaceImporter, const Teuchos::ArrayRCP<LocalOrdinal>& regionInterfaceLIDs ) const
407 {
408 XPETRA_MONITOR("TpetraCrsMatrix::apply(region)");
409 RCP<const Map< LocalOrdinal, GlobalOrdinal, Node >> regionInterfaceMap = regionInterfaceImporter->getTargetMap();
410 mtx_->localApply(toTpetra(X), toTpetra(Y), mode, alpha, beta);
411 if (sumInterfaceValues)
412 {
413 // preform communication to propagate local interface
414 // values to all the processor that share interfaces.
416 matvecInterfaceTmp->doImport(Y, *regionInterfaceImporter, INSERT);
417
418 // sum all contributions to interface values
419 // on all ranks
421 ArrayRCP<Scalar> interfaceData = matvecInterfaceTmp->getDataNonConst(0);
422 for(LocalOrdinal interfaceIdx = 0; interfaceIdx < static_cast<LocalOrdinal>(interfaceData.size()); ++interfaceIdx) {
423 YData[regionInterfaceLIDs[interfaceIdx]] += interfaceData[interfaceIdx];
424 }
425 }
426 }
427
428 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
430
431 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
433
434 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
435 std::string TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::description() const { XPETRA_MONITOR("TpetraCrsMatrix::description"); return mtx_->description(); }
436
437 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
438 void TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const { XPETRA_MONITOR("TpetraCrsMatrix::describe"); mtx_->describe(out, verbLevel); }
439
440 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
442 XPETRA_MONITOR("TpetraCrsMatrix::setObjectLabel");
444 mtx_->setObjectLabel(objectLabel);
445 }
446
447
448
449 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
451 : mtx_(Teuchos::rcp(new Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>(*(matrix.mtx_),Teuchos::Copy))) {}
452
453 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
455 XPETRA_MONITOR("TpetraCrsMatrix::getLocalDiagCopy");
456 XPETRA_DYNAMIC_CAST(TpetraVectorClass, diag, tDiag, "Xpetra::TpetraCrsMatrix.getLocalDiagCopy() only accept Xpetra::TpetraVector as input arguments.");
457 mtx_->getLocalDiagCopy(*tDiag.getTpetra_Vector());
458 }
459
460 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
462 XPETRA_MONITOR("TpetraCrsMatrix::getLocalDiagOffsets");
463 mtx_->getLocalDiagOffsets(offsets);
464 }
465
466 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
468 XPETRA_MONITOR("TpetraCrsMatrix::getLocalDiagCopy");
469 mtx_->getLocalDiagCopy(*(toTpetra(diag)), offsets);
470 }
471
472 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
474 XPETRA_MONITOR("TpetraCrsMatrix::replaceDiag");
475 Tpetra::replaceDiagonalCrsMatrix(*mtx_, *(toTpetra(diag)));
476 }
477
478 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
480 XPETRA_MONITOR("TpetraCrsMatrix::leftScale");
481 mtx_->leftScale(*(toTpetra(x)));
482 }
483
484 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
486 XPETRA_MONITOR("TpetraCrsMatrix::rightScale");
487 mtx_->rightScale(*(toTpetra(x)));
488 }
489
490 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
492
493 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
496 XPETRA_MONITOR("TpetraCrsMatrix::doImport");
497
498 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, source, tSource, "Xpetra::TpetraCrsMatrix::doImport only accept Xpetra::TpetraCrsMatrix as input arguments.");//TODO: remove and use toTpetra()
500 //mtx_->doImport(toTpetraCrsMatrix(source), *tImporter.getTpetra_Import(), toTpetra(CM));
501 mtx_->doImport(*v, toTpetra(importer), toTpetra(CM));
502 }
503
505 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
508 XPETRA_MONITOR("TpetraCrsMatrix::doExport");
509
510 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, dest, tDest, "Xpetra::TpetraCrsMatrix::doImport only accept Xpetra::TpetraCrsMatrix as input arguments.");//TODO: remove and use toTpetra()
512 mtx_->doExport(*v, toTpetra(importer), toTpetra(CM));
513
514 }
515
516 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
519 XPETRA_MONITOR("TpetraCrsMatrix::doImport");
520
521 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, source, tSource, "Xpetra::TpetraCrsMatrix::doImport only accept Xpetra::TpetraCrsMatrix as input arguments.");//TODO: remove and use toTpetra()
523 mtx_->doImport(*v, toTpetra(exporter), toTpetra(CM));
524
525 }
526
527 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
530 XPETRA_MONITOR("TpetraCrsMatrix::doExport");
531
532 XPETRA_DYNAMIC_CAST(const TpetraCrsMatrixClass, dest, tDest, "Xpetra::TpetraCrsMatrix::doImport only accept Xpetra::TpetraCrsMatrix as input arguments.");//TODO: remove and use toTpetra()
534 mtx_->doExport(*v, toTpetra(exporter), toTpetra(CM));
535
536 }
537
538 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
540 XPETRA_MONITOR("TpetraCrsMatrix::removeEmptyProcessesInPlace");
541 mtx_->removeEmptyProcessesInPlace(toTpetra(newMap));
542 }
543
544 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
546 return ! mtx_.is_null ();
547 }
548
549 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
550 TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::TpetraCrsMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) : mtx_(mtx) { }
551
552 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
554
556 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
558
560 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
562 const MultiVector & B,
563 MultiVector & R) const {
564 Tpetra::Details::residual(*mtx_,toTpetra(X),toTpetra(B),toTpetra(R));
565 }
566
567
570// End of TpetrCrsMatrix class definition //
573
574} // Xpetra namespace
575
576#endif // XPETRA_TPETRACRSMATRIX_DEF_HPP
#define XPETRA_MONITOR(funcName)
#define XPETRA_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
size_type size() const
void resize(const size_type n, const T &val=T())
size_type size() const
virtual void setObjectLabel(const std::string &objectLabel)
bool is_null() const
T * get() const
static Teuchos::RCP< MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Build(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &map, size_t NumVectors, bool zeroOut=true)
Constructor specifying the number of non-zeros for all rows.
virtual Teuchos::ArrayRCP< Scalar > getDataNonConst(size_t j)=0
View of the local values in a particular vector of this multivector.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
Get the underlying Tpetra matrix.
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
TpetraCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying fixed number of entries for each row.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage.
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_matrix_type local_matrix_type
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
void getLocalDiagCopy(Vector &diag) const
Get a copy of the diagonal entries owned by this node, with local row indices.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
void rightScale(const Vector &x)
Right scale operator with given vector values.
void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
void residual(const MultiVector &X, const MultiVector &B, MultiVector &R) const
Compute a residual R = B - (*this) * X.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
void resumeFill(const RCP< ParameterList > &params=null)
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
void leftScale(const Vector &x)
Left scale operator with given vector values.
void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
bool hasMatrix() const
Does this have an underlying matrix.
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
std::string description() const
A simple one-line description of this object.
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
Get the underlying Tpetra matrix.
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
bool isFillActive() const
Returns true if the matrix is in edit mode.
void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)
Replaces the current domainMap and importer with the user-specified objects.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph.
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
void expertStaticFillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &importer=Teuchos::null, const RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > &exporter=Teuchos::null, const RCP< ParameterList > &params=Teuchos::null)
Expert static fill complete.
void replaceDiag(const Vector &diag)
Replace the diagonal entries of the matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
void setObjectLabel(const std::string &objectLabel)
void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Map associated with the range of this operator, which must be compatible with Y....
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Copy
Xpetra namespace
size_t global_size_t
Global size_t object.
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)
CombineMode
Xpetra::Combine Mode enumerable type.