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::BlockCrsMatrix
178 <scalar_type,local_ordinal_type,global_ordinal_type,node_type> block_crs_matrix_type;
186 typedef typename Container<MatrixType>::row_matrix_type row_matrix_type;
188 static_assert (std::is_same<MatrixType, row_matrix_type>::value,
189 "Ifpack2::BlockTriDiContainer: MatrixType must be a Tpetra::RowMatrix specialization.");
211 const Teuchos::Array<Teuchos::Array<local_ordinal_type> >& partitions,
212 const Teuchos::RCP<const import_type>& importer,
230 const Teuchos::Array<Teuchos::Array<local_ordinal_type> >& partitions,
231 bool overlapCommAndComp =
false,
bool useSequentialMethod =
false);
236 struct ComputeParameters {
243 magnitude_type addRadiallyToDiagonal = Kokkos::ArithTraits<magnitude_type>::zero();
247 struct ApplyParameters {
250 bool zeroStartingSolution =
false;
252 scalar_type dampingFactor = Kokkos::ArithTraits<scalar_type>::one();
255 int maxNumSweeps = 1;
265 magnitude_type tolerance = Kokkos::ArithTraits<magnitude_type>::zero();
273 int checkToleranceEvery = 1;
281 void setParameters(
const Teuchos::ParameterList& List)
override;
283 void clearBlocks()
override;
290 void initialize ()
override;
293 void compute ()
override;
296 void applyInverseJacobi (
const mv_type& X, mv_type& Y,
297 scalar_type dampingFactor,
298 bool zeroStartingSolution =
false,
299 int numSweeps = 1)
const override;
302 ComputeParameters createDefaultComputeParameters ()
const;
315 void compute (
const ComputeParameters& input);
318 ApplyParameters createDefaultApplyParameters ()
const;
326 int applyInverseJacobi (
const mv_type& X, mv_type& Y,
327 const ApplyParameters& input)
const;
332 const magnitude_type getNorms0 ()
const;
336 const magnitude_type getNormsFinal ()
const;
341 apply (const_host_view_type X,
344 Teuchos::ETransp mode = Teuchos::NO_TRANS,
345 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
346 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const override;
351 weightedApply (const_host_view_type X,
353 const_host_view_type W,
355 Teuchos::ETransp mode = Teuchos::NO_TRANS,
356 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
357 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const override;
366 std::ostream& print (std::ostream& os)
const override;
373 std::string description ()
const override;
377 describe (Teuchos::FancyOStream &out,
378 const Teuchos::EVerbosityLevel verbLevel =
379 Teuchos::Describable::verbLevel_default)
const override;
384 static std::string getName();
391 Teuchos::RCP<BlockTriDiContainerDetails::ImplObject<MatrixType> > impl_;
394 void initInternal (
const Teuchos::RCP<const row_matrix_type>& matrix,
395 const Teuchos::Array<Teuchos::Array<local_ordinal_type> >& partitions,
396 const Teuchos::RCP<const import_type> &importer,
397 const bool overlapCommAndComp,
398 const bool useSeqMethod);
400 void clearInternal();
414 typedef typename MatrixType::scalar_type scalar_type;
415 typedef typename Kokkos::ArithTraits<scalar_type>::magnitudeType magnitude_type;
416 typedef typename Container<MatrixType>::local_ordinal_type local_ordinal_type;
417 typedef typename Container<MatrixType>::global_ordinal_type global_ordinal_type;
419 typedef typename Container<MatrixType>::mv_type mv_type;
420 typedef typename Container<MatrixType>::import_type import_type;
423 typedef typename Container<MatrixType>::ConstHostView const_host_view_type;
424 typedef typename Container<MatrixType>::row_matrix_type row_matrix_type;
426 static_assert (std::is_same<MatrixType, row_matrix_type>::value,
427 "Ifpack2::BlockTriDiContainer: MatrixType must be a Tpetra::RowMatrix specialization.");
431 const Teuchos::Array<Teuchos::Array<local_ordinal_type> >& partitions,
432 const Teuchos::RCP<const import_type>& importer,
435 TEUCHOS_TEST_FOR_EXCEPT_MSG(
true,
"Error: BlockTriDiContainer is not available for this scalar_type");
439 void clearBlocks()
override {}
444 scalar_type dampingFactor,
445 bool zeroStartingSolution =
false,
446 int numSweeps = 1)
const override {}
452 Teuchos::ETransp mode = Teuchos::NO_TRANS,
453 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
454 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const override {}
459 const_host_view_type W,
461 Teuchos::ETransp mode = Teuchos::NO_TRANS,
462 scalar_type alpha = Teuchos::ScalarTraits<scalar_type>::one(),
463 scalar_type beta = Teuchos::ScalarTraits<scalar_type>::zero())
const override {}
465 std::ostream&
print (std::ostream& os)
const override {
466 return os <<
"Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
469 std::string description ()
const override {
470 return "Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
474 describe (Teuchos::FancyOStream &out,
475 const Teuchos::EVerbosityLevel verbLevel =
476 Teuchos::Describable::verbLevel_default)
const override {
477 out <<
"Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
480 static std::string getName() {
481 return "Ifpack2::BlockTriDiContainer::ImplNotAvailTag";
void weightedApply(const_host_view_type X, host_view_type Y, const_host_view_type W, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const override
Compute Y := alpha * diag(D) * M^{-1} (diag(D) * X) + beta*Y.
Definition Ifpack2_BlockTriDiContainer_decl.hpp:457
void apply(const_host_view_type X, host_view_type Y, int blockIndex, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const override
Compute Y := alpha * M^{-1} X + beta*Y.
Definition Ifpack2_BlockTriDiContainer_decl.hpp:449