Xpetra  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Xpetra_CrsMatrixWrap_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 
47 // WARNING: This code is experimental. Backwards compatibility should not be expected.
48 
49 #ifndef XPETRA_CRSMATRIXWRAP_DEF_HPP
50 #define XPETRA_CRSMATRIXWRAP_DEF_HPP
51 
52 #include <Tpetra_KokkosCompat_DefaultNode.hpp>
53 
54 #include <Teuchos_SerialDenseMatrix.hpp>
55 #include <Teuchos_Hashtable.hpp>
56 
57 #include "Xpetra_ConfigDefs.hpp"
58 #include "Xpetra_Exceptions.hpp"
59 
60 #include "Xpetra_MultiVector.hpp"
61 #include "Xpetra_CrsGraph.hpp"
62 #include "Xpetra_CrsMatrix.hpp"
64 
65 #include "Xpetra_Matrix.hpp"
66 
68 
69 namespace Xpetra {
70 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
72  : finalDefaultView_(false) {
75 }
76 
77 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
79  size_t maxNumEntriesPerRow)
80  : finalDefaultView_(false) {
81  matrixData_ = CrsMatrixFactory::Build(rowMap, maxNumEntriesPerRow);
83 }
84 
85 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
87  const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc)
88  : finalDefaultView_(false) {
89  matrixData_ = CrsMatrixFactory::Build(rowMap, NumEntriesPerRowToAlloc);
91 }
92 
93 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
94 CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::CrsMatrixWrap(const RCP<const Map> &rowMap, const RCP<const Map> &colMap, size_t maxNumEntriesPerRow)
95  : finalDefaultView_(false) {
96  // Set matrix data
97  matrixData_ = CrsMatrixFactory::Build(rowMap, colMap, maxNumEntriesPerRow);
98 
99  // Default view
101 }
102 
103 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
104 CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::CrsMatrixWrap(const RCP<const Map> &rowMap, const RCP<const Map> &colMap, const ArrayRCP<const size_t> &NumEntriesPerRowToAlloc)
105  : finalDefaultView_(false) {
106  // Set matrix data
107  matrixData_ = CrsMatrixFactory::Build(rowMap, colMap, NumEntriesPerRowToAlloc);
108 
109  // Default view
111 }
112 
113 #ifdef HAVE_XPETRA_TPETRA
114 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
115 CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::CrsMatrixWrap(const RCP<const Map> &rowMap, const RCP<const Map> &colMap, const local_matrix_type &lclMatrix, const Teuchos::RCP<Teuchos::ParameterList> &params)
116  : finalDefaultView_(false) {
117  // Set matrix data
118  matrixData_ = CrsMatrixFactory::Build(rowMap, colMap, lclMatrix, params);
119 
120  // Default view
122 }
123 
124 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
125 CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::CrsMatrixWrap(const local_matrix_type &lclMatrix, const RCP<const Map> &rowMap, const RCP<const Map> &colMap,
126  const RCP<const Map> &domainMap, const RCP<const Map> &rangeMap,
127  const Teuchos::RCP<Teuchos::ParameterList> &params)
128  : finalDefaultView_(false) {
129  // Set matrix data
130  matrixData_ = CrsMatrixFactory::Build(lclMatrix, rowMap, colMap, domainMap, rangeMap, params);
131 
132  // Default view
134 }
135 #else
136 #ifdef __GNUC__
137 #warning "Xpetra Kokkos interface for CrsMatrix is enabled (HAVE_XPETRA_KOKKOS_REFACTOR) but Tpetra is disabled. The Kokkos interface needs Tpetra to be enabled, too."
138 #endif
139 #endif
140 
141 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
143  : finalDefaultView_(matrix->isFillComplete()) {
144  // Set matrix data
145  matrixData_ = matrix;
146 
147  // Default view
149 }
150 
151 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
152 CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::CrsMatrixWrap(const RCP<const CrsGraph> &graph, const RCP<ParameterList> &paramList)
153  : finalDefaultView_(false) {
154  // Set matrix data
155  matrixData_ = CrsMatrixFactory::Build(graph, paramList);
156 
157  // Default view
159 }
160 
161 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
164 
165  const RCP<ParameterList> &paramList)
166  : finalDefaultView_(false) {
167  // Set matrix data
168  matrixData_ = CrsMatrixFactory::Build(graph, values, paramList);
169 
170  // Default view
172 }
173 
174 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
176 
177 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
178 void CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::insertGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals) {
179  matrixData_->insertGlobalValues(globalRow, cols, vals);
180 }
181 
182 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
183 void CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::insertLocalValues(LocalOrdinal localRow, const ArrayView<const LocalOrdinal> &cols, const ArrayView<const Scalar> &vals) {
184  matrixData_->insertLocalValues(localRow, cols, vals);
185 }
186 
187 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
189  const ArrayView<const GlobalOrdinal> &cols,
190  const ArrayView<const Scalar> &vals) { matrixData_->replaceGlobalValues(globalRow, cols, vals); }
191 
192 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
194  const ArrayView<const LocalOrdinal> &cols,
195  const ArrayView<const Scalar> &vals) { matrixData_->replaceLocalValues(localRow, cols, vals); }
196 
197 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
198 void CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::setAllToScalar(const Scalar &alpha) { matrixData_->setAllToScalar(alpha); }
199 
200 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
202  matrixData_->scale(alpha);
203 }
204 
205 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
207  matrixData_->resumeFill(params);
208 }
209 
210 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
211 void CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::fillComplete(const RCP<const Map> &domainMap, const RCP<const Map> &rangeMap, const RCP<Teuchos::ParameterList> &params) {
212  matrixData_->fillComplete(domainMap, rangeMap, params);
213 
214  // Update default view with the colMap because colMap can be <tt>null</tt> until fillComplete() is called.
215  updateDefaultView();
216 }
217 
218 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
220  matrixData_->fillComplete(params);
221 
222  // Update default view with the colMap because colMap can be <tt>null</tt> until fillComplete() is called.
223  updateDefaultView();
224 }
225 
226 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
228  return matrixData_->getGlobalNumRows();
229 }
230 
231 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
233  return matrixData_->getGlobalNumCols();
234 }
235 
236 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
238  return matrixData_->getLocalNumRows();
239 }
240 
241 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
243  return matrixData_->getGlobalNumEntries();
244 }
245 
246 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
248  return matrixData_->getLocalNumEntries();
249 }
250 
251 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
253  return matrixData_->getNumEntriesInLocalRow(localRow);
254 }
255 
256 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
258  return matrixData_->getNumEntriesInGlobalRow(globalRow);
259 }
260 
261 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
263  return matrixData_->getGlobalMaxNumRowEntries();
264 }
265 
266 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
268  return matrixData_->getLocalMaxNumRowEntries();
269 }
270 
271 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
273  return matrixData_->isLocallyIndexed();
274 }
275 
276 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
278  return matrixData_->isGloballyIndexed();
279 }
280 
281 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
283  return matrixData_->isFillComplete();
284 }
285 
286 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
288  const ArrayView<LocalOrdinal> &Indices,
289  const ArrayView<Scalar> &Values,
290  size_t &NumEntries) const {
291  matrixData_->getLocalRowCopy(LocalRow, Indices, Values, NumEntries);
292 }
293 
294 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
295 void CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values) const {
296  matrixData_->getGlobalRowView(GlobalRow, indices, values);
297 }
298 
299 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
300 void CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values) const {
301  matrixData_->getLocalRowView(LocalRow, indices, values);
302 }
303 
304 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
306  matrixData_->getLocalDiagCopy(diag);
307 }
308 
309 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
311  matrixData_->getLocalDiagOffsets(offsets);
312 }
313 
314 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
316  matrixData_->getLocalDiagCopy(diag, offsets);
317 }
318 
319 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
320 typename ScalarTraits<Scalar>::magnitudeType CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getFrobeniusNorm() const {
321  return matrixData_->getFrobeniusNorm();
322 }
323 
324 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
326  matrixData_->leftScale(x);
327 }
328 
329 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
331  matrixData_->rightScale(x);
332 }
333 
334 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
336  return matrixData_->haveGlobalConstants();
337 }
338 
339 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
342  Teuchos::ETransp mode,
343  Scalar alpha,
344  Scalar beta) const {
345  matrixData_->apply(X, Y, mode, alpha, beta);
346 }
347 
348 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
351  Teuchos::ETransp mode,
352  Scalar alpha,
353  Scalar beta,
354  bool sumInterfaceValues,
355  const RCP<Import<LocalOrdinal, GlobalOrdinal, Node>> &regionInterfaceImporter,
356  const Teuchos::ArrayRCP<LocalOrdinal> &regionInterfaceLIDs) const {
357  matrixData_->apply(X, Y, mode, alpha, beta, sumInterfaceValues, regionInterfaceImporter, regionInterfaceLIDs);
358 }
359 
360 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
361 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getDomainMap() const {
362  return matrixData_->getDomainMap();
363 }
364 
365 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
366 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getRangeMap() const {
367  return matrixData_->getRangeMap();
368 }
369 
370 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
371 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> &CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getColMap() const { return getColMap(Matrix::GetCurrentViewLabel()); }
372 
373 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
374 const RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> &CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getColMap(viewLabel_t viewLabel) const {
375  TEUCHOS_TEST_FOR_EXCEPTION(Matrix::operatorViewTable_.containsKey(viewLabel) == false, Xpetra::Exceptions::RuntimeError, "Xpetra::Matrix.GetColMap(): view '" + viewLabel + "' does not exist.");
376  updateDefaultView(); // If CrsMatrix::fillComplete() have been used instead of CrsMatrixWrap::fillComplete(), the default view is updated.
377  return Matrix::operatorViewTable_.get(viewLabel)->GetColMap();
378 }
379 
380 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
382  matrixData_->removeEmptyProcessesInPlace(newMap);
383  this->operatorViewTable_.get(this->GetCurrentViewLabel())->SetRowMap(matrixData_->getRowMap());
384  this->operatorViewTable_.get(this->GetCurrentViewLabel())->SetColMap(matrixData_->getColMap());
385 }
386 
387 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
388 const Teuchos::RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>> CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getMap() const {
389  return matrixData_->getMap();
390 }
391 
392 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
395  const CrsMatrixWrap &sourceWrp = dynamic_cast<const CrsMatrixWrap &>(source);
396  matrixData_->doImport(*sourceWrp.getCrsMatrix(), importer, CM);
397 }
398 
399 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
402  const CrsMatrixWrap &destWrp = dynamic_cast<const CrsMatrixWrap &>(dest);
403  matrixData_->doExport(*destWrp.getCrsMatrix(), importer, CM);
404 }
405 
406 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
409  const CrsMatrixWrap &sourceWrp = dynamic_cast<const CrsMatrixWrap &>(source);
410  matrixData_->doImport(*sourceWrp.getCrsMatrix(), exporter, CM);
411 }
412 
413 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
416  const CrsMatrixWrap &destWrp = dynamic_cast<const CrsMatrixWrap &>(dest);
417  matrixData_->doExport(*destWrp.getCrsMatrix(), exporter, CM);
418 }
419 
420 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
422  return "Xpetra::CrsMatrixWrap";
423 }
424 
425 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
426 void CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const {
427  // Teuchos::EVerbosityLevel vl = verbLevel;
428  // if (vl == VERB_DEFAULT) vl = VERB_LOW;
429  // RCP<const Comm<int> > comm = this->getComm();
430  // const int myImageID = comm->getRank(),
431  // numImages = comm->getSize();
432 
433  // if (myImageID == 0) out << this->description() << std::endl;
434 
435  matrixData_->describe(out, verbLevel);
436 
437  // Teuchos::OSTab tab(out);
438 }
439 
440 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
442  Teuchos::LabeledObject::setObjectLabel(objectLabel);
443  matrixData_->setObjectLabel(objectLabel);
444 }
445 
446 #ifdef HAVE_XPETRA_TPETRA
447 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
450  return matrixData_->getLocalMatrixHost();
451 }
452 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
455  return matrixData_->getLocalMatrixDevice();
456 }
457 #else
458 #ifdef __GNUC__
459 #warning "Xpetra Kokkos interface for CrsMatrix is enabled (HAVE_XPETRA_KOKKOS_REFACTOR) but Tpetra is disabled. The Kokkos interface needs Tpetra to be enabled, too."
460 #endif
461 #endif
462 
463 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
465 
466 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
467 RCP<const Xpetra::CrsGraph<LocalOrdinal, GlobalOrdinal, Node>> CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getCrsGraph() const { return matrixData_->getCrsGraph(); }
468 
469 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
470 RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>> CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>::getCrsMatrix() const { return matrixData_; }
471 
472 // Default view is created after fillComplete()
473 // Because ColMap might not be available before fillComplete().
474 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
476  // Create default view
477  this->defaultViewLabel_ = "point";
478  this->CreateView(this->GetDefaultViewLabel(), matrixData_->getRowMap(), matrixData_->getColMap());
479 
480  // Set current view
481  this->currentViewLabel_ = this->GetDefaultViewLabel();
482 }
483 
484 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
486  if ((finalDefaultView_ == false) && matrixData_->isFillComplete()) {
487  // Update default view with the colMap
488  Matrix::operatorViewTable_.get(Matrix::GetDefaultViewLabel())->SetColMap(matrixData_->getColMap());
489  finalDefaultView_ = true;
490  }
491 }
492 
493 // Expert only
494 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
496  // Clear the old view table
497  Teuchos::Hashtable<viewLabel_t, RCP<MatrixView>> dummy_table;
498  Matrix::operatorViewTable_ = dummy_table;
499 
500  finalDefaultView_ = M->isFillComplete();
501  // Set matrix data
502  matrixData_ = M;
503 
504  // Default view
505  CreateDefaultView();
506 }
507 
508 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
510  return matrixData_->GetStorageBlockSize();
511 }
512 
513 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
518  matrixData_->residual(X, B, R);
519 }
520 
521 } // namespace Xpetra
522 
523 #endif // ifndef XPETRA_CRSMATRIXWRAP_DEF_HPP
void setObjectLabel(const std::string &objectLabel)
RCP< CrsMatrix > getCrsMatrix() const
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale matrix using the given vector entries.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs.
virtual local_matrix_type getLocalMatrixDevice() const
global_size_t getGlobalNumRows() const
Returns the number of global rows in this matrix.
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
const Teuchos::RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
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...
const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
bool hasCrsGraph() const
Supports the getCrsGraph() call.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this.
LocalOrdinal GetStorageBlockSize() const
Returns the block size of the storage mechanism, which is usually 1, except for Tpetra::BlockCrsMatri...
Exception throws to report errors in the internal logical of the program.
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map > &newMap)
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.
void doExport(const Matrix &dest, const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
Teuchos::Hashtable< viewLabel_t, RCP< MatrixView > > operatorViewTable_
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Get Frobenius norm of the matrix.
void doImport(const Matrix &source, const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
std::string description() const
Return a simple one-line description of this object.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
const viewLabel_t & GetCurrentViewLabel() const
CrsMatrix::local_matrix_type local_matrix_type
void getLocalDiagCopy(Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
CrsMatrixWrap(const RCP< const Map > &rowMap)
Constructor for a dynamic profile matrix (Epetra only)
void residual(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &R) const
Compute a residual R = B - (*this) * X.
RCP< const CrsGraph > getCrsGraph() const
Returns the CrsGraph associated with 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 getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
std::string viewLabel_t
size_t global_size_t
Global size_t object.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs.
void fillComplete(const RCP< const Map > &domainMap, const RCP< const Map > &rangeMap, const RCP< Teuchos::ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
virtual ~CrsMatrixWrap()
Destructor.
const RCP< const Map > & getColMap() const
Returns the Map that describes the column distribution in this matrix. This might be null until fillC...
void resumeFill(const RCP< ParameterList > &params=null)
void replaceCrsMatrix(RCP< CrsMatrix > &M)
Expert only.
bool isFillComplete() const
Returns true if fillComplete() has been called and the matrix is in compute mode. ...
virtual void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalar.
const RCP< const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Concrete implementation of Xpetra::Matrix.
global_size_t getGlobalNumCols() const
Returns the number of global columns in the matrix.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs.
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
CombineMode
Xpetra::Combine Mode enumerable type.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale matrix using the given vector entries.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
const viewLabel_t & GetDefaultViewLabel() const
virtual void apply(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
KokkosSparse::CrsMatrix< impl_scalar_type, LocalOrdinal, execution_space, void, typename local_graph_type::size_type > local_matrix_type
The specialization of Kokkos::CrsMatrix that represents the part of the sparse matrix on each MPI pro...
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the specified global row.
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.
virtual local_matrix_type::HostMirror getLocalMatrixHost() const
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
Xpetra-specific matrix class.
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...