43 #ifndef IFPACK2_BLOCKTRIDICONTAINER_DECL_HPP
44 #define IFPACK2_BLOCKTRIDICONTAINER_DECL_HPP
49 #include "Ifpack2_config.h"
50 #include "Ifpack2_Container.hpp"
51 #include "Tpetra_MultiVector.hpp"
52 #include "Tpetra_Map.hpp"
53 #include "Tpetra_RowMatrix.hpp"
54 #include "Tpetra_BlockCrsMatrix_decl.hpp"
55 #include <type_traits>
111 namespace BlockTriDiContainerDetails {
116 struct ImplSimdTag {};
117 struct ImplSacadoTag {};
119 template<
typename T>
struct ImplTag {
typedef ImplNotAvailTag type; };
120 template<>
struct ImplTag<float> {
typedef ImplSimdTag type; };
121 template<>
struct ImplTag<double> {
typedef ImplSimdTag type; };
122 template<>
struct ImplTag<std::complex<float> > {
typedef ImplSimdTag type; };
123 template<>
struct ImplTag<std::complex<double> > {
typedef ImplSimdTag type; };
132 template <
typename MatrixType,
141 template <
typename MatrixType>
152 typedef MatrixType matrix_type;
155 typedef typename MatrixType::scalar_type scalar_type;
157 typedef typename Kokkos::ArithTraits<scalar_type>::magnitudeType magnitude_type;
159 typedef typename Container<MatrixType>::local_ordinal_type local_ordinal_type;
161 typedef typename Container<MatrixType>::global_ordinal_type global_ordinal_type;
163 typedef typename Container<MatrixType>::node_type node_type;
165 typedef typename Container<MatrixType>::mv_type mv_type;
166 typedef typename Container<MatrixType>::map_type map_type;
167 typedef typename Container<MatrixType>::vector_type vector_type;
168 typedef typename Container<MatrixType>::import_type import_type;
171 typedef typename Container<MatrixType>::ConstHostView const_host_view_type;
172 typedef host_view_type HostView;
173 typedef const_host_view_type ConstHostView;
177 typedef Tpetra::CrsMatrix
178 <scalar_type,local_ordinal_type,global_ordinal_type,node_type> crs_matrix_type;
179 typedef Tpetra::BlockCrsMatrix
180 <scalar_type,local_ordinal_type,global_ordinal_type,node_type> block_crs_matrix_type;
190 typedef typename Container<MatrixType>::row_matrix_type row_matrix_type;
192 static_assert (std::is_same<MatrixType, row_matrix_type>::value,
193 "Ifpack2::BlockTriDiContainer: MatrixType must be a Tpetra::RowMatrix specialization.");
235 const int n_subparts_per_part = 1,
236 bool overlapCommAndComp =
false,
237 bool useSequentialMethod =
false,
238 const int block_size = -1);
243 struct ComputeParameters {
250 magnitude_type addRadiallyToDiagonal = Kokkos::ArithTraits<magnitude_type>::zero();
254 struct ApplyParameters {
257 bool zeroStartingSolution =
false;
259 scalar_type dampingFactor = Kokkos::ArithTraits<scalar_type>::one();
262 int maxNumSweeps = 1;
272 magnitude_type tolerance = Kokkos::ArithTraits<magnitude_type>::zero();
280 int checkToleranceEvery = 1;
290 void clearBlocks()
override;
297 void initialize ()
override;
300 void compute ()
override;
303 void applyInverseJacobi (
const mv_type& X, mv_type& Y,
304 scalar_type dampingFactor,
305 bool zeroStartingSolution =
false,
306 int numSweeps = 1)
const override;
309 ComputeParameters createDefaultComputeParameters ()
const;
322 void compute (
const ComputeParameters& input);
325 ApplyParameters createDefaultApplyParameters ()
const;
333 int applyInverseJacobi (
const mv_type& X, mv_type& Y,
334 const ApplyParameters& input)
const;
339 const magnitude_type getNorms0 ()
const;
343 const magnitude_type getNormsFinal ()
const;
348 apply (const_host_view_type X,
358 weightedApply (const_host_view_type X,
360 const_host_view_type W,
373 std::ostream& print (std::ostream& os)
const override;
380 std::string description ()
const override;
391 static std::string getName();
399 int n_subparts_per_part_;
404 const bool overlapCommAndComp,
405 const bool useSeqMethod,
406 const int block_size = -1);
408 void clearInternal();
418 template <
typename MatrixType>
422 typedef typename MatrixType::scalar_type scalar_type;
423 typedef typename Kokkos::ArithTraits<scalar_type>::magnitudeType magnitude_type;
424 typedef typename Container<MatrixType>::local_ordinal_type local_ordinal_type;
425 typedef typename Container<MatrixType>::global_ordinal_type global_ordinal_type;
427 typedef typename Container<MatrixType>::mv_type mv_type;
428 typedef typename Container<MatrixType>::import_type import_type;
431 typedef typename Container<MatrixType>::ConstHostView const_host_view_type;
432 typedef typename Container<MatrixType>::row_matrix_type row_matrix_type;
434 static_assert (std::is_same<MatrixType, row_matrix_type>::value,
435 "Ifpack2::BlockTriDiContainer: MatrixType must be a Tpetra::RowMatrix specialization.");
447 void clearBlocks()
override {}
451 void applyInverseJacobi (
const mv_type& X, mv_type& Y,
452 scalar_type dampingFactor,
453 bool zeroStartingSolution =
false,
454 int numSweeps = 1)
const override {}
457 apply (const_host_view_type X,
465 weightedApply (const_host_view_type X,
467 const_host_view_type W,
473 std::ostream&
print (std::ostream& os)
const override {
474 return os <<
"Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
477 std::string description ()
const override {
478 return "Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
485 out <<
"Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
488 static std::string getName() {
489 return "Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
496 #endif // IFPACK2_BLOCKTRIDICONTAINER_DECL_HPP
Definition: Ifpack2_BlockTriDiContainer_decl.hpp:115
void compute() override
Extract the local diagonal blocks and prepare the solver.
Definition: Ifpack2_BlockTriDiContainer_decl.hpp:450
Store and solve local block tridiagonal linear problems.
Definition: Ifpack2_BlockTriDiContainer_decl.hpp:134
void initialize() override
Do all set-up operations that only require matrix structure.
Definition: Ifpack2_BlockTriDiContainer_decl.hpp:449
#define TEUCHOS_TEST_FOR_EXCEPT_MSG(throw_exception_test, msg)
forward declaration
Definition: Ifpack2_BlockTriDiContainer_decl.hpp:126
typename mv_type::dual_view_type::t_host HostView
Definition: Ifpack2_Container_decl.hpp:139
std::ostream & print(std::ostream &os) const override
Print basic information about the container to os.
Definition: Ifpack2_BlockTriDiContainer_decl.hpp:473
void setParameters(const Teuchos::ParameterList &List) override
Set parameters, if any.
Definition: Ifpack2_BlockTriDiContainer_decl.hpp:446
static const EVerbosityLevel verbLevel_default
Interface for creating and solving a set of local linear problems.
Definition: Ifpack2_Container_decl.hpp:112