Ifpack2 Templated Preconditioning Package Version 1.0
Loading...
Searching...
No Matches
Ifpack2_BlockRelaxation_def.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_BLOCKRELAXATION_DEF_HPP
44#define IFPACK2_BLOCKRELAXATION_DEF_HPP
45
47#include "Ifpack2_LinearPartitioner.hpp"
48#include "Ifpack2_LinePartitioner.hpp"
50#include "Ifpack2_Details_UserPartitioner_def.hpp"
51#include <Ifpack2_Parameters.hpp>
52#include "Teuchos_TimeMonitor.hpp"
53
54namespace Ifpack2 {
55
56template<class MatrixType,class ContainerType>
58setMatrix (const Teuchos::RCP<const row_matrix_type>& A)
59{
60 if (A.getRawPtr () != A_.getRawPtr ()) { // it's a different matrix
61 IsInitialized_ = false;
62 IsComputed_ = false;
63 Partitioner_ = Teuchos::null;
64 W_ = Teuchos::null;
65
66 if (! A.is_null ()) {
67 IsParallel_ = (A->getRowMap ()->getComm ()->getSize () != 1);
68 Teuchos::RCP<const block_crs_matrix_type> A_bcrs =
69 Teuchos::rcp_dynamic_cast<const block_crs_matrix_type> (A);
70 hasBlockCrsMatrix_ = !A_bcrs.is_null();
71 }
72 if (! Container_.is_null ()) {
73 //This just frees up the container's memory.
74 //The container will be fully re-constructed during initialize().
75 Container_->clearBlocks();
76 }
77 NumLocalBlocks_ = 0;
78
79 A_ = A;
80 computeImporter();
81 }
82}
83
84template<class MatrixType,class ContainerType>
86BlockRelaxation (const Teuchos::RCP<const row_matrix_type>& A)
87:
88 Container_ (Teuchos::null),
89 Partitioner_ (Teuchos::null),
90 PartitionerType_ ("linear"),
91 NumSweeps_ (1),
92 NumLocalBlocks_(0),
93 containerType_ ("TriDi"),
94 PrecType_ (Ifpack2::Details::JACOBI),
95 ZeroStartingSolution_ (true),
96 hasBlockCrsMatrix_ (false),
97 DoBackwardGS_ (false),
98 OverlapLevel_ (0),
99 nonsymCombine_(0),
100 schwarzCombineMode_("ZERO"),
101 DampingFactor_ (STS::one ()),
102 IsInitialized_ (false),
103 IsComputed_ (false),
104 NumInitialize_ (0),
105 NumCompute_ (0),
106 NumApply_ (0),
107 InitializeTime_ (0.0),
108 ComputeTime_ (0.0),
109 NumLocalRows_ (0),
110 NumGlobalRows_ (0),
111 NumGlobalNonzeros_ (0),
112 W_ (Teuchos::null),
113 Importer_ (Teuchos::null)
114{
115 setMatrix(A);
116}
117
118template<class MatrixType,class ContainerType>
121{}
122
123template<class MatrixType,class ContainerType>
124Teuchos::RCP<const Teuchos::ParameterList>
126getValidParameters () const
127{
128 Teuchos::RCP<Teuchos::ParameterList> validParams = Teuchos::parameterList ("Ifpack2::BlockRelaxation");
129
130 validParams->set("relaxation: container", "TriDi");
131 validParams->set("relaxation: backward mode",false);
132 validParams->set("relaxation: type", "Jacobi");
133 validParams->set("relaxation: sweeps", 1);
134 validParams->set("relaxation: damping factor", STS::one());
135 validParams->set("relaxation: zero starting solution", true);
136 validParams->set("block relaxation: decouple dofs", false);
137 validParams->set("schwarz: compute condest", false); // mfh 24 Mar 2015: for backwards compatibility ONLY
138 validParams->set("schwarz: combine mode", "ZERO"); // use string mode for this
139 validParams->set("schwarz: use reordering", true);
140 validParams->set("schwarz: filter singletons", false);
141 validParams->set("schwarz: overlap level", 0);
142 validParams->set("partitioner: type", "greedy");
143 validParams->set("partitioner: local parts", 1);
144 validParams->set("partitioner: overlap", 0);
145 validParams->set("partitioner: combine mode", "ZERO"); // use string mode for this
146 Teuchos::Array<Teuchos::ArrayRCP<int>> tmp0;
147 validParams->set("partitioner: parts", tmp0);
148 Teuchos::Array<Teuchos::ArrayRCP<typename MatrixType::global_ordinal_type> > tmp1;
149 validParams->set("partitioner: global ID parts", tmp1);
150 validParams->set("partitioner: nonsymmetric overlap combine", false);
151 validParams->set("partitioner: maintain sparsity", false);
152 validParams->set("fact: ilut level-of-fill", 1.0);
153 validParams->set("fact: absolute threshold", 0.0);
154 validParams->set("fact: relative threshold", 1.0);
155 validParams->set("fact: relax value", 0.0);
156
157 Teuchos::ParameterList dummyList;
158 validParams->set("Amesos2",dummyList);
159 validParams->sublist("Amesos2").disableRecursiveValidation();
160 validParams->set("Amesos2 solver name", "KLU2");
161
162 Teuchos::ArrayRCP<int> tmp;
163 validParams->set("partitioner: map", tmp);
164
165 validParams->set("partitioner: line detection threshold", 0.0);
166 validParams->set("partitioner: PDE equations", 1);
167 Teuchos::RCP<Tpetra::MultiVector<typename Teuchos::ScalarTraits<typename MatrixType::scalar_type>::magnitudeType,
168 typename MatrixType::local_ordinal_type,
169 typename MatrixType::global_ordinal_type,
170 typename MatrixType::node_type> > dummy;
171 validParams->set("partitioner: coordinates",dummy);
172
173 return validParams;
174}
175
176template<class MatrixType,class ContainerType>
177void
179setParameters (const Teuchos::ParameterList& pl)
180{
181 // CAG: Copied form Relaxation
182 // FIXME (aprokop 18 Oct 2013) Casting away const is bad here.
183 // but otherwise, we will get [unused] in pl
184 this->setParametersImpl(const_cast<Teuchos::ParameterList&>(pl));
185}
186
187template<class MatrixType,class ContainerType>
188void
190setParametersImpl (Teuchos::ParameterList& List)
191{
192 if (List.isType<double>("relaxation: damping factor")) {
193 // Make sure that ST=complex can run with a damping factor that is
194 // a double.
195 scalar_type df = List.get<double>("relaxation: damping factor");
196 List.remove("relaxation: damping factor");
197 List.set("relaxation: damping factor",df);
198 }
199
200 // Note that the validation process does not change List.
201 Teuchos::RCP<const Teuchos::ParameterList> validparams;
202 validparams = this->getValidParameters();
203 List.validateParameters (*validparams);
204
205 if (List.isParameter ("relaxation: container")) {
206 // If the container type isn't a string, this will throw, but it
207 // rightfully should.
208
209 // If its value does not match the currently registered Container types,
210 // the ContainerFactory will throw with an informative message.
211 containerType_ = List.get<std::string> ("relaxation: container");
212 // Intercept more human-readable aliases for the *TriDi containers
213 if(containerType_ == "Tridiagonal") {
214 containerType_ = "TriDi";
215 }
216 if(containerType_ == "Block Tridiagonal") {
217 containerType_ = "BlockTriDi";
218 }
219 }
220
221 if (List.isParameter ("relaxation: type")) {
222 const std::string relaxType = List.get<std::string> ("relaxation: type");
223
224 if (relaxType == "Jacobi") {
225 PrecType_ = Ifpack2::Details::JACOBI;
226 }
227 else if (relaxType == "MT Split Jacobi") {
228 PrecType_ = Ifpack2::Details::MTSPLITJACOBI;
229 }
230 else if (relaxType == "Gauss-Seidel") {
231 PrecType_ = Ifpack2::Details::GS;
232 }
233 else if (relaxType == "Symmetric Gauss-Seidel") {
234 PrecType_ = Ifpack2::Details::SGS;
235 }
236 else {
237 TEUCHOS_TEST_FOR_EXCEPTION
238 (true, std::invalid_argument, "Ifpack2::BlockRelaxation::"
239 "setParameters: Invalid parameter value \"" << relaxType
240 << "\" for parameter \"relaxation: type\".");
241 }
242 }
243
244 if (List.isParameter ("relaxation: sweeps")) {
245 NumSweeps_ = List.get<int> ("relaxation: sweeps");
246 }
247
248 // Users may specify this as a floating-point literal. In that
249 // case, it may have type double, even if scalar_type is something
250 // else. We take care to try the various types that make sense.
251 if (List.isParameter ("relaxation: damping factor")) {
252 if (List.isType<double> ("relaxation: damping factor")) {
253 const double dampingFactor =
254 List.get<double> ("relaxation: damping factor");
255 DampingFactor_ = static_cast<scalar_type> (dampingFactor);
256 }
257 else if (List.isType<scalar_type> ("relaxation: damping factor")) {
258 DampingFactor_ = List.get<scalar_type> ("relaxation: damping factor");
259 }
260 else if (List.isType<magnitude_type> ("relaxation: damping factor")) {
261 const magnitude_type dampingFactor =
262 List.get<magnitude_type> ("relaxation: damping factor");
263 DampingFactor_ = static_cast<scalar_type> (dampingFactor);
264 }
265 else {
266 TEUCHOS_TEST_FOR_EXCEPTION
267 (true, std::invalid_argument, "Ifpack2::BlockRelaxation::"
268 "setParameters: Parameter \"relaxation: damping factor\" "
269 "has the wrong type.");
270 }
271 }
272
273 if (List.isParameter ("relaxation: zero starting solution")) {
274 ZeroStartingSolution_ = List.get<bool> ("relaxation: zero starting solution");
275 }
276
277 if (List.isParameter ("relaxation: backward mode")) {
278 DoBackwardGS_ = List.get<bool> ("relaxation: backward mode");
279 }
280
281 if (List.isParameter ("partitioner: type")) {
282 PartitionerType_ = List.get<std::string> ("partitioner: type");
283 }
284
285 // Users may specify this as an int literal, so we need to allow
286 // both int and local_ordinal_type (not necessarily same as int).
287 if (List.isParameter ("partitioner: local parts")) {
288 if (List.isType<local_ordinal_type> ("partitioner: local parts")) {
289 NumLocalBlocks_ = List.get<local_ordinal_type> ("partitioner: local parts");
290 }
291 else if (! std::is_same<int, local_ordinal_type>::value &&
292 List.isType<int> ("partitioner: local parts")) {
293 NumLocalBlocks_ = List.get<int> ("partitioner: local parts");
294 }
295 else {
296 TEUCHOS_TEST_FOR_EXCEPTION
297 (true, std::invalid_argument, "Ifpack2::BlockRelaxation::"
298 "setParameters: Parameter \"partitioner: local parts\" "
299 "has the wrong type.");
300 }
301 }
302
303 if (List.isParameter ("partitioner: overlap level")) {
304 if (List.isType<int> ("partitioner: overlap level")) {
305 OverlapLevel_ = List.get<int> ("partitioner: overlap level");
306 }
307 else if (! std::is_same<int, local_ordinal_type>::value &&
308 List.isType<local_ordinal_type> ("partitioner: overlap level")) {
309 OverlapLevel_ = List.get<local_ordinal_type> ("partitioner: overlap level");
310 }
311 else {
312 TEUCHOS_TEST_FOR_EXCEPTION
313 (true, std::invalid_argument, "Ifpack2::BlockRelaxation::"
314 "setParameters: Parameter \"partitioner: overlap level\" "
315 "has the wrong type.");
316 }
317 }
318 // when using global ID parts, assume that some blocks overlap even if
319 // user did not explicitly set the overlap level in the input file.
320 if ( ( List.isParameter("partitioner: global ID parts")) && (OverlapLevel_ < 1)) OverlapLevel_ = 1;
321
322 if (List.isParameter ("partitioner: nonsymmetric overlap combine"))
323 nonsymCombine_ = List.get<bool> ("partitioner: nonsymmetric overlap combine");
324
325 if (List.isParameter ("partitioner: combine mode"))
326 schwarzCombineMode_ = List.get<std::string> ("partitioner: combine mode");
327
328 std::string defaultContainer = "TriDi";
329 if(std::is_same<ContainerType, Container<MatrixType> >::value)
330 {
331 //Generic container template parameter, container type specified in List
332 Ifpack2::getParameter(List, "relaxation: container", defaultContainer);
333 }
334 // check parameters
335 if (PrecType_ != Ifpack2::Details::JACOBI) {
336 OverlapLevel_ = 0;
337 }
338 if (NumLocalBlocks_ < static_cast<local_ordinal_type> (0)) {
339 NumLocalBlocks_ = A_->getLocalNumRows() / (-NumLocalBlocks_);
340 }
341
342 decouple_ = false;
343 if(List.isParameter("block relaxation: decouple dofs"))
344 decouple_ = List.get<bool>("block relaxation: decouple dofs");
345
346 // other checks are performed in Partitioner_
347
348 // NTS: Sanity check to be removed at a later date when Backward mode is enabled
349 TEUCHOS_TEST_FOR_EXCEPTION(
350 DoBackwardGS_, std::runtime_error,
351 "Ifpack2::BlockRelaxation:setParameters: Setting the \"relaxation: "
352 "backward mode\" parameter to true is not yet supported.");
353
354 // copy the list as each subblock's constructor will
355 // require it later
356 List_ = List;
357}
358
359template<class MatrixType,class ContainerType>
360Teuchos::RCP<const Teuchos::Comm<int> >
362{
363 TEUCHOS_TEST_FOR_EXCEPTION
364 (A_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::getComm: "
365 "The matrix is null. You must call setMatrix() with a nonnull matrix "
366 "before you may call this method.");
367 return A_->getComm ();
368}
369
370template<class MatrixType,class ContainerType>
371Teuchos::RCP<const Tpetra::RowMatrix<typename MatrixType::scalar_type,
372 typename MatrixType::local_ordinal_type,
373 typename MatrixType::global_ordinal_type,
374 typename MatrixType::node_type> >
376 return A_;
377}
378
379template<class MatrixType,class ContainerType>
380Teuchos::RCP<const Tpetra::Map<typename MatrixType::local_ordinal_type,
381 typename MatrixType::global_ordinal_type,
382 typename MatrixType::node_type> >
384getDomainMap () const
385{
386 TEUCHOS_TEST_FOR_EXCEPTION
387 (A_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::"
388 "getDomainMap: The matrix is null. You must call setMatrix() with a "
389 "nonnull matrix before you may call this method.");
390 return A_->getDomainMap ();
391}
392
393template<class MatrixType,class ContainerType>
394Teuchos::RCP<const Tpetra::Map<typename MatrixType::local_ordinal_type,
395 typename MatrixType::global_ordinal_type,
396 typename MatrixType::node_type> >
398getRangeMap () const
399{
400 TEUCHOS_TEST_FOR_EXCEPTION
401 (A_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::"
402 "getRangeMap: The matrix is null. You must call setMatrix() with a "
403 "nonnull matrix before you may call this method.");
404 return A_->getRangeMap ();
405}
406
407template<class MatrixType,class ContainerType>
408bool
410hasTransposeApply () const {
411 return true;
412}
413
414template<class MatrixType,class ContainerType>
415int
417getNumInitialize () const {
418 return NumInitialize_;
419}
420
421template<class MatrixType,class ContainerType>
422int
424getNumCompute () const
425{
426 return NumCompute_;
427}
428
429template<class MatrixType,class ContainerType>
430int
432getNumApply () const
433{
434 return NumApply_;
435}
436
437template<class MatrixType,class ContainerType>
438double
440getInitializeTime () const
441{
442 return InitializeTime_;
443}
444
445template<class MatrixType,class ContainerType>
446double
448getComputeTime () const
449{
450 return ComputeTime_;
451}
452
453template<class MatrixType,class ContainerType>
454double
456getApplyTime () const
457{
458 return ApplyTime_;
459}
460
461
462template<class MatrixType,class ContainerType>
464 TEUCHOS_TEST_FOR_EXCEPTION(
465 A_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::getNodeSmootherComplexity: "
466 "The input matrix A is null. Please call setMatrix() with a nonnull "
467 "input matrix, then call compute(), before calling this method.");
468 // Relaxation methods cost roughly one apply + one block-diagonal inverse per iteration
469 // NOTE: This approximates all blocks as dense, which may overstate the cost if you have a sparse (or banded) container.
470 size_t block_nnz = 0;
471 for (local_ordinal_type i = 0; i < NumLocalBlocks_; ++i)
472 block_nnz += Partitioner_->numRowsInPart(i) *Partitioner_->numRowsInPart(i);
473
474 return block_nnz + A_->getLocalNumEntries();
475}
476
477template<class MatrixType,class ContainerType>
478void
480apply (const Tpetra::MultiVector<typename MatrixType::scalar_type,
481 typename MatrixType::local_ordinal_type,
482 typename MatrixType::global_ordinal_type,
483 typename MatrixType::node_type>& X,
484 Tpetra::MultiVector<typename MatrixType::scalar_type,
485 typename MatrixType::local_ordinal_type,
486 typename MatrixType::global_ordinal_type,
487 typename MatrixType::node_type>& Y,
488 Teuchos::ETransp mode,
489 scalar_type alpha,
490 scalar_type beta) const
491{
492 TEUCHOS_TEST_FOR_EXCEPTION
493 (A_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::apply: "
494 "The matrix is null. You must call setMatrix() with a nonnull matrix, "
495 "then call initialize() and compute() (in that order), before you may "
496 "call this method.");
497 TEUCHOS_TEST_FOR_EXCEPTION(
498 ! isComputed (), std::runtime_error, "Ifpack2::BlockRelaxation::apply: "
499 "isComputed() must be true prior to calling apply.");
500 TEUCHOS_TEST_FOR_EXCEPTION(
501 X.getNumVectors () != Y.getNumVectors (), std::invalid_argument,
502 "Ifpack2::BlockRelaxation::apply: X.getNumVectors() = "
503 << X.getNumVectors () << " != Y.getNumVectors() = "
504 << Y.getNumVectors () << ".");
505 TEUCHOS_TEST_FOR_EXCEPTION(
506 mode != Teuchos::NO_TRANS, std::logic_error, "Ifpack2::BlockRelaxation::"
507 "apply: This method currently only implements the case mode == Teuchos::"
508 "NO_TRANS. Transposed modes are not currently supported.");
509 TEUCHOS_TEST_FOR_EXCEPTION(
510 alpha != Teuchos::ScalarTraits<scalar_type>::one (), std::logic_error,
511 "Ifpack2::BlockRelaxation::apply: This method currently only implements "
512 "the case alpha == 1. You specified alpha = " << alpha << ".");
513 TEUCHOS_TEST_FOR_EXCEPTION(
514 beta != Teuchos::ScalarTraits<scalar_type>::zero (), std::logic_error,
515 "Ifpack2::BlockRelaxation::apply: This method currently only implements "
516 "the case beta == 0. You specified beta = " << beta << ".");
517
518 const std::string timerName ("Ifpack2::BlockRelaxation::apply");
519 Teuchos::RCP<Teuchos::Time> timer = Teuchos::TimeMonitor::lookupCounter (timerName);
520 if (timer.is_null ()) {
521 timer = Teuchos::TimeMonitor::getNewCounter (timerName);
522 }
523
524 double startTime = timer->wallTime();
525
526 {
527 Teuchos::TimeMonitor timeMon (*timer);
528
529 // If X and Y are pointing to the same memory location,
530 // we need to create an auxiliary vector, Xcopy
531 Teuchos::RCP<const MV> X_copy;
532 {
533 if (X.aliases(Y)) {
534 X_copy = rcp (new MV (X, Teuchos::Copy));
535 } else {
536 X_copy = rcpFromRef (X);
537 }
538 }
539
540 if (ZeroStartingSolution_) {
541 Y.putScalar (STS::zero ());
542 }
543
544 // Flops are updated in each of the following.
545 switch (PrecType_) {
546 case Ifpack2::Details::JACOBI:
547 ApplyInverseJacobi(*X_copy,Y);
548 break;
549 case Ifpack2::Details::GS:
550 ApplyInverseGS(*X_copy,Y);
551 break;
552 case Ifpack2::Details::SGS:
553 ApplyInverseSGS(*X_copy,Y);
554 break;
555 case Ifpack2::Details::MTSPLITJACOBI:
556 //note: for this method, the container is always BlockTriDi
557 Container_->applyInverseJacobi(*X_copy, Y, DampingFactor_, ZeroStartingSolution_, NumSweeps_);
558 break;
559 default:
560 TEUCHOS_TEST_FOR_EXCEPTION
561 (true, std::logic_error, "Ifpack2::BlockRelaxation::apply: Invalid "
562 "PrecType_ enum value " << PrecType_ << ". Valid values are Ifpack2::"
563 "Details::JACOBI = " << Ifpack2::Details::JACOBI << ", Ifpack2::Details"
564 "::GS = " << Ifpack2::Details::GS << ", and Ifpack2::Details::SGS = "
565 << Ifpack2::Details::SGS << ". Please report this bug to the Ifpack2 "
566 "developers.");
567 }
568 }
569
570 ApplyTime_ += (timer->wallTime() - startTime);
571 ++NumApply_;
572}
573
574template<class MatrixType,class ContainerType>
575void
577applyMat (const Tpetra::MultiVector<typename MatrixType::scalar_type,
578 typename MatrixType::local_ordinal_type,
579 typename MatrixType::global_ordinal_type,
580 typename MatrixType::node_type>& X,
581 Tpetra::MultiVector<typename MatrixType::scalar_type,
582 typename MatrixType::local_ordinal_type,
583 typename MatrixType::global_ordinal_type,
584 typename MatrixType::node_type>& Y,
585 Teuchos::ETransp mode) const
586{
587 A_->apply (X, Y, mode);
588}
589
590template<class MatrixType,class ContainerType>
591void
593initialize ()
594{
595 using Teuchos::rcp;
596 typedef Tpetra::RowGraph<local_ordinal_type, global_ordinal_type, node_type>
597 row_graph_type;
598
599 TEUCHOS_TEST_FOR_EXCEPTION
600 (A_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::initialize: "
601 "The matrix is null. You must call setMatrix() with a nonnull matrix "
602 "before you may call this method.");
603
604 Teuchos::RCP<Teuchos::Time> timer =
605 Teuchos::TimeMonitor::getNewCounter ("Ifpack2::BlockRelaxation::initialize");
606 double startTime = timer->wallTime();
607
608 { // Timing of initialize starts here
609 Teuchos::TimeMonitor timeMon (*timer);
610 IsInitialized_ = false;
611
612 // Check whether we have a BlockCrsMatrix
613 Teuchos::RCP<const block_crs_matrix_type> A_bcrs =
614 Teuchos::rcp_dynamic_cast<const block_crs_matrix_type> (A_);
615 hasBlockCrsMatrix_ = !A_bcrs.is_null();
616 if (A_bcrs.is_null ()) {
617 hasBlockCrsMatrix_ = false;
618 }
619 else {
620 hasBlockCrsMatrix_ = true;
621 }
622
623 NumLocalRows_ = A_->getLocalNumRows ();
624 NumGlobalRows_ = A_->getGlobalNumRows ();
625 NumGlobalNonzeros_ = A_->getGlobalNumEntries ();
626
627 // NTS: Will need to add support for Zoltan2 partitions later Also,
628 // will need a better way of handling the Graph typing issue.
629 // Especially with ordinal types w.r.t the container.
630 Partitioner_ = Teuchos::null;
631
632 if (PartitionerType_ == "linear") {
633 Partitioner_ =
634 rcp (new Ifpack2::LinearPartitioner<row_graph_type> (A_->getGraph ()));
635 } else if (PartitionerType_ == "line") {
636 Partitioner_ =
638 } else if (PartitionerType_ == "user") {
639 Partitioner_ =
640 rcp (new Ifpack2::Details::UserPartitioner<row_graph_type> (A_->getGraph () ) );
641 } else {
642 // We should have checked for this in setParameters(), so it's a
643 // logic_error, not an invalid_argument or runtime_error.
644 TEUCHOS_TEST_FOR_EXCEPTION
645 (true, std::logic_error, "Ifpack2::BlockRelaxation::initialize: Unknown "
646 "partitioner type " << PartitionerType_ << ". Valid values are "
647 "\"linear\", \"line\", and \"user\".");
648 }
649
650 // need to partition the graph of A
651 Partitioner_->setParameters (List_);
652 Partitioner_->compute ();
653
654 // get actual number of partitions
655 NumLocalBlocks_ = Partitioner_->numLocalParts ();
656
657 // Note: Unlike Ifpack, we'll punt on computing W_ until compute(), which is where
658 // we assume that the type of relaxation has been chosen.
659
660 if (A_->getComm()->getSize() != 1) {
661 IsParallel_ = true;
662 } else {
663 IsParallel_ = false;
664 }
665
666 // We should have checked for this in setParameters(), so it's a
667 // logic_error, not an invalid_argument or runtime_error.
668 TEUCHOS_TEST_FOR_EXCEPTION
669 (NumSweeps_ < 0, std::logic_error, "Ifpack2::BlockRelaxation::initialize: "
670 "NumSweeps_ = " << NumSweeps_ << " < 0.");
671
672 // Extract the submatrices
673 ExtractSubmatricesStructure ();
674
675 // Compute the weight vector if we're doing overlapped Jacobi (and
676 // only if we're doing overlapped Jacobi).
677 if (PrecType_ == Ifpack2::Details::JACOBI && OverlapLevel_ > 0) {
678 TEUCHOS_TEST_FOR_EXCEPTION
679 (hasBlockCrsMatrix_, std::runtime_error,
680 "Ifpack2::BlockRelaxation::initialize: "
681 "We do not support overlapped Jacobi yet for Tpetra::BlockCrsMatrix. Sorry!");
682
683 // weight of each vertex
684 W_ = rcp (new vector_type (A_->getRowMap ()));
685 W_->putScalar (STS::zero ());
686 {
687 Teuchos::ArrayRCP<scalar_type > w_ptr = W_->getDataNonConst(0);
688
689 for (local_ordinal_type i = 0 ; i < NumLocalBlocks_ ; ++i) {
690 for (size_t j = 0 ; j < Partitioner_->numRowsInPart(i) ; ++j) {
691 local_ordinal_type LID = (*Partitioner_)(i,j);
692 w_ptr[LID] += STS::one();
693 }
694 }
695 }
696 // communicate to sum together W_[k]'s (# of blocks/patches) that update
697 // kth dof) and have this information in overlapped/extended vector.
698 // only needed when Schwarz combine mode is ADD as opposed to ZERO (which is RAS)
699
700 if (schwarzCombineMode_ == "ADD") {
701 typedef Tpetra::MultiVector< typename MatrixType::scalar_type, typename MatrixType::local_ordinal_type, typename MatrixType::global_ordinal_type,typename MatrixType::node_type> scMV;
702 Teuchos::RCP<const import_type> theImport = A_->getGraph()->getImporter();
703 if (!theImport.is_null()) {
704 scMV nonOverLapW(theImport->getSourceMap(), 1, false);
705 Teuchos::ArrayRCP<scalar_type > w_ptr = W_->getDataNonConst(0);
706 Teuchos::ArrayRCP<scalar_type> nonOverLapWArray = nonOverLapW.getDataNonConst(0);
707 nonOverLapW.putScalar(STS::zero ());
708 for (int ii = 0; ii < (int) theImport->getSourceMap()->getLocalNumElements(); ii++) nonOverLapWArray[ii] = w_ptr[ii];
709 nonOverLapWArray = Teuchos::null;
710 w_ptr = Teuchos::null;
711 nonOverLapW.doExport (*W_, *theImport, Tpetra::ADD);
712 W_->doImport( nonOverLapW, *theImport, Tpetra::INSERT);
713 }
714
715 }
716 W_->reciprocal (*W_);
717 }
718 } // timing of initialize stops here
719
720 InitializeTime_ += (timer->wallTime() - startTime);
721 ++NumInitialize_;
722 IsInitialized_ = true;
723}
724
725
726template<class MatrixType,class ContainerType>
727void
729compute ()
730{
731 using Teuchos::rcp;
732
733 TEUCHOS_TEST_FOR_EXCEPTION
734 (A_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::compute: "
735 "The matrix is null. You must call setMatrix() with a nonnull matrix, "
736 "then call initialize(), before you may call this method.");
737
738 if (! isInitialized ()) {
739 initialize ();
740 }
741
742 Teuchos::RCP<Teuchos::Time> timer =
743 Teuchos::TimeMonitor::getNewCounter ("Ifpack2::BlockRelaxation::compute");
744
745 double startTime = timer->wallTime();
746
747 {
748 Teuchos::TimeMonitor timeMon (*timer);
749
750 // reset values
751 IsComputed_ = false;
752
753 Container_->compute(); // compute each block matrix
754 }
755
756 ComputeTime_ += (timer->wallTime() - startTime);
757 ++NumCompute_;
758 IsComputed_ = true;
759}
760
761template<class MatrixType, class ContainerType>
762void
765{
766 TEUCHOS_TEST_FOR_EXCEPTION
767 (Partitioner_.is_null (), std::runtime_error, "Ifpack2::BlockRelaxation::"
768 "ExtractSubmatricesStructure: Partitioner object is null.");
769
770 std::string containerType = ContainerType::getName ();
771 if (containerType == "Generic") {
772 // ContainerType is Container<row_matrix_type>. Thus, we need to
773 // get the container name from the parameter list. We use "TriDi"
774 // as the default container type.
775 containerType = containerType_;
776 }
777 //Whether the Container will define blocks (partitions)
778 //in terms of individual DOFs, and not by nodes (blocks).
779 bool pointIndexed = decouple_ && hasBlockCrsMatrix_;
780 Teuchos::Array<Teuchos::Array<local_ordinal_type> > blockRows;
781 if(decouple_)
782 {
783 //dofs [per node] is how many blocks each partition will be split into
784 auto A_bcrs = Teuchos::rcp_dynamic_cast<const block_crs_matrix_type>(A_);
785 local_ordinal_type dofs = hasBlockCrsMatrix_ ?
786 A_bcrs->getBlockSize() : List_.get<int>("partitioner: PDE equations");
787 blockRows.resize(NumLocalBlocks_ * dofs);
788 if(hasBlockCrsMatrix_)
789 {
790 for(local_ordinal_type i = 0; i < NumLocalBlocks_; i++)
791 {
792 size_t blockSize = Partitioner_->numRowsInPart(i);
793 //block i will be split into j different blocks,
794 //each corresponding to a different dof
795 for(local_ordinal_type j = 0; j < dofs; j++)
796 {
797 local_ordinal_type blockIndex = i * dofs + j;
798 blockRows[blockIndex].resize(blockSize);
799 for(size_t k = 0; k < blockSize; k++)
800 {
801 //here, the row and dof are combined to get the point index
802 //(what the row would be if A were a CrsMatrix)
803 blockRows[blockIndex][k] = (*Partitioner_)(i, k) * dofs + j;
804 }
805 }
806 }
807 }
808 else
809 {
810 //Rows in each partition are distributed round-robin to the blocks -
811 //that's how MueLu stores DOFs in a non-block matrix
812 for(local_ordinal_type i = 0; i < NumLocalBlocks_; i++)
813 {
814 //#ifdef HAVE_IFPACK2_DEBUG
815 TEUCHOS_TEST_FOR_EXCEPTION(Partitioner_->numRowsInPart(i) % dofs != 0, std::logic_error,
816 "Expected size of all blocks (partitions) to be divisible by MueLu dofs/node.");
817 size_t blockSize = Partitioner_->numRowsInPart(i) / dofs;
818 //#endif
819 //block i will be split into j different blocks,
820 //each corresponding to a different dof
821 for(local_ordinal_type j = 0; j < dofs; j++)
822 {
823 local_ordinal_type blockIndex = i * dofs + j;
824 blockRows[blockIndex].resize(blockSize);
825 for(size_t k = 0; k < blockSize; k++)
826 {
827 blockRows[blockIndex][k] = (*Partitioner_)(i, k * dofs + j);
828 }
829 }
830 }
831 }
832 }
833 else
834 {
835 //No decoupling - just get the rows directly from Partitioner
836 blockRows.resize(NumLocalBlocks_);
837 for(local_ordinal_type i = 0; i < NumLocalBlocks_; ++i)
838 {
839 const size_t numRows = Partitioner_->numRowsInPart (i);
840 blockRows[i].resize(numRows);
841 // Extract a list of the indices of each partitioner row.
842 for (size_t j = 0; j < numRows; ++j)
843 {
844 blockRows[i][j] = (*Partitioner_) (i,j);
845 }
846 }
847 }
848 //right before constructing the
849 Container_ = ContainerFactory<MatrixType>::build(containerType, A_, blockRows, Importer_, pointIndexed);
850 Container_->setParameters(List_);
851 Container_->initialize();
852}
853
854template<class MatrixType,class ContainerType>
855void
856BlockRelaxation<MatrixType,ContainerType>::
857ApplyInverseJacobi (const MV& X, MV& Y) const
858{
859 const size_t NumVectors = X.getNumVectors ();
860
861 MV AY (Y.getMap (), NumVectors);
862
863 // Initial matvec not needed
864 int starting_iteration = 0;
865 if (OverlapLevel_ > 0)
866 {
867 //Overlapping jacobi, with view of W_
868 auto WView = W_->getLocalViewHost (Tpetra::Access::ReadOnly);
869 if(ZeroStartingSolution_) {
870 auto XView = X.getLocalViewHost (Tpetra::Access::ReadOnly);
871 auto YView = Y.getLocalViewHost (Tpetra::Access::ReadWrite);
872 Container_->DoOverlappingJacobi(XView, YView, WView, DampingFactor_, nonsymCombine_);
873 starting_iteration = 1;
874 }
875 const scalar_type ONE = STS::one();
876 for(int j = starting_iteration; j < NumSweeps_; j++)
877 {
878 applyMat (Y, AY);
879 AY.update (ONE, X, -ONE);
880 {
881 auto AYView = AY.getLocalViewHost (Tpetra::Access::ReadOnly);
882 auto YView = Y.getLocalViewHost (Tpetra::Access::ReadWrite);
883 Container_->DoOverlappingJacobi (AYView, YView, WView, DampingFactor_, nonsymCombine_);
884 }
885 }
886 }
887 else
888 {
889 //Non-overlapping
890 if(ZeroStartingSolution_)
891 {
892 auto XView = X.getLocalViewHost (Tpetra::Access::ReadOnly);
893 auto YView = Y.getLocalViewHost (Tpetra::Access::ReadWrite);
894 Container_->DoJacobi (XView, YView, DampingFactor_);
895 starting_iteration = 1;
896 }
897 const scalar_type ONE = STS::one();
898 for(int j = starting_iteration; j < NumSweeps_; j++)
899 {
900 applyMat (Y, AY);
901 AY.update (ONE, X, -ONE);
902 {
903 auto AYView = AY.getLocalViewHost (Tpetra::Access::ReadOnly);
904 auto YView = Y.getLocalViewHost (Tpetra::Access::ReadWrite);
905 Container_->DoJacobi (AYView, YView, DampingFactor_);
906 }
907 }
908 }
909}
910
911template<class MatrixType,class ContainerType>
912void
913BlockRelaxation<MatrixType,ContainerType>::
914ApplyInverseGS (const MV& X, MV& Y) const
915{
916 using Teuchos::Ptr;
917 using Teuchos::ptr;
918 size_t numVecs = X.getNumVectors();
919 //Get view of X (is never modified in this function)
920 auto XView = X.getLocalViewHost (Tpetra::Access::ReadOnly);
921 auto YView = Y.getLocalViewHost (Tpetra::Access::ReadWrite);
922 //Pre-import Y, if parallel
923 Ptr<MV> Y2;
924 bool deleteY2 = false;
925 if(IsParallel_)
926 {
927 Y2 = ptr(new MV(Importer_->getTargetMap(), numVecs));
928 deleteY2 = true;
929 }
930 else
931 Y2 = ptr(&Y);
932 if(IsParallel_)
933 {
934 for(int j = 0; j < NumSweeps_; ++j)
935 {
936 //do import once per sweep
937 Y2->doImport(Y, *Importer_, Tpetra::INSERT);
938 auto Y2View = Y2->getLocalViewHost (Tpetra::Access::ReadWrite);
939 Container_->DoGaussSeidel(XView, YView, Y2View, DampingFactor_);
940 }
941 }
942 else
943 {
944 auto Y2View = Y2->getLocalViewHost (Tpetra::Access::ReadWrite);
945 for(int j = 0; j < NumSweeps_; ++j)
946 {
947 Container_->DoGaussSeidel(XView, YView, Y2View, DampingFactor_);
948 }
949 }
950 if(deleteY2)
951 delete Y2.get();
952}
953
954template<class MatrixType,class ContainerType>
955void
956BlockRelaxation<MatrixType,ContainerType>::
957ApplyInverseSGS (const MV& X, MV& Y) const
958{
959 using Teuchos::Ptr;
960 using Teuchos::ptr;
961 //Get view of X (is never modified in this function)
962 auto XView = X.getLocalViewHost (Tpetra::Access::ReadOnly);
963 auto YView = Y.getLocalViewHost (Tpetra::Access::ReadWrite);
964 //Pre-import Y, if parallel
965 Ptr<MV> Y2;
966 bool deleteY2 = false;
967 if(IsParallel_)
968 {
969 Y2 = ptr(new MV(Importer_->getTargetMap(), X.getNumVectors()));
970 deleteY2 = true;
971 }
972 else
973 Y2 = ptr(&Y);
974 if(IsParallel_)
975 {
976 for(int j = 0; j < NumSweeps_; ++j)
977 {
978 //do import once per sweep
979 Y2->doImport(Y, *Importer_, Tpetra::INSERT);
980 auto Y2View = Y2->getLocalViewHost (Tpetra::Access::ReadWrite);
981 Container_->DoSGS(XView, YView, Y2View, DampingFactor_);
982 }
983 }
984 else
985 {
986 auto Y2View = Y2->getLocalViewHost (Tpetra::Access::ReadWrite);
987 for(int j = 0; j < NumSweeps_; ++j)
988 {
989 Container_->DoSGS(XView, YView, Y2View, DampingFactor_);
990 }
991 }
992 if(deleteY2)
993 delete Y2.get();
994}
995
996template<class MatrixType,class ContainerType>
997void BlockRelaxation<MatrixType,ContainerType>::computeImporter () const
998{
999 using Teuchos::RCP;
1000 using Teuchos::rcp;
1001 using Teuchos::Ptr;
1002 using Teuchos::ArrayView;
1003 using Teuchos::Array;
1004 using Teuchos::Comm;
1005 using Teuchos::rcp_dynamic_cast;
1006 if(IsParallel_)
1007 {
1008 if(hasBlockCrsMatrix_)
1009 {
1010 const RCP<const block_crs_matrix_type> bcrs =
1011 rcp_dynamic_cast<const block_crs_matrix_type>(A_);
1012 int bs_ = bcrs->getBlockSize();
1013 RCP<const map_type> oldDomainMap = A_->getDomainMap();
1014 RCP<const map_type> oldColMap = A_->getColMap();
1015 // Because A is a block CRS matrix, import will not do what you think it does
1016 // We have to construct the correct maps for it
1017 global_size_t numGlobalElements = oldColMap->getGlobalNumElements() * bs_;
1018 global_ordinal_type indexBase = oldColMap->getIndexBase();
1019 RCP<const Comm<int>> comm = oldColMap->getComm();
1020 ArrayView<const global_ordinal_type> oldColElements = oldColMap->getLocalElementList();
1021 Array<global_ordinal_type> newColElements(bs_ * oldColElements.size());
1022 for(int i = 0; i < oldColElements.size(); i++)
1023 {
1024 for(int j = 0; j < bs_; j++)
1025 newColElements[i * bs_ + j] = oldColElements[i] * bs_ + j;
1026 }
1027 RCP<map_type> colMap(new map_type(numGlobalElements, newColElements, indexBase, comm));
1028 // Create the importer
1029 Importer_ = rcp(new import_type(oldDomainMap, colMap));
1030 }
1031 else if(!A_.is_null())
1032 {
1033 Importer_ = A_->getGraph()->getImporter();
1034 if(Importer_.is_null())
1035 Importer_ = rcp(new import_type(A_->getDomainMap(), A_->getColMap()));
1036 }
1037 }
1038 //otherwise, Importer_ is not needed and is left NULL
1039}
1040
1041template<class MatrixType, class ContainerType>
1042std::string
1044description () const
1045{
1046 std::ostringstream out;
1047
1048 // Output is a valid YAML dictionary in flow style. If you don't
1049 // like everything on a single line, you should call describe()
1050 // instead.
1051 out << "\"Ifpack2::BlockRelaxation\": {";
1052 if (this->getObjectLabel () != "") {
1053 out << "Label: \"" << this->getObjectLabel () << "\", ";
1054 }
1055 // out << "Initialized: " << (isInitialized () ? "true" : "false") << ", ";
1056 // out << "Computed: " << (isComputed () ? "true" : "false") << ", ";
1057 if (A_.is_null ()) {
1058 out << "Matrix: null, ";
1059 }
1060 else {
1061 // out << "Matrix: not null"
1062 // << ", Global matrix dimensions: ["
1063 out << "Global matrix dimensions: ["
1064 << A_->getGlobalNumRows () << ", " << A_->getGlobalNumCols () << "], ";
1065 }
1066
1067 // It's useful to print this instance's relaxation method. If you
1068 // want more info than that, call describe() instead.
1069 out << "\"relaxation: type\": ";
1070 if (PrecType_ == Ifpack2::Details::JACOBI) {
1071 out << "Block Jacobi";
1072 } else if (PrecType_ == Ifpack2::Details::GS) {
1073 out << "Block Gauss-Seidel";
1074 } else if (PrecType_ == Ifpack2::Details::SGS) {
1075 out << "Block Symmetric Gauss-Seidel";
1076 } else if (PrecType_ == Ifpack2::Details::MTSPLITJACOBI) {
1077 out << "MT Split Jacobi";
1078 } else {
1079 out << "INVALID";
1080 }
1081
1082 // BlockCrs if we have that
1083 if(hasBlockCrsMatrix_)
1084 out<<", BlockCrs";
1085
1086 // Print the approximate # rows per part
1087 int approx_rows_per_part = A_->getLocalNumRows()/Partitioner_->numLocalParts();
1088 out <<", blocksize: "<<approx_rows_per_part;
1089
1090 out << ", overlap: " << OverlapLevel_;
1091
1092 out << ", " << "sweeps: " << NumSweeps_ << ", "
1093 << "damping factor: " << DampingFactor_ << ", ";
1094
1095 std::string containerType = ContainerType::getName();
1096 out << "block container: " << ((containerType == "Generic") ? containerType_ : containerType);
1097 if(List_.isParameter("partitioner: PDE equations"))
1098 out << ", dofs/node: "<<List_.get<int>("partitioner: PDE equations");
1099
1100
1101 out << "}";
1102 return out.str();
1103}
1104
1105template<class MatrixType,class ContainerType>
1106void
1108describe (Teuchos::FancyOStream& out,
1109 const Teuchos::EVerbosityLevel verbLevel) const
1110{
1111 using std::endl;
1112 using std::setw;
1113 using Teuchos::TypeNameTraits;
1114 using Teuchos::VERB_DEFAULT;
1115 using Teuchos::VERB_NONE;
1116 using Teuchos::VERB_LOW;
1117 using Teuchos::VERB_MEDIUM;
1118 using Teuchos::VERB_HIGH;
1119 using Teuchos::VERB_EXTREME;
1120
1121 Teuchos::EVerbosityLevel vl = verbLevel;
1122 if (vl == VERB_DEFAULT) vl = VERB_LOW;
1123 const int myImageID = A_->getComm()->getRank();
1124
1125 // Convention requires that describe() begin with a tab.
1126 Teuchos::OSTab tab (out);
1127
1128 // none: print nothing
1129 // low: print O(1) info from node 0
1130 // medium:
1131 // high:
1132 // extreme:
1133 if (vl != VERB_NONE && myImageID == 0) {
1134 out << "Ifpack2::BlockRelaxation<"
1135 << TypeNameTraits<MatrixType>::name () << ", "
1136 << TypeNameTraits<ContainerType>::name () << " >:";
1137 Teuchos::OSTab tab1 (out);
1138
1139 if (this->getObjectLabel () != "") {
1140 out << "label: \"" << this->getObjectLabel () << "\"" << endl;
1141 }
1142 out << "initialized: " << (isInitialized () ? "true" : "false") << endl
1143 << "computed: " << (isComputed () ? "true" : "false") << endl;
1144
1145 std::string precType;
1146 if (PrecType_ == Ifpack2::Details::JACOBI) {
1147 precType = "Block Jacobi";
1148 } else if (PrecType_ == Ifpack2::Details::GS) {
1149 precType = "Block Gauss-Seidel";
1150 } else if (PrecType_ == Ifpack2::Details::SGS) {
1151 precType = "Block Symmetric Gauss-Seidel";
1152 }
1153 out << "type: " << precType << endl
1154 << "global number of rows: " << A_->getGlobalNumRows () << endl
1155 << "global number of columns: " << A_->getGlobalNumCols () << endl
1156 << "number of sweeps: " << NumSweeps_ << endl
1157 << "damping factor: " << DampingFactor_ << endl
1158 << "nonsymmetric overlap combine" << nonsymCombine_ << endl
1159 << "backwards mode: "
1160 << ((PrecType_ == Ifpack2::Details::GS && DoBackwardGS_) ? "true" : "false")
1161 << endl
1162 << "zero starting solution: "
1163 << (ZeroStartingSolution_ ? "true" : "false") << endl;
1164 }
1165}
1166
1167}//namespace Ifpack2
1168
1169
1170// Macro that does explicit template instantiation (ETI) for
1171// Ifpack2::BlockRelaxation. S, LO, GO, N correspond to the four
1172// template parameters of Ifpack2::Preconditioner and
1173// Tpetra::RowMatrix.
1174//
1175// We only instantiate for MatrixType = Tpetra::RowMatrix. There's no
1176// need to instantiate for Tpetra::CrsMatrix too. All Ifpack2
1177// preconditioners can and should do dynamic casts if they need a type
1178// more specific than RowMatrix. This keeps build time short and
1179// library and executable sizes small.
1180
1181#ifdef HAVE_IFPACK2_EXPLICIT_INSTANTIATION
1182
1183#define IFPACK2_BLOCKRELAXATION_INSTANT(S,LO,GO,N) \
1184 template \
1185 class Ifpack2::BlockRelaxation< \
1186 Tpetra::RowMatrix<S, LO, GO, N>, \
1187 Ifpack2::Container<Tpetra::RowMatrix<S, LO, GO, N> > >;
1188
1189#endif // HAVE_IFPACK2_EXPLICIT_INSTANTIATION
1190
1191#endif // IFPACK2_BLOCKRELAXATION_DEF_HPP
Ifpack2::BlockRelaxation class declaration.
Declaration of a user-defined partitioner in which the user can define a partition of the graph....
Block relaxation preconditioners (or smoothers) for Tpetra::RowMatrix and Tpetra::CrsMatrix sparse ma...
Definition: Ifpack2_BlockRelaxation_decl.hpp:91
Teuchos::RCP< const row_matrix_type > getMatrix() const
The input matrix of this preconditioner's constructor.
Definition: Ifpack2_BlockRelaxation_def.hpp:375
Teuchos::RCP< const map_type > getDomainMap() const
Returns the Tpetra::Map object associated with the domain of this operator.
Definition: Ifpack2_BlockRelaxation_def.hpp:384
int getNumApply() const
Returns the number of calls to apply().
Definition: Ifpack2_BlockRelaxation_def.hpp:432
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
Return a list of all the parameters that this class accepts.
Definition: Ifpack2_BlockRelaxation_def.hpp:126
BlockRelaxation(const Teuchos::RCP< const row_matrix_type > &Matrix)
Constructor.
Definition: Ifpack2_BlockRelaxation_def.hpp:86
void compute()
compute the preconditioner for the specified matrix, diagonal perturbation thresholds and relaxation ...
Definition: Ifpack2_BlockRelaxation_def.hpp:729
size_t getNodeSmootherComplexity() const
Get a rough estimate of cost per iteration.
Definition: Ifpack2_BlockRelaxation_def.hpp:463
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
Definition: Ifpack2_BlockRelaxation_def.hpp:1108
MatrixType::local_ordinal_type local_ordinal_type
The type of local indices in the input MatrixType.
Definition: Ifpack2_BlockRelaxation_decl.hpp:100
virtual ~BlockRelaxation()
Destructor.
Definition: Ifpack2_BlockRelaxation_def.hpp:120
int getNumCompute() const
Returns the number of calls to compute().
Definition: Ifpack2_BlockRelaxation_def.hpp:424
double getApplyTime() const
Returns the time spent in apply().
Definition: Ifpack2_BlockRelaxation_def.hpp:456
int getNumInitialize() const
Returns the number of calls to initialize().
Definition: Ifpack2_BlockRelaxation_def.hpp:417
void initialize()
Initialize.
Definition: Ifpack2_BlockRelaxation_def.hpp:593
Teuchos::RCP< const map_type > getRangeMap() const
Returns the Tpetra::Map object associated with the range of this operator.
Definition: Ifpack2_BlockRelaxation_def.hpp:398
void applyMat(const MV &X, MV &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS) const
Applies the matrix to a Tpetra::MultiVector.
Definition: Ifpack2_BlockRelaxation_def.hpp:577
double getComputeTime() const
Returns the time spent in compute().
Definition: Ifpack2_BlockRelaxation_def.hpp:448
void setParameters(const Teuchos::ParameterList &params)
Sets all the parameters for the preconditioner.
Definition: Ifpack2_BlockRelaxation_def.hpp:179
double getInitializeTime() const
Returns the time spent in initialize().
Definition: Ifpack2_BlockRelaxation_def.hpp:440
void apply(const MV &X, MV &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, scalar_type alpha=Teuchos::ScalarTraits< scalar_type >::one(), scalar_type beta=Teuchos::ScalarTraits< scalar_type >::zero()) const
Applies the preconditioner to X, returns the result in Y.
Definition: Ifpack2_BlockRelaxation_def.hpp:480
std::string description() const
A one-line description of this object.
Definition: Ifpack2_BlockRelaxation_def.hpp:1044
Teuchos::RCP< const Teuchos::Comm< int > > getComm() const
The communicator over which the input matrix is distributed.
Definition: Ifpack2_BlockRelaxation_def.hpp:361
virtual void setMatrix(const Teuchos::RCP< const row_matrix_type > &A)
Change the matrix to be preconditioned.
Definition: Ifpack2_BlockRelaxation_def.hpp:58
MatrixType::scalar_type scalar_type
The type of the entries of the input MatrixType.
Definition: Ifpack2_BlockRelaxation_decl.hpp:97
Partition in which the user can define a nonoverlapping partition of the graph in any way they choose...
Definition: Ifpack2_Details_UserPartitioner_decl.hpp:69
Ifpack2::LinePartitioner: A class to define partitions into a set of lines.
Definition: Ifpack2_LinePartitioner_decl.hpp:78
A class to define linear partitions.
Definition: Ifpack2_LinearPartitioner_decl.hpp:60
Ifpack2 implementation details.
Preconditioners and smoothers for Tpetra sparse matrices.
Definition: Ifpack2_AdditiveSchwarz_decl.hpp:74
void getValidParameters(Teuchos::ParameterList &params)
Fills a list which contains all the parameters possibly used by Ifpack2.
Definition: Ifpack2_Parameters.cpp:51
void getParameter(const Teuchos::ParameterList &params, const std::string &name, T &value)
Set a value from a ParameterList if a parameter with the specified name exists.
Definition: Ifpack2_Parameters.hpp:59
static Teuchos::RCP< BaseContainer > build(std::string containerType, const Teuchos::RCP< const MatrixType > &A, const Teuchos::Array< Teuchos::Array< local_ordinal_type > > &partitions, const Teuchos::RCP< const import_type > importer, bool pointIndexed)
Build a specialization of Ifpack2::Container given a key that has been registered.
Definition: Ifpack2_ContainerFactory_def.hpp:89