Ifpack2 Templated Preconditioning Package Version 1.0
Loading...
Searching...
No Matches
Ifpack2_Container_decl.hpp
1/*@HEADER
2// ***********************************************************************
3//
4// Ifpack2: Templated Object-Oriented Algebraic Preconditioner Package
5// Copyright (2009) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
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
43#ifndef IFPACK2_CONTAINER_DECL_HPP
44#define IFPACK2_CONTAINER_DECL_HPP
45
48
49#include "Ifpack2_ConfigDefs.hpp"
50#include "Tpetra_RowMatrix.hpp"
51#include "Teuchos_Describable.hpp"
52#include <Tpetra_Map.hpp>
53#include <Tpetra_BlockCrsMatrix.hpp>
54#include <Teuchos_ParameterList.hpp>
55#include <iostream>
56
57#ifndef DOXYGEN_SHOULD_SKIP_THIS
58namespace Teuchos {
59 // Forward declaration to avoid include.
60 class ParameterList;
61} // namespace Teuchos
62#endif // DOXYGEN_SHOULD_SKIP_THIS
63
64namespace Ifpack2 {
65
111template<class MatrixType>
112class Container : public Teuchos::Describable
113{
114protected:
115 using scalar_type = typename MatrixType::scalar_type;
116 using local_ordinal_type = typename MatrixType::local_ordinal_type;
117 using global_ordinal_type = typename MatrixType::global_ordinal_type;
118 using node_type = typename MatrixType::node_type;
119 using SC = scalar_type;
120 using LO = local_ordinal_type;
121 using GO = global_ordinal_type;
122 using NO = node_type;
123 using import_type = Tpetra::Import<LO, GO, NO>;
124 using mv_type = Tpetra::MultiVector<SC, LO, GO, NO>;
125 using vector_type = Tpetra::Vector<SC, LO, GO, NO>;
126 using map_type = Tpetra::Map<LO, GO, NO>;
127 using crs_matrix_type = Tpetra::CrsMatrix<SC, LO, GO, NO>;
128 using block_crs_matrix_type = Tpetra::BlockCrsMatrix<SC, LO, GO, NO>;
129 using row_matrix_type = Tpetra::RowMatrix<SC, LO, GO, NO>;
130
131 static_assert(std::is_same<MatrixType, row_matrix_type>::value,
132 "Ifpack2::Container: Please use MatrixType = Tpetra::RowMatrix.");
133
135 using ISC = typename Kokkos::Details::ArithTraits<SC>::val_type;
136
139 using HostView = typename mv_type::dual_view_type::t_host;
140 using ConstHostView = typename HostView::const_type;
141
142public:
153 Container (const Teuchos::RCP<const row_matrix_type>& matrix,
154 const Teuchos::Array<Teuchos::Array<LO> >& partitions,
155 bool pointIndexed);
156
158 virtual ~Container();
159
175 Teuchos::ArrayView<const LO> getBlockRows(int blockIndex) const;
176
185 virtual void initialize () = 0;
186
188 void setBlockSizes(const Teuchos::Array<Teuchos::Array<LO> >& partitions);
189
190 void getMatDiag() const;
191
193 bool isInitialized () const;
194
196 bool isComputed () const;
197
206 virtual void compute () = 0;
207
208 void DoJacobi(ConstHostView X, HostView Y, SC dampingFactor) const;
209 void DoOverlappingJacobi(ConstHostView X, HostView Y, ConstHostView W, SC dampingFactor, bool nonsymScaling) const;
210 void DoGaussSeidel(ConstHostView X, HostView Y, HostView Y2, SC dampingFactor) const;
211 void DoSGS(ConstHostView X, HostView Y, HostView Y2, SC dampingFactor) const;
212
214 virtual void setParameters (const Teuchos::ParameterList& List) = 0;
215
228 virtual void
229 apply(ConstHostView X,
230 HostView Y,
231 int blockIndex,
232 Teuchos::ETransp mode = Teuchos::NO_TRANS,
233 SC alpha = Teuchos::ScalarTraits<SC>::one(),
234 SC beta = Teuchos::ScalarTraits<SC>::zero()) const = 0;
235
237 virtual void
238 weightedApply(ConstHostView X,
239 HostView Y,
240 ConstHostView D,
241 int blockIndex,
242 Teuchos::ETransp mode = Teuchos::NO_TRANS,
243 SC alpha = Teuchos::ScalarTraits<SC>::one(),
244 SC beta = Teuchos::ScalarTraits<SC>::zero()) const = 0;
245
248 // The underlying container implements the splitting <tt>A = D + R</tt>. Only
249 // it can have efficient access to D and R, as these are constructed in the
250 // symbolic and numeric phases.
251 //
252 // This is the first performance-portable implementation of a block
253 // relaxation, and it is supported currently only by BlockTriDiContainer.
254 virtual void applyInverseJacobi (const mv_type& /* X */, mv_type& /* Y */,
255 SC dampingFactor,
256 bool /* zeroStartingSolution = false */,
257 int /* numSweeps = 1 */) const = 0;
258
260 virtual void applyMV (const mv_type& X, mv_type& Y) const;
261
263 virtual void weightedApplyMV (const mv_type& X,
264 mv_type& Y,
265 vector_type& W) const;
266
267 virtual void clearBlocks();
268
270 virtual std::ostream& print (std::ostream& os) const = 0;
271
274 static std::string getName();
275
276protected:
277
279 virtual void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid,
280 SC dampingFactor, LO i) const;
281
283 Teuchos::RCP<const row_matrix_type> inputMatrix_;
284
286 Teuchos::RCP<const crs_matrix_type> inputCrsMatrix_;
287
289 Teuchos::RCP<const block_crs_matrix_type> inputBlockMatrix_;
290
294 Teuchos::Array<LO> blockRows_; //size: total # of local rows (in all local blocks)
296 Teuchos::Array<LO> blockSizes_; //size: # of blocks
298 Teuchos::Array<LO> blockOffsets_;
300 mutable Teuchos::RCP<vector_type> Diag_;
318
323
324 LO maxBlockSize_;
325};
326
327namespace Details
328{
329 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
330 struct StridedRowView;
331}
332
341
342template<class MatrixType, class LocalScalarType>
343class ContainerImpl : public Container<MatrixType>
344{
346
347protected:
348 using local_scalar_type = LocalScalarType;
349 using SC = typename Container<MatrixType>::scalar_type;
350 using LO = typename Container<MatrixType>::local_ordinal_type;
351 using GO = typename Container<MatrixType>::global_ordinal_type;
352 using NO = typename Container<MatrixType>::node_type;
354 using typename Container<MatrixType>::import_type;
355 using typename Container<MatrixType>::row_matrix_type;
356 using typename Container<MatrixType>::crs_matrix_type;
357 using typename Container<MatrixType>::block_crs_matrix_type;
358 using typename Container<MatrixType>::mv_type;
359 using typename Container<MatrixType>::vector_type;
360 using typename Container<MatrixType>::map_type;
361 using typename Container<MatrixType>::ISC;
363 using LSC = LocalScalarType;
364 using LISC = typename Kokkos::Details::ArithTraits<LSC>::val_type;
365
366 using local_mv_type = Tpetra::MultiVector<LSC, LO, GO, NO>;
367
368 using typename Container<MatrixType>::HostView;
369 using typename Container<MatrixType>::ConstHostView;
370 using HostViewLocal = typename local_mv_type::dual_view_type::t_host;
371 using HostSubviewLocal = Kokkos::View<LISC**, Kokkos::LayoutStride, typename HostViewLocal::memory_space>;
372 using ConstHostSubviewLocal = Kokkos::View<const LISC**, Kokkos::LayoutStride, typename HostViewLocal::memory_space>;
373
374 static_assert(std::is_same<MatrixType, row_matrix_type>::value,
375 "Ifpack2::Container: Please use MatrixType = Tpetra::RowMatrix.");
377 //
378
385 LO translateRowToCol(LO row);
386
389
394 virtual void extract() = 0;
395
396public:
397
398 ContainerImpl (const Teuchos::RCP<const row_matrix_type>& matrix,
399 const Teuchos::Array<Teuchos::Array<LO> >& partitions,
400 bool pointIndexed);
401
403 virtual ~ContainerImpl();
404
413 virtual void initialize () = 0;
414
423 virtual void compute () = 0;
424
426 virtual void setParameters (const Teuchos::ParameterList& List);
427
440 virtual void
441 apply(ConstHostView X,
442 HostView Y,
443 int blockIndex,
444 Teuchos::ETransp mode = Teuchos::NO_TRANS,
445 SC alpha = Teuchos::ScalarTraits<SC>::one(),
446 SC beta = Teuchos::ScalarTraits<SC>::zero()) const;
447
449 virtual void
450 weightedApply(ConstHostView X,
451 HostView Y,
452 ConstHostView D,
453 int blockIndex,
454 Teuchos::ETransp mode = Teuchos::NO_TRANS,
455 SC alpha = Teuchos::ScalarTraits<SC>::one(),
456 SC beta = Teuchos::ScalarTraits<SC>::zero()) const;
457
460 // The underlying container implements the splitting <tt>A = D + R</tt>. Only
461 // it can have efficient access to D and R, as these are constructed in the
462 // symbolic and numeric phases.
463 //
464 // This is the first performance-portable implementation of a block
465 // relaxation, and it is supported currently only by BlockTriDiContainer.
466 virtual void applyInverseJacobi (const mv_type& /* X */, mv_type& /* Y */,
467 SC dampingFactor,
468 bool /* zeroStartingSolution = false */,
469 int /* numSweeps = 1 */) const;
470
472 void applyMV (const mv_type& X, mv_type& Y) const;
473
475 void weightedApplyMV (const mv_type& X,
476 mv_type& Y,
477 vector_type& W) const;
478
479 virtual void clearBlocks();
480
482 virtual std::ostream& print (std::ostream& os) const = 0;
483
486 static std::string getName();
487
488protected:
489 //Do Gauss-Seidel on only block i (this is used by DoGaussSeidel and DoSGS)
490 void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid,
491 SC dampingFactor, LO i) const;
492
498 virtual void
499 solveBlock(ConstHostSubviewLocal X,
500 HostSubviewLocal Y,
501 int blockIndex,
502 Teuchos::ETransp mode,
503 const LSC alpha,
504 const LSC beta) const;
505
507 mutable HostViewLocal X_local_; //length: blockRows_.size()
508 mutable HostViewLocal Y_local_; //length: blockRows_.size()
509
520 mutable HostViewLocal weightedApplyScratch_;
521
523 mutable std::vector<HostSubviewLocal> X_localBlocks_;
524
526 mutable std::vector<HostSubviewLocal> Y_localBlocks_;
527};
528
529namespace Details {
535 template<typename Scalar, typename LocalOrdinal, typename GlobalOrdinal, typename Node>
537 {
538 using SC = Scalar;
539 using LO = LocalOrdinal;
540
541 using block_crs_matrix_type = Tpetra::BlockCrsMatrix<SC, LO, GlobalOrdinal, Node>;
542
543 using h_inds_type = typename block_crs_matrix_type::local_inds_host_view_type;
544 using h_vals_type = typename block_crs_matrix_type::values_host_view_type;
546 StridedRowView(h_vals_type vals_, h_inds_type inds_, int blockSize_, size_t nnz_);
547
549 // StridedRowView(const SC* vals_, const LO* inds_, int blockSize_, size_t nnz_);
550
552 StridedRowView(Teuchos::Array<SC>& vals_, Teuchos::Array<LO>& inds_);
553
554 SC val(size_t i) const;
555 LO ind(size_t i) const;
556
557 size_t size() const;
558
559 private:
560 h_vals_type vals;
561 h_inds_type inds;
562 int blockSize;
563 size_t nnz;
564 //These arrays are only used if the inputMatrix_ doesn't support row views.
565 Teuchos::Array<SC> valsCopy;
566 Teuchos::Array<LO> indsCopy;
567 };
568} // namespace Details
569
570} // namespace Ifpack2
571
573template <class MatrixType>
574std::ostream& operator<<(std::ostream& os, const Ifpack2::Container<MatrixType>& obj);
575
576namespace Teuchos {
577
582template<class MatrixType>
583class TEUCHOSCORE_LIB_DLL_EXPORT TypeNameTraits< ::Ifpack2::Container<MatrixType> >
584{
585 public:
586 static std::string name () {
587 return std::string ("Ifpack2::Container<") +
588 TypeNameTraits<MatrixType>::name () + ">";
589 }
590
591 static std::string concreteName (const ::Ifpack2::Container<MatrixType>&) {
592 return name ();
593 }
594};
595
596} // namespace Teuchos
597
598#endif // IFPACK2_CONTAINER_HPP
The implementation of the numerical features of Container (Jacobi, Gauss-Seidel, SGS)....
Definition: Ifpack2_Container_decl.hpp:344
std::vector< HostSubviewLocal > X_localBlocks_
Views for holding pieces of X corresponding to each block.
Definition: Ifpack2_Container_decl.hpp:523
HostViewLocal X_local_
Scratch vectors used in apply().
Definition: Ifpack2_Container_decl.hpp:507
virtual std::ostream & print(std::ostream &os) const =0
Print basic information about the container to os.
virtual void setParameters(const Teuchos::ParameterList &List)
Set parameters, if any.
Definition: Ifpack2_Container_def.hpp:530
Details::StridedRowView< SC, LO, GO, NO > getInputRowView(LO row) const
View a row of the input matrix.
Definition: Ifpack2_Container_def.hpp:898
virtual void extract()=0
Extract the submatrices identified by the local indices set by the constructor. BlockMatrix may be an...
std::vector< HostSubviewLocal > Y_localBlocks_
Views for holding pieces of Y corresponding to each block.
Definition: Ifpack2_Container_decl.hpp:526
virtual void weightedApply(ConstHostView X, HostView Y, ConstHostView D, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const
Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y.
Definition: Ifpack2_Container_def.hpp:739
virtual void compute()=0
Extract the local diagonal blocks and prepare the solver.
void applyMV(const mv_type &X, mv_type &Y) const
Wrapper for apply with MVs, used in unit tests (never called by BlockRelaxation)
Definition: Ifpack2_Container_def.hpp:544
virtual void apply(ConstHostView X, HostView Y, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const
Compute Y := alpha * M^{-1} X + beta*Y.
Definition: Ifpack2_Container_def.hpp:628
LocalScalarType LSC
The internal representation of LocalScalarType in Kokkos::View.
Definition: Ifpack2_Container_decl.hpp:363
void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid, SC dampingFactor, LO i) const
Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS)
Definition: Ifpack2_Container_def.hpp:286
virtual void solveBlock(ConstHostSubviewLocal X, HostSubviewLocal Y, int blockIndex, Teuchos::ETransp mode, const LSC alpha, const LSC beta) const
Definition: Ifpack2_Container_def.hpp:572
HostViewLocal weightedApplyScratch_
Definition: Ifpack2_Container_decl.hpp:520
static std::string getName()
Definition: Ifpack2_Container_def.hpp:565
LO translateRowToCol(LO row)
Definition: Ifpack2_Container_def.hpp:585
virtual void initialize()=0
Do all set-up operations that only require matrix structure.
virtual ~ContainerImpl()
Destructor.
Definition: Ifpack2_Container_def.hpp:526
void weightedApplyMV(const mv_type &X, mv_type &Y, vector_type &W) const
Wrapper for weightedApply with MVs, used in unit tests (never called by BlockRelaxation)
Definition: Ifpack2_Container_def.hpp:553
virtual void applyInverseJacobi(const mv_type &, mv_type &, SC dampingFactor, bool, int) const
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y).
Definition: Ifpack2_Container_def.hpp:534
Interface for creating and solving a set of local linear problems.
Definition: Ifpack2_Container_decl.hpp:113
virtual ~Container()
Destructor.
Definition: Ifpack2_Container_def.hpp:114
bool hasBlockCrs_
Whether the input matrix is a BlockCRS matrix.
Definition: Ifpack2_Container_decl.hpp:310
virtual void apply(ConstHostView X, HostView Y, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const =0
Compute Y := alpha * M^{-1} X + beta*Y.
Teuchos::RCP< const crs_matrix_type > inputCrsMatrix_
The input matrix, dynamic cast to CrsMatrix. May be null.
Definition: Ifpack2_Container_decl.hpp:286
bool IsParallel_
Whether the problem is distributed across multiple MPI processes.
Definition: Ifpack2_Container_decl.hpp:302
int numBlocks_
The number of blocks (partitions) in the container.
Definition: Ifpack2_Container_decl.hpp:292
Teuchos::ArrayView< const LO > getBlockRows(int blockIndex) const
Local indices of the rows of the input matrix that belong to this block.
Definition: Ifpack2_Container_def.hpp:118
static std::string getName()
Definition: Ifpack2_Container_def.hpp:190
virtual void applyInverseJacobi(const mv_type &, mv_type &, SC dampingFactor, bool, int) const =0
Compute Y := (1 - a) Y + a D^{-1} (X - R*Y).
Teuchos::RCP< vector_type > Diag_
Diagonal elements.
Definition: Ifpack2_Container_decl.hpp:300
virtual void weightedApply(ConstHostView X, HostView Y, ConstHostView D, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, SC alpha=Teuchos::ScalarTraits< SC >::one(), SC beta=Teuchos::ScalarTraits< SC >::zero()) const =0
Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y.
int bcrsBlockSize_
If hasBlockCrs_, the number of DOFs per vertex. Otherwise 1.
Definition: Ifpack2_Container_decl.hpp:312
Teuchos::RCP< const block_crs_matrix_type > inputBlockMatrix_
The input matrix, dynamic cast to BlockCrsMatrix. May be null.
Definition: Ifpack2_Container_decl.hpp:289
virtual void setParameters(const Teuchos::ParameterList &List)=0
Set parameters, if any.
typename Kokkos::Details::ArithTraits< SC >::val_type ISC
Internal representation of Scalar in Kokkos::View.
Definition: Ifpack2_Container_decl.hpp:135
virtual void DoGSBlock(ConstHostView X, HostView Y, HostView Y2, HostView Resid, SC dampingFactor, LO i) const
Do one step of Gauss-Seidel on block i (used by DoGaussSeidel and DoSGS)
Definition: Ifpack2_Container_def.hpp:196
virtual void applyMV(const mv_type &X, mv_type &Y) const
Wrapper for apply with MultiVector.
Definition: Ifpack2_Container_def.hpp:176
LO NumLocalRows_
Number of local rows in input matrix.
Definition: Ifpack2_Container_decl.hpp:304
bool isInitialized() const
Whether the container has been successfully initialized.
Definition: Ifpack2_Container_def.hpp:165
void setBlockSizes(const Teuchos::Array< Teuchos::Array< LO > > &partitions)
Initialize arrays with information about block sizes.
Definition: Ifpack2_Container_def.hpp:125
LO scalarsPerRow_
Definition: Ifpack2_Container_decl.hpp:317
Teuchos::Array< LO > blockSizes_
Number of rows in each block.
Definition: Ifpack2_Container_decl.hpp:296
GO NumGlobalNonzeros_
Number of nonzeros in input matrix.
Definition: Ifpack2_Container_decl.hpp:308
virtual void initialize()=0
Do all set-up operations that only require matrix structure.
virtual std::ostream & print(std::ostream &os) const =0
Print basic information about the container to os.
virtual void weightedApplyMV(const mv_type &X, mv_type &Y, vector_type &W) const
Wrapper for weightedApply with MultiVector.
Definition: Ifpack2_Container_def.hpp:183
Teuchos::RCP< const row_matrix_type > inputMatrix_
The input matrix to the constructor.
Definition: Ifpack2_Container_decl.hpp:283
bool pointIndexed_
(If hasBlockCrs_) Whether the blocks are described using sub-block row indices instead of full block ...
Definition: Ifpack2_Container_decl.hpp:314
bool IsInitialized_
If true, the container has been successfully initialized.
Definition: Ifpack2_Container_decl.hpp:320
Teuchos::Array< LO > blockRows_
Local indices of the rows of the input matrix that belong to this block.
Definition: Ifpack2_Container_decl.hpp:294
GO NumGlobalRows_
Number of global rows in input matrix.
Definition: Ifpack2_Container_decl.hpp:306
bool isComputed() const
Whether the container has been successfully computed.
Definition: Ifpack2_Container_def.hpp:170
typename mv_type::dual_view_type::t_host HostView
Definition: Ifpack2_Container_decl.hpp:139
Teuchos::Array< LO > blockOffsets_
Starting index in blockRows_ of local row indices for each block.
Definition: Ifpack2_Container_decl.hpp:298
bool IsComputed_
If true, the container has been successfully computed.
Definition: Ifpack2_Container_decl.hpp:322
virtual void compute()=0
Extract the local diagonal blocks and prepare the solver.
Ifpack2 implementation details.
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:74
Structure for read-only views of general matrix rows.
Definition: Ifpack2_Container_decl.hpp:537