42 #ifndef IFPACK2_BANDEDCONTAINER_DECL_HPP 43 #define IFPACK2_BANDEDCONTAINER_DECL_HPP 50 #include "Ifpack2_Details_LapackSupportsScalar.hpp" 51 #include "Tpetra_MultiVector.hpp" 52 #include "Tpetra_Map.hpp" 53 #include "Tpetra_RowMatrix.hpp" 54 #include "Teuchos_SerialDenseVector.hpp" 55 #include "Teuchos_SerialBandDenseMatrix.hpp" 104 template<
class MatrixType,
105 class LocalScalarType,
109 template<
class MatrixType,
class LocalScalarType>
121 typedef MatrixType matrix_type;
123 typedef LocalScalarType local_scalar_type;
125 typedef typename Kokkos::Details::ArithTraits<local_scalar_type>::val_type local_impl_scalar_type;
128 typedef typename Container<MatrixType>::scalar_type scalar_type;
130 typedef typename Container<MatrixType>::local_ordinal_type local_ordinal_type;
132 typedef typename Container<MatrixType>::global_ordinal_type global_ordinal_type;
134 typedef typename Container<MatrixType>::node_type node_type;
136 typedef typename Container<MatrixType>::mv_type mv_type;
137 typedef typename Container<MatrixType>::map_type map_type;
138 typedef Tpetra::MultiVector<local_scalar_type, local_ordinal_type, global_ordinal_type, node_type> local_mv_type;
139 typedef typename Container<MatrixType>::vector_type vector_type;
141 typedef typename Container<MatrixType>::import_type import_type;
143 typedef typename Container<MatrixType>::HostView HostView;
144 typedef typename local_mv_type::dual_view_type::t_host HostViewLocal;
146 static_assert(std::is_same<MatrixType,
147 Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> >::value,
148 "Ifpack2::BandedContainer: Please use MatrixType = Tpetra::RowMatrix.");
158 typedef typename Container<MatrixType>::row_matrix_type row_matrix_type;
182 const Teuchos::Array<Teuchos::Array<local_ordinal_type> >& partitions,
183 const Teuchos::RCP<const import_type>& importer,
185 scalar_type DampingFactor);
188 const Teuchos::Array<local_ordinal_type>& localRows);
198 virtual bool isInitialized ()
const {
199 return IsInitialized_;
203 virtual bool isComputed ()
const {
208 virtual void setParameters (
const Teuchos::ParameterList& List);
215 virtual void initialize ();
218 virtual void compute ();
228 Teuchos::ETransp mode = Teuchos::NO_TRANS,
229 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
230 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const;
234 weightedApply (HostView& X,
239 Teuchos::ETransp mode = Teuchos::NO_TRANS,
240 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
241 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const;
250 virtual std::ostream& print (std::ostream& os)
const;
257 virtual std::string description ()
const;
261 describe (Teuchos::FancyOStream &out,
262 const Teuchos::EVerbosityLevel verbLevel =
263 Teuchos::Describable::verbLevel_default)
const;
268 static std::string getName();
272 BandedContainer (
const BandedContainer<MatrixType, LocalScalarType>& rhs);
291 applyImpl (HostViewLocal& X,
295 Teuchos::ETransp mode,
296 const local_scalar_type alpha,
297 const local_scalar_type beta)
const;
300 std::vector<Teuchos::SerialBandDenseMatrix<int, local_scalar_type> > diagBlocks_;
303 mutable std::vector<HostViewLocal> X_local;
306 mutable std::vector<HostViewLocal> Y_local;
309 Teuchos::Array<int> ipiv_;
317 Teuchos::Array<local_ordinal_type> kl_;
318 Teuchos::Array<local_ordinal_type> ku_;
321 local_scalar_type* scalars_;
324 Teuchos::Array<local_ordinal_type> scalarOffsets_;
327 template<
class MatrixType,
class LocalScalarType>
328 class BandedContainer<MatrixType, LocalScalarType, false> :
329 public Container<MatrixType> {
339 typedef MatrixType matrix_type;
341 typedef LocalScalarType local_scalar_type;
343 typedef typename Kokkos::Details::ArithTraits<local_scalar_type>::val_type local_impl_scalar_type;
346 typedef typename Container<MatrixType>::scalar_type scalar_type;
348 typedef typename Container<MatrixType>::local_ordinal_type local_ordinal_type;
350 typedef typename Container<MatrixType>::global_ordinal_type global_ordinal_type;
352 typedef typename Container<MatrixType>::node_type node_type;
354 typedef typename Container<MatrixType>::mv_type mv_type;
355 typedef typename Container<MatrixType>::map_type map_type;
356 typedef Tpetra::MultiVector<local_scalar_type, local_ordinal_type, global_ordinal_type, node_type> local_mv_type;
357 typedef typename Container<MatrixType>::vector_type vector_type;
358 typedef typename Container<MatrixType>::partitioner_type partitioner_type;
359 typedef typename Container<MatrixType>::import_type import_type;
361 typedef typename Container<MatrixType>::HostView HostView;
362 typedef typename local_mv_type::dual_view_type::t_host HostViewLocal;
364 static_assert(std::is_same<MatrixType,
365 Tpetra::RowMatrix<scalar_type, local_ordinal_type, global_ordinal_type, node_type> >::value,
366 "Ifpack2::BandedContainer: Please use MatrixType = Tpetra::RowMatrix.");
376 typedef typename Container<MatrixType>::row_matrix_type row_matrix_type;
399 BandedContainer (
const Teuchos::RCP<const row_matrix_type>& matrix,
400 const Teuchos::Array<Teuchos::Array<local_ordinal_type> >& partitions,
401 const Teuchos::RCP<const import_type>& importer,
403 scalar_type DampingFactor);
405 BandedContainer (
const Teuchos::RCP<const row_matrix_type>& matrix,
406 const Teuchos::Array<local_ordinal_type>& localRows);
409 virtual ~BandedContainer ();
416 virtual bool isInitialized ()
const {
417 return IsInitialized_;
421 virtual bool isComputed ()
const {
426 virtual void setParameters (
const Teuchos::ParameterList& List);
433 virtual void initialize ();
436 virtual void compute ();
446 Teuchos::ETransp mode = Teuchos::NO_TRANS,
447 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
448 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const;
452 weightedApply (HostView& X,
457 Teuchos::ETransp mode = Teuchos::NO_TRANS,
458 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
459 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const;
468 virtual std::ostream& print (std::ostream& os)
const;
475 virtual std::string description ()
const;
479 describe (Teuchos::FancyOStream &out,
480 const Teuchos::EVerbosityLevel verbLevel =
481 Teuchos::Describable::verbLevel_default)
const;
486 static std::string getName();
490 BandedContainer (
const BandedContainer<MatrixType, LocalScalarType>& rhs);
509 applyImpl (HostViewLocal& X,
513 Teuchos::ETransp mode,
514 const local_scalar_type alpha,
515 const local_scalar_type beta)
const;
518 std::vector<Teuchos::SerialBandDenseMatrix<int, local_scalar_type> > diagBlocks_;
521 mutable std::vector<HostViewLocal> X_local;
524 mutable std::vector<HostViewLocal> Y_local;
527 Teuchos::Array<int> ipiv_;
535 Teuchos::Array<local_ordinal_type> kl_;
536 Teuchos::Array<local_ordinal_type> ku_;
539 local_scalar_type* scalars_;
542 Teuchos::Array<local_ordinal_type> scalarOffsets_;
547 #endif // IFPACK2_BANDEDCONTAINER_DECL_HPP Ifpack2::Container class declaration.
Declaration and definition of the Ifpack2::Details::MultiVectorLocalGatherScatter class...
Ifpack2::Partitioner:
Definition: Ifpack2_Partitioner.hpp:179
Interface for creating and solving a local linear problem.
Definition: Ifpack2_Container.hpp:114
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:72
Type traits class that says whether Teuchos::LAPACK has a valid implementation for the given ScalarTy...
Definition: Ifpack2_Details_LapackSupportsScalar.hpp:17
Store and solve a local Banded linear problem.
Definition: Ifpack2_BandedContainer_decl.hpp:107