Tpetra parallel linear algebra  Version of the Day
TpetraExt_MatrixMatrix_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Tpetra: Templated Linear Algebra Services Package
5 // Copyright (2008) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 // @HEADER
41 
42 #ifndef TPETRA_MATRIXMATRIX_DECL_HPP
43 #define TPETRA_MATRIXMATRIX_DECL_HPP
44 
45 #include <string>
46 #include <Teuchos_RCP.hpp>
47 #include <Teuchos_Array.hpp>
48 #include "Tpetra_ConfigDefs.hpp"
49 #include "Tpetra_CrsMatrix.hpp"
50 #include "Tpetra_Vector.hpp"
51 #include "TpetraExt_MMHelpers.hpp"
52 
53 
59 namespace Tpetra {
60 
61 namespace MatrixMatrix {
62 
93 template <class Scalar,
94  class LocalOrdinal,
95  class GlobalOrdinal,
96  class Node>
97 void Multiply(
98  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
99  bool transposeA,
100  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
101  bool transposeB,
102  CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
103  bool call_FillComplete_on_result=true,
104  const std::string & label = std::string());
105 
119 template <class Scalar,
120  class LocalOrdinal,
121  class GlobalOrdinal,
122  class Node>
123 void Add(
124  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
125  bool transposeA,
126  Scalar scalarA,
127  CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
128  Scalar scalarB );
129 
130 
160 template <class Scalar,
161  class LocalOrdinal,
162  class GlobalOrdinal,
163  class Node>
164 Teuchos::RCP<CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
165 add (const Scalar& alpha,
166  const bool transposeA,
167  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
168  const Scalar& beta,
169  const bool transposeB,
170  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
171  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >& domainMap=Teuchos::null,
172  const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >& rangeMap=Teuchos::null,
173  const Teuchos::RCP<Teuchos::ParameterList>& params=Teuchos::null);
174 
204 template <class Scalar,
205  class LocalOrdinal,
206  class GlobalOrdinal,
207  class Node>
208 void Add(
209  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
210  bool transposeA,
211  Scalar scalarA,
212  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
213  bool transposeB,
214  Scalar scalarB,
215  RCP<CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > C);
216 
217 
239  template <class Scalar,
240  class LocalOrdinal,
241  class GlobalOrdinal,
242  class Node>
243  void Jacobi(Scalar omega,
244  const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Dinv,
245  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
246  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
247  CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
248  bool call_FillComplete_on_result=true,
249  const std::string & label = std::string());
250 
251 } // namespace MatrixMatrix
252 
253 namespace MMdetails{
254 
255 template<class Scalar,
256  class LocalOrdinal,
257  class GlobalOrdinal,
258  class Node>
259 void mult_AT_B_newmatrix(
260  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& A,
261  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& B,
262  CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
263  const std::string & label = std::string());
264 
265 
266 template<class Scalar,
267  class LocalOrdinal,
268  class GlobalOrdinal,
269  class Node>
270 void mult_A_B(
271  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Aview,
272  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Bview,
273  CrsWrapper<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
274  const std::string & label = std::string());
275 
276 template<class Scalar,
277  class LocalOrdinal,
278  class GlobalOrdinal,
279  class Node>
280 void mult_A_B_newmatrix(
281  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Aview,
282  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Bview,
283  CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
284  const std::string & label = std::string());
285 
286 
287 template<class Scalar,
288  class LocalOrdinal,
289  class GlobalOrdinal,
290  class Node>
291 void jacobi_A_B_newmatrix(
292  Scalar omega,
293  const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Dinv,
294  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Aview,
295  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Bview,
296  CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& C,
297  const std::string & label = std::string());
298 
299 
300 template<class Scalar,
301  class LocalOrdinal,
302  class GlobalOrdinal,
303  class Node>
304 void import_and_extract_views(
305  const CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& M,
306  RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> > targetMap,
307  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Mview,
308  RCP<const Import<LocalOrdinal,GlobalOrdinal, Node> > prototypeImporter = Teuchos::null,
309  bool userAssertsThereAreNoRemotes=false,
310  const std::string & label = std::string());
311 
312 template<class Scalar,
313  class LocalOrdinal,
314  class GlobalOrdinal,
315  class Node>
316 void setMaxNumEntriesPerRow(
317  CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Mview);
318 
319 }//end namespace MMdetails
320 
321 } // end of Tpetra namespace
322 
323 #endif // TPETRA_MATRIXMATRIX_DECL_HPP
324 
Namespace Tpetra contains the class and methods constituting the Tpetra library.
void Add(const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool transposeA, Scalar scalarA, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, Scalar scalarB)
void Jacobi(Scalar omega, const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Dinv, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &C, bool call_FillComplete_on_result=true, const std::string &label=std::string())
void Multiply(const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool transposeA, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, bool transposeB, CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &C, bool call_FillComplete_on_result=true, const std::string &label=std::string())
Sparse matrix-matrix multiply.
Teuchos::RCP< CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > add(const Scalar &alpha, const bool transposeA, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta, const bool transposeB, const CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Compute the sparse matrix sum C = scalarA * Op(A) + scalarB * Op(B), where Op(X) is either X or its t...