Ifpack2 Templated Preconditioning Package Version 1.0
Loading...
Searching...
No Matches
Ifpack2_Details_ScaledDampedResidual_decl.hpp
Go to the documentation of this file.
1/*
2//@HEADER
3// ***********************************************************************
4//
5// Ifpack2: Templated Object-Oriented Algebraic Preconditioner Package
6// Copyright (2009) Sandia Corporation
7//
8// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
9// license for use of this work by or on behalf of the U.S. Government.
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// ***********************************************************************
39//@HEADER
40*/
41
47
48#ifndef IFPACK2_DETAILS_SCALEDDAMPEDRESIDUAL_DECL_HPP
49#define IFPACK2_DETAILS_SCALEDDAMPEDRESIDUAL_DECL_HPP
50
51#include "Ifpack2_config.h"
52#include "Tpetra_CrsMatrix_fwd.hpp"
53#include "Tpetra_MultiVector_fwd.hpp"
54#include "Tpetra_Operator_fwd.hpp"
55#include "Tpetra_Vector_fwd.hpp"
56#include "Tpetra_Export_fwd.hpp"
57#include "Tpetra_Import_fwd.hpp"
58#include "Teuchos_RCP.hpp"
59#include <memory>
60
61namespace Ifpack2 {
62namespace Details {
63
76template<class TpetraOperatorType>
78private:
79 using SC = typename TpetraOperatorType::scalar_type;
80 using LO = typename TpetraOperatorType::local_ordinal_type;
81 using GO = typename TpetraOperatorType::global_ordinal_type;
82 using NT = typename TpetraOperatorType::node_type;
83
84 using crs_matrix_type = Tpetra::CrsMatrix<SC, LO, GO, NT>;
85 using multivector_type = Tpetra::MultiVector<SC, LO, GO, NT>;
86 using operator_type = Tpetra::Operator<SC, LO, GO, NT>;
87 using vector_type = Tpetra::Vector<SC, LO, GO, NT>;
88
89public:
90 ScaledDampedResidual (const Teuchos::RCP<const operator_type>& A);
91
92 void
93 setMatrix (const Teuchos::RCP<const operator_type>& A);
94
95 void
96 compute (multivector_type& W,
97 const SC& alpha,
98 vector_type& D_inv,
99 multivector_type& B,
100 multivector_type& X,
101 const SC& beta);
102
103private:
104 using import_type = Tpetra::Import<LO, GO, NT>;
105 using export_type = Tpetra::Export<LO, GO, NT>;
106
107 Teuchos::RCP<const operator_type> A_op_;
108 Teuchos::RCP<const crs_matrix_type> A_crs_;
109 Teuchos::RCP<const import_type> imp_;
110 Teuchos::RCP<const export_type> exp_;
111 std::unique_ptr<vector_type> X_colMap_;
112 std::unique_ptr<multivector_type> V1_;
113
114 Teuchos::RCP<vector_type> W_vec_, B_vec_, X_vec_;
115
116 // Do the Import, if needed, and return the column Map version of X.
117 vector_type&
118 importVector (vector_type& X_domMap);
119
120 bool canFuse (const multivector_type& B) const;
121
122 void
123 unfusedCase (multivector_type& W,
124 const SC& alpha,
125 vector_type& D_inv,
126 multivector_type& B,
127 const operator_type& A,
128 multivector_type& X,
129 const SC& beta);
130
131 void
132 fusedCase (vector_type& W,
133 const SC& alpha,
134 vector_type& D_inv,
135 vector_type& B,
136 const crs_matrix_type& A,
137 vector_type& X,
138 const SC& beta);
139};
140
141} // namespace Details
142} // namespace Ifpack2
143
144#endif // IFPACK2_DETAILS_SCALEDDAMPEDRESIDUAL_DECL_HPP
Compute scaled damped residual for Chebyshev.
Definition: Ifpack2_Details_ScaledDampedResidual_decl.hpp:77
Ifpack2 implementation details.
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:74