Xpetra_CrsMatrix.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_CRSMATRIX_HPP
47 #define XPETRA_CRSMATRIX_HPP
48 
49 /* this file is automatically generated - do not edit (see script/interfaces.py) */
50 
51 #include <Kokkos_DefaultNode.hpp>
52 #include "Xpetra_ConfigDefs.hpp"
53 #include "Xpetra_RowMatrix.hpp"
54 #include "Xpetra_DistObject.hpp"
55 #include "Xpetra_CrsGraph.hpp"
56 #include "Xpetra_Vector.hpp"
57 
58 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
59 #include <Kokkos_StaticCrsGraph.hpp>
60 #include <Kokkos_CrsMatrix.hpp>
61 #endif
62 
63 namespace Xpetra {
64 
65  template <class Scalar = double,
66  class LocalOrdinal = CrsGraph<>::local_ordinal_type,
67  class GlobalOrdinal =
69  class Node =
71  class CrsMatrix
72  : public RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>,
73  public DistObject<char, LocalOrdinal,GlobalOrdinal,Node>
74  {
75  public:
76  typedef Scalar scalar_type;
77  typedef LocalOrdinal local_ordinal_type;
78  typedef GlobalOrdinal global_ordinal_type;
79  typedef Node node_type;
80 
82 
83 
85  virtual ~CrsMatrix() { }
86 
88 
90 
92  virtual void
93  insertGlobalValues (GlobalOrdinal globalRow,
95  const ArrayView<const Scalar>& vals) = 0;
96 
98  virtual void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)= 0;
99 
101  virtual void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)= 0;
102 
104  virtual void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)= 0;
105 
107  virtual void setAllToScalar(const Scalar &alpha)= 0;
108 
110  virtual void scale(const Scalar &alpha)= 0;
111 
113  //** \warning This is an expert-only routine and should not be called from user code. */
114  virtual void allocateAllValues(size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values)=0;
115 
117  virtual void setAllValues(const ArrayRCP<size_t> & rowptr, const ArrayRCP<LocalOrdinal> & colind, const ArrayRCP<Scalar> & values)=0;
118 
120  virtual void getAllValues(ArrayRCP<const size_t>& rowptr, ArrayRCP<const LocalOrdinal>& colind, ArrayRCP<const Scalar>& values) const = 0;
121 
123 
125 
126 
128  virtual void resumeFill(const RCP< ParameterList > &params=null)= 0;
129 
131  virtual void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)= 0;
132 
134  virtual void fillComplete(const RCP< ParameterList > &params=null)= 0;
135 
138 
140  virtual void expertStaticFillComplete(const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & domainMap,
141  const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & rangeMap,
142  const RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > &importer=Teuchos::null,
143  const RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > &exporter=Teuchos::null,
144  const RCP<ParameterList> &params=Teuchos::null) = 0;
146 
148 
149 
152 
155 
158 
160  virtual global_size_t getGlobalNumRows() const = 0;
161 
163  virtual global_size_t getGlobalNumCols() const = 0;
164 
166  virtual size_t getNodeNumRows() const = 0;
167 
169  virtual global_size_t getGlobalNumEntries() const = 0;
170 
172  virtual size_t getNodeNumEntries() const = 0;
173 
175  virtual size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const = 0;
176 
178  virtual global_size_t getGlobalNumDiags() const = 0;
179 
181  virtual size_t getNodeNumDiags() const = 0;
182 
184  virtual size_t getGlobalMaxNumRowEntries() const = 0;
185 
187  virtual size_t getNodeMaxNumRowEntries() const = 0;
188 
190  virtual bool isLocallyIndexed() const = 0;
191 
193  virtual bool isGloballyIndexed() const = 0;
194 
196  virtual bool isFillComplete() const = 0;
197 
199  virtual bool isFillActive() const = 0;
200 
202  virtual typename ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const = 0;
203 
205  virtual bool supportsRowViews() const = 0;
206 
208  virtual void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const = 0;
209 
211  virtual void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView<GlobalOrdinal> &indices, const ArrayView<Scalar> &values, size_t &numEntries) const = 0;
212 
214  virtual void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const = 0;
215 
218 
220  virtual void getLocalDiagOffsets(Teuchos::ArrayRCP<size_t> &offsets) const = 0;
221 
224 
225  virtual void removeEmptyProcessesInPlace(const RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >& newMap) = 0;
226 
228 
230 
231 
234 
237 
240 
242 
244 
245 
247  virtual std::string description() const = 0;
248 
251 
253 
255 
256 #ifdef HAVE_XPETRA_KOKKOS_REFACTOR
257  typedef typename Kokkos::Details::ArithTraits<Scalar>::val_type impl_scalar_type;
258  typedef typename node_type::execution_space execution_space;
259 
260  // that is the local_graph_type in Tpetra::CrsGraph...
261  typedef Kokkos::StaticCrsGraph<LocalOrdinal,
262  Kokkos::LayoutLeft,
263  execution_space> local_graph_type;
267  typedef Kokkos::CrsMatrix<impl_scalar_type, LocalOrdinal, execution_space,void,
268  typename local_graph_type::size_type> local_matrix_type;
269 
271  virtual local_matrix_type getLocalMatrix () const = 0;
272 #endif
273 
275 
276  // Adding these functions by hand, as they're in the skip list.
277 
279  virtual size_t getNodeNumCols() const = 0;
280 
282  virtual void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const = 0;
283 
285  virtual bool hasMatrix() const = 0;
286 
287 
288  }; // CrsMatrix class
289 
290 } // Xpetra namespace
291 
292 #define XPETRA_CRSMATRIX_SHORT
293 #endif // XPETRA_CRSMATRIX_HPP
virtual 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)=0
Expert static fill complete.
virtual global_size_t getGlobalNumRows() const =0
Number of global elements in the row map of this matrix.
virtual void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const =0
Computes the sparse matrix-multivector multiplication.
virtual size_t getNodeMaxNumRowEntries() const =0
Returns the maximum number of entries across all rows/columns on this node.
virtual void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Replace matrix entries, using local IDs.
virtual void removeEmptyProcessesInPlace(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)=0
virtual const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const =0
Returns the Map that describes the column distribution in this matrix.
GlobalOrdinal global_ordinal_type
virtual void scale(const Scalar &alpha)=0
Scale the current values of a matrix, this = alpha*this.
virtual bool hasMatrix() const =0
Does this have an underlying matrix.
virtual void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)=0
Signal that data entry is complete, specifying domain and range maps.
virtual void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Insert matrix entries, using global IDs.
virtual const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const =0
Returns the Map associated with the range of this operator, which must be compatible with Y...
virtual const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const =0
Returns the Map that describes the row distribution in this matrix.
Xpetra namespace
virtual size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const =0
Returns the current number of entries on this node in the specified local row.
virtual bool isLocallyIndexed() const =0
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
virtual void setAllToScalar(const Scalar &alpha)=0
Set all matrix entries equal to scalarThis.
LocalOrdinal local_ordinal_type
virtual void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)=0
Sets the 1D pointer arrays of the graph.
virtual const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const =0
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
virtual size_t getNodeNumEntries() const =0
Returns the local number of entries in this matrix.
virtual global_size_t getGlobalNumDiags() const =0
Returns the number of global diagonal entries, based on global row/column index comparisons.
virtual void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const =0
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage...
virtual std::string description() const =0
A simple one-line description of this object.
LocalOrdinal local_ordinal_type
virtual void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)=0
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
virtual void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const =0
Extract a const, non-persisting view of local indices in a specified row of the matrix.
virtual global_size_t getGlobalNumCols() const =0
Number of global columns in the matrix.
GlobalOrdinal global_ordinal_type
virtual bool supportsRowViews() const =0
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
virtual void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Replace matrix entries, using global IDs.
virtual void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const =0
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
virtual void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const =0
Print the object with some verbosity level to an FancyOStream object.
virtual void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)=0
Replaces the current domainMap and importer with the user-specified objects.
virtual void resumeFill(const RCP< ParameterList > &params=null)=0
virtual size_t getNodeNumDiags() const =0
Returns the number of local diagonal entries, based on global row/column index comparisons.
virtual void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)=0
Insert matrix entries, using local IDs.
size_t global_size_t
Global size_t object.
static const EVerbosityLevel verbLevel_default
virtual global_size_t getGlobalNumEntries() const =0
Returns the global number of entries in this matrix.
virtual size_t getGlobalMaxNumRowEntries() const =0
Returns the maximum number of entries across all rows/columns on all nodes.
virtual void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const =0
Get a copy of the diagonal entries owned by this node, with local row indices.
virtual ~CrsMatrix()
Destructor.
virtual ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const =0
Returns the Frobenius norm of the matrix.
virtual bool isFillComplete() const =0
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
virtual void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const =0
Get offsets of the diagonal entries in the matrix.
virtual size_t getNodeNumCols() const =0
Returns the number of matrix columns owned on the calling node.
virtual RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const =0
Returns the CrsGraph associated with this matrix.
virtual bool isGloballyIndexed() const =0
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
virtual size_t getNodeNumRows() const =0
Returns the number of matrix rows owned on the calling node.
virtual bool isFillActive() const =0
Returns true if the matrix is in edit mode.
virtual void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const =0
Gets the 1D pointer arrays of the graph.
virtual void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const =0
Extract a const, non-persisting view of global indices in a specified row of the matrix.